基于高通量測序技術分析青藏高原牦牛和犏牛乳中微生物多樣性的研究

馬 靜，張琳琳，柴沙駝，王 迅，崔占鴻，孫 璐, ，劉書杰,

(1.青海大學畜牧獸醫科學院，青海西寧 810016；2.青海高原牦牛研究中心，青海西寧 810016；3.青海省高原放牧家畜動物營養與飼料科學重點實驗室，青海西寧 810016)

青藏高原及其毗鄰地區的牦牛是當地牧民生產生活及畜牧經濟主要的畜種之一，能很好地適應高海拔缺氧的環境并為人類提供乳、肉和運輸等基本資源[1]。牦牛乳被稱為天然濃縮乳，所含的營養物質如蛋白質、脂肪、必需氨基酸、礦物質(磷除外)大多均高于其他牛乳[2-4]。牦牛與奶(肉)牛的F1后代為犏牛[5]，具有雜種優勢，適應性更強，且母犏牛比母牦牛的乳產量高。同時有研究發現雖然牦牛乳和犏牛乳都具有豐富的脂肪酸，犏牛乳具有高端牛奶的資源優勢，乳質優良，氣味芬芳，營養豐富[6-8]。牛乳作為人類飲食生活的重要組成成分，幾乎適合所有年齡段的人食用，乳中豐富的微生物區系與人類健康密切相關，一些微生物可能會降低炎性腸病，腹瀉，便秘，食物過敏甚至結腸癌的發病率[9]。牛乳除了含有多種對機體健康有益的微生物之外，還存在一些對人類健康有威脅的致病菌如：李斯特菌、沙門氏菌、大腸桿菌等[10]。

隨著科學技術的發展，牛乳微生物組成的神秘面紗被逐漸揭開。早期研究中，采用傳統的美蘭實驗法對牦牛乳微生物的含量有了初步的了解[11]。此外，肜豪峰等[12]采用稀釋平板涂布法分離計數對牦牛乳及其制品中微生物區系進行了分析。微生物在極端環境下經過長期的自然選擇，具備相對特殊的結構、機能以及遺傳特性[13-14]。其群落結構組成不僅受環境差異影響還和來源動物有關，宿主動物的來源不同直接影響著乳微生物群落結構組成[15]。已知犏牛乳品質的雜種優勢，但其微生物組成及種類還不清楚。因此，開展牦牛乳和犏牛乳中的微生物組的研究非常必要。

近年來越來越多的證據表明，高通量測序技術可以克服基于微生物培養的檢測方法的局限性，因此被用來快速鑒定微生物和發現新微生物[16-17]，廣泛的用于各個領域。如16S rDNA來檢測食品中的微生物群落結構及鑒定綿羊糞便微生物群的多樣性和功能[18-19]。綜上所述，牦牛乳和犏牛乳作為青藏高原特色乳其營養價值較高，但目前對兩種原乳微生物結構組成差異及優勢菌群了解較少，因此本試驗通過提取牦牛乳和犏牛乳中微生物的DNA進行高通量測序，更加全面的展現了牦牛和犏牛乳微生物組成差異及多樣性，以期為牦牛乳及犏牛乳制品的開發利用及優勢菌群的發掘等提供基礎依據。

1 材料與方法

1.1 材料與儀器

牛乳樣品采集信息 采樣時間：2019年7月至8月。采樣地點：在青海省貴南縣塔秀鄉達隆村夏曲駿周太畜牧業專業合作社。選取自然放牧狀態下全天候隨群放牧的未進行任何補飼的體況一致的泌乳中期4～5胎次的健康牦牛(yak)和犏牛(dzo)各5頭，每頭各采3份，共采生鮮牦牛乳15份，生鮮犏牛乳15份，兩種生鮮乳各五組，一組因污染剔除，牦牛乳和犏牛乳各剩四組(yak5.2,yak5.6,yak5.7,yak5.10;dzo16.2,dzo16.6,dzo16.7,dzo16.10)，每組乳樣取50 mL，分為兩部分：取45 mL至離心管(康寧/Corning公司)中，干冰速凍后于-20 °C冰箱儲存，用于乳常規營養成分分析。取5 mL于無菌無酶凍存管(康寧/Corning公司)中，液氮速凍后運回實驗室-80 °C冰箱儲存，用于測定微生物多樣性及組成。

FOSS多功能乳制品分析儀(MilkoScan FT1) 上海瑞玢國際貿易有限公司；DCC體細胞檢測儀 利拉伐天津有限公司；超低溫冰箱(-80 ℃) 美國Thermo公司；Fresco21型高速冷凍離心機 美國Thermo公司；Tprofessiona PCR儀 德國Biometral公司；Gel DOC XR凝膠成像系統、PowerPacUniversal水平電泳儀、SUBCELL GT電泳槽(20 cm×25 cm)美國Bio-Rad公司。

1.2 實驗方法

1.2.1 牦牛乳和犏牛乳常規營養成分分析 采用FOSS FT-1乳品分析儀測定牦牛乳和犏牛乳常規營養成分。

1.2.2 基因組總DNA提取及PCR擴增、測序 采用SDS法提取基因組總DNA，通過瓊脂糖凝膠電泳檢測DNA完整性，并用無菌水將其稀釋至1 ng/μL作為PCR模板，使用帶 Barcode 的特異性引物515F(5′-GTGCCAGCMGCCGCGG-3′)和806R(5′-GGA CTACHVGGGTWTCTAAT-3′)，New England Biolabs

公司的Phusion? High-Fidelity PCR Master Mix with GC Buffer及高效高保真酶對特定區域進行PCR擴增。擴增條件：95 ℃預變性2 min后，95 ℃變性30 s，55 ℃復性30 s，72 ℃延伸1 min，共30個循環，72 ℃延伸5 min。

所得PCR產物用2%瓊脂糖凝膠電泳檢測，對目的條帶使用Qiagen公司提供的膠回收試劑盒(MinElute Gel Extraction Kit)回收產物。然后使用TruSeq? DNA PCR-Free Sample Preparation Kit建庫試劑盒進行文庫構建，構建好的文庫經過Qubit和Q-PCR定量，使用NovaSeq6000進行上機測序。PCR擴增、PCR產物的混樣、純化，文庫的構建和上機測序流程均由天津諾禾致源生物信息科技有限公司完成。

1.3 數據分析

1.3.1 測序數據處理 根據Barcode序列和PCR擴增引物序列，截去Barcode和引物序列后使用FLASH(Version1.2.http://ccb.jhu.edu/software/FLAS H/)對每個樣本的Reads進行拼接，得到的拼接序列為原始Tags數據；拼接得到的原始Tags數據，經過嚴格的過濾處理得到高質量的Tags數據。參照Qiime(Version1.9.1.http://qiime.org/scripts/split_libr aries_fastq.html)的Tags質量控制流程，進行Tags截取及長度過濾之后進行去除嵌合體序列的處理，Tags序列通過(https://github.com/torognes/vsearch/)與物種注釋數據庫進行比對檢測嵌合體序列，并最終去除其中的嵌合體序列，得到最終的有效數據。

1.3.2 菌群多樣性分析 通過Uparse軟件(Uparse version7.0.1001.http://www.drive5.com/uparse/)對所有樣本進行聚類，用Mothur方法與SILVA132(http://www.arb-silva.de/)的SSUrRNA數據庫進行物種注釋分析(設定閾值為0.8～1)，獲得在各個分類水平分類學信息。通過Qiime軟件(Version 1.9.1)計算Chao1、Shannon、Simpson指數，使用R軟件(Version 2.15.3)繪制稀釋曲線及Alpha多樣性指數組間差異分析。

基于Qiime軟件(Version 1.9.1)計算Unifrac，并通過R軟件進行Beta多樣性指數組間差異分析。然后對菌群組成進行多級物種差異判別分析(LEfSe)，默認設置LDA Score的篩選值為4，且利用R軟件分析組間差異顯著的物種。

2 結果與分析

2.1 牦牛乳和犏牛乳常規營養成分對比

牦牛乳和犏牛乳的9個常規營養成分檢測，并利用SAS 9.4在0.05水平下對其進行顯著性單因素方差分析，結果如表1所示。牦牛乳中各項指標的檢測結果均高于犏牛乳，且二者蛋白、非脂固形物、總乳固形物及密度的差異顯著(P<0.05)。

表1 牦牛和犏牛乳常規營養成分Table 1 Nutritional component contents of yak and cattle-yak milk

2.2 微生物多樣性分析

2.2.1 樣品測序結果及多樣性指數 經PCR-free文庫構建、測序、Reads拼接后平均每樣品測得93918條序列，經質控過濾平均得到75922條有效數據，質控有效數據量達63698，質控有效率達68.28%。基于97%一致性對OTUs進行注釋，共得到5516個OTUs。如圖1所示，稀釋曲線逐漸趨于平緩，表明測序結果合理，且間接反映了物種的豐富度，即樣本的當前測序深度足夠反映該群落樣本所包含的微生物多樣性。如圖2所示，Rank Abundance曲線直觀地反映樣品中物種的豐富度和均勻度，水平方向的曲線跨度體現了物種豐富度，曲線在橫軸上的范圍較大，表明物種豐度較高；垂直方向的曲線漸進平緩，表明物種分布均勻。

圖1 Alpha多樣性指數稀釋曲線圖Fig.1 Rarefaction curve for Alpha diversity index

圖2 Rank Abundance曲線Fig.2 Rank Abundance curve

如表2所示，采用Shannon和Simpson指數對牦牛乳和犏牛乳中微生物的多樣性進行評價，采用Observed species和Chao1對二者微生物菌群的豐度進行評價，結果表明犏牛乳中所含物種數較牦牛乳高，二者微生物菌群豐度差異顯著，但菌群多樣性相似。

表2 基于T檢驗的Alpha多樣性指數分析Table 2 Alpha diversity index analysis based on T test

2.2.2 微生物群落組成分析

2.2.2.1 物種韋恩圖分析 由圖3可知，共測得5516個OTUs，共有OTUs 2506個，牦牛乳組(yak)和犏牛乳組(dzo)分別特有1054和1956個OTUs，犏牛乳特有的OUTs相對牦牛乳較高，表明犏牛乳中的微生物群落較牦牛乳豐富。

圖3 韋恩圖Fig.3 Venn graph

2.2.2.2 微生物門分類水平的比較 牦牛乳和犏牛乳在門水平上前10種微生物相對豐度比較結果如圖4顯示，牦牛乳和犏牛乳中變形菌門(Proteobacteria)的相對豐度分別為29.80%、45.36%；厚壁菌門(Firmicutes)的相對豐度分別為35.99%、25.79%；放線菌門(Actinobacteria)的相對豐度分別為8.41%、7.39%。從菌群門水平排名前十的熱圖(圖5)可知，犏牛乳中芽單胞菌門(Gemmatimonadetes)、酸桿菌門(Acidobacteria)及變形菌門等較為常見，而牦牛乳則以厚壁菌門、擬桿菌門(Bacteroidetes)及柔膜菌門(Tenericutes)為主。

圖4 門水平下微生物的相對豐度Fig.4 Relative abundance of microorganisms at phylum level

圖5 微生物門水平熱圖Fig.5 The heat map of microorganisms at phylum level

2.2.2.3 微生物屬分類水平的比較 牦牛乳和犏牛乳微生物前10屬水平比較如圖6所示，結果表明牦牛乳中的優勢菌屬為未分類的藍藻細菌屬(unidentified-Cyanobacteria)；犏牛乳中優勢菌屬為慢生根瘤菌屬(Bradyrhizobium)，其次是耶爾森菌屬(Yersinia)、乳桿菌屬(Lactobacillus)及貪銅菌屬(Cupriavidus)。兩種牛乳中都含有這幾種微生物，但所占比例差異較大。耶爾森菌屬和貪銅菌屬在犏牛乳中所占比例明顯高于牦牛乳，而乳桿菌屬的含量在犏牛乳中較牦牛乳相比更低。

圖6 屬水平下微生物的相對豐度Fig.6 Relative abundance of microorganisms at genus level

2.2.2.4 微生物種分類水平的比較 牦牛乳和犏牛乳微生物前10種水平比較表明，埃爾坎尼中慢生根瘤菌種(Bradyrhizobium elkanii)和Kosakonia oryzae為牦牛乳和犏牛乳中的優勢菌種，其中埃爾坎尼中慢生根瘤菌種在犏牛乳中的相對豐度較高，達11.53%，比牦牛乳中所占比例高；菌種Kosakonia oryzae在犏牛乳(5.53%)的占比同樣比牦牛乳(3.36%)高(圖7)。

圖7 種水平下微生物的相對豐度Fig.7 Relative abundance of microorganisms at species level

2.2.2.5 物種主坐標分析 基于unweightedUnifrac距離進行物種主坐標分析(圖8)，結果顯示，第一主成分對樣品差異的解釋率為25.57%，第二主成分對樣品差異的解釋率為18.39%。牦牛乳和犏牛乳分別聚集在一起，兩組樣本距離較近，由此可知，牦牛乳和犏牛乳微生物群落差異輕微。

圖8 基于unweightedUniFrac距離的PCoA主坐標分析Fi g.8 PCoA analysis based on unweightedUniFrac distance

2.2.3 組間差異分析 對牦牛、犏牛乳的LEfSe分析，LDA線性判別分析結果見圖9A，LDA score>4的物種共有八種，其中莫拉菌科(Moraxellaceae)、假單 胞 菌 目(Pseudomonadales)、莫 拉 氏 菌 屬(Moraxella)及肉桿菌科(Carnobacteriaceae)對牦牛乳的菌群結構具有顯著影響；δ-變形菌綱(Deltaproteobacteria)、unidentified_Gammaproteobacteria、伯克氏菌科(Burkholderiaceae)及α-變形菌綱(Alphaproteobacteria)則與犏牛乳微生物群落結構相關。

圖9 LEfSe組間差異分析Fig.9 LEfSe analysis of differences between groups

2.3 乳常規和微生物相關性分析

圖10為優勢物種Spearman相關系數分析熱圖，揭示了環境因子與物種之間的關系。縱向為環境因子信息(脂肪、蛋白質、乳糖、非脂固形物和總乳固形物)，橫向為物種信息(微生物)，中間熱圖對應的值為Spearman相關系數r，介于-1和1之間，r<0為負相關，r>0為正相關。由圖可知，牦牛乳和犏牛乳中的共性優勢菌門中的變形菌門與乳脂肪、蛋白質、乳糖、非脂固形物及總乳固形物均呈負相關；放線菌門與乳糖呈負相關，與其他的指標呈正相關；牦牛乳中主要的厚壁菌門與乳脂肪、蛋白質、乳糖、非脂固形物及總乳固形物均呈正相關；異常球菌-棲熱菌門(Deinococcus-Thermus) 與牦牛乳和犏牛乳中具有顯著差異的蛋白質呈顯著正相關(P<0.05)；犏牛乳中常見的芽單胞菌門與脂肪、乳糖及總乳固形物呈正相關，與蛋白質、非脂固形物呈負相關。

圖10 Spearman相關性分析熱圖Fig.10 Spearman correlation analysis heat map

3 討論

牦牛乳作為藏區牧民營養攝取的主要來源之一，對多種疾病有較好的預防作用，有助于降低代謝綜合征發病率[20]。本試驗通過對牦牛乳和犏牛乳常規營養成分的測定，得出牦牛乳和犏牛乳的蛋白質、非脂固形物、總乳固形物及密度差異顯著(P<0.05)。引起乳營養成分含量差異的原因有很多，如牛群遺傳特性、泌乳生理差異及日糧營養結構等[21]。而該試驗在同一飼養環境中，且年齡、泌乳期、胎次、健康狀況基本相同的前提下進行的，則出現這種差異的主要原因可能取決于品種及個體遺傳因素。

隨著科學技術的進步，人們已經了解到家畜的原乳有自己獨特的內原微生物群落，對人類和家畜后代的健康都有影響[22]。本試驗結果顯示犏牛乳的細菌豐富度較牦牛乳高，而兩種乳的細菌多樣性相似。犏牛作為牦牛和黃牛的雜交品種，其雜種優勢明顯，由韋恩圖分析結果可知牦牛乳和犏牛乳二者所共有的微生物物種很多，但犏牛乳中特有的微生物與物種高于牦牛乳，推測犏牛雜種優勢及遺傳特性可能是其乳中細菌豐富度較高的原因[23-24]。在微生物門分類水平上，與張敏[25]的研究結果相似。此外，本研究結果牦牛乳中擬桿菌門的含量也較多。擬桿菌門的擬桿菌通過抑制脂肪細胞中脂蛋白脂酶活性，從而對脂質代謝產生積極影響，且該菌與家畜腸道免疫有關[26]。當機體脂質代謝異常時則會引起許多常見的疾病，對機體健康造成威脅。因此，牦牛乳可能對家畜脂質代謝及提高食用者及有犢牛的腸道免疫有積極影響，更有利于脂肪在小腸的消化吸收為機體提供所需要的能量，但這需要進一步試驗來證明。本試驗在屬水平上主要的優勢屬為慢生根瘤菌屬、耶爾森菌屬及未分類的藍藻細菌屬，與發酵乳制品優勢菌屬相比完全不同，鮮乳通過不同的加工方式加工后其細菌的多樣性及種屬結構可能會發生改變。耶爾森菌屬中的部分菌種會通過破壞細胞死亡途徑，擾亂炎癥過程并利用免疫細胞促進疾病，對人類具有致病性[27]。而結果中含有該細菌，可能是因為牛乳在擠奶過程中被外界環境所污染，有研究表明，耶爾森菌屬多存在于糞便中，由此推斷，有害菌群的產生可能是因為牦牛飼養管理不當，采樣時牛乳可能被糞便污染所致[28]。乳桿菌屬是一類能夠利用碳水化合物產生乳酸的桿狀細菌，長期以來被用于制造乳制品，該菌屬有25個菌種部分具有產生抗氧化劑并且具有益生菌潛力，通過抗炎反應參與機體的免疫調節[29]。犏牛乳中乳桿菌屬的比例更高，可能是更理想的原料奶。

通過 Spearman相關性分析比較乳常規和微生物的關系，牦牛乳和犏牛乳中的變形菌門與乳脂肪、蛋白質、乳糖、非脂固形物及總乳固形物呈負相關，而厚壁菌門與乳常規營養成分呈正相關。有研究報道，厚壁菌門與小鼠肥胖密切相關，促進腸道對單糖的吸收，誘導脂肪形成[30-31]。細菌門水平分類比較結果表明，牦牛乳中的厚壁菌門含量較犏牛乳高，而乳常規營養成分分析結果同樣也表明牦牛乳中脂肪和乳糖含量比犏牛乳高。牦牛乳中脂肪和乳糖含量較高的原因可能與此有關。

4 結論

營養成分分析：牦牛乳和犏牛乳蛋白質、總乳固形物、非脂固形物及密度差異顯著(P<0.05)，其他乳常規營養成分差異不顯著。微生物多樣性分析：牦牛乳和犏牛乳微生物豐度差異明顯，且犏牛乳的微生物豐度及多樣性較高。通過Spearman相關性分析發現兩種乳營養成分含量的差異與其微生物組成之間具有相關性。這為進一步研究青藏高原牦牛乳和犏牛乳的營養調控提供了理論依據。