Phylogenetic Conservation and Environmental Threshold of Soil Microbes in Response to Environmental Changes

仪器设备

1. 普通Windows系统个人电脑，内存8 G，需求内存空间约500 M

软件

1. R (v3.6.3)

注：本教程是基于已经在个人电脑上安装好的相关软件和依赖包进行的。

实验步骤

1. 数据准备
   本分析用到两个数据，即土壤微生物群落数据，环境因子数据 (与采样点一一对应)。
2. 数据准备
   导入已抽平的微生物群落数据及环境因子数据 (图1和图2)
   Micro<-read.csv("otu_data.csv", header = T,row.names=1)
   Env<-read.csv("env.csv",header=T,row.names=1)
   
   
   图1. 微生物群落数据集
   
   
   图2. 环境因子数据集
   
   
3. 系统发育信号的计算
   
   3.1

   微生物群落特征数据集的获取
   为了获得潜在性状信息，通过微生物类群与环境变量之间Spearman相关来定义每个OTU的生态偏好。例如，OTU与pH的正负相关被定义为acid-preferred或alkaline-preferred，相应的在性状数据集中：正相关为1，负相关为-1，没有显著相关的为0。
         相应的计算过程如下：
   library(Hmisc) #所需依赖包：lattice、survival、ggplot2和Formula
   Micro<-Micro/11730 #11730为每个样本的最小抽平序列数
   Micro<-t(Micro)
   cor<-rcorr(as.matrix(cbind(Micro,Env)),type="spear")
   r<-cor$r
   q<-cor$P
   q[q<0.05]<-1
   q[q>0.05&q<1]<-0
   g<-r*q
   heat<-g[1:ncol(Micro),-(1:ncol(Micro))]
   traits<-replace(heat, heat<0,-1)
   traits<-replace(traits, traits>0,1)
   write.csv(traits,"traits.csv") #写出相应的性状数据集（图3）
   
   
   图3. 微生物群落性状数据集
   
   
   3.2

   Blomberg’s K的计算
   library(picante) #所需依赖包：vegan、permute和nlme
   phy<-read.tree("Micro.tre") #读取微生物类群对应的系统发育树
   phy<-multi2di(phy) #无根树转为有根树
   traits <- traits[phy$tip.label, ]
   traits<-as.matrix(traits)
   Micro_Blombergk<-multiPhylosignal(traits,phy) #计算Blomberg’s K
   write.csv(Micro_Blombergk," Micro_Blombergk.csv")
   写出计算结果 (图4)
   
   
   图4. 布朗K值数据集
   
         K值描述了一个分类单元与系统发育的关系，正如布朗运动的系统发育信号强度度量所期望的那样。K值接近于0表示进化是随机的或收敛的，而K值大于1则表示强烈的系统发育信号和性状保守性。
   3.3

   绘制系统发育信号柱状图 (图5)
   data<-read.csv("Micro_Blombergk.csv",header=T,row.names=1)
   data<-t(data)
   data<-data[1,]
   data<-as.matrix(data)
   data<-t(data)
   barplot(data,axes=T,width=2,axis.lty=1,ylab="Blomberg’s K ",
   legend=rownames(data),beside=TRUE,col=c("#366CA8"))
   
   
   图5. 系统发育信号柱状图
   
   
4. 环境阈值的计算 (图6和7)
   
   4.1

   对输入的微生物群落数据集Micro进行筛选，标准为：每个OTU至少在3个样本中存在，获得的数据集为Micro1
   Micro1<- Micro[rowSums(Micro >= 1)>= 3,]
   Micro1<- Micro1/11730 #11730为每个样本的最小抽平序列数
   Micro1<-t(Micro1)
   library(TITAN2)
   CEC_titan<-titan(Env[,3], Micro1,numPerm=250,
   boot=TRUE,nBoot=500,imax=FALSE,ivTot=FALSE,
   pur.cut=0.95,rel.cut=0.95,memory=TRUE)
   plot_sumz(CEC_titan,filter=TRUE,xlab=expression("CEC"),pch1=20,pch2=20, col1="#fdbb2d",col2="#1E9600")
   
   
   图6. 微生物群落负响应种 (Z—) 和正响应种 (Z+) 指示总分沿CEC梯度突变点的响应曲线
   
   plot_taxa(CEC_titan, z.med=F,leg.x=.8, leg.y=5, xlab=expression("CEC"), cex.taxa=.9,cex.axis=1.35,cex=1,fil1="white",fil2="white",pch1=20, pch2=20,
   col1="#fdbb2d",col2="#1E9600",pur.cut=0.95, pval.cut=0.05, rel.cut=0.95)
   
   
   图7. 微生物群落沿CEC梯度负响应种和正响应种分布图
   
   
   4.2

   重复第4.1条步骤，分别探究微生物群落沿其他环境梯度的环境阈值。

致谢

本方案主要来源于课题组先前发表的相关文章 (Jiao et al., 2020; Zhang et al., 2020)。相关研究得到了国家自然科学基金 (41807030)、国家重点研发计划项目(2016YFD0200306) 和中国博士后科学基金 (2018M630041) 等项目的资助。

参考文献

1. Baker, M. E. and King, R. S. (2010). A new method for detecting and interpreting biodiversity and ecological community thresholds. Methods Ecol Evol 1: 25-37.
2. Isobe, K., Allison, S. D., Khalili, B., Martiny, A. C. and Martiny, J. B. H. (2019). Phylogenetic conservation of bacterial responses to soil nitrogen addition across continents. Nat Commun 10(1): 2499.
3. Jiao, S. and Lu, Y. (2020). Abundant fungi adapt to broader environmental gradients than rare fungi in agricultural fields. Glob Chang Biol 26(8).
4. Oliverio, A. M., Bradford, M. A. and Fierer, N. (2016). Identifying the microbial taxa that consistently respond tosoil warming across time and space. Glob Chang Biol 23(5)
5. Zhang, Z., Zhang, J. and Jiao, S. (2020). Fungi show broader environmental thresholds in wet than dry agricultural soils with distinct biogeographic patterns. Sci Total Environ 750: 141761.