用 SPIEC-EASI 进行微生物网络分析

生信菜鸟团 · 公众号 · 生物 · 2020-10-18 22:00

正文

在之前的文章中我们提到，由于微生物组数据之间的非独立性，往往会导致传统相关性计算方法（Pearson）出现偏差。除了 SparCC 之外， SPIEC-EASI （ SP arse I nvers E C ovariance Estimation for E cological A ssociation I nference，读作 speakeasy ）也是一种常用的统计方法。这种方法将针对组成数据开发的数据转换方法与稀疏图形模型推理框架相结合，使用稀疏邻域和逆协方差选择算法构建微生物组网络，该流程封装于 SpiecEasi R 包中。

Sparse and compositionally robust inference of microbial ecological networks. PLoS Comput. Biol. 2015; 11: e1004226

GitHub 地址：https://github.com/zdk123/SpiecEasi

SPIEC-EASI 分析流程

SPIEC-EASI 流程主要分三个阶段进行：首先，对数据进行预处理，进行中心对数比转换（CLR），确保组成数据的鲁棒性。第二步，选择图形模型推断方法：1）邻域选择（MB 方法）或 2）逆协方差选择（glasso 方法，默认值）。SPIEC-EASI 网络推论假设基础网络稀疏。通过一种叫 Stability Approach to Regularization Selection （StARS）的 bootstrap 方法，对数据集进行随机子采样，在所选边集中找到可靠网络，以推断正确的模型稀疏度。SPIEC-EASI 输出包括网络（来自逆协方差网络的非零项）和一个逆协方差矩阵。

安装

library(devtools)install_github("zdk123/SpiecEasi")library(SpiecEasi)

基本用法

SpiecEasi 一行命令即可完成整个流程，这里我们使用基于 American Gut 的示例数据，同时比较 MB 方法和 glasso 方法以及 SparCC 结果之间的区别。

data(amgut1.filt)se.mb.amgut  spiec.easi(amgut1.filt, method='mb', lambda.min.ratio=1e-2,                          nlambda=20, pulsar.params=list(rep.num=50))se.gl.amgut  spiec




    
.easi(amgut1.filt, method='glasso', lambda.min.ratio=1e-2,                          nlambda=20, pulsar.params=list(rep.num=50))sparcc.amgut  sparcc(amgut1.filt)## 定义 SparCC 阈值sparcc.graph  abs(sparcc.amgut$Cor) >= 0.3diag(sparcc.graph)  0library(Matrix)sparcc.graph  Matrix(sparcc.graph, sparse=TRUE)## 创建 igraph 对象ig.mb      adj2igraph(getRefit(se.mb.amgut))ig.gl      adj2igraph(getRefit(se.gl.amgut))ig.sparcc  adj2igraph(sparcc.graph)

用 igraph 可视化：

library(igraph)## set size of vertex proportional to clr-meanvsize     rowMeans(clr(amgut1.filt, 1))+6am.coord  layout.fruchterman.reingold(ig.mb)par(mfrow=c(1,3))plot(ig.mb, layout=am.coord, vertex.size=vsize, vertex.label=NA, main="MB")plot(ig.gl, layout=am.coord, vertex.size=vsize, vertex.label=NA, main="glasso")plot(ig.sparcc, layout=am.coord, vertex.size=vsize, vertex.label=NA, main="sparcc")

由于 SPIEC-EASI 基于惩罚估计量，因此边缘权重无法直接与 SparCC（或Pearson/Spearman相关系数）相比较，需要进行调整。

library(




    
Matrix)secor   cov2cor(getOptCov(se.gl.amgut))sebeta  symBeta(getOptBeta(se.mb.amgut), mode='maxabs')elist.gl      summary(triu(secor*getRefit(se.gl.amgut), k=1))elist.mb      summary(sebeta)elist.sparcc  summary(sparcc.graph*sparcc.amgut$Cor)hist(elist.sparcc[,3], main='', xlab='edge weights')hist(elist.mb[,3], add=TRUE, col='forestgreen')hist(elist.gl[,3], add=TRUE, col='red')

比较不同方法推断出的网络性质。

dd.gl      degree.distribution(ig.gl)dd.mb      degree.distribution(ig.mb)dd.sparcc  degree.distribution(ig.sparcc)plot(0:(length(dd.sparcc)-1), dd.sparcc, ylim=c(0,.35), type='b',      ylab="Frequency", xlab="Degree", main="Degree Distributions")points(0:(length(dd.gl)-1), dd.gl, col="red" , type='b')points(0:(length(dd.mb)-1), dd.mb, col="forestgreen", type='b')legend("topright", c("MB", "glasso", "sparcc"),        col=c("forestgreen", "red", "black"), pch=1, lty=1)

用 phyloseq 对象构建网络

phyloseq 对象可直接导入 SpiecEasi 流程进行相关性分析。

用 SPIEC-EASI 进行微生物网络分析

正文

相关性分析和条件独立性

SPIEC-EASI 分析流程

安装

基本用法

用 phyloseq 对象构建网络

请到「今天看啥」查看全文