LEfSe(LDA Effect Size)是一种用于组间显著差异物种检测的工具,能够在不同组间寻找统计学差异显著的Biomarker。通过将统计显著性检验与生物学一致性、效应相关性检验相结合,确定最有可能解释类别之间的差异特征(基因、功能、进化枝、分类单元、生物体)。
原理:
- 采用非参数因组Kruskal-Wallis秩和检验检测不同分组间丰度差异显著的物种;
- 微生物多样性分析中,样本物种丰度分布不确定,而Kruskal-Wallis秩和检验假设各样本服从概率分布具有相同的中位数,不考虑样本总体分布类型,数据是否符合方差齐性、是否符合正态分布。但此检验并未识别出这些差异发生在那些样本之间及差异的大小。
- 再利用Wilcoxon秩和检验对显著差异物种类中的所有子类成对比较检验生物学一致性;
- Wilcoxon秩和检验应用于不同类的子类的比较,把观测值和零假设的中心位置之差的绝对值的秩分别按照不同符号相加作为其检验统计量,检验成对观测数据之差是否来自均值为0的总体(即产生数据的总体是否具有相同的均值)。如果一个特征出现在两个类样本中,每个类有3个子类,那么不同子类之间的所有9个比较都必须违反原假设,并且中位数之间差异的所有符合必须一致,才可被认为生物marker。
- 最后用线性判别分析(LDA)对数据进行降维并评估差异显著物种的影响力(即LDA score),并进行可视化展示。
- LDA线性判别分析是一种有监督数据降维方法,在机器学习、数据挖掘领域比较热门的一个算法。主要思想是将高维空间中的数据投影到低维的空间中,尽量使得投影后同一类样本聚集在一起,不同类别间的样本距离越远越好,即投影后类内方差小而类间方差最大。建立LDA模型,其中类作为因变量,特征值、子类、主题值作为自变量,用于评估它们的影响大小。
分析结果:
LEfSe分析结果包括LDA值分布柱状图、进化分支图、差异物种在不同分组间的样本相对丰度比较图。
- LDA值分布柱状图:展示LDA Score大于预设值(默认2)的显著差异物种,即具有统计学差异的Biomarker。柱状图颜色代表不同组别,长短代表LDA Score,不同组间显著差异物种的影响程度。
- 进化分支图:由内至外辐射的圆圈代表由门至属的分类级别。在不同类别上的每一个小圆圈代表该水平下的一个分类,小圆圈直径大小与相对丰度大小呈正比。着色原则:无显著差异的物种统一着色为黄色,差异显著的物种(Biomarker)跟随组进行着色,红色节点表示在红色组别中起到重要作用的微生物,绿色节点表示在绿色组别中起到重要作用的微生物。若图中某一组缺失,说明此组无显著差异的物种。未在图中显示的Biomarker对应的字母编号与物种名称绘制在右侧。
安装
conda
conda install -c bioconda lefse
conda activate base
#lefse默认安装在base下
docker
docker pull biobakery/lefse
docker run -ti biobakery/lefse bash
使用
格式化输入文件
lefse_format_input.py hmp_aerobiosis_small.txt hmp_aerobiosis_small.in -c 1 -s 2 -u 3 -o 1000000
#参数:-f {c,r}:设置特征值在行或列,默认r(行)
# -c [1..n_feats]:设置那个特征值作为分类,默认1,第一行/列作为分类;
# -s [1..n_feats]: 设置那个特征值作为子分类,默认-1,没有子分类;
# -o float : 设置归一化值,默认-1.0,没有归一化;
# -u [1..n_feats]: 设置那个特征值作为主题,默认-1,没有主题;
# -m {f,s}: 设置缺失值采用的策略,f移除缺失值的特征,s移除缺失值的样本,默认f;
# -n int: 设置每个子类的最小基数(基数低的子类会被归为一组,如果基数仍然很低,则不与他们进行成对比较)
# -biom_c BIOM_CLASS:对于输入文件是biom格式,设置那个特征值用作分类;
# -biom_s BIOM_SUBCLASS:对于输入文件是biom格式,设置那个特征值作为子类;
生成LEfSe分析结果
lefse_run.py hmp_aerobiosis_small.in hmp_aerobiosis_small.res
#参数:-o str: 设置用于导出结果的文件(只有简明的文字形式)
# -a float: 设置Anova测试的alpha值,默认0.05;
# -w float: 设置Wilcoxon的alpha值,默认0.05;
# -l float: 设置LDA分数绝对值的阈值,默认2.0;
# --nlogs int max log ingluence of LDA coeff
# --wilc int: 是否执行Wicoxon步骤,默认1;
# -r str: select LDA or SVM for effect size,默认LDA
# --svm_norm int: whether to normalize the data in [0,1] for SVM feature waiting (default strongly suggested)
# -b int:设置LDA的迭代次数,默认30;
# -e int: 设置是否只在同名子类中进行wilcoxon测试,默认0;
# -c int: set whether perform the wilcoxon test ing the Curtis's approach [BETA VERSION] (default 0);
# -f float: set the subsampling fraction value for each bootstrap iteration,默认0.66666;
# -s {0,1,2}: set the multiple testing correction options. 0 no correction (more strict, default), 1 correction for independent comparisons, 2 correction for dependent comparison;
# --min_c int: 执行wilcoxon测试的每个子类的最小样本数,默认10;
# -t str: 设置分析的标题,默认不带扩展名的输入文件名;
# -y {0,1} (for multiclass tasks) set whether the test is performed in a one-against-one ( 1 - more strict!) or in a one- against-all setting ( 0 - less strict) (default 0);
结果可视化
lefse_plot_res.py hmp_aerobiosis_small.res hmp_aerobiosis_small.png
#参数:--format {png,svg,pdf}:输出图像格式;
# --dpi DPI
# --title TITLE
# --title_font_size TITLE_FONT_SIZE
# --class_legend_font_size CLASS_LEGEND_FONT_SIZE
# --width WIDTH
# --height HEIGHT :仅使用垂直的直方图;
# --left_space LS
# --right_space RS
# --orientation {h,v}
# --autoscale {0,1}
# --background_color {k,w} set the color of the background
# --subclades N_SCL 要显示的标签级别数(-1从叶子开始, 1表示所有级别,默认为 1);
# --max_feature_len MAX_FEATURE_LEN:特征字符串的最大长度,默认60;
# --otu_only:Plot only species resolved OTUs (as opposed to all levels)
# --report_features Report important features to STDOUT
lefse_plot_cladogram.py hmp_aerobiosis_small.res hmp_aerobiosis_small.cladogram.png --format png
#参数:--width WIDTH
# --height HEIGHT
# --top TOP set maximum y limit (-1.0 means automatic limit)
# --bot BOT set minimum y limit (default 0.0, -1.0 means automatic limit)
# --title_font_size TITLE_FONT_SIZE
# --class_font_size CLASS_FONT_SIZE
# --class_label_pos {up,down}
# --subcl_mean {y,n}
# --subcl_median {y,n}
# --font_size FONT_SIZE
# -n flt unused
# --format {png,pdf,svg} the format for the output file
# -f {all,diff,one}: 是否绘制所有特征(all),根据LEfSe绘制差异丰富的特征(diff),仅绘制一个(使用--feature_name给出的特征);
# --feature_name FEATURE_NAME:要绘制的特征名(由.分隔级别)
# --feature_num FEATURE_NUM:要绘制的特征数量;
# --archive {zip,none}
# --background_color {k,w}: 设置背景颜色
# --dpi DPI
输入文件示例
制表符分隔的文本文件,由数字特征表、类、子类和主题构成。特征值可以是reads数,也可以是相对丰度。菌群分类名称用“|”或“.”分隔。
bodysite mucosal mucosal mucosal mucosal mucosal non_mucosal non_mucosal non_mucosal non_mucosal non_mucosal
subsite oral gut oral oral gut skin nasal skin ear nasal
id 1023 1023 1672 1876 1672 159005010 1023 1023 1023 1672
Bacteria 0.99999 0.99999 0.999993 0.999989 0.999997 0.999927 0.999977 0.999987 0.999997 0.999993
Bacteria|Actinobacteria 0.311037 0.000864363 0.00446132 0.0312045 0.000773642 0.359354 0.761108 0.603002 0.95913 0.753688
Bacteria|Bacteroidetes 0.0689602 0.804293 0.00983343 0.0303561 0.859838 0.0195298 0.0212741 0.145729 0.0115617 0.0114511
Bacteria|Firmicutes 0.494223 0.173411 0.715345 0.813046 0.124552 0.177961 0.189178 0.188964 0.0226835 0.192665
Bacteria|Proteobacteria 0.0914284 0.0180378 0.265664 0.109549 0.00941215 0.430869 0.0225884 0.0532684 0.00512034 0.0365453
Bacteria|Firmicutes|Clostridia 0.090041 0.170246 0.004831
喜欢的小伙伴可关注下方公众号