孤獨癥譜系障礙兒童腸道菌群多樣性研究及功能預測分析

劉智程 吳德 屈愛娜 王璐璐

（安徽醫科大學第一附屬醫院小兒神經康復中心，安徽 合肥 230022）

孤獨癥譜系障礙（autism spectrum disorder，ASD）是一種以社交和溝通障礙、重復刻板行為及狹窄興趣為核心癥狀的神經發育障礙性疾病。近年來，全球范圍內ASD患病率快速上升。據統計，2018年美國8歲兒童中每44人就有1人患ASD［1］。2020年，我國6～12歲兒童ASD患病率為7‰［2］。目前，典型ASD尚無有效的治療方法，給患兒家庭和社會帶來了巨大負擔。

ASD的病因尚不清楚，臨床實踐觀察到近一半ASD兒童共患便秘、腹瀉等胃腸道問題，因此，有研究者認為微生物-腸-腦軸（microbiota-gutbrain axis，MGBA）紊亂可能與ASD發病相關。已有為數不多的研究發現，ASD兒童腸道菌群的多樣性或結構與正常發育（typically developing，TD）兒童不同，但具體研究結果并不統一。Lou等［3］及Kang等［4］研究均表明ASD兒童腸道菌群的多樣性降低；Finegold等［5］研究則顯示ASD兒童腸道菌群的多樣性升高；此外，也有研究報道ASD兒童腸道菌群的α多樣性無明顯變化［6］。這些研究均基于操作分類單元（operational taxonomic units，OTUs）聚類的方法進行物種注釋，其物種注釋和多樣性分析的準確性均低于基于擴增子序列變體（amplicon sequence variants，ASVs）的分析方法。本研究采用基于ASVs進行物種注釋的16S rDNA高通量測序技術，前瞻性地分析ASD兒童腸道菌群的分布及多樣性特征，并預測分析2組菌群代謝功能的差異，為進一步探索兒童ASD與腸道菌群的相關性，以及ASD的微生態療法提供科學理論依據。

1 資料與方法

1.1 研究對象

選取安徽醫科大學第一附屬醫院小兒神經康復中心門診確診為ASD的30例兒童納入ASD組，同時招募20例TD兒童作為對照（TD組）。ASD組需符合美國精神障礙診斷與統計手冊第5版中ASD的 診 斷 標 準［7］，且 兒 童 孤 獨 癥 評 定 量 表（Childhood Autism Rating Scale，CARS）評 分≥30分。此外，2組需符合以下共同的納排標準。納入標準：（1）年齡≤12歲；（2）監護人同意此研究并簽署知情同意書。排除標準：（1）早產及過期產兒童；（2）有神經系統或遺傳代謝性疾病史；（3）有消化道畸形或消化道手術史；（4）合并胃腸道急慢性感染；（5）有糞便采樣前1個月內抗生素、微生態制劑、質子泵抑制劑或激素類藥物使用史。本研究經安徽醫科大學第一附屬醫院倫理委員會批準（PJ2021-14-10）。

1.2 基本信息

受試者主要監護人填寫問卷，包括受試者年齡、性別、出生方式和喂養方式、居住環境及嬰兒期抗生素使用史等基本信息。

1.3 樣本采集與保存

收集2組兒童上午10時之前的新鮮糞便2～3 g于無菌采樣管中，立即低溫轉運至-80℃冰箱保存。

1.4 DNA提取、文庫制備及測序

采用QIAamp Fast DNA Stool Mini Kit（QIAGEN）提取樣本DNA，采用NanoDrop ND-1000（Thermo Fisher Scientific，Waltham，MA，USA）分光光度計測定DNA數量和質量，瓊脂糖凝膠電泳檢測DNA完整性。利用引物515F（5'-GTGCCAGCMGCCGCGGTAA-3'）和806R（5'-GGACTACHVGGGTWTCTAAT-3'）PCR擴 增16S rDNA V4區。用AMPure XP Beads（Beckman Coulter，Indianapolis，IN）試劑盒純化PCR產物，使用Qubit（Thermofisher，Qubit 4.0）進行文庫定量。采用Illumina NovaSeq6000平臺進行雙端測序。使用Vsearch 2.4.4組裝拼接雙端序列，質控并過濾低質量數據、接頭序列或PCR錯誤等，再經過Vsearch 2.15.0去重、降噪、去嵌合體等步驟得到ASVs，樣本中含量低于總序列0.001%的ASVs將被去除。通過Qiime2基于Silva138數據庫對ASVs的代表序列進行物種注釋。

1.5 生物信息學分析

主要使用Qiime2和R包3.2.0對測序數據進行分析。使用Qiime2計算α多樣性指數（Chao1、Shannon和Simpson）評估物種的豐富度和均勻度。繪制稀釋性曲線評價測序深度是否合理。采用Anosim分析檢驗組間差異是否顯著大于組內，從而判斷分組是否有意義。采用基于非加權UniFrac距離的主坐標分析（principal component analysis，PCoA）評估組間物種組成差異。采用Wilcoxon秩和檢驗分析選取2組樣本在各分類水平上的差異物種。采用線性判別影響大小（linear discriminant analysis effect size，LEfSe）分析方法篩選組間生物標志物（即差異顯著物種）。構建隨機森林模型對樣本進行分類預測，采用受試者工作特征（receiver operating characteristic，ROC）曲線驗證模型的可靠性。采用Bugbase軟件預測菌群表型特征。基于Picrust2 2.3.0-b軟件預測菌群代謝功能。

1.6 統計學分析

使用SPSS 24.0軟件進行統計學分析。非正態分布的計量資料以中位數（四分位數間距）［M（P25，P75）］表示，組間比較采用Wilcoxon秩和檢驗；計數資料采用例數和百分率（%）表示，組間比較使用卡方檢驗。P＜0.05為差異有統計學意義。

2 結果

2.1 2組兒童基本信息比較

ASD組和TD組兒童入組年齡、性別、出生方式、喂養方式、嬰兒期抗生素使用史、鄉村居住史及寵物喂養史方面比較差異均無統計學意義（P＞0.05），見表1。

表1 2組兒童基本信息比較

2.2 測序數據統計

測序讀長范圍245～292 bp。50例樣本測序得到總有效序列數為6 428 113。測得樣本序列數最小為107 491，均數為128 562。稀釋性曲線顯示測序深度合理，見圖1。得到ASVs總數為1 354，其中146個為ASD組特有，40個為TD組特有。

圖1 ASVs稀釋性曲線 曲線趨于平緩，提示測序深度已基本覆蓋樣本中所有物種。

2.3 物種多樣性分析

α多樣性指數中Chao1指數反映物種豐富度，Shannon和Simpson指數均反映物種的豐富度和均勻 度。ASD組 和TD組 樣 本Chao1、Shannon及Simpson指數比較差異均無統計學意義（P＞0.05），見表2。

表2 2組兒童腸道菌群α多樣性分析［M（P25，P75）］

β多樣性分析利用不同樣本的豐度變化或進化關系來計算樣本間距離，從而反映樣本或組間微生物組成的相似性。Anosim分析可判斷分組是否有意義，R＞0說明組間差異大于組內，P值說明了差異的顯著性。PCoA是一種非約束性數據降維分析方法，可用來研究樣本群落組成的相似性。基于非加權UniFrac距離的Anosim及PCoA分析顯示，2組兒童腸道菌群的分布不同，見圖2。

圖2 非加權UniFrac距離Anosim分析（A）和PCoA組間差異分析（B） A：Between代表組間差異，ASD、TD分別表示各組內差異；圖中各箱線圖的凹槽互不重疊，說明各組中位數差異顯著。B：每個點代表一個樣本，點間距反映樣本相似性；圖中ASD組和TD組散點部分分離，提示2組間存在一定差異。

2.4 物種分布及差異分析

ASVs經物種注釋共得到8門、12綱、26目、44科、107屬。ASD組和TD組樣本中相對豐度位列前四的菌門相同，依次為厚壁菌門、擬桿菌門、放線菌門和變形菌門，這4種菌門的總豐度在各組中占比均超過99%，其中厚壁菌門和擬桿菌門總豐度在各組中占比均超過89%。2組樣本中相對豐度位列前三的菌屬均為擬桿菌屬、糞桿菌屬和雙歧桿菌屬。在門和綱水平，2組樣本的菌群分布差異無統計學意義（P＞0.05）。在屬水平，ASD組中有6個菌屬的相對豐度大于TD組，分別為巨單胞菌屬（P＜0.001）、巴恩斯氏菌屬（P=0.004）、小桿菌屬（P=0.006）、巨球菌屬（P=0.020）、瘤胃球菌屬扭鏈群（P=0.028）及梭桿菌屬（P=0.040）。2組樣本屬水平主要菌群構成見圖3。

圖3 2組樣本屬水平物種相對豐度分布圖 右側圖例列出了屬水平相對豐度位列前十的物種名稱，其余物種歸為“其他”。

2.5 LEfSe篩選組間生物標志物

應用LEfSe分析方法選取各組優勢物種，即組間生物標志物。結果顯示，ASD組的優勢物種為巨單胞菌屬、小桿菌屬、巨球菌屬、消化鏈球菌科、巴恩斯氏菌屬、瘤胃球菌屬扭鏈群及鏈型桿菌屬等；TD組優勢物種為毛螺菌科、消化球菌科、Saccharimonadaceae、TM7x、阿里松氏菌屬、帕拉普氏菌屬、乳桿菌目_g_P5D1_392及Saccharimonadales等。見圖4。

圖4 LEfSe分析篩選各組生物標志物 柱狀圖的長度代表差異顯著物種的影響大小（LDA score），其絕對值大于預設值（默認為2）表示該物種為差異有統計學意義的標志物。

2.6 隨機森林組間預測

構建隨機森林模型對ASD組和TD組樣本進行鑒別，結果顯示以巨單胞菌屬、罕見小球菌屬及小桿菌屬等屬水平細菌為特征可以較準確地區分2組樣本（錯誤率為30.83%），其中巨單胞菌屬對分類預測的效果影響最大，見圖5。繪制ROC曲線評估隨機森林分類模型的優劣，曲線下面積（area under curve，AUC）為0.82，說明該模型分組預測具有一定的準確性，見圖6。

圖5 隨機森林物種重要性點圖 橫坐標（平均準確度下降）代表重要性水平，縱坐標為模型中對分類效果起主要作用的菌屬，按作用從大到小排列。

圖6 評估隨機森林分類模型的ROC曲線 AUC越接近于1，分類結果越準確；AUC在0.7～0.9時，分類結果有一定準確性。縱坐標（靈敏度）越高代表分類的準確率越高，橫坐標（1-特異度）越低代表誤判率越低，曲線越靠近左上角，說明模型的分類效果越好。

2.7 菌群表型特征預測

ASD組和TD組兒童腸道菌群的表型，如革蘭陽性、革蘭陰性、生物膜形成、致病性、移動元件含量、氧需求及氧化脅迫耐受等差異均無統計學意義（P＞0.05），見表3。

表3 2組兒童腸道菌群表型特征預測分析［M（P25，P75）］

2.8 菌群功能預測分析

ASD組菌群的細胞內吞、鞘糖脂生物合成-乳糖和新乳糖系列代謝途徑的豐度較TD組高（P值分別為0.031、0.018）；而ASD組非同源末端連接代謝途徑的豐度較TD組低（P=0.009）。

為了進一步發掘樣本中與腸-腦軸模塊相關的代謝信息，將測序數據基于KEGG數據庫進行基因注釋，并將注釋信息轉化為腸道代謝模塊和腸-腦軸模塊的分析結果。腸道代謝模塊分析顯示，ASD組菌群的谷氨酸降解Ⅰ（MF0030）、色氨酸降解（MF0025）、異化性硝酸鹽還原（MF0101）、丙酮酸脫氫酶復合體（MF0072）、甘油降解Ⅱ（MF0061）及丁酸鹽生成Ⅱ（MF0089）等代謝功能的豐度較TD組低，差異有統計學意義（P值分別為0.002、0.008、0.029、0.036、0.042、0.046）。腸-腦軸模塊分析顯示，ASD組菌群的谷氨酸降解Ⅰ（MGB050）、色氨酸降解（MGB049）、多巴胺降解（MGB023）及丁酸鹽生成Ⅱ（MGB053）等代謝功能的豐度較TD組低（P值分別為0.002、0.008、0.022、0.046），γ-氨 基 丁 酸 降 解（MGB019）功能的豐度較TD組高（P=0.006）。見圖7～8。

圖7 腸道代謝模塊菌群功能預測 MF0030：谷氨酸降解Ⅰ；MF0025：色氨酸降解；MF0101：異化性硝酸鹽還原；MF0072：丙酮酸脫氫酶復合體；MF0061：甘油降解Ⅱ；MF0089：丁酸鹽生成Ⅱ。ASD組和TD組各指標的比較均P＜0.05。

圖8 腸-腦軸模塊菌群功能預測 MGB050：谷氨酸降解Ⅰ；MGB019：γ-氨基丁酸降解；MGB049：色氨酸降解；MGB023：多巴胺降解；MGB053：丁酸鹽生成Ⅱ。ASD組和TD組各指標的比較均P＜0.05。

3 討論

腸道菌群是人體腸道微生態系統的主要成分，數量約為人體細胞總數的10倍，被稱為人體的“第二基因組”［8］。腸道菌群參與人體營養物質代謝、維生素合成、腸道黏膜生物屏障形成、免疫系統成熟及免疫應答調節等多種生理活動，進而影響人體多系統功能［9］。研究表明，腸道菌群可能對腦功能有重要影響，并且與ASD密切相關［10］。

多數相關研究表明，ASD兒童腸道菌群的α多樣性較健康兒童低［4］。而本研究則顯示，2組樣本的α多樣性指數Chao1、Shannon及Simpson比較差異均無統計學意義，可能與樣本數量少有關。本研究發現2組樣本菌群結構存在差異，屬水平共發現6個差異顯著物種，分別為巨單胞菌屬、巴恩斯氏菌屬、小桿菌屬、巨球菌屬、瘤胃球菌屬扭鏈群及梭桿菌屬。

關于ASD兒童腸道巨單胞菌屬含量的研究報道并不一致，部分研究結果相反［11］。Zou等［12］發表ASD兒童腸道巨單胞菌屬的含量增多，與本研究相同。巨單胞菌屬在腸道可以將葡萄糖發酵成乙酸和丙酸，用于合成脂質和膽固醇，也是宿主的能量來源之一［13］。此外，該菌屬在肥胖者腸道中含量增多［14］。趙銳豪等［15］發現巴恩斯氏菌屬是ASD兒童腸道的優勢菌群之一，與本研究結果相似。另一項研究表明，相較于TD兒童，合并便秘的ASD兒童腸道巴恩斯氏菌屬的豐度升高［16］。也有研究顯示，ASD兒童腸道巴恩斯氏菌屬的平均豐度降低［17］。巴恩斯氏菌屬可競爭性抑制人體腸道致病菌，并與免疫調節及炎癥性疾病相關［18-20］。但該菌屬與ASD的具體相關性仍不清楚，需進一步探索。研究顯示ASD患者腸道中小桿菌屬的含量降低［11］。而本研究結果卻與之相反，這可能與采用的統計方法不同或本研究樣本量小有關。研究表明，對食物過敏的兒童腸道小桿菌屬的含量降低［21］。因此，研究ASD與小桿菌屬之間關系時尚需考慮過敏因素。研究顯示ASD兒童腸道巨球菌屬的含量增多［11］，與本研究結果相同。而一項病例對照研究則顯示，ASD男孩腸道中巨球菌屬的豐度降低［22］。ASD兒童腸道巨球菌屬的豐度是否受性別影響需進一步驗證。巨球菌屬可利用乳酸產生丁酸等短鏈脂肪酸［23］。此外，肥胖者腸道中該菌屬的含量增多［24］。然而，巨球菌屬是否參與ASD的病理過程仍未可知。本研究中，與TD兒童相比，ASD兒童腸道瘤胃球菌屬扭鏈群的相對豐度升高，這與Ahmed等［25］研究結果相似。一項前瞻性研究顯示，早期兒童腸道扭鏈瘤胃球菌的相對豐度與ASD兒童社交反應量表-2得分呈負相關［26］。可見，扭鏈瘤胃球菌與ASD具有相關性。本研究顯示，ASD兒童腸道梭桿菌屬的相對豐度升高。有研究表明，相較于TD兒童，有便秘癥狀的ASD兒童腸道梭桿菌屬含量增多［15］；母體免疫激活子代大鼠的腸道梭桿菌屬豐度升高，并與重復行為和焦慮樣癥狀相關［27］。由此可見，ASD可能與腸道梭桿菌屬之間存在相關性。

盡管多數相關研究均認為ASD與腸道菌群密切相關，但進一步探索二者之間因果關系的研究卻很少。相關研究多傾向于認為腸道菌群異常可能參與了ASD的病理過程，雖然具體機制尚不明確，但目前較統一的認識是基于MGBA理論。這一理論認為，微生物及其代謝產物可通過神經內分泌、神經免疫及神經通路等途徑影響腦功能［28］。結合本研究，我們不僅發現ASD兒童腸道菌群結構上的異常，而且通過功能預測分析還發現，ASD兒童腸道菌群的代謝功能與TD兒童相比存在顯著差異，如丁酸鹽生成，色氨酸、谷氨酸及γ-氨基丁酸等活性物質的降解功能等，而這些物質均對腦功能有重要影響。因此，腸道菌群異常可能通過MGBA介導神經炎癥、免疫異常及神經遞質紊亂等病理過程，進而參與ASD的發生與發展。

本研究采用了高通量16S rDNA測序技術，創新之處在于使用了去噪方法獲得ASVs，基于ASVs的物種注釋相較于以往的OTUs聚類更加準確，且與其他研究的ASVs數據結果進行比較時，更具可比性和一致性。本研究也存在不足之處：首先，樣本量相對較小，針對相關因素進行多組分析較困難；其次，未對血液及糞便代謝組學進行檢測分析，以驗證基于擴增子測序對腸道菌群代謝功能預測的結果。期望以后的研究采用大樣本、考慮多因素并使用代謝組學等更完善的研究方法更加深入地探究腸道菌群與ASD的關系。

綜上所述，ASD兒童和TD兒童腸道菌群的構成存在顯著差異，臨床實踐中通過調節腸道菌群結構以改善微生態可能有助于ASD兒童獲得更好的預后；預測2組兒童腸道菌群的代謝功能存在顯著不同，其中腸道菌群的丁酸鹽生成、谷氨酸及γ-氨基丁酸降解等代謝功能可能是未來相關研究需關注的重點。

利益沖突聲明：所有作者均聲明不存在利益沖突。