微生物群落多样性测序与功能分析Word文件下载.docx
《微生物群落多样性测序与功能分析Word文件下载.docx》由会员分享,可在线阅读,更多相关《微生物群落多样性测序与功能分析Word文件下载.docx(24页珍藏版)》请在冰豆网上搜索。
微生物菌群多样性测序受DNA提取和扩增影响很大,不同的扩增区段和扩增引物甚至PCR循环数的差异都会对结果有所影响。
因而建议同一项目不同样品的都采用相同的条件和测序方法,这样相互之间才存在可比性。
∙完成PCR之后的产物一般可以直接上测序仪测序,在上机测序前我们需要对所有样本进行定量和均一化,通常要进行荧光定量PCR。
完成定量的样品混合后就可以上机测序。
∙16srDNA测序目前可以采用多种不同的测序仪进行测序,包括罗氏的454,Illumina的MiSeq,Life的PGM或Pacbio的RSII三代测序仪。
不同的仪器各有优缺点,目前最主流的是Illumina公司的MiSeq,因为其在通量、长度和价格三者之间最为平衡。
MiSeq测序仪可以产生2x300bp的测序读长,一次可以产生15Gb的测序数据远远大于其他测序仪的测序通量。
方法/步骤
1.1
16srDNA分析基本流程:
2.2
原始数据处理:
原始测序数据需要去除接头序列,并将双端测序序列进行拼接成单条序列。
根据测序barcode序列区分不同的样本序列。
过滤低质量序列和无法比对到16srDNA数据库的序列。
3.3
OTU分类和统计:
OTU(operationaltaxonomicunits)是在系统发生学研究或群体遗传学研究中,为了便于进行分析,人为给某一个分类单元(品系,种,属,分组等)设置的同一标志。
通常按照97%的相似性阈值将序列划分为不同的OTU,每一个OTU通常被视为一个微生物物种。
相似性小于97%就可以认为属于不同的种,相似性小于93%-95%,可以认为属于不同的属。
样品中的微生物多样性和不同微生物的丰度都是基于对OTU的分析。
使用QIIME(version1.8.0)工具包进行统计注释。
4.4
我们还可以对这些种属的构成进行柱状图显示:
横坐标中每一个条形图代表一个样本,纵坐标代表该分类层级的序列数目或比例。
同一种颜色代表相同的分类级别。
图中的每根柱子中的颜色表示该样本在不同级别(门、纲、目等)的序列数目,序列数目只计算级别最低的分类,例如在属中计算过了,则在科中则不重复计算。
Q:
为什么要选择V3-V4区的测序长度?
为什么有些文献是V6区,有什么区别?
A:
16SrRNA总长约1540bp,包含9个可变区。
由于高通量测序的测序长度的限制,不可能将16SrRNA的9个可变区全部测序,所以在PCR扩增时往往只能选择1-3个可变区作为扩增片段。
Kozich等评估了Miseq测序仪分析的不同16SrRNA可变区的准确性发现,测定V4区效果最佳。
根据我们的测序长度,v3-v4区是最佳选择。
5.5
我们还需要对样本之间或分组之间的OTU进行比较获得韦恩图:
注意,韦恩图目前一般最多只能显示5个样本或分组,过多的样本无法无法进行韦恩图绘制
6.6
样品构成丰度:
稀释曲线
微生物多样性分析中需要验证测序数据量是否足以反映样品中的物种多样性,稀释曲线(丰富度曲线)可以用来检验这一指标。
稀释曲线是用来评价测序量是否足以覆盖所有类群,并间接反映样品中物种的丰富程度。
稀释曲线是利用已测得16SrDNA序列中已知的各种OTU的相对比例,来计算抽取n个(n小于测得reads序列总数)reads时出现OTU数量的期望值,然后根据一组n值(一般为一组小于总序列数的等差数列)与其相对应的OTU数量的期望值做出曲线来。
当曲线趋于平缓或者达到平台期时也就可以认为测序深度已经基本覆盖到样品中所有的物种;
反之,则表示样品中物种多样性较高,还存在较多未被测序检测到的物种。
下图中的稀释曲线
横坐标代表随机抽取的序列数量;
纵坐标代表观测到的OTU数量。
样本曲线的延伸终点的横坐标位置为该样本的测序数量,如果曲线趋于平坦表明测序已趋于饱和,增加测序数据无法再找到更多的OTU;
反之表明不饱和,增加数据量可以发现更多OTU。
7.7
Shannon-Winner曲线
Shannon-Wiener曲线,是利用shannon指数来进行绘制的,反映样品中微生物多样性的指数,利用各样品的测序量在不同测序深度时的微生物多样性指数构建曲线,以此反映各样本在不同测序数量时的微生物多样性。
当曲线趋向平坦时,说明测序数据量足够大,可以反映样品中绝大多数的微生物物种信息。
与上图一样,横坐标代表随机抽取的序列数量;
纵坐标代表的是反映物种多样性的Shannon指数。
其中曲线的最高点也就是该样本的Shannon指数,指数越高表明样品的物种多样性越高。
Shannon指数怎么算的?
Shannon指数公式:
其中,Sobs=
实际测量出的OTU数目;
ni=
含有i条序列的OTU数目;
N
=
所有的序列数。
8.8
Rank-Abundance曲线
用于同时解释样品多样性的两个方面,即样品所含物种的丰富程度和均匀程度。
物种的丰富程度由曲线在横轴上的长度来反映,曲线越宽,表示物种的组成越丰富;
物种组成的均匀程度由曲线的形状来反映,曲线越平坦,表示物种组成的均匀程度越高。
一般超过20个样本图就会变得非常复杂而且不美观,所以一般20个样本以下会做该图,图片保存为结果目录中rank.pdf。
横坐标代表物种排序的数量;
纵坐标代表观测到的相对丰度。
样本曲线的延伸终点的横坐标位置为该样本的物种数量,如果曲线越平滑下降表明样本的物种多样性越高,而曲线快速陡然下降表明样本中的优势菌群所占比例很高,多样性较低。
9.9
Alpha多样性(样本内多样性)
Alpha多样性是指一个特定区域或者生态系统内的多样性,常用的度量指标有Chao1丰富度估计量(Chao1richnessestimator)、香农-威纳多样性指数(Shannon-wienerdiversityindex)、辛普森多样性指数(Simpsondiversityindex)等。
计算菌群丰度:
Chao、ace;
计算菌群多样性:
Shannon、Simpson。
Simpson指数值越大,说明群落多样性越高;
Shannon指数越大,说明群落多样性越高。
表中显示前10个样本,如果样本大于10个,详见结果目录中的alpha_div.txt。
能不能解释下每个指数(如chao1、shannon)?
Chao1:
是用chao1算法估计群落中含OTU数目的指数,chao1在生态学中常用来估计物种总数,由Chao(1984)最早提出。
Chao1值越大代表物种总数越多。
Schao1=Sobs+n1(n1-1)/2(n2+1)
其中Schao1为估计的OTU数,Sobs为观测到的OTU数,n1为只有一条序列的OTU数目,n2为只有两条序列的OTU数目。
Shannon:
用来估算样品中微生物的多样性指数之一。
它与Simpson多样性指数均为常用的反映alpha多样性的指数。
Shannon值越大,说明群落多样性越高。
Ace:
用来估计群落中含有OTU数目的指数,由Chao提出,是生态学中估计物种总数的常用指数之一,与Chao1的算法不同。
Simpson:
用来估算样品中微生物的多样性指数之一,由EdwardHughSimpson(1949)提出,在生态学中常用来定量的描述一个区域的生物多样性。
Simpson指数值越大,说明群落多样性越高。
辛普森多样性指数=随机取样的两个个体属于不同种的概率
=1-随机取样的两个个体属于同种的概率
10.10
Beta多样性分析(样品间差异分析)
Beta多样性度量时空尺度上物种组成的变化,
是生物多样性的重要组成部分,
与许多生态学和进化生物学问题密切相关,
因此在最近10年间成为生物多样性研究的热点问题之一。
PCoA分析
PCoA(principalco-ordinatesanalysis)是一种研究数据相似性或差异性的可视化方法,通过一系列的特征值和特征向量进行排序后,选择主要排在前几位的特征值,PCoA可以找到距离矩阵中最主要的坐标,结果是数据矩阵的一个旋转,它没有改变样品点之间的相互位置关系,只是改变了坐标系统。
通过PCoA可以观察个体或群体间的差异。
每一个点代表一个样本,相同颜色的点来自同一个分组,两点之间距离越近表明两者的群落构成差异越小。
PCoA有多张图,分别代表的PCoA1-2,2-3,3-1。
11.11
NMDS分析(非度量多维尺度分析)
NMDS(NonmetricMultidimensionalScaling)常用于比对样本组之间的差异,可以基于进化关系或数量距离矩阵。
横轴和纵轴:
表示基于进化或者数量距离矩阵的数值在二维表中成图。
与PCA分析的主要差异在于考量了进化上的信息。
12.12
PCA分析
主成分分析PCA(Principalcomponentanalysis)是一种研究数据相似性或差异性的可视化方法,通过一系列的特征值和特征向量进行排序后,选择主要的前几位特征值,采取降维的思想,PCA可以找到距离矩阵中最主要的坐标,结果是数据矩阵的一个旋转,它没有改变样品点之间的相互位置关系,只是改变了坐标系统。
详细关于主成分分析的解释推荐大家看一篇文章,。
通过PCA可以观察个体或群体间的差异。
以上三个图可能遇到的问题:
1:
PCA,PcoA,NMDS分析分别是基于什么数据画的?
回答:
PCA,PcoA,NMDS分析均是基于OTU分类taxon数据所画,用的是R语言Vegan包中的相关函数画成,其中PcoA与NMDS还要基于样本之间的距离矩阵才能画成。
2:
PCA分析如果图中大部分点集中在一起,少数点在很远的外围,是什么原因造成的?
是因为样本OTU分类时候,少数样本某些菌含量特别高所造成,导致这些样本偏离正常范围,建议单独拿出这些样本观察,看是否是实验错误。
3:
PCA分析时,不是有PC1,PC2,PC3三个坐标吗?
是给出三张图吗?
还是三维立体图?
PCA作图时,会得出PC1,PC2,PC3三个坐标,可以根据PC12,PC13,PC23分别作图,一般给出的是PC12的图,当PC12图质量不好,看不出明显的样本分类效果时,可以看PC13或PC23的图分类是否清晰,也可以用R语言rgl包做出PC123三维图。
QIIME本身结果中有提供PCA的三维图结果,可以通过网页打开。
13.13
LDA差异贡献分析
PCA和LDA的差别在于,PCA,它所作的只是将整组数据整体映射到最方便表示这组数据的坐标轴上,映射时没有利用任何数据内部的分类信息,是无监督的,而LDA是由监督的,增加了种属之间的信息关系后,结合显著性差异标准测试(克鲁斯卡尔-沃利斯检验和两两Wilcoxon测试)和线性判别分析的方法进行特征选择。
除了可以检测重要特征,他还可以根据效应值进行功能特性排序,这些功能特性可以解释顶部的大部分生物学差异。
详细说明可以参考这篇文章。
不同颜色代表不同样本或组之间的显著差异物种。
使用LefSe软件分析获得,其中显著差异的logarithmicLDAscore设为2。
问题:
LDA分析有什么用?
组间差异显著物种又可以称作生物标记物(biomarkers),该分析主要是想找到组间在丰度上有显著差异的物种。
14.物种进化树的样本群落分布图
是将不同样本的群落构成及分布以物种分类树的形式在一个环图中展示。
数据经过分析后,将物种分类树和分类丰度信息通过软件GraPhlAn(http:
//huttenhower.sph.harvard.edu/GraPhlAn)进行绘制。
其目的是将物种之间的进化关系以及不同样本的物种分布丰度和最高分布样本的信息在一个视觉集中的环图中一次展示,其提供的信息量较其他图最为丰富。
中间为物种进化分类树,不同颜色的分支代表不同的纲(具体的代表颜色见右上角的图例),之后外圈的灰色标示字母的环表示的是本次研究中比例最高的15个科(字母代表的科参见左上角的图例)。
之后的外圈提供的是热力图,如果样本数<
=10个则绘制样本,如果样本数超过10个则按照分组绘制,每一环为一个样本,根据其丰度绘制的热力图。
最外圈为柱状图,绘制的是该属所占比例最高的样本的丰度和样本颜色(样本颜色见环最下方的样本名字的颜色)。
其中热力图和柱状图取值均为原比例值x10000后进行log2转换后的值
参考文献:
1.Vazquez-BaezaY,PirrungM,GonzalezA,KnightR.2013.Emperor:
Atoolforvisualizinghigh-throughputmicrobialcommunitydata.Gigascience2
(1):
16.
2.Legendre,P.andLegendre,L.1998.NumericalEcology.SecondEnglishEdition.DevelopmentsinEnvironmentalModelling20.Elsevier,Amsterdam.
3.SegataN,IzardJ,WaldronL,etal.Metagenomicbiomarkerdiscoveryandexplanation[J].GenomeBiol,2011,12(6):
R60.
4.LangilleMGI,ZaneveldJ,CaporasoJG,McDonaldD,KnightsD,ReyesJAetal.(2013).Predictivefunctionalprofilingofmicrobialcommunitiesusing16SrRNAmarkergenesequences.NatBiotechnol31:
814–821.
15.物种相关性分析
根据各个物种在各个样品中的丰度以及变化情况,计算物种之间的相关性,包括正相关和负相关。
相关性分析使用CCREPE算法,首先对原始16s测序数据的种属数量进行标准化,然后进行Spearman和Pearson秩相关分析并进行统计检验,计算出各个物种之间的相关性,之后在所有物种中根据simscore绝对值的大小,挑选出相关性最高的前100组数据,基于Cytoscap绘制共表达分析网络图,网络图采用两种不同的形式表现出来。
物种相关性网络图A:
图中每一个点代表一个物种,存在相关性的物种用连线连接,其中,红色的连线代表负相关,绿色的先代表正相关,连线颜色的深浅代表相关性的高低。
物种相关性网络图B:
图中每一个点代表一个物种,点的大小表示与其他物种的关联关系的多少,其中与之有相关性的物种数越多,点的半径和字体越大,连线的粗细代表两物种之间相关性的大小,连线越粗,相关性越高。
SchwagerE,WeingartG,BielskiC,etal.CCREPE:
CompositionalityCorrectedbyPermutationandRenormalization[J].2014.
16.聚类分析
根据OUT数据进行标准化处理(1wlog10)之后,选取数目最多的前60个物种,基于Rheatmap进行作图,热图中的每一个色块代表一个样品的一个属的丰度,样品横向排列,属纵向排列,两个热图,差异是是否对样品进行聚类,从聚类中可以了解样品之间的相似性以及属水平上的群落构成相似性。
如果聚类结果中出现大面积的白或黑是因为大量的菌含量非常低,导致都没有数值,可以在绘制之前进行标准化操作,对每一类菌单独自身进行Z标准化。
17.群落功能差异分析
通过对已有测序微生物基因组的基因功能的构成进行分析后,我们可以通过16s测序获得的物种构成推测样本中的功能基因的构成,从而分析不同样本和分组之间在功能上的差异(PICRUStNatureBiotechnology,1-10.82013)。
通过对宏基因组测序数据功能分析和对应16s预测功能分析结果的比较发现,此方法的准确性在84%-95%,对肠道微生物菌群和土壤菌群的功能分析接近95%,能非常好的反映样品中的功能基因构成。
为了能够通过16s测序数据来准确的预测出功能构成,首先需要对原始16s测序数据的种属数量进行标准化,因为不同的种属菌包含的16s拷贝数不相同。
然后将16s的种属构成信息通过构建好的已测序基因组的种属功能基因构成表映射获得预测的功能结果。
(根据属这个水平,对不同样本间的物种丰度进行显著性差异两两检验,我们这里的检验方法使用STAMP中的two-sample中T-TEST方法,Pvalue值过滤为0.05,作Extenterrorbar图。
)
此处提供COG,KO基因预测以及KEGG代谢途径预测。
用户也可自行使用我们提供的文件和软件(STAMP)对不同层级以及不同分组之间进行统计分析和制图,以及选择不同的统计方法和显著性水平。
DonovanH.Parks1,GeneW.Tyson,STAMP:
statisticalanalysisoftaxonomicandfunctionalprofiles,Bioinformatics(2014)30(21):
3123-3124.doi:
10.1093
18.COG构成差异分析图
图中不同颜色代表不同的分组,列出了COG构成在组间存在显著差异的功能分类以及在各组的比例,此外右侧还给出了差异的比例和置信区间以及P-value。
19.KEGG代谢途径差异分析图
通过KEGG代谢途径的预测差异分析,我们可以了解到不同分组的样品之间在微生物群落的功能基因在代谢途径上的差异,以及变化的高低。
为我们了解群落样本的环境适应变化的代谢过程提供一种简便快捷的方法。
图解读:
图中不同颜色代表不同的分组,列出了在第三层级的构成在组间存在显著差异的KEGG代谢途径第三层分类以及在各组的比例,此外右侧还给出了差异的比例和置信区间以及P-value。
本例图所显示的是第三层级的KEGG代谢途径的差异分析,也可以针对第二或第一层的分级进行分析。
20.基因的差异分析图
除了能对大的基因功能分类和代谢途径进行预测外,我们还能提供精细的功能基因的数量和构成的预测,以及进行样本间以及组间的差异分析,并给出具有统计意义和置信区间的分析结果。
这一分析将我们对于样本群落的差异进一步深入到了每一类基因的层面。
图中不同颜色代表不同的分组,列出了在组间/样本间存在显著差异的每一个功能基因(酶)以及在各组的比例,此外右侧还给出了差异的比例和置信区间以及P-value。
21.在获得标准报告后如果希望单独修改分组或对某些组之间进行显著性差异分析,可以使用STAMP软件在自己的电脑上进行数据分析。
STAMP提供了丰富的统计检验方法和图形化结果的输出。
在使用STAMP之前需要首先准备需要的spf格式文件和样品分组信息表。
在我们的报告中已经将KEGG和KO以及COG的结果文件后经过转换生成了适用于STAMP软件打开的spf格式文件,还有对应的分组信息表文件groupfile.txt。
以下是使用STAMP时的一些相关问题,详细的STAMP使用教程可以参考我们提供的STAMP使用教程。
1、stamp作图用的原始数据的来源?
STAMP可以直接使用来自QIIME的biom文件和PICUST的KEGG和ko文件,groupfile.txt文件的格式为tab-saperatedvalue(tab键隔开的数据)
2、分组问题:
导入数据之后,viewà
grouplegend,在窗口右侧会出现分组栏,根据需要进行分组。
3、Unclassiffied选项中,remainUncla