宏基因组学课件

宏基因组学

Microorganismsandenvironment微生物是重要的分解者，参与营养物质的循环，如碳循环，氮循环，了解微生物的种群组成，对其在环境中发挥的功能具有重要的意义。CommunityFingerprintingmethods基于不同的微生物具有不同的基因分型的可视化方法T-RFLPDGGEARISA共同的缺陷，只能“定性”，而不能定量，而且这个“定性”，并不具有序列信息这样的唯一性。灵敏度不够，丰度低的菌群无法获得很好反映。应用方法AS(16SrDNA)法vsDSS法Flowdiagramoftwonext-generationsequencingapproachesformicrobiomestudies.(A)SequencingofDNAampliconstargetingspecific16SrRNAgenefragments(hypervariableregions).TheobtainedsequencedataarecomparedwithsequencesinsmallsubunitribosomalRNAgenedatabasesandareusedfortaxonomicprofilinganddiversityanalyses.(B)DirectsequencingofrandomDNAfragments,alsocalledmetagenomicshotgunsequencing.Asequencecontigisacontiguous,overlappingsequenceresultingfromthere-assemblyofthesmallDNAfragments.Theobtainedsequencedataarecomparedwithfull-genomereferencedatabasesandareusedtodescribethepredominantfunctionsofthemicrobialcommunities,aswellastoidentifythemicrobialtaxa.Stepsindicatedingraydifferbetweenthetwomethods数据分析流程AS法ClusteringintoOTUsTaxonomyassignmentAlphadiversityBetadiversityOrdinationanalysisClassificationHypothesisgenerationHypothesistestingHypothesisvalidationDSS法FragmentassemblyintocontigBuildinggenecatalogGeneannotationandfunctionpredictionTaxonomyassignmentAlphadiversityBetadiversityDimensionreductionHypothesisgenerationHypothesistestingHypothesisvalidation16SrDNA法Why16s？普遍存在、长度适中、进化缓慢、相对保守。WhyMiseq？Miseq测序平台读长达到2x300PE通量高，V3试剂可产出15Gb的数据，从而降低了单条序列的测序成本，提高了对低丰度菌群的捕获能力无需进行细菌培养，简化了实验操作16srDNA保守区和高变区测序区域选择多样性测序可以同时对不同微生物的16s/18s、ITS等不同区域及功能基因进行测序。根据已发引物或引物比对文献，给出如下推荐：细菌16s区域：真菌18s和ITS区域：MiSeq引物[1]PNAS.2011,108:4680-4687.[2]EnvironmentalMicrobiology.2012,14(9):2457-2466[3]Res.2013,73:5905-5913.[4]TheISMEJournal.2013,7:1391-1401.[5]AcademicPress.1990,315-322.[6]TheISMEJournal.2010,4:1340-1351.24建议测序量[1]Hepatology.2011,54(2):562-572.[2]TheISMEJournal.2012,6:1858-1868.[3]TheISMEJournal,2013,1-13.[4]GenomeRes.2013,23(12):2103-14.[5]ApplEnvironMicrobiol.2012,78(20):7173-84.[6]TheISMEJournal.2013,7:210-220.[7]MarinePollutionBulletin.72(2013)181-187.测序前的二三事目标（Objective）直接点的目标，就是要了解样本中有哪些微生物，它们的组成和丰度了解你的样品土壤水样动物粪便组织Makeadecision（你想获得哪部分微生物的信息）细菌/古菌真菌Important！！这些选择关系到实验的设计样本采集建议水体样本水样的取样深度和范围，可根据研究目的进行确定。水样或略浑浊的水样：取水样≥1L，低温环境下运输至实验室，0.22µm过滤后将滤膜放置冻存管中，干冰运输寄送；水样泥样需提供约5-10g，有必要时需借助采样器，进行采样。参考文献：AntonievanLeeuwenhoek.2013,103:1309-1319.粪便样本取人粪便时，为了避免尿液的污染，可先排尽尿液。用无菌勺子取粪便内部（注：粪便表面与空气接触容易产生变化，不建议取表面样品），一式三份，每份5g以上。取完后将样品存于专门的粪盒，取样之后存于-80℃，干冰运输。鼠粪便需要提供1g以上。

参考文献：PLoSONE.2013,8(10).肠道内容物用无菌解剖刀，在无菌状态下取出腹腔肠道并用PBS冲洗，放置-80℃保存。参考文献：TheISMEJournal.2012,6:2091-2106.18样本采集建议土壤样本土壤样本取样深度和范围，可根据研究目的进行确定。可多点采样法进行土壤样品收集，需-80℃保存，土壤每个样品需要提供约5-10g。参考文献：TheISMEJournal.2013,7:652-659.物体表面微生物样本将物体置于无菌容器内，加入适量的PBS浸没物体，利用摇床等仪器旋转震荡，使表面微生物与物体脱离，收集水样，低温高速离心，收集沉淀。参考文献：PLoSONE.2012,7(5).口腔样本口腔表面样本可采用拭子形式或漱口水形式提供，唾液可直接进行收集后送样，干冰送样。参考文献：Microbiology.2010,156:3282-3287.TerminologyAmplicon：扩增子，通过设计标准引物对16s（或者其他标志基因）的一个或者多个区进行扩增Clustering：聚类，按照一定的相似性阈值将序列分成不同的组OTU（OperationalTaxonomicUnit）：分类操作单元，相当于我们要分析的“个体”，但是对细菌来说OTU和物种并不一致，是基于序列相似性97%的分类。基本策略通过目标基因或者目标基因区域PCR扩增产物（amplicons）构建文库进行测序细菌和古菌适用16SrRNA的一个区或者几个区真菌适合用internaltranscribedspacer（ITS）真核生物适用18SrRNA相关的数据库Greengenes()(16S)RibosomalDatabaseProject()(16S)Silva(http://www.arb-silva.de)16S+18SUnite(http://unite.ut.ee)ITS流程图OTUpicking通过reads序列相似性，将reads聚类形成OTU，一般我们将阈值设置为97%OTU1OTU2OTU3OTU4OTU5AssigntaxonomicinformationOTU1OTU2OTU3OTU4OTU5通过数据库可以为OTU的代表性序列提供分类学注释RDPSILVAGreengenes…genus1genus2genus3genus4genus5单样本物种组成单样品中不同分类水平下的物种组成分布。通过饼图可以直观的展示优势菌群。ClassGenusOTU分布韦恩图

Α-多样性分析在获得OTU数据之后，可以进行α-多样性分析α-多样性是和β-多样性相对而生的概念α-多样性反映的是样本内的物种丰度以及组成。Β-多样性则反映的是样本间的物种丰度和组成的差异比较。α-多样性衡量方法包括Chao指数ACE指数CoverageShannon曲线Richnessrarefaction曲线Α-多样性衡量指标测序深度指数：

Coverage-theGood’scoverage(/wiki/Coverage)计算菌群多样性(Communitydiversity)：Shannon-theShannonindex(/wiki/Shannon)Simpson-theSimpsonindex(/wiki/Simpson)计算菌群丰度(Communityrichness)：Chao-theChao1estimator(/wiki/Chao)Ace-theACEestimator(/wiki/Ace)Α-多样性衡量指标Rarefactioncurve(稀释性曲线)

用途：用来比较测序数据量不同的样本中物种的丰富度，也可以用来说明样本的测序数据量是否合理。原理：从样本中随机抽取一定数量的个体，统计这些个体所代表的物种数目，并以个体数与物种数来构建曲线。采用对序列进行随机抽样的方法，以抽到的序列数与它们所能代表OTU的数目做图。解读：当曲线趋向平坦时，说明测序数据量合理，更多的数据量只会产生少量新的OTU，反之则表明继续测序还可能产生较多新的OTU。因此，通过作稀释性曲线，可得出样本的测序深度情况。绘图方法：使用97%相似度构建OTU，利用qiime或mothur做rarefaction分析，利用R语言工具制作曲线图。RarefactioncurveShannon-WienerCurve用途：反映样本中微生物多样性利用各样本的测序量在不同测序深度时的微生物多样性指数构建曲线，以此反映各样本在不同测序数量时的微生物多样性。曲线解读：当曲线趋向平坦时，说明测序数据量足够大，可以反映样本中绝大多数的微生物信息。计算方法：使用97%相似度构建OTU，利用mothur计算不同随机抽样下的shannon值，利用R语言工具制作曲线图。Shannon-WienerCurveRank-AbundanceCurve用途：分析多样性构建原理：统计单一样本中，每一个OTU所含的序列数，将OTUs按丰度（所含有的序列条数）由大到小等级排序，再以OTU等级为横坐标，以每个OTU中所含的序列数（也可用OTU中序列数的相对百分含量）为纵坐标做图。结果解读：用来解释多样性的两个方面，即物种丰度（abundance）和物种均匀度（eveness）。在水平方向，物种的丰度由曲线的宽度来反映，物种的丰度越高，曲线在横轴上的范围越大；曲线的形状（平滑程度）反映了样本中物种的均度，曲线越平缓，物种分布越均匀。绘图方法：利用R语言工具制作曲线图。Rank-AbundanceCurve注：横坐标表示OTU等级，“500”代表样本中按照丰度排列第500位的OTU；纵坐标表示该等级OTU中序列数的相对百分含量，即属于该OTU的序列数除以总序列数，纵坐标轴上数字，例如“100”代表相对丰度为100%，“10”代表相对丰度为10%，依次类推。beta多样性分析-群落组成最简单的比较是在各个分类水平下菌群的丰度的比较。多个样本菌群组成以及丰度的变化可以通过条形图展示。样品间物种丰度热图各纵列为样本，行为该分类水平下物种物种丰度热图，用物种丰度矩阵绘制，图中每一列代表一个样本，行代表菌群，颜色块代表相对物种丰度值，即Z值，颜色越绿表示相对丰度越高，颜色越红反之。另外热图对样本及菌群均做了聚类，样本菌群分布越类似则样本距离越近，在图上方聚类树中的位置越靠近；菌群聚类亦如此，在不同样本中具有相似的分布模式菌群被聚在一起，这些菌群可能就具有相关性。排序分析（Ordination）每个样本都具有自己的alpha多样性的特征,如菌群和丰度组成。那样本之间可以比较吗？如何比较？手段是排序（Ordination）分析。Ordination是通过数学方法，将各个样本的菌群和丰度转换为距离矩阵，并将生态距离映射到二维坐标轴上点之间的距离。Uncontrainedmethods依赖于排序的方法，大都依赖于降维。为什么要降维？假设两个样本都只有一种物种，丰度不一样假设每个样本有2个…3个不同的物种假设每个样本有1000个物种所幸的是，我们的数学家给我们带来了降维的方法PCAPCoANMDSNonmetricMultidimensionalScalingPCA方法主成分分析是一种对数据进行简化分析的技术：--去除噪音和冗余--将原有的复杂数据降维--揭示隐藏在复杂数据背后的简单结构简单且无参数限制。通过分析不同样本OTU（97%相似性）组成可以反映样本间的差异和距离PCoA方法主坐标分析和主成分分析类似，也是一种非约束性的数据降维分析方法，可用来研究样本群落组成的相似性或差异性。它首先假定有多个梯度（gradient），但是重要性依次递减，计算距离的方式是ecologically-relevantdistances如Bray-Curtisdistance等，而如果用欧氏距离就会变成PCA。在计算获得了距离矩阵的坐标后，我们就能在2d或者3d图上通过点之间的距离直观的反映样本之间的距离。NMDS方法目标：用最低的维度可靠的代表样本之间的距离.依赖于距离的最稳健的方法非线性模型维度可以指定最大限度保持了点之间的距离，但是可视化的时候点比较铺开，视觉上聚类效果不明显软件：Qiime计算beta多样性距离矩阵，R语言vegan软件包作NMDS分析和作图。Unifrac方法基于系统发育树的样品分类。分析中首先将OTU种子序列比对到Greengenes核心16SrDNA序列中，根据多序列队列构建以种子序列为节点的系统发育树，而后通过进化树计算Unifracmetric，Unifracmetric是通过计算进化树中不同环境样本间Uniquebranch长度总和获得。Unifracmetric值在0‐1之间，值越小说明样本间相似度越高。计算Unifrac值时有两种方式，一种为未加权重方式（只考虑物种有或无），另外一种为加权重计算方式（考虑物种丰度）。HAONOHAOPLANTYANNOYANPLANTHAONO00.3565040.1678750.551072HAOPLANT0.35650400.3501960.406638YANNO0.1678750.35019600.531851YANPLANT0.5510720.4066380.5318510OTU_IDHAONOHAOPLANTYANNOYANPLANTTotalOtu110001Otu210001Otu312803380466Otu410078248481311Otu52015017Otu62800028Otu710001Contrainedmethods非限制性排序仅仅利用了物种（species）的数据，而限制性排序方法同时利用了物种数据以及环境参数矩阵（回归分析）。限制性排序的目的是尝试解释物种组成差异和环境因子之间的关联。其方法主要包括RDACCARDA方法RedundancyAnalysis(RDA)和PCA类似，基于线性模型物种和环境因子以箭头表示RDA适用于shortgradient,反应short-term实验RDA可以反映物种的正相关关系坐标轴能解释的方差依次递减（和PCA类似）CCA方法CanonicalCorrespondenceAnalysis(CCA),基于单峰模型CCA聚焦于样本的组成，也就是相对丰度CCA无法检测物种间的正相关性CCA是CA和多元回归的结合使得物种得分和样本得分的关联最大化样本得分受到解释变量的线性组合限制特征值比CA小RDA/CCA用途：多元直接梯度分析，反映菌群与环境因子之间关系。排序方法：将对应分析与多元回归分析相结合，每一步计算均与环境因子进行回归，RDA基于线性模型，CCA则基于单峰模型。分析可以检测环境因子、样本、菌群三者之间的关系或者两两之间的关系。RDA或CCA模型的选择：用species-sample数据（97%相似性的样本OTU表）做DCA分析，看分析结果中Lengthsofgradient的第

一轴的大小，如果大于4.0，就应该选CCA，如果3.0-4.0之间，选RDA和CCA均可，如果小于3.0，RDA的结果要好于CCA。要点：通过bioenv函数判断环境因子与样本群落分布差异的最大Pearson相关系数，通过最大相关系数得到环境因子子集；将样本

物种分布表与环境因子或环境因子子集分别做CCA或者RDA分析；

通过类似于ANOVA的permutest分析判断CCA或者RDA分析显著性。软件：R语言vegan包中rda或者cca分析和作图。RDA/CCAClassificationBuildcluster:RpackagehclustManipulateanddisplaytree+dragram:Rpackagerect.hcclust,veganfuctionsordihull,ordispiderandordiellipse.HypothesisTestadonis:MultivariateANOVAbasedondissimilarities--similartoANOVA--testifgroup“means”aresignificantlydifferentbetadisper--testifthegroup“variances”areequalmantel--testifgroupsaresignificantlycorrespondingtoenvironmentalfactorsIntegrativeAnalysisEnvironmentsQIIMEa.Principalcoordinatesanalysisplotsformice,H-miceandtwins.Colorscorrespondtoseparatesamplesbyspeciesandtimepoint,andareconsistentthroughoutthepanels.UnweightedUniFracdistancehistogramsbetweenthedataforfecalmicrobiotaofhumantwins;humandonorsfortheH-micestudy;day56post-transplantH-miceonalow-fat(LF)andplantpolysaccharide–rich(PP)diet;day1H-mice(LFandPPdiet);andday0H-mice.Taxonomicclassificationsarepresentedattheclasslevel.b.c.AlphadiversityrarefactionplotsofphylogeneticdiversityfortheH-micesamples.d.OTUnetworkconnectivityofH-micetimeseriesdata.CONV-D,conventionalizedmice;CONV-R,conventionallyraisedmice;andGF,germ-freemice.QIIMEDSS法测序、组装物种丰度计算构建基因目录使用MetaGene对contig进行ORF预测；选择长度>100的DNA序列翻译成氨基酸序列；将预测的基因序列用CD-HIT进行聚类(参数：95%identity,90%coverage)在每个类群中选择最长的基因作为代表序列，即构成non-redundantgenecatalog(非冗余基因集)。基因注释和功能预测NR:KEGG:http://www.genome.jp/kegg/eggNOG:http://eggnogdb.embl.de/

FGR：单样品(phylumlevel)分类谱多样品(Phylumlevel)分类谱单样本COG功能分类统计图KEGG代谢途径注：矩形框表示催化反应进行的酶（红色边框：当前数据比对上的基因对应的KO相关的酶），框内编号即EC编号；圆圈表示代谢产物（酶促反应的反应物或产物）；实线箭头的方向表示酶促反应进行的方向；虚线箭头表示此产物可再通过中间产物与其它代谢途径发生关系；圆角矩形框代表其它代谢途径。QiM,etal.(2011).PLoSONE宏转录组物种组成Top10PhymmBL%NR%KEGG%Firmicutes4945623.894794637351.244852381245.22268Proteobacteria7074234.179171749819.33631109921.07872Bacteroidetes179778.6856321730819.12634822115.61295Actinobacteria126496.11139625452.81237225654.871332Euryarchaeota109085.27022722822.52174219963.790713Spirochaetes50122.4215615901.75704210131.923844Fibrobacteres29431.4219185360.5923119531.809895Chloroflexi4590.221767510.0563583490.662805Fusobacteria7210.3483532100.2320623250.617225Ascomycota155447.5101223750.4143971920.364638Sum18641190.064948876898.093775052595.9548宏转录组物种组成vs功能特征宏转录组代谢途径ECOrg#read（THL）Actinobacteria3Alphaproteobacteria10Gammaproteobacteria7Betaproteobacteria6Spirochaetes7Euryarchaeota3Firmicutes3857（BHBD）Actinobacteria3Gammaproteobacteria1Betaproteobacteria5Fusobacteria4Euryarchaeota2Spirochaetes1Firmicutes225（CRO）Fusobacteria2Deltaproteobacteria1Spirochaetes18Firmicutes9（BCD）Betaproteobacteria1Bacteroidetes4Fusobacteria2Spirochaetes2

Firmicutes42宏转录组代谢途径QiM,etal.(2011).PLoSONE宏转录组聚类分析HumanMicrobiomeHowhumanarewe?SimilaritiesbetweenIndividualsSimilaritiesbetweenIndividualsEachbodysitecontainsadistinctMicrobiomeACatalogofReferenceGenomesfromtheHumanMicrobiomeFig.1.Phylogenetictreeof16SrDNAsequences.Thetreewascreatedusing~150016SrDNAsrepresentingsinglespecies.OrganismssequencedaspartoftheHMPprojectarehighlightedinblue.Additionalcoloringindicatesseparationbyphylum:yellow,Actinobacteria;darkgreen,Bacteroidetes;lightgreen,Cyanobacteria;red,Firmicutes;cyan,Fusobacteria;darkred,Planctomycetes;gray,Proteobacteria;magenta,Spirochaetes;lightpink,TM7;tan,Tenericutes.Thepurposeofthisanalysisisnotthedetailsofthebranchingstructure(whichincludeminorknownartifacts),buttheoveralldistributionoftheHMPstrains(inblue)aroundthetreeoflife./content/sci/328/5981/994.full.pdfMetaHIT127AhumangutmicrobialgenecatalogueestablishedbymetagenomicsequencingVol464|4March2010|doi:10.1038/nature08821128PredictedORFsinthehumangutmicrobiomea,Numberofuniquegenesasafunctionoftheextentofsequencing.ThegeneaccumulationcurvecorrespondstotheSobs(MaoTau)values(numberofobservedgenes),calculatedusingEstimateS21(version8.2.0)onrandomlychosen100samplesb,Coverageofgenesfrom89frequentgutmicrobialspeciesc,Numberoffunctionscapturedbynumberofsamplesinvestigated--known(wellcharacterized)orthologousgroups(OGs;bottom),--knownplusunknownorthologousgroups(including,forexample,putative,predicted,conservedhypotheticalfunctions;middle)--orthologousgroupsplusnovelgenefamilies(.20proteins)recoveredfromthemetagenome(top).Boxesdenotetheinterquartilerange(IQR)betweenthefirstandthirdquartiles(25thand75thpercentiles,respectively)andthelineinsidedenotesthemedian.Whiskersdenotethelowestandhighestvalueswithin1.5timesIQRfromthefirstandthirdquartiles,respectively.Circlesdenoteoutliersbeyondthewhiskers.Relativeabundanceof57frequentmicrobialgenomesamongindividualsofthecohort.BacterialspeciesabundancedifferentiatesIBDpatientsandhealthyindividualsPrincipalcomponentanalysiswithhealthstatusasinstrumentalvariables,basedontheabundanceof155specieswith$1%genomecoveragebytheIlluminareadsinatleast1individualofthecohort,wascarriedoutwith14healthyindividualsand25IBDpatients(21ulcerativecolitisand4Crohn’sdisease)fromSpain(SupplementaryTable1).Twofirstcomponents(PC1andPC2)wereplottedandrepresented7.3%ofwholeinertia.Individuals(representedbypoints)wereclusteredandcentreofgravitycomputedforeachclass;P-valueofthelinkbetweenhealthstatusandspeciesabundancewasassessedusingaMonte-Carlotest(999replicates).Characterizationoftheminimalgutgenomeandmetagenomea,ProjectionoftheminimalgutgenomeontheKEGGpathwaysusingtheiPathtool38.b,Functionalcompositionofthem