微生物群落多樣性測(cè)序與功能分析

《微生物群落多樣性測(cè)序與功能分析》由會(huì)員分享，可在線閱讀，更多相關(guān)《微生物群落多樣性測(cè)序與功能分析（27頁珍藏版）》請(qǐng)?jiān)谘b配圖網(wǎng)上搜索。

微生物群落多樣性測(cè)序與功能分析 微生物群落測(cè)序是指對(duì)微生物群體進(jìn)行高通量測(cè)序，通過分析測(cè)序序列的構(gòu)成分析特定環(huán)境中微生物群體的構(gòu)成情況或基因的組成以及功能。借助不同環(huán)境下微生物群落的構(gòu)成差異分析我們可以分析微生物與環(huán)境因素或宿主之間的關(guān)系，尋找標(biāo)志性菌群或特定功能的基因。對(duì)微生物群落進(jìn)行測(cè)序包括兩類，一類是通過16s rDNA，18s rDNA，ITS區(qū)域進(jìn)行擴(kuò)增測(cè)序分析微生物的群體構(gòu)成和多樣性；還有一類是宏基因組測(cè)序，是不經(jīng)過分離培養(yǎng)微生物，而對(duì)所有微生物DNA進(jìn)行測(cè)序，從而分析微生物群落構(gòu)成，基因構(gòu)成，挖掘有應(yīng)用價(jià)值的基因資源。 以16s rDNA擴(kuò)增進(jìn)行測(cè)序分析主要用于微生物群落多樣性和構(gòu)成的分析，目前的生物信息學(xué)分析也可以基于16s rDNA的測(cè)序?qū)ξ⑸锶郝涞幕驑?gòu)成和代謝途徑進(jìn)行預(yù)測(cè)分析，大大拓展了我們對(duì)于環(huán)境微生物的微生態(tài)認(rèn)知。 目前我們根據(jù)16s的測(cè)序數(shù)據(jù)可以將微生物群落分類到種（species）（一般只能對(duì)部分菌進(jìn)行種的鑒定），甚至對(duì)亞種級(jí)別進(jìn)行分析， 幾個(gè)概念： 16S rDNA（或16S rRNA）：16S rRNA基因是編碼原核生物核糖體小亞基的基因，長(zhǎng)度約為1542bp，其分子大小適中，突變率小，是細(xì)菌系統(tǒng)分類學(xué)研究中最常用和最有用的標(biāo)志。16S rRNA基因序列包括9個(gè)可變區(qū)和10個(gè)保守區(qū)，保守區(qū)序列反映了物種間的親緣關(guān)系，而可變區(qū)序列則能體現(xiàn)物種間的差異。16S rRNA基因測(cè)序以細(xì)菌16S rRNA基因測(cè)序?yàn)橹?核心是研究樣品中的物種分類、物種豐度以及系統(tǒng)進(jìn)化。 OTU：operational taxonomic units (OTUs)在微生物的免培養(yǎng)分析中經(jīng)常用到，通過提取樣品的總基因組DNA，利用16S rRNA或ITS的通用引物進(jìn)行PCR擴(kuò)增，通過測(cè)序以后就可以分析樣品中的微生物多樣性，那怎么區(qū)分這些不同的序列呢，這個(gè)時(shí)候就需要引入operational taxonomic units，一般情況下，如果序列之間，比如不同的 16S rRNA序列的相似性高于97%就可以把它定義為一個(gè)OTU，每個(gè)OTU對(duì)應(yīng)于一個(gè)不同的16S rRNA序列，也就是每個(gè)OTU對(duì)應(yīng)于一個(gè)不同的細(xì)菌（微生物）種。通過OTU分析，就可以知道樣品中的微生物多樣性和不同微生物的豐度。 測(cè)序區(qū)段：由于16s rDNA較長(zhǎng)（1.5kb），我們只能對(duì)其中經(jīng)常變化的區(qū)域也就是可變區(qū)進(jìn)行測(cè)序。16s rDNA包含有9個(gè)可變區(qū)，分別是v1-v9。一般我們對(duì)v3-v4雙可變區(qū)域進(jìn)行擴(kuò)增和測(cè)序，也有對(duì)v1-v3區(qū)進(jìn)行擴(kuò)增測(cè)序。 工具/原料 16s rDNA測(cè)序首先需要提取環(huán)境樣品的DNA，這些DNA可以來自土壤、糞便、空氣或水體等任何來源。 提取DNA后需要經(jīng)過質(zhì)檢和純化，一般16s rDNA測(cè)序擴(kuò)增對(duì)DNA的總量要求并不高，總量大于100ng，濃度大于10ng/ul一般都可以滿足要求。如果是來自和寄主共生的環(huán)境如昆蟲的腸道微生物，提取時(shí)可能包括了寄主本身的大量DNA，對(duì)DNA的總量要求會(huì)提高。微生物菌群多樣性測(cè)序受DNA提取和擴(kuò)增影響很大，不同的擴(kuò)增區(qū)段和擴(kuò)增引物甚至PCR循環(huán)數(shù)的差異都會(huì)對(duì)結(jié)果有所影響。因而建議同一項(xiàng)目不同樣品的都采用相同的條件和測(cè)序方法，這樣相互之間才存在可比性。 完成PCR之后的產(chǎn)物一般可以直接上測(cè)序儀測(cè)序，在上機(jī)測(cè)序前我們需要對(duì)所有樣本進(jìn)行定量和均一化，通常要進(jìn)行熒光定量PCR。完成定量的樣品混合后就可以上機(jī)測(cè)序。 16s rDNA測(cè)序目前可以采用多種不同的測(cè)序儀進(jìn)行測(cè)序，包括羅氏的454，Illumina的MiSeq，Life的PGM或Pacbio的RSII三代測(cè)序儀。不同的儀器各有優(yōu)缺點(diǎn)，目前最主流的是Illumina公司的MiSeq，因?yàn)槠湓谕?、長(zhǎng)度和價(jià)格三者之間最為平衡。MiSeq測(cè)序儀可以產(chǎn)生2x300bp的測(cè)序讀長(zhǎng)，一次可以產(chǎn)生15Gb的測(cè)序數(shù)據(jù)遠(yuǎn)遠(yuǎn)大于其他測(cè)序儀的測(cè)序通量。 方法/步驟 1. 1 16s rDNA分析基本流程： 2. 2 原始數(shù)據(jù)處理： 原始測(cè)序數(shù)據(jù)需要去除接頭序列，并將雙端測(cè)序序列進(jìn)行拼接成單條序列。 根據(jù)測(cè)序barcode序列區(qū)分不同的樣本序列。 過濾低質(zhì)量序列和無法比對(duì)到16s rDNA數(shù)據(jù)庫的序列。 3. 3 OTU分類和統(tǒng)計(jì)： OTU(operational taxonomic units) 是在系統(tǒng)發(fā)生學(xué)研究或群體遺傳學(xué)研究中，為了便于進(jìn)行分析，人為給某一個(gè)分類單元（品系，種，屬，分組等）設(shè)置的同一標(biāo)志。通常按照 97% 的相似性閾值將序列劃分為不同的 OTU，每一個(gè) OTU 通常被視為一個(gè)微生物物種。相似性小于97%就可以認(rèn)為屬于不同的種，相似性小于93%-95%，可以認(rèn)為屬于不同的屬。樣品中的微生物多樣性和不同微生物的豐度都是基于對(duì)OTU的分析。 使用QIIME（version 1.8.0）工具包進(jìn)行統(tǒng)計(jì)注釋。 使用QIIME（version 1.9.0, http://bio.cug.edu.cn/qiime/）的ucluster方法根據(jù)97%的序列相似度將所有序列進(jìn)行同源比對(duì)并聚類成operational taxonomic units (OTUs)。然后與數(shù)據(jù)庫GreenGenes（version gg_13_8, http://greengenes.lbl.gov/cgi-bin/JD_Tutorial/nph-16S.cgi）進(jìn)行比對(duì)，比對(duì)方法uclust，identity 0.9 。 然后對(duì)每個(gè)OTUs進(jìn)行reads數(shù)目統(tǒng)計(jì)。 下面的2個(gè)表，其中一個(gè)表是對(duì)每個(gè)樣本的測(cè)序數(shù)量和OTU數(shù)目進(jìn)行統(tǒng)計(jì)，并且在表栺中列出了測(cè)序覆蓋的完整度（顯示前10個(gè)樣本）。 另一個(gè)表是對(duì)每個(gè)樣本在分類字水平上的數(shù)量進(jìn)行統(tǒng)計(jì)，并且在表栺中列出了在每個(gè)分類字水平上的物種數(shù)目（顯示前10個(gè)樣本）。 可以看到絕大部分的OTU都分類到了屬（Genus），也有很多分類到了種（Species）。但是仍然有很多無法完全分類到種一級(jí)，這是由于環(huán)境微生物本身存在非常豐富的多樣性，還有大量的菌仍然沒有被測(cè)序和發(fā)現(xiàn)。 測(cè)序數(shù)目統(tǒng)計(jì)表主要是對(duì)每個(gè)樣本的測(cè)序數(shù)量和OTU數(shù)目進(jìn)行統(tǒng)計(jì)，并且在表格中列出了測(cè)序覆蓋的完整度（顯示前10個(gè)樣本，如果樣本超過10個(gè)，請(qǐng)查看結(jié)果中otu_stat.txt文件） 其中 SampleName表示樣本名稱；SampleSize表示樣本序列總數(shù)；OTUsNumber表示注釋上的OTU數(shù)目；OTUsSeq表示注釋上OTU的樣本序列總數(shù)。 Coverage是指各樣品文庫的覆蓋率，其數(shù)值越高，則樣本中序列沒有被測(cè)出的概率越低。該指數(shù)實(shí)際反映了本次測(cè)序結(jié)果是否代表樣本的真實(shí)情況。 計(jì)算公式為：C=1-n1/N 其中n1 = 只含有一條序列的OTU的數(shù)目； N = 抽樣中出現(xiàn)的總的序列數(shù)目。 分類水平統(tǒng)計(jì)表主要是對(duì)每個(gè)樣本在分類學(xué)水平上的數(shù)量進(jìn)行統(tǒng)計(jì)，并且在表格中列出了在每個(gè)分類學(xué)水平上的物種數(shù)目（只顯示前10個(gè)樣本，如果樣本超過10個(gè)，請(qǐng)查看結(jié)果中taxon_all.txt文件） 其中SampleName表示樣本名稱；Phylum表示分類到門的OTU數(shù)量；Class表示分類到綱的OTU數(shù)量；Order表示分類到目的OTU數(shù)量；Family表示分類到科的OTU數(shù)量；Genus表示分類到屬的OTU數(shù)量；Species表示分類到種的OTU數(shù)量。 4. 4 我們還可以對(duì)這些種屬的構(gòu)成進(jìn)行柱狀圖顯示： 橫坐標(biāo)中每一個(gè)條形圖代表一個(gè)樣本，縱坐標(biāo)代表該分類層級(jí)的序列數(shù)目或比例。同一種顏色代表相同的分類級(jí)別。圖中的每根柱子中的顏色表示該樣本在不同級(jí)別（門、綱、目等）的序列數(shù)目，序列數(shù)目只計(jì)算級(jí)別最低的分類，例如在屬中計(jì)算過了，則在科中則不重復(fù)計(jì)算。 Q: 為什么要選擇V3-V4區(qū)的測(cè)序長(zhǎng)度？為什么有些文獻(xiàn)是V6區(qū)，有什么區(qū)別？ A: 16S rRNA總長(zhǎng)約1540 bp，包含 9個(gè)可變區(qū)。由于高通量測(cè)序的測(cè)序長(zhǎng)度的限制，不可能將16S rRNA的9個(gè)可變區(qū)全部測(cè)序，所以在PCR擴(kuò)增時(shí)往往只能選擇1-3個(gè)可變區(qū)作為擴(kuò)增片段。Kozich 等評(píng)估了Miseq測(cè)序儀分析的不同16S rRNA可變區(qū)的準(zhǔn)確性發(fā)現(xiàn)，測(cè)定 V4 區(qū)效果最佳。根據(jù)我們的測(cè)序長(zhǎng)度，v3-v4區(qū)是最佳選擇。 5. 5 我們還需要對(duì)樣本之間或分組之間的OTU進(jìn)行比較獲得韋恩圖： 注意，韋恩圖目前一般最多只能顯示5個(gè)樣本或分組，過多的樣本無法無法進(jìn)行韋恩圖繪制 6. 6 樣品構(gòu)成豐度： 稀釋曲線 微生物多樣性分析中需要驗(yàn)證測(cè)序數(shù)據(jù)量是否足以反映樣品中的物種多樣性，稀釋曲線（豐富度曲線）可以用來檢驗(yàn)這一指標(biāo)。 稀釋曲線是用來評(píng)價(jià)測(cè)序量是否足以覆蓋所有類群，并間接反映樣品中物種的豐富程度。稀釋曲線是利用已測(cè)得16S rDNA序列中已知的各種OTU的相對(duì)比例，來計(jì)算抽取n個(gè)（n小于測(cè)得reads序列總數(shù)）reads時(shí)出現(xiàn)OTU數(shù)量的期望值，然后根據(jù)一組n值（一般為一組小于總序列數(shù)的等差數(shù)列）與其相對(duì)應(yīng)的OTU數(shù)量的期望值做出曲線來。當(dāng)曲線趨于平緩或者達(dá)到平臺(tái)期時(shí)也就可以認(rèn)為測(cè)序深度已經(jīng)基本覆蓋到樣品中所有的物種；反之，則表示樣品中物種多樣性較高，還存在較多未被測(cè)序檢測(cè)到的物種。 下圖中的稀釋曲線 橫坐標(biāo)代表隨機(jī)抽取的序列數(shù)量；縱坐標(biāo)代表觀測(cè)到的OTU數(shù)量。樣本曲線的延伸終點(diǎn)的橫坐標(biāo)位置為該樣本的測(cè)序數(shù)量，如果曲線趨于平坦表明測(cè)序已趨于飽和，增加測(cè)序數(shù)據(jù)無法再找到更多的OTU；反之表明不飽和，增加數(shù)據(jù)量可以發(fā)現(xiàn)更多OTU。 7. 7 Shannon-Winner曲線 Shannon-Wiener 曲線，是利用shannon指數(shù)來進(jìn)行繪制的，反映樣品中微生物多樣性的指數(shù)，利用各樣品的測(cè)序量在不同測(cè)序深度時(shí)的微生物多樣性指數(shù)構(gòu)建曲線，以此反映各樣本在不同測(cè)序數(shù)量時(shí)的微生物多樣性。 當(dāng)曲線趨向平坦時(shí)，說明測(cè)序數(shù)據(jù)量足夠大，可以反映樣品中絕大多數(shù)的微生物物種信息。 與上圖一樣，橫坐標(biāo)代表隨機(jī)抽取的序列數(shù)量；縱坐標(biāo)代表的是反映物種多樣性的Shannon指數(shù)。 樣本曲線的延伸終點(diǎn)的橫坐標(biāo)位置為該樣本的測(cè)序數(shù)量，如果曲線趨于平坦表明測(cè)序已趨于飽和，增加測(cè)序數(shù)據(jù)無法再找到更多的OTU；反之表明不飽和，增加數(shù)據(jù)量可以發(fā)現(xiàn)更多OTU。 其中曲線的最高點(diǎn)也就是該樣本的Shannon指數(shù)，指數(shù)越高表明樣品的物種多樣性越高。 Q: Shannon指數(shù)怎么算的？ A: Shannon指數(shù)公式： 其中，Sobs=實(shí)際測(cè)量出的OTU數(shù)目；ni=含有i 條序列的OTU數(shù)目；N=所有的序列數(shù)。 8. 8 Rank-Abundance曲線 用于同時(shí)解釋樣品多樣性的兩個(gè)方面，即樣品所含物種的豐富程度和均勻程度。 物種的豐富程度由曲線在橫軸上的長(zhǎng)度來反映，曲線越寬，表示物種的組成越豐富； 物種組成的均勻程度由曲線的形狀來反映，曲線越平坦，表示物種組成的均勻程度越高。 一般超過20個(gè)樣本圖就會(huì)變得非常復(fù)雜而且不美觀，所以一般20個(gè)樣本以下會(huì)做該圖，圖片保存為結(jié)果目錄中rank.pdf。 橫坐標(biāo)代表物種排序的數(shù)量；縱坐標(biāo)代表觀測(cè)到的相對(duì)豐度。 樣本曲線的延伸終點(diǎn)的橫坐標(biāo)位置為該樣本的物種數(shù)量，如果曲線越平滑下降表明樣本的物種多樣性越高，而曲線快速陡然下降表明樣本中的優(yōu)勢(shì)菌群所占比例很高，多樣性較低。 9. 9 Alpha多樣性（樣本內(nèi)多樣性） Alpha多樣性是指一個(gè)特定區(qū)域或者生態(tài)系統(tǒng)內(nèi)的多樣性，常用的度量指標(biāo)有Chao1 豐富度估計(jì)量（Chao1 richness estimator） 、香農(nóng) - 威納多樣性指數(shù)（Shannon-wiener diversity index）、辛普森多樣性指數(shù)（Simpson diversity index）等。 計(jì)算菌群豐度：Chao、ace； 計(jì)算菌群多樣性：Shannon、Simpson。 Simpson指數(shù)值越大，說明群落多樣性越高；Shannon指數(shù)越大，說明群落多樣性越高。表中顯示前10個(gè)樣本，如果樣本大于10個(gè)，詳見結(jié)果目錄中的alpha_div.txt。 Q: 能不能解釋下每個(gè)指數(shù)（如chao1、shannon）？ A: Chao1：是用chao1 算法估計(jì)群落中含OTU 數(shù)目的指數(shù)，chao1 在生態(tài)學(xué)中常用來估計(jì)物種總數(shù)，由Chao (1984) 最早提出。Chao1值越大代表物種總數(shù)越多。 Schao1=Sobs+n1(n1-1)/2(n2+1) 其中Schao1為估計(jì)的OTU數(shù)，Sobs為觀測(cè)到的OTU數(shù)，n1為只有一條序列的OTU數(shù)目，n2為只有兩條序列的OTU數(shù)目。 Shannon：用來估算樣品中微生物的多樣性指數(shù)之一。它與 Simpson 多樣性指數(shù)均為常用的反映 alpha 多樣性的指數(shù)。Shannon值越大，說明群落多樣性越高。 Ace：用來估計(jì)群落中含有OTU 數(shù)目的指數(shù)，由Chao 提出，是生態(tài)學(xué)中估計(jì)物種總數(shù)的常用指數(shù)之一，與Chao1 的算法不同。 Simpson：用來估算樣品中微生物的多樣性指數(shù)之一，由Edward Hugh Simpson ( 1949) 提出，在生態(tài)學(xué)中常用來定量的描述一個(gè)區(qū)域的生物多樣性。Simpson 指數(shù)值越大，說明群落多樣性越高。 辛普森多樣性指數(shù)=隨機(jī)取樣的兩個(gè)個(gè)體屬于不同種的概率 =1-隨機(jī)取樣的兩個(gè)個(gè)體屬于同種的概率 10. 10 Beta多樣性分析（樣品間差異分析） Beta多樣性度量時(shí)空尺度上物種組成的變化,是生物多樣性的重要組成部分,與許多生態(tài)學(xué)和進(jìn)化生物學(xué)問題密切相關(guān),因此在最近10年間成為生物多樣性研究的熱點(diǎn)問題之一。 PCoA分析 PCoA（principal co-ordinates analysis）是一種研究數(shù)據(jù)相似性或差異性的可視化方法，通過一系列的特征值和特征向量進(jìn)行排序后，選擇主要排在前幾位的特征值，PCoA 可以找到距離矩陣中最主要的坐標(biāo)，結(jié)果是數(shù)據(jù)矩陣的一個(gè)旋轉(zhuǎn)，它沒有改變樣品點(diǎn)之間的相互位置關(guān)系，只是改變了坐標(biāo)系統(tǒng)。通過PCoA 可以觀察個(gè)體或群體間的差異。 每一個(gè)點(diǎn)代表一個(gè)樣本，相同顏色的點(diǎn)來自同一個(gè)分組，兩點(diǎn)之間距離越近表明兩者的群落構(gòu)成差異越小。PCoA有多張圖，分別代表的PCoA1-2,2-3,3-1。 11. 11 NMDS分析（非度量多維尺度分析） NMDS（Nonmetric Multidimensional Scaling）常用于比對(duì)樣本組之間的差異，可以基于進(jìn)化關(guān)系或數(shù)量距離矩陣。 橫軸和縱軸：表示基于進(jìn)化或者數(shù)量距離矩陣的數(shù)值 在二維表中成圖。 與PCA分析的主要差異在于考量了進(jìn)化上的信息。 每一個(gè)點(diǎn)代表一個(gè)樣本，相同顏色的點(diǎn)來自同一個(gè)分組，兩點(diǎn)之間距離越近表明兩者的群落構(gòu)成差異越小。 12. 12 PCA分析 主成分分析PCA（Principal component analysis）是一種研究數(shù)據(jù)相似性或差異性的可視化方法，通過一系列的特征值和特征向量進(jìn)行排序后，選擇主要的前幾位特征值，采取降維的思想，PCA 可以找到距離矩陣中最主要的坐標(biāo)，結(jié)果是數(shù)據(jù)矩陣的一個(gè)旋轉(zhuǎn)，它沒有改變樣品點(diǎn)之間的相互位置關(guān)系，只是改變了坐標(biāo)系統(tǒng)。詳細(xì)關(guān)于主成分分析的解釋推薦大家看一篇文章，http://blog.csdn.net/aywhehe/article/details/5736659 。通過PCA 可以觀察個(gè)體或群體間的差異。 每一個(gè)點(diǎn)代表一個(gè)樣本，相同顏色的點(diǎn)來自同一個(gè)分組，兩點(diǎn)之間距離越近表明兩者的群落構(gòu)成差異越小。 以上三個(gè)圖可能遇到的問題： 1：PCA，PcoA，NMDS分析分別是基于什么數(shù)據(jù)畫的？ 回答：PCA，PcoA，NMDS分析均是基于OTU分類taxon數(shù)據(jù)所畫，用的是R語言Vegan包中的相關(guān)函數(shù)畫成，其中PcoA與NMDS還要基于樣本之間的距離矩陣才能畫成。 2：PCA分析如果圖中大部分點(diǎn)集中在一起，少數(shù)點(diǎn)在很遠(yuǎn)的外圍，是什么原因造成的？ 回答：是因?yàn)闃颖綩TU分類時(shí)候，少數(shù)樣本某些菌含量特別高所造成，導(dǎo)致這些樣本偏離正常范圍，建議單獨(dú)拿出這些樣本觀察，看是否是實(shí)驗(yàn)錯(cuò)誤。 3：PCA分析時(shí)，不是有PC1，PC2，PC3三個(gè)坐標(biāo)嗎？是給出三張圖嗎？還是三維立體圖？ 回答：PCA作圖時(shí)，會(huì)得出PC1，PC2，PC3三個(gè)坐標(biāo)，可以根據(jù)PC12,PC13,PC23分別作圖，一般給出的是PC12的圖，當(dāng)PC12圖質(zhì)量不好，看不出明顯的樣本分類效果時(shí)，可以看PC13或PC23的圖分類是否清晰，也可以用R語言rgl包做出PC123三維圖。 QIIME本身結(jié)果中有提供PCA的三維圖結(jié)果，可以通過網(wǎng)頁打開。 13. 13 LDA差異貢獻(xiàn)分析 PCA和LDA的差別在于，PCA，它所作的只是將整組數(shù)據(jù)整體映射到最方便表示這組數(shù)據(jù)的坐標(biāo)軸上，映射時(shí)沒有利用任何數(shù)據(jù)內(nèi)部的分類信息，是無監(jiān)督的，而LDA是由監(jiān)督的，增加了種屬之間的信息關(guān)系后，結(jié)合顯著性差異標(biāo)準(zhǔn)測(cè)試(克魯斯卡爾-沃利斯檢驗(yàn)和兩兩Wilcoxon測(cè)試)和線性判別分析的方法進(jìn)行特征選擇。除了可以檢測(cè)重要特征，他還可以根據(jù)效應(yīng)值進(jìn)行功能特性排序，這些功能特性可以解釋頂部的大部分生物學(xué)差異。詳細(xì)說明可以參考這篇文章http://blog.csdn.net/sunmenggmail/article/details/8071502 。 不同顏色代表不同樣本或組之間的顯著差異物種。使用LefSe軟件分析獲得，其中顯著差異的logarithmic LDA score設(shè)為2。 問題：LDA分析有什么用？ 回答：組間差異顯著物種又可以稱作生物標(biāo)記物（biomarkers），該分析主要是想找到組間在豐度上有顯著差異的物種。 14. 物種進(jìn)化樹的樣本群落分布圖 是將不同樣本的群落構(gòu)成及分布以物種分類樹的形式在一個(gè)環(huán)圖中展示。數(shù)據(jù)經(jīng)過分析后，將物種分類樹和分類豐度信息通過軟件GraPhlAn(http://huttenhower.sph.harvard.edu/GraPhlAn )進(jìn)行繪制。其目的是將物種之間的進(jìn)化關(guān)系以及不同樣本的物種分布豐度和最高分布樣本的信息在一個(gè)視覺集中的環(huán)圖中一次展示，其提供的信息量較其他圖最為豐富。 中間為物種進(jìn)化分類樹，不同顏色的分支代表不同的綱（具體的代表顏色見右上角的圖例），之后外圈的灰色標(biāo)示字母的環(huán)表示的是本次研究中比例最高的15個(gè)科（字母代表的科參見左上角的圖例）。之后的外圈提供的是熱力圖，如果樣本數(shù)<=10個(gè)則繪制樣本，如果樣本數(shù)超過10個(gè)則按照分組繪制，每一環(huán)為一個(gè)樣本，根據(jù)其豐度繪制的熱力圖。最外圈為柱狀圖，繪制的是該屬所占比例最高的樣本的豐度和樣本顏色（樣本顏色見環(huán)最下方的樣本名字的顏色）。其中熱力圖和柱狀圖取值均為原比例值x10000后進(jìn)行l(wèi)og2轉(zhuǎn)換后的值 參考文獻(xiàn)： 1. Vazquez-Baeza Y, Pirrung M, Gonzalez A, Knight R. 2013. Emperor: A tool for visualizing high-throughput microbial community data. Gigascience 2(1):16. 2. Legendre, P. and Legendre, L. 1998. Numerical Ecology. Second English Edition. Developments in Environmental Modelling 20. Elsevier, Amsterdam. 3. Segata N, Izard J, Waldron L, et al. Metagenomic biomarker discovery and explanation[J]. Genome Biol, 2011, 12(6): R60. 4. Langille MGI, Zaneveld J, Caporaso JG, McDonald D, Knights D, Reyes JA et al. (2013). Predictive functional profiling of microbial communities using 16S rRNA marker gene sequences. Nat Biotechnol 31: 814–821. 15. 物種相關(guān)性分析 根據(jù)各個(gè)物種在各個(gè)樣品中的豐度以及變化情況，計(jì)算物種之間的相關(guān)性，包括正相關(guān)和負(fù)相關(guān)。 相關(guān)性分析使用CCREPE算法，首先對(duì)原始16s測(cè)序數(shù)據(jù)的種屬數(shù)量進(jìn)行標(biāo)準(zhǔn)化，然后進(jìn)行Spearman和Pearson秩相關(guān)分析并進(jìn)行統(tǒng)計(jì)檢驗(yàn)，計(jì)算出各個(gè)物種之間的相關(guān)性，之后在所有物種中根據(jù)simscore絕對(duì)值的大小，挑選出相關(guān)性最高的前100組數(shù)據(jù)，基于Cytoscap繪制共表達(dá)分析網(wǎng)絡(luò)圖，網(wǎng)絡(luò)圖采用兩種不同的形式表現(xiàn)出來。 物種相關(guān)性網(wǎng)絡(luò)圖A：圖中每一個(gè)點(diǎn)代表一個(gè)物種，存在相關(guān)性的物種用連線連接，其中，紅色的連線代表負(fù)相關(guān)，綠色的先代表正相關(guān)，連線顏色的深淺代表相關(guān)性的高低。 物種相關(guān)性網(wǎng)絡(luò)圖B：圖中每一個(gè)點(diǎn)代表一個(gè)物種，點(diǎn)的大小表示與其他物種的關(guān)聯(lián)關(guān)系的多少，其中與之有相關(guān)性的物種數(shù)越多，點(diǎn)的半徑和字體越大，連線的粗細(xì)代表兩物種之間相關(guān)性的大小，連線越粗，相關(guān)性越高。 參考文獻(xiàn)： Schwager E, Weingart G, Bielski C, et al. CCREPE: Compositionality Corrected by Permutation and Renormalization[J]. 2014. 16. 聚類分析 根據(jù)OUT數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化處理（1wlog10）之后，選取數(shù)目最多的前60個(gè)物種，基于R heatmap進(jìn)行作圖，熱圖中的每一個(gè)色塊代表一個(gè)樣品的一個(gè)屬的豐度，樣品橫向排列，屬縱向排列，兩個(gè)熱圖，差異是是否對(duì)樣品進(jìn)行聚類，從聚類中可以了解樣品之間的相似性以及屬水平上的群落構(gòu)成相似性。 如果聚類結(jié)果中出現(xiàn)大面積的白或黑是因?yàn)榇罅康木糠浅５?，?dǎo)致都沒有數(shù)值，可以在繪制之前進(jìn)行標(biāo)準(zhǔn)化操作，對(duì)每一類菌單獨(dú)自身進(jìn)行Z標(biāo)準(zhǔn)化。 17. 群落功能差異分析 通過對(duì)已有測(cè)序微生物基因組的基因功能的構(gòu)成進(jìn)行分析后，我們可以通過16s測(cè)序獲得的物種構(gòu)成推測(cè)樣本中的功能基因的構(gòu)成，從而分析不同樣本和分組之間在功能上的差異（PICRUSt Nature Biotechnology, 1-10. 8 2013）。 通過對(duì)宏基因組測(cè)序數(shù)據(jù)功能分析和對(duì)應(yīng)16s預(yù)測(cè)功能分析結(jié)果的比較發(fā)現(xiàn)，此方法的準(zhǔn)確性在84%-95%，對(duì)腸道微生物菌群和土壤菌群的功能分析接近95%，能非常好的反映樣品中的功能基因構(gòu)成。 為了能夠通過16s測(cè)序數(shù)據(jù)來準(zhǔn)確的預(yù)測(cè)出功能構(gòu)成，首先需要對(duì)原始16s測(cè)序數(shù)據(jù)的種屬數(shù)量進(jìn)行標(biāo)準(zhǔn)化，因?yàn)椴煌姆N屬菌包含的16s拷貝數(shù)不相同。然后將16s的種屬構(gòu)成信息通過構(gòu)建好的已測(cè)序基因組的種屬功能基因構(gòu)成表映射獲得預(yù)測(cè)的功能結(jié)果。（根據(jù)屬這個(gè)水平，對(duì)不同樣本間的物種豐度進(jìn)行顯著性差異兩兩檢驗(yàn)，我們這里的檢驗(yàn)方法使用STAMP中的two-sample中T-TEST方法，Pvalue值過濾為0.05，作Extent error bar圖。） 此處提供COG，KO基因預(yù)測(cè)以及KEGG代謝途徑預(yù)測(cè)。用戶也可自行使用我們提供的文件和軟件（STAMP）對(duì)不同層級(jí)以及不同分組之間進(jìn)行統(tǒng)計(jì)分析和制圖，以及選擇不同的統(tǒng)計(jì)方法和顯著性水平。 參考文獻(xiàn)： Donovan H. Parks1 ,Gene W. Tyson,STAMP: statistical analysis of taxonomic and functional profiles, Bioinformatics(2014)30(21):3123-3124.doi:10.1093 18. COG構(gòu)成差異分析圖 圖中不同顏色代表不同的分組，列出了COG構(gòu)成在組間存在顯著差異的功能分類以及在各組的比例，此外右側(cè)還給出了差異的比例和置信區(qū)間以及P-value。 19. KEGG代謝途徑差異分析圖 通過KEGG代謝途徑的預(yù)測(cè)差異分析，我們可以了解到不同分組的樣品之間在微生物群落的功能基因在代謝途徑上的差異，以及變化的高低。為我們了解群落樣本的環(huán)境適應(yīng)變化的代謝過程提供一種簡(jiǎn)便快捷的方法。 圖解讀：圖中不同顏色代表不同的分組，列出了在第三層級(jí)的構(gòu)成在組間存在顯著差異的KEGG代謝途徑第三層分類以及在各組的比例，此外右側(cè)還給出了差異的比例和置信區(qū)間以及P-value。 本例圖所顯示的是第三層級(jí)的KEGG代謝途徑的差異分析，也可以針對(duì)第二或第一層的分級(jí)進(jìn)行分析。 20. 基因的差異分析圖 除了能對(duì)大的基因功能分類和代謝途徑進(jìn)行預(yù)測(cè)外，我們還能提供精細(xì)的功能基因的數(shù)量和構(gòu)成的預(yù)測(cè)，以及進(jìn)行樣本間以及組間的差異分析，并給出具有統(tǒng)計(jì)意義和置信區(qū)間的分析結(jié)果。 這一分析將我們對(duì)于樣本群落的差異進(jìn)一步深入到了每一類基因的層面。 圖解讀：圖中不同顏色代表不同的分組，列出了在組間/樣本間存在顯著差異的每一個(gè)功能基因（酶）以及在各組的比例，此外右側(cè)還給出了差異的比例和置信區(qū)間以及P-value。 21. 在獲得標(biāo)準(zhǔn)報(bào)告后如果希望單獨(dú)修改分組或?qū)δ承┙M之間進(jìn)行顯著性差異分析，可以使用STAMP軟件在自己的電腦上進(jìn)行數(shù)據(jù)分析。STAMP提供了豐富的統(tǒng)計(jì)檢驗(yàn)方法和圖形化結(jié)果的輸出。 在使用STAMP之前需要首先準(zhǔn)備需要的spf格式文件和樣品分組信息表。在我們的報(bào)告中已經(jīng)將KEGG和KO以及COG的結(jié)果文件后經(jīng)過轉(zhuǎn)換生成了適用于STAMP軟件打開的spf格式文件，還有對(duì)應(yīng)的分組信息表文件groupfile.txt。 以下是使用STAMP時(shí)的一些相關(guān)問題，詳細(xì)的STAMP使用教程可以參考我們提供的STAMP使用教程。 1、stamp作圖用的原始數(shù)據(jù)的來源？ STAMP 可以直接使用來自QIIME的biom文件和PICUST的KEGG和ko 文件，groupfile.txt文件的格式為tab-saperated value (tab鍵隔開的數(shù)據(jù)) 2、分組問題：導(dǎo)入數(shù)據(jù)之后，viewgroup legend ,在窗口右側(cè)會(huì)出現(xiàn)分組欄，根據(jù)需要進(jìn)行分組。 3、Unclassiffied選項(xiàng)中，remain Unclassiffied reads、remove Unclassiffied reads、和use only for calculating frequency profiles 方法的區(qū)別？ remain Unclassiffied reads和use only for calculating frequency profiles方法會(huì)保留所有的數(shù)據(jù)，而remove Unclassiffied reads僅僅保留有確定分組信息的數(shù)據(jù)。 4、Statistical test 中，Welch’s t-test、t-test、white’s non-parametric t-test的區(qū)別，各自優(yōu)缺點(diǎn)？ 為了確保統(tǒng)計(jì)學(xué)意義和準(zhǔn)確度和精確性，需要足夠多的樣本數(shù)目，t-test檢驗(yàn)可以在最少樣本數(shù)為4的時(shí)候確保高的準(zhǔn)確度和精確性。 當(dāng)兩個(gè)樣本之間具有相同方差的時(shí)候，用t-test更為準(zhǔn)確，當(dāng)兩個(gè)樣本沒有相同方差，Welch’s t-test更為準(zhǔn)確。 當(dāng)樣本數(shù)目少于8的時(shí)候，可以使用white’s non-parametric t-test，該計(jì)算時(shí)間較長(zhǎng)，當(dāng)樣本數(shù)目過多的時(shí)候不宜使用該方法。 5、Two-group 中type: one side 和 two side 的區(qū)別？ One side 只會(huì)顯示前一個(gè)group與后一個(gè)group差異的比例，而two side 兩者之間的比例均會(huì)顯示。 6、 STAMP在使用時(shí)首先打開了一個(gè)分析文件，如果新打開一個(gè)可能會(huì)導(dǎo)致顯示錯(cuò)誤？ 目前版本的STAMP存在一些小問題，一次分析只能使用一個(gè)數(shù)據(jù)文件，如果要打開新的需要關(guān)閉軟件后再打開。 22. 環(huán)境因子分析 RDA分析 CCA/RDA分析 基于對(duì)應(yīng)分析發(fā)展的一種排序方法，將對(duì)應(yīng)分析與多元回歸分析相結(jié)合，每一步計(jì)算均與環(huán)境因子進(jìn)行回歸，又稱多元直接梯度分析。主要用來反映菌群與環(huán)境因子之間的關(guān)系。 RDA 是基于線性模型，CCA是基于單峰模型。分析可以檢測(cè)環(huán)境因子、樣品、菌群三者之間的關(guān)系或者兩兩之間的關(guān)系。 橫軸和縱軸：RDA 和CCA 分析，模型不同，橫縱坐標(biāo)上的刻度為每個(gè)樣品或者物種在與環(huán)境因子進(jìn)行回歸分析計(jì)算時(shí)產(chǎn)生的值，可以繪制于二維圖形中。 圖解讀：冗余分析可以基于所有樣品的OTU作圖，也可以基于樣品中優(yōu)勢(shì)物種作圖； 箭頭射線：箭頭分別代表不同的環(huán)境因子； 夾角：環(huán)境因子之間的夾角為銳角時(shí)表示兩個(gè)環(huán)境因子之間呈正相關(guān)關(guān)系，鈍角時(shí)呈負(fù)相關(guān)關(guān)系。環(huán)境因子的射線越長(zhǎng)，說明該影響因子的影響程度越大； 不同顏色的點(diǎn)表示不同組別的樣品或者同一組別不同時(shí)期的樣品，圖中的拉丁文代表物種名稱，可以將關(guān)注的優(yōu)勢(shì)物種也納入圖中； 環(huán)境因子數(shù)量要少于樣本數(shù)量，同時(shí)在分析時(shí)，需要提供環(huán)境因子的數(shù)據(jù)，比如 pH值，測(cè)定的溫度值等。 23. 有其他問題可以聯(lián)系谷禾信息。 END 注意事項(xiàng) 假設(shè)您要對(duì)人體腸道微生物菌群進(jìn)行測(cè)序，一般可以選擇對(duì)糞便進(jìn)行取樣。一般建議采取新鮮樣品，因?yàn)槲⑸锉旧硎腔畹娜后w，樣品長(zhǎng)期保存或不在原環(huán)境下保存會(huì)改變?cè)芯旱臉?gòu)成，最終導(dǎo)致我們得到的菌群構(gòu)成發(fā)生偏差。比如糞便樣品如果4度低溫保存了一段時(shí)間，則其中部分耐低溫的菌可能仍然在持續(xù)繁殖，這樣最終樣品的菌群會(huì)發(fā)生偏差。如果是送往公司測(cè)序也建議首先對(duì)樣品進(jìn)行DNA提取后再寄送，因?yàn)樵紭悠返募乃瓦^程也可能導(dǎo)致菌群變化。