LEfSe学习[通俗易懂]

#!/bin/sh
# in this script we show how to perform the biomarker discovery operation
# using LEfSe. The scripts require LEfSe to be installed and in the system path

# convert the sample names in the table of abundance into classes (i.e. the two bodysites)
sed 's/\([A-Z][A-Z]\)_\w*/\1/g' output/merged_abundance_table.txt > tmp/merged_abundance_table.4lefse.txt

# first LEfSe step: format the input specifying that the class info is in the first row
format_input.py tmp/merged_abundance_table.4lefse.txt tmp/merged_abundance_table.lefse -c 1 -o 1000000

# run the LEfSe biomarker discovery tool with default options apart for the 
# threshold on the LDA effect size which is increaset to 4
run_lefse.py tmp/merged_abundance_table.lefse tmp/merged_abundance_table.lefse.out -l 4

# Plot the resulting list of biomarkers with the corresponsing effect size
plot_res.py --dpi 300 tmp/merged_abundance_table.lefse.out output_images/lefse_biomarkers.png

# Plot the biomarkers on the underlying cladogram
plot_cladogram.py --dpi 300 --format png tmp/merged_abundance_table.lefse.out output_images/lefse_biomarkers_cladogram.png 

# Plot one features specifically (Firmicutes in this case) 
plot_features.py -f one --feature_name "k__Bacteria.p__Firmicutes" tmp/merged_abundance_table.lefse tmp/merged_abundance_table.lefse.out Firmicutes.png

# Plot all biomarkes saving the images in one zip archive ("-f diff" is for plotting biomarkers only, with "-f all" one can plot all input features)
plot_features.py -f diff --archive zip tmp/merged_abundance_table.lefse tmp/merged_abundance_table.lefse.out biomarkers.zip







## 格式化输入表
format_input.py merged_abundance_table.4lefse.txt merged_abundance_table.lefse -c 1 -o 1000000

## LEfSe生物标记物发现工具可以与默认统计选项一起使用（-l 阈值的设置）
run_lefse.py merged_abundance_table.lefse merged_abundance_table.lefse.out -l 4


## 绘制具有相应效应大小的生成标记
plot_res.py --dpi 300 merged_abundance_table.lefse.out lefse_biomarkers.png


## 基础分类树上的生物标记
plot_cladogram.py --dpi 300 --format png merged_abundance_table.lefse.out lefse_biomarkers_cladogram.png 


## 
plot_features.py -f diff --archive zip merged_abundance_table.lefse merged_abundance_table.lefse.out biomarkers.zip


PS：
### merged_abundance_table.4lefse.txt
ID    group1    group1    group1    group1    group2    group2    group2    group2        
k__Archaea      0.01107 0.0     0.01438 0.0     0.00591 0.00909 0.0     0.00264 0.06548 0.06238
k__Archaea|p__Euryarchaeota     0.01107 0.0     0.01438 0.0     0.00591 0.00909 0.0     0.00264 0.06548 0.06238
k__Archaea|p__Euryarchaeota|c__Methanobacteria  0.01107 0.0     0.01438 0.0     0.00591 0.00909 0.0     0.00264 0.06548 0.06238
k__Archaea|p__Euryarchaeota|c__Methanobacteria|o__Methanobacteriales    0.01107 0.0     0.01438 0.0     0.00591 0.00909 0.0     0.00264 0.06548 0.06238
k__Archaea|p__Euryarchaeota|c__Methanobacteria|o__Methanobacteriales|f__Methanobacteriaceae     0.01107 0.0     0.01438 0.0     0.00591 0.00909 0.0     0.00264 0.06548 0.06238
k__Archaea|p__Euryarchaeota|c__Methanobacteria|o__Methanobacteriales|f__Methanobacteriaceae|g__Methanobrevibacter       0.01107 0.0     0.01438 0.0     0.00591 0.00909 0.0     0.00264 0.06548 0.06238