人類微生物組研究設計、樣本採集和生物信息分析指南

人類微生物組研究設計、樣本採集和生物信息分析指南

A guide to human microbiome research: study design, sample collection, and bioinformatics analysis

Chinese Medical Journal [IF: 1.585]

DOI: https://doi.org/10.1097/CM9.0000000000000871

Review: 2020-6-26

錢旭波1, 陳同2, 徐益萍1, 陳雷3, 孫馥香4, 盧美萍1, 劉永鑫5,6

浙江大學醫學院附屬兒童醫院風溼、免疫和變態反應科

中國中醫科學院中藥資源中心

首都醫科大學附屬復興醫院

易漢博基因科技（北京）有限公司

中國科學院遺傳與發育生物學研究所

中國科學院大學，生物互作卓越創新中心

錢旭波和陳同為共同第一作者

通訊作者：盧美萍，浙江大學醫學院附屬兒童醫院風溼、免疫和變態反應科，中國浙江杭州竹竿巷57號，郵編：310003，郵箱：[e-mail]

摘要

這篇綜述的目的是為醫學研究人員，特別是那些沒有生物信息學背景的研究者提供簡單易懂的微生物組學知識，包括研究中常用的概念、技術和分析方法等。首先，我們介紹了基本概念，例如微生物群（microbiota）、微生物組（microbiome）和宏基因組（metagenome）等。然後，我們討論了研究設計方案、樣本量計算方法以及提高研究可靠性的方法。我們特別強調了陽性和陰性對照的重要性。接下來，我們討論了微生物組研究中常用的統計分析方法，重點關注多重比較的問題以及組間β多樣性分析的方法。最後，我們介紹了生物信息學分析的具體流程。總之，嚴謹的研究設計是獲得有意義結果的關鍵步驟，而適當的統計方法對於準確解釋微生物組數據很重要。通過閱讀這篇文章，研究者能獲得研究設計、樣本採集和生物信息分析等全方位的微生物組學知識。

關鍵詞：微生物組、研究設計、統計分析、樣本量、生物信息分析、分析流程

1. 前言

隨著測序技術和數據分析方法的發展，近幾年醫學微生物組研究領域出現了一些令人矚目的成果[1-3]，比如微生物組與代謝性疾病[4-6]、消化系統疾病[7-10]和心血管系統疾病[11]之間的關係日益明確。這些發展和發現增加了醫生在微生物組研究方面的興趣，進而也湧現出了大量有價值的論文[12]。另外，隨著QIIME 2[13]和多組學方法[1, 9]等先進技術和分析流程的出現，微生物組分析方法也不斷進步。然而，理解和掌握這些技術和分析流程並非易事，特別對於醫生來說更是如此。

本文的目的是為研究者，特別是那些沒有生物信息學背景的醫生提供易懂的微生物組學知識，這些知識包括詳細的微生物組學基本概念、科研設計方法、樣本採集和保存方法、統計分析方法以及生物信息分析方法。我們希望醫生們通過閱讀此文能夠快速掌握以上知識和方法，進而有效地挖掘數據背後的生物學意義。

2. 基本概念2.1 Microbiota、Microbiome等術語

Microbiota(微生物群/微生物組)是指定植在人體特定部位的微生物，包括細菌、古菌、病毒、真菌和原生動物[14, 15]。在醫學研究中，如果測序技術採用的是16S rRNA基因（又稱為rDNA），則microbiota是指細菌和古菌。Microbiome是指整個微生境，包括微生物、基因組和周圍環境[14, 15]。不過，microbiota和microbiome有時存在混用情況。我們建議，如果你的研究僅涉及微生物本身，則應該使用microbiota，否則應該使用microbiome（圖1）。例如，如果研究者想探索腸道短鏈脂肪酸與微生物的關係，使用microbiome更合適。宏基因組（metagenome）是指微生物基因組的集合[14]，一般用鳥槍法宏基因組測序獲得，宏基因組學則是研究宏基因組的學科[12, 14]。病毒組（virome）指人體內或表面的病毒集合，包括內源性逆轉錄病毒、真核生物病毒和噬菌體[16]。研究病毒組的學科就是病毒組學。作者註：Microbiota國內有些學者翻譯為「微生物群」，microbiome翻譯為「微生物組」。不過中文文獻用「微生物組」或「××菌群」即可，多數情況下不需要區分是microbiota或microbiome。

圖1：微生物組、微生物群、宏基因組和16S rDNA的概念。

（A）微生物組（microbiome）的概念不僅涵蓋微生物，而且涵蓋周圍的環境條件。微生物群（microbiota）僅指微生物本身。（B）宏基因組是指微生物的所有基因組，而16S rDNA僅涵蓋基因組的一部分。（C）α多樣性衡量樣本中的多樣性，而β多樣性比較樣本之間的物種差異。

2.2 細菌層級分類

細菌分類最常用的層級為門、綱、目、科、屬、種、株。例如，臨床上十分常見的大腸埃希菌的層級分類見表1。

表 1: 大腸埃希菌細菌層級分類分類層級分類名稱門變形菌門綱丙型變形菌綱目腸桿菌目科腸桿菌科屬埃希氏桿菌屬種埃希氏菌株EIEC112ac株2.3 操作分類單元和擴增子序列變異

操作分類單元（operational taxonomic units，OTUs）的構建對於標記基因（擴增子）數據分析非常重要[17]。OTU是指一組高度相似的序列，通常將具有97%相似性的一組序列歸為一個OTU[18, 19]。不過，這種OTU的方法有顯著的缺點，它人為地設置一個相似性閾值，漏掉了細微的和真正的生物學序列差異[20]。最近開發的擴增子序列變異（amplicon sequence variants, ASVs）方法可以解決這些問題，它使用序列變異信息將序列數據解析為準確的序列特徵。ASV具有單核苷酸解析度，並且具有比OTU相似或更好的敏感性和特異性[20]。注意，OTU或ASV不等於物種，一個OTU / ASV可能包括多個物種，反之亦然[21]。

2.4 α-多樣性

α-多樣性是指樣本內的多樣性，常見的樣本有糞便，唾液或支氣管肺泡灌洗液等[15]。醫學研究中經常使用3種α多樣性指數：Chao 1指數、香農指數和辛普森指數。Chao 1指數主要反映物種數量（richness），它計算時考慮以下三個因素：物種數量、單條序列數量和雙條序列數量[22]。這意味著它不能反映微生物組的豐度（abundance）。香農指數結合了豐度和均勻度信息[23]，它賦予稀有物種更多的權重[22]，這意味著當稀有物種的數量增加時，它的值會更大。香農指數的值通常不超過5.0；它的值越高，α多樣性就越豐富[22]。辛普森指數也整合了豐度和均勻度，不過與香農指數比較，計算時它對常見物種有更大權重。它的值介於0-1之間，這個值越大，α多樣性越豐富[22]。在以上指數中，richness是指一個樣本中物種的數量[17, 24]，而abundance（豐度）指物種的原始序列讀數[24]。如果原始序列讀數被轉換成百分比後，它就稱為相對豐度。

2.5 β-多樣性

β-多樣性是指樣本或組間的微生物組差異，通常用於了解兩組微生物組組成的差異是否顯著。在這裡，我們關注兩個常用的β多樣性指數：Bray-Curtis相異性和UniFrac距離。Bray-Curtis相異性是一種用於量化兩個樣本或組間的物種組成差異的指標，其值的範圍是0到1，其中0表示兩個樣本或組間具有相同物種，而1則表示它們不共享任何物種[25]。此外，它在計算時給予常見物種更大的權重[23]。請注意，Bray-Curtis相異性不是真正的距離度量指標，因此用「Bray-Curtis相異性」的叫法比「Bray-Curtis距離」更恰當[22]。

UniFrac距離可以不加權，也可以加權，它基於系統發育距離估算微生物組樣本或組間的差異[26]。未加權的UniFrac距離只考慮了物種是否存在，它對於檢測稀有物種的數量變化很敏感，但是在計算中忽略了豐度信息[27]。加權UniFrac距離計算時納入了豐度信息[28]，並減少了稀有物種的權重[29]。

2.6 排序

排序用於探索數據結構，由降維後的正交軸圖形表示。排序圖是可視化β多樣性的有效方法。排序可以分為2大類：非約束排序和約束排序[30-32]。如果圖形上的點不受環境因素（樣本元數據）的約束，這種排序叫做非約束排序，否則叫約束排序[32]。常用的非約束標準包括主成分分析（principal component analysis, PCA）、對應分析（correspondence analysis, CA）、主坐標分析（principal coordinate analysis, PCoA）和非度量多維標度（non-metric multidimensional scaling, NMDS）[30, 32]。常用的約束排序有冗餘分析（redundancy analysis, RDA）和典範對應分析（canonical correspondence analysis, CCA）[31, 32]。

微生物組信息是高維數據。PCA通過將數據以幾何方式投影到較少的維度上來簡化複雜性，它在計算中使用歐幾裡得（Euclidean）距離[30]。通常情況下它並不適用於物種豐度數據的分析，因為PCA分析的數據必須是線性的[30]。但是如果物種數據經過Hellinger轉換，則PCA可以用於物種數據分析[30]。相反，CA適合於物種豐度數據分析，而且無需預先轉換數據。在CA分析中，所有樣本均使用Pearson卡方距離進行排序[30]。但是請注意，稀有物種可能會對CA分析產生過大影響[33]。如果研究人員希望基於相異性指標來對樣本或特徵進行排序，那麼PCoA是一個不錯的選擇。在微生物組研究中，PCoA分析最常使用Bray-Curtis相異性和UniFrac距離。NMDS用於表示排序圖中樣本的相對位置。與PCoA相似，NMDS分析可以使用任何距離或相異矩陣。參考文獻[30]詳細介紹了PCoA和NMDS之間的差異，在大多數情況下PCoA比較常用。

RDA是一種結合了PCA和回歸的約束排序，它的響應矩陣是微生物組數據，解釋矩陣是臨床指標（樣本元數據）。RDA對於顯示微生物組數據是否受臨床指標影響很有用。但是請注意，由於PCA計算過程要求響應矩陣的數據結構必須是線性的，因此可能需要對數據進行預轉換。最後，CCA其實就是CA的約束版本，它具有CA的基本特性和缺點[31]。

3. 研究設計3.1 研究設計方案

嚴謹的研究設計對於獲得準確而有意義的結果很重要。醫學微生物組研究中最常使用的研究方法包括橫斷面研究、病例對照研究、縱向研究和隨機對照試驗（randomized controlled trial, RCT）。前3種是不應用幹預因素的觀察性研究，而最後一個是典型的實驗性研究。

橫斷面研究分為描述性橫斷面研究和分析性橫斷面研究[34]。前者僅是描述性的，主要用於調查一個或多個人群中的微生物組成，而後者則用於探討微生物組與健康結果之間的關聯。但是，微生物組與健康結果之間的關聯可能源於混雜因素，例如性別[35]、年齡[36]、體重指數（body mass index, BMI）[37]、飲食[5, 38]、季節[39]和藥物治療[40, 41]。此外，橫斷面研究時，微生物組和結果是同時測量的，因此很難確定它們之間的因果關係。通常，橫斷面研究僅用於探索微生物組的基本特徵，並且可以作為後續研究的初步實驗。

在大多數情況下，微生物組被視為暴露（exposure），疾病被視為結局（outcome）。在這些假設下，傳統的病例對照研究很少用於微生物組研究，因為以前的暴露（微生物組）信息很難獲得。但是，如果暴露和結局對調，則可以使用病例對照研究設計方案。

同樣，在上述假設下進行前瞻性隊列研究也很困難，因為很難知道哪些微生物是潛在的暴露。而且，定義可用作暴露或非暴露因素的特定微生物組並非易事，因此難以將研究對象確定為暴露或非暴露個體。在實踐中，有或沒有疾病的個體通常歸入研究組或對照組，然後在不同時間點前瞻性地收集含有微生物組的樣本[17]。也就是說，前瞻性隊列研究中的研究對象通常根據臨床結局而不是特定的微生物組模式進行分組。

RCT或其他實驗研究的目的是評估幹預措施的有效性。幹預措施可以是藥物或微生物組。例如，糞菌移植研究中的幹預措施是微生物群[42, 43]。

值得注意的是，對照組的選擇應恰當。以上這些研究設計中應注意匹配混雜因素，這部分內容將在下面討論。有時對照的選擇很困難，尤其是在臨床研究中幹預措施是微生物群本身的情況下。在這種情況下，如果其他研究設計不合適，那麼進行有對照組的前後自身對照試驗（controlled before-after trial）或歷史對照試驗將是一個不錯的選擇[44]。

3.2 定義納入和排除標準

定義確切的納入標準和排除標準可以使組間更好地匹配，並且有利於控制混淆因素，比如年齡[36, 45]、性別[35]、BMI[46]、飲食[47]、季節因素[39]、藥物治療[40, 41]、種族[48]、地理區域[45]和共存疾病等[7]。年齡可顯著影響微生物組，對於那些小於16歲的人更是如此[36, 45]。因此，對於涉及兒童的研究，年齡必須很好地匹配。飲食是另一個對微生物組改變有影響的因素，所以也要進行匹配[47]。為了增加組間的可比性，地理區域因素在研究設計時也需要考慮在內[45]。由於藥物治療對於微生物組有顯著影響，所以入組前數月內接受過藥物治療的患者應該排除在外[41, 49]，這裡講的數月通常指入組前3~6個月[49]。

3.3 微生物組研究的樣本量和檢驗效能計算

在進行研究設計時估計樣本量大小非常重要。適當的樣本量可使微生物組研究識別出組間的差異，並節省資源和時間。但是，樣本量和檢驗效能計算對於研究者來說仍然是一個挑戰[50]。微生物組研究中最常用的樣本量和效能計算方法可以用t檢驗、方差分析、χ2檢驗和Dirichlet多項式模型[51]。以t檢驗為例，分3個步驟確定樣本大小和效能計算。首先，通過初步實驗獲得少量擴增子數據。其次，使用R包vegan計算出每個樣品的香農指數[52]。最後一步是使用R軟體包pwr中的power.t.test()函數計算樣本量和效能。當研究者僅關注兩組之間物種多樣性的差異時，可使用t檢驗計算樣本量和效能。在參考文獻[51]中有樣本量和效能計算的詳細介紹。

3.4 陰性和陽性對照的重要性

微生物組研究的結果可能會受到多種因素的影響，例如DNA提取試劑盒、採樣方法、汙染和測序方法等[53]，不過可以通過使用陰性和陽性對照來減少這些影響。不幸的是，以前的研究中只有30％報告使用了陰性對照，只有10％報告使用了陽性對照[53]。使用對照對於準確認識微生物組非常重要，尤其是當樣本的微生物含量較低時。以前的研究發現，過去被認為是無菌的標本（例如胎盤和關節液）可能會被微生物定植[54]。但是，這些陽性的結果可能是由其他因素導致的，例如汙染。有趣的是，這些低生物含量標本在採用陰性和／或陽性對照後已被證明是無菌的[55]。因此，我們建議當樣本為低生物含量樣本（例如血液、羊水、腦脊液、關節液和胎盤等）時，應考慮使用陰性和陽性對照。值得注意的是，陰性對照和陽性對照在病毒學研究中也很重要，因為病毒和細菌通常是同時進行檢測的[16]。此外，R包decontam可用於鑑定和去除擴增子和宏基因組學數據中的汙染物序列[56]。

3.5 測序方法的選擇

微生物組研究中使用的測序方法包括擴增子測序、宏基因組測序和宏轉錄組測序。擴增子測序包括適用於細菌和古菌的16S rDNA測序以及適用於真菌的內部轉錄間隔區（internal transcribed spacer, ITS）測序。每種測序方法的優缺點在這兩篇參考文獻中有詳細討論[17, 57]。簡而言之，擴增子測序很便宜，可應用於受宿主DNA汙染的低生物含量標本，但一般僅能注釋到「屬」層級，並且易受某些固有偏倚來源的影響，例如PCR循環數[58]。宏基因組測序方法對樣品中存在的所有DNA進行測序，包括細菌、病毒、真核生物和宿主的DNA。它不僅將其分類學解析度擴展到「種」或「株」的水平，而且還提供了潛在功能信息[17]。但是，擴增子和宏基因組測序方法都無法區分死微生物或活微生物[17]。轉錄組測序僅產生群落的活躍功能信息。鑑於這些測序方法的優缺點不同，建議將多種測序方法整合在一起以優化研究設計。簡而言之，測序方法的選擇主要取決於實驗成本和樣本質量。擴增子測序通常用於獲得微生物群落的概況[59]，並且通常適用於大規模研究[6, 60]。如果您有足夠的項目資金，並且想要獲得菌株水平的解析度和潛在功能，甚至想要恢復整個基因組，宏基因組測序是一種首選方法[61-65]。

3.6 提高研究可靠性的方法

簡單的橫斷面研究在微生物組研究中的意義有限。在本小節中我們討論了提高研究可靠性的方法。首先，首選縱向研究或RCT研究，而不是橫斷面研究或病例對照研究[17, 66]。其次，應計算樣本量[51]。第三，混淆因素應匹配，元數據（即各種臨床指標等信息）應仔細收集。第四，應詳細定義納入和排除標準。例如，幼年特發性關節炎有幾種亞型，每種亞型可能代表不同的疾病[67]。研究者應確定患者組中是否包括所有亞型。第五，最好考慮使用陰性和/或陽性對照[68]。第六，整合其他組學方法，例如代謝組學、轉錄組學和蛋白質組學，這對於全面了解微生物群落的結構和功能至關重要[17]。因此，應考慮獲取微生物群落代謝物概況和／或其他多組學數據。目前，僅探索微生物群落結構的研究不被視為論證效率強的研究設計[17]。最後，建議在動物模型中驗證從臨床試驗獲得的初步結果。

表2列出了設計臨床微生物組研究需要考慮的因素，圖2展示了典型的工作流程。實驗研究需要考慮的因素見參考文獻[49]。

表 2: 臨床微生物組研究設計需要考慮的要素核對表需要考慮的要素核對詳情研究設計類型□橫斷面研究  □病例對照研究  □隊列研究  □RCT  □其他：性別□已匹配  □未匹配  □其他：年齡□已匹配  □未匹配  □其他：BMI□已匹配  □未匹配  □其他：種族□已匹配  □未匹配  □其他：地理區域□已匹配  □未匹配  □其他：飲食□組間已經均衡並已記錄：列出詳細信息；□未記錄季節因素□樣本收集自相同季節 □樣本收集自不同季節藥物治療入組前使用了哪些藥物？使用了多久？納入標準□已定義好  □定義不清晰排除標準□已定義好  □定義不清晰樣本量□已計算  □未計算測序方法□擴增子  □宏基因組  □其他陰性和／或陽性對照□有陰性對照  □無陰性對照 □有陽性對照  □無陽性對照多組學方法□代謝組  □轉錄組  □蛋白組樣本類型□糞便  □結腸灌洗液  □腔內刷  □組織鉗出物  □黏膜下組織  □關節液  □尿液  □牙菌斑  □唾液  □皮膚  □其他：動物模型驗證□結果將在動物模型中驗證 □結果不將在動物模型中驗證

RCT：隨機對照試驗

4. 樣本類型、保存和儲藏4.1 樣本類型

人類微生物組研究的樣本類型包括糞便、結腸灌洗液和腔內刷等（表2）。樣本類型的選擇取決於感興趣的研究假設。例如，糞便樣本易於收集，可用於大規模和縱向研究。另一方面，活檢樣本對於探索微生物群與宿主之間的相互作用更有用[69]。注意，在一項研究中應該固定採樣位置，因為人體的不同部位定植著不同的微生物群[70, 71]。

圖2：人類微生物組研究的典型流程。

4.2 保存和儲藏

樣品保存和儲藏的方法應適合實驗方法和樣品類型。最通用的方法是直接冷凍樣品，它可用於各種測序和實驗方法，例如擴增子、宏基因組、轉錄組測序和代謝組學測定。建議將樣品收集後15分鐘內保存在-20℃下[72, 73]，然後在收集24小時內用乾冰轉移到-80℃冰箱中儲藏。不過樣本通常是在家裡而不是在醫院收集的，在這種情況下可以使用保存液。保存液中保存的樣本可以在環境溫度下保存一周以上[74]。請注意，樣品的保存和儲藏方法應一致，以最大程度地減少潛在的混淆因素幹擾。

5. 微生物組研究中的統計分析方法

醫學研究者通常熟悉單變量統計方法，例如t檢驗、方差分析、χ2檢驗和秩和檢驗。因此，我們在這裡僅討論與多重比較和其他多元統計方法有關的問題。我們首先討論多重比較會遇到的問題及其解決方案，包括P值調整和使用錯誤發現率（FDR）。然後，我們討論其他多元統計方法，例如置換多元方差分析（permutational multivariate analysis of variance, PERMANOVA）和Mantel檢驗。

5.1 多重比較的問題及解決方法

由於微生物組數據是高維的，因此多重比較經常在微生物組研究中使用。例如，特徵表（feature table）具有成百上千個OTU或ASV，並且每個OTU或ASV都可以進行多次比較。醫學研究者經常遇到的另一個例子可能更容易理解。假設一項研究分為3組，例如A組、B組和C組，而研究者想比較這3組之間的差異。在這種情況下就應調整P值，因為每個組都進行了2次比較，即A組與B組，A組與C組，B組與C組。如果有任何組或變量需要進行多次比較則必須進行P值調整，以便減少假陽性率[75]。

調整P值的經典方法是控制family-wise錯誤率，即Ⅰ類錯誤或α水平。Bonferroni是校正α水平最常用的方法。校正P值的計算非常容易：單個檢驗的α值除以檢驗次數。因此，對於上述具有3個檢驗次數的例子，調整後的P值為0.05 / 3 = 0.017，即只有P <0.017的檢驗結果才被認為是有意義的[75]。請注意，Bonferroni校正僅適用於多重比較次數較少的假設檢驗，否則會導致較高的假陰性率（圖3）[75]。

解決多重比較問題的另一種方法是控制錯誤發現率（false discovery rate, FDR），它是I類錯誤或假陽性的數量與所有被拒絕的無效假設的預期比例。例如，如果100個陽性假設檢驗結果中有5個是錯誤發現，則FDR為5％。在微生物組研究中，通常使用「Benjamini-Hochberg（BH）校正的P值」而不是原始P值。校正後的P = 原始P * m/i，其中m是檢驗次數，i是每個P值從小到大排序的序號[75]。如果校正後的P值小於你選擇的所選FDR，則認為該檢驗是有統計學意義的。與Bonferroni方法相比，BH方法不那麼保守（即校正強度不是很大），BH法通常用於微生物組特徵的多重比較。Bonferroni和BH是最常用的P值校正方法[76]，這兩種P值校正方法的校正強度見圖3所示。

圖3：不同P值校正方法的校正強度

該圖顯示，Benjamini-Hochberg校正強度小於Bonferroni。隨著原始P值的增加，Bonferroni校正法生成的校正後P值快速接近1.0。

5.2 PERMANOVA檢驗

有幾種統計方法或模型可以用於組間β多樣性比較，比如PERMANOVA、Mantel檢驗、相似性分析（ANOSIM）和多響應置換程序（multi-response permutation procedures, MRPP）。PERMANOVA最常用，並且被認為是以上檢驗方法中檢驗效能最大的一種[77]，它可通過R包vegan中的函數adonis()實現[52]。vegan包可計算4種常用相異性或距離度量：Bray-Curtis相異性、Jaccard距離以及加權和未加權UniFrac距離[29]。如果PERMANOVA檢驗的P值小於0.05，則表明不同組間的β多樣性差異具有統計學意義；該檢驗的另一個輸出結果是R2，它表示總方差可以用分組因素來解釋的比例[29]。

5.3 Mantel檢驗

Mantel檢驗通常用於分析元數據矩陣和微生物組矩陣之間的關聯[77]，它可使用R包vegan中的mantel()函數實現[52, 77]。該檢驗的輸出至少2個主要統計量：P值和r。與其他類型的相關係數類似，r的值範圍是-1 ~ +1[29]。例如，假設研究人員想知道元數據種的分組因素（例如吸菸狀態）是否對腸道微生物組的組成產生影響。如果P＜0.05並且r＞0，這表明吸菸組和不吸菸組之間腸道微生物組的組成不同，元數據矩陣和微生物組矩陣呈正相關。

6. 生物信息分析6.1 擴增子數據分析：從原始數據到物種分類表

有幾種流行的軟體或分析流程（pipeline）可用於擴增子數據分析，例如QIIME 2[13]、USEARCH[78]、VSEARCH[79]和mothur[80]。前兩者具有許多優點，並已被許多研究者使用和推薦。每種軟體或分析流程的優缺點已在我們先前的論文中詳細描述[81]遺傳：微生物組數據分析方法與應用 和 Protein Cell：擴增子和宏基因組數據分析實用指南。擴增子分析的主要步驟見圖4A。我們通常從fastq格式的原始雙端Illumina數據開始，最終輸出是一個特徵表，也稱為OTU表或ASV表。

第一步是從原始數據中恢復純淨的擴增子序列，因為原始數據包括人造序列產物（artifact），例如引物和標籤（barcode）。它包括3個主要過程：合併雙端序列，通過標籤拆分序列和去除引物。由於原始數據沒有統一的標準格式，因此我們需要設計適合上述過程的分析流程。另外，我們也可以使用基因測序公司提供的純淨擴增子數據。圖4B顯示了用於恢復純淨擴增子序列的典型分析流程。

第二步是濾除低質量序列，以便減少背景「噪音」。

第三步是識別非冗餘序列並且計數。高質量序列仍然有許多人造序列產物，例如錯誤序列和嵌合體。非冗餘序列的計數是找出可靠序列的關鍵信息。

第四步是選擇代表性序列（特徵）。此步驟基於唯一序列，並通過將序列聚類成OTU或降噪生成ASV來實現[18, 82]。此步驟還包括de novo檢測和去除嵌合體。

第五步是有參嵌合體檢測，這是可選項[83]。通過將序列比對到資料庫中，例如rRNA資料庫SILVA[84]，可以進一步過濾特徵序列。應當指出的是，該步驟可以降低假陽性率但易於導致假陰性結果。

最後，通過將純淨的擴增子數據與特徵序列進行比較來生成特徵表（圖4A）。然後使用基於RDP[85]、SILVA或Greengenes[86]資料庫的分類器實現特徵序列的物種分類。此外，基於16S rRNA基因譜，使用PICRUSt[87, 88]、FAPROTAX[87, 89]和BugBase[90]等工具可實現功能預測。

6.2 宏基因組分析：從原始數據到物種和功能分類表

擴增子測序僅能獲得微生物組部分的分類學信息，而且PCR過程很容易產生偏倚和嵌合體[83]。鳥槍宏基因組測序比擴增子測序提供更詳細的基因組信息和更高的分類學解析度[66]。與擴增子方法相比，宏基因組學分析更為複雜，但是它提供了更準確的物種分類、多維度的功能信息，甚至是末培養微生物的基因組草圖。宏基因組分析流程如圖4C所示。

第一步是預處理原始序列數據。原始數據包含低質量的汙染序列以及與宿主相關序列。我們可以使用FastQC軟體（http://www.bioinformatics.babraham.ac.uk/ projects/fastqc/）進行數據質量檢查，然後使用KneadData流程進行質量控制[91]並去除宿主DNA[92]。有關更多KneadData的信息，請訪問 http://huttenhower.sph.harvard. edu/kneaddata 。

第二步是使用基於序列的方法分析物種分類和功能代謝特徵。人類微生物組具有高質量的基因集（gene catalog）和基因組[64, 65]，因此我們建議使用HUMAnN2[93]工具並採用基於序列的方法進行物種分類和代謝通路分析，該方法高效且易於操作。但是，這種方法只使用一小部分序列信息，而且分析結果受到已知資料庫的限制[66]。

如果需要發現新物種或基因功能，則需要進行第三步。有幾個好的軟體工具可以用於將純淨序列組裝為重疊群（contigs），例如MEGAHIT[94]和metaSPAdes[95]。然後通過MetaProdigal[96]或Prokka[97]從長序列中預測基因。另外，其他軟體工具也可以用於從短序列中預測編碼基因，例如MetaGeneAnnotator[98]、MetaGeneMark[99]、Glimmer-MG[100]、MetaGUN[101]、FragGeneScan[102]和Orphelia[103]。為了減少重複基因，在分析多個樣品或批次時需要使用CD-HIT構建非冗餘基因集[104]。通過採用Bowtie 2[92]或Salmon[105]工具進行比對的方法可以計算基因豐度。目前至少有20個軟體工具可用於宏基因組數據物種分類[106]。我們建議使用超快速分類器Kraken 2，它可以提供快速、準確和「種」級別的分類結果[107]。至於功能注釋，許多研究人員都推薦使用DIAMOND[108]，它是一種快速、敏感的蛋白質比對工具[108]。每個資料庫都提供了獨特的功能視角，例如，京都基因與基因組百科全書（Kyoto Encyclopedia of Genes and Genomes, KEGG）[109]、EggNOG（一個提供直系同源關係、功能注釋和基因進化歷史的資料庫）[110]、碳水化合物活性酶資料庫（Carbohydrate-Active enZYmes Database, CAZy）[111]、致病菌的毒力因子（Virulence Factors of Pathogenic Bacteria, VFDF）[112]和綜合抗生素抗性資料庫（Comprehensive Antibiotic Resistance Database, CARD）[113]。宏基因組通常包含100~1000個物種[64]，很難釐清彼此關係。分箱算法可以恢復無法培養的高豐度菌的基因組草圖，並重建系統發育和代謝通路。

最後一步是使用metaWRAP[114]或DASTool[115]執行分箱流程（圖4C）。這些軟體工具有逐步操作教程，並且在其網站上提供了有關人類微生物組的一些樣本數據集[81]。另外，幾個集成的分析流程，例如MOCAT 2[116]、bioBakery[98]、IMP[117]和微生物組助手（Microbiome Helper）[118]，可以執行上述部分或全部分析步驟。你可以在微信公眾號「宏基因組」中找到一些受歡迎軟體的中文教程。

現在你已經獲得了物種分類和功能信息文件。通過STAMP或LEfSe可以輕鬆找到你感興趣的生物標記[119, 120]。使用R語言或ImageGP（http://www.ehbio.com/ ImageGP）可以將所有結果可視化。

圖4：人類微生物組研究的生物信息學分析流程

（A）擴增子數據分析的主要步驟。（B）擴增子數據預處理的典型流程圖：從原始的雙端序列到純淨的擴增子。（C）宏基因組測序數據的分析流程。（a）預處理。它涉及刪除低質量序列、接頭和宿主序列。輸出文件是純淨序列。（b）基於序列的分析。它將序列與資料庫比對來推斷物種分類和代謝特徵。（c）基於組裝的分析。它將短序列組裝為長序列，預測基因，構建非冗餘基因集，並與資料庫比對進行物種分類和功能注釋。（d）分箱。它涉及恢復未培養微生物的基因組草圖，並重建系統發育和代謝通路。KEGG：京都基因與基因組百科全書（Kyoto encyclopedia of genes and genomes）；eggNOG：基因進化譜系：非監督直系同源群（Evolutionary genealogy of genes: non-supervised orthologous groups）；CAZy：碳水化合物活性酶資料庫（Carbohydrate-active enzymes database）；CARD：抗性基因綜合資料庫（Comprehensive antibiotic resistance database）；VFDB：毒力因子資料庫（Virulence factor database）。

7. 病毒組在人類疾病中的作用

近年來病毒組在人類疾病中的作用吸引了醫學研究者的關注[121]。使用病毒組學的方法已發現了許多令人信服的研究成果[122]，其中一些技術已經用於臨床[123]。在微生物組研究中，病毒組學與其他多組學方法整合後顯示出廣闊的應用前景。但是，病毒組學研究仍然面臨一些挑戰。例如，至少40％的病毒序列無法注釋[124]。此外，病毒的測序結果容易受到背景噪音的影響[17]。最後，很難獲得用於病毒組研究的商業化陽性對照，即病毒模擬群落[16]。

8. 總結和結論

本文討論了用於微生物組研究的研究設計、樣本收集、統計方法和生物信息學分析方法。在「研究設計」部分，我們強調了研究設計的重要性，特別是設計方案、樣本量計算以及用於提高研究可靠性的多種措施。研究設計非常重要，因為不好的研究設計可能會產生無意義的數據。在「統計分析」部分，我們介紹了詳細的多重比較P值校正方法。選擇合適的統計方法對於準確解釋微生物組數據很重要。最後，「生物信息學分析」部分介紹了用於分析微生物組數據分析的方法。本文圖中使用的腳本可從 https://github.com/YongxinLiu/Qian2020CMJ 獲得。

綜上所述，對於微生物組研究而言，嚴謹的研究設計在獲得有意義的結果方面具有舉足輕重的作用，而適當的統計方法對於準確解釋微生物組數據非常重要。循序漸進的分析流程為研究者掌握最新生物信息學分析方法提供了幫助。

參考文獻

略，詳見原文

Xu-Bo Qian, Tong Chen, Yi-Ping Xu, Lei Chen, Fu-Xiang Sun, Mei-Ping Lu & Yong-Xin Liu. (2020). A guide to human microbiome research: study design, sample collection, and bioinformatics analysis. Chinese Medical Journal Publish Ahead of Print, doi: https://doi.org/10.1097/cm9.0000000000000871

相關文章猜你喜歡

10000+：菌群分析 寶寶與貓狗 梅毒狂想曲 提DNA發Nature Cell專刊 腸道指揮大腦

系列教程：微生物組入門 Biostar 微生物組  宏基因組

專業技能：學術圖表 高分文章 生信寶典 不可或缺的人

一文讀懂：宏基因組 寄生蟲益處 進化樹

必備技能：提問 搜索  Endnote

文獻閱讀 熱心腸 SemanticScholar Geenmedical

擴增子分析：圖表解讀 分析流程 統計繪圖

16S功能預測   PICRUSt  FAPROTAX  Bugbase Tax4Fun

在線工具：16S預測培養基 生信繪圖

科研經驗：雲筆記  雲協作 公眾號

編程模板: Shell  R Perl

生物科普:  腸道細菌 人體上的生命 生命大躍進  細胞暗戰 人體奧秘  

寫在後面

為鼓勵讀者交流、快速解決科研困難，我們建立了「宏基因組」專業討論群，目前己有國內外5000+ 一線科研人員加入。參與討論，獲得專業解答，歡迎分享此文至朋友圈，並掃碼加主編好友帶你入群，務必備註「姓名-單位-研究方向-職稱/年級」。PI請明示身份，另有海內外微生物相關PI群供大佬合作交流。技術問題尋求幫助，首先閱讀《如何優雅的提問》學習解決問題思路，仍未解決群內討論，問題不私聊，幫助同行。

學習16S擴增子、宏基因組科研思路和分析實戰，關注「宏基因組」

點擊閱讀原文