基于高通量測序分析年糕菌群結構和優勢菌屬

蔡懷依 雷麗萍 婁永江 李勇勇

(寧波大學海洋學院，寧波 315211)

年糕(rice pastry)是一種具有悠久歷史的半方便性的米制品，有濃郁的地方傳統和特色[1]。其營養豐富，水分含量高，水分活度大，且整個加工過程中年糕原輔料半成品及成品都直接與水或空氣接觸，這也直接為微生物提供了適宜的生長條件[2]。通過分析年糕儲藏過程中的菌群結構，從而有目的的選用合適的抑菌劑，能為年糕品質劣化控制提供基礎資料。

目前，分析微生物群體的多樣性及群落結構的經典方法是分離、培養以及鑒定，需要進行一系列繁雜的形態特征和生理生化實驗，且能夠培養分離出的微生物僅占樣品的1%～10%，無法解析微生物的組成及豐度情況[3]。隨著宏基因組學概念的提出和高通量測序技術的發展，使樣品中不可培養的微生物、低豐度的微生物均能被檢測出，更完整地反映樣品中微生物的群落特征，使菌群的分析更準確和快速[4-7]。本研究采用宏基因組結合高通量測序技術分析年糕菌群結構和優勢菌屬，為進一步探討微生物與年糕品質的相關性提供參考。

1 材料與方法

1.1 材料與試劑

年糕：寧波江北五橋糧油食品有限責任公司提供；E.Z.N.A.Soil DNA Kit D5625-01：OMEGA；Qubit2.0 DNA檢測試劑盒Q10212：Life；Taq DNA Polymerase Ep0406：Thermo；Agencourt AMPure XP A63882：Beckman；實驗用水：經Milli-Q Integral 5超純水系統純化。

1.2 儀器與設備

Pico-21臺式離心機：Thermo Fisher；GL-88B漩渦混合器：海門市其林貝爾儀器制造有限公司；TND03-H-H混勻型干式恒溫器：深圳拓能達科技有限公司；DYY-6C型電泳儀、DYCZ-21電泳槽：北京市六一儀器廠；凝膠成像系統：美國UVP；Q32866 Qubit? 2.0熒光計：Invitrogen；T100TM Thermal Cyeler PCR儀：BIO-RAD；Research plus移液器： Eppendorf。

1.3 實驗方法

1.3.1 樣品預處理

新鮮年糕采取液氮研磨的方法將樣品磨成粉末狀，用于基因組DNA的提取。

1.3.2 宏基因組DNA的提取

采用E.Z.N.ATM Mag-Bind Soil DNA Kit的試劑盒(OMEGA公司)提取總DNA，提取的總DNA經1%瓊脂糖凝膠電泳檢測完整性，超微量分光光度計(Thermo NanoDrop 2000)檢測濃度。檢測合格后于-20 ℃保存，用于后續實驗。

1.3.3 PCR擴增及高通量測序

利用Qubit2.0 DNA檢測試劑盒對基因組DNA精確定量，以確定PCR反應加入的DNA量。PCR所用的引物已經融合了Miseq測序平臺的V4-V5通用引物，其中515F引物：CCCTACACGACGCTCTTCCGAT CTN (barcode) GTGCCAGCMGCCGCGGTAA，909R引物AGTTCCTTGGCACCCGAGAATTCCACCCCGY CAATTCMTTTRAGT。

通過兩輪PCR擴增并完成接頭序列的連接。第一輪PCR反應體系為30 μL反應液，含2×Taq master Mix 15 μL，Bar-PCR primer F(10uM) 1 μL，Primer R (10uM) 1 μL，Genomic DNA 10 ng， H2O補足至30 μL 。PCR反應條件為：94 ℃預變性 3 min；94 ℃30 s，45 ℃20 s，65 ℃ 30 s進行5個循環；之后94 ℃20 s，55 ℃ 20 s，72 ℃30 s進行20個循環，最后72 ℃延伸5 min。第二輪擴增引入Illumina橋式PCR兼容引物， PCR 反應條件為：95 ℃預變性30 s，然后95 ℃15 s，55 ℃15 s，72 ℃30 s進行5個循環，最后72 ℃延伸5 min。PCR 結束后，利用瓊脂糖電泳進行鑒定，之后進行DNA純化回收。在25 μL PCR產物中加入體積0.8倍的磁珠(Agencourt AMPure XP)，震蕩充分懸浮后放在磁力架上吸附5 min，小心的用移液槍吸出上清。加入30 μL0.8倍的磁珠洗滌液，震蕩充分懸浮后放在磁力架上吸附5 min，小心吸出上清。加入90 μL WashBuffer，反向放置在磁力架上，使磁珠吸附到PCR管的另外一面，充分吸附后吸出上清。將PCR管或8聯管放在55 ℃烘箱5 min，使里面的酒精完全揮發。加入30 μL Elution Buffer洗脫。將PCR管放在吸附架上5 min，充分吸附，移出上清到干凈的1.5 mL離心管中，定量備用。利用Qubit 2.0 DNA檢測試劑盒對回收的DNA精確定量，以方便按照1∶1的等量混合后測序。等量混合時，每個樣品DNA量取10 ng，最終上機測序濃度為20 pmol。

1.3.4 數據分析

Illumina MiseqTM得到的原始圖像數據文件經CASAVA堿基識別 (Base Calling) 分析轉化為原始測序序列(Sequenced Reads) ，其中包含測序序列 (reads) 的序列信息以及其對應的測序質量信息。Miseq測序序列中含有barcode序列，以及測序時加入的引物和接頭序列。首先需要去除引物接頭序列，再根據PE reads之間的overlap關系，將成對的reads 拼接 (merge) 成一條序列，然后按照barcode標簽序列識別并區分樣品得到各樣本數據，最后對各樣本數據的質量進行質控過濾，得到各樣本有效數據。去除3′端測序引物接頭，Read1 3′端測序接頭為TGGAATTCTCGGGTGCCAAGGAACTC，根據PE reads之間的overlap關系將成對reads拼接(merge) 成一條序列，拼接序列的overlap區域允許的最大錯配比率是0.1，根據各樣本barcode序列從融合后數據中分割出各樣本數據，去除各樣本中reads尾部質量值在20以下的堿基。設置10 bp的端口，如果窗口內的平均質量值低于20，從窗口開始去除后端的堿基。切除reads中含N部分序列，并去除數據中的短序列，長度閾值200 bp，隨后再對低復雜度的序列進行過濾。將多條序列按其序列間的距離進行聚類，對相似性在 0.97 以上的序列進行歸并，生成操作分類單元(OTU)[8-10]。根據聚類分析結果，計算 ACE、Chao1、Shannon、Simpson 指數進行 alpha 多樣性分析，其中 ACE、Chao1 指數是對菌群豐度進行評估，Shannon、Simpson 指數是對菌群多樣性進行評估。采用 RDP classifier 軟件對序列進行物種分類，根據分類學分析比對結果，在門、屬等水平上對樣品的菌群結構進行種類和豐度分析[11]。

2 結果與分析

2.1 序列長度分布

在年糕樣品總DNA的16S rDNA V4-V5區中測得的原始reads數目34 778條，原始序列平均長度為414.32(圖1)。質量控制之后剩余reads數目33 090條，序列平均長度為374.47(圖2)。從序列長度的分布來看， 與16S rDNA- V4-V5 區序列長度大致吻合。去除嵌合體與靶區域之前序列總數33 090條，比對到細胞器組織序列數目6 763條，非靶區域序列數目0條，嵌合體數目159條，處理后剩余序列26 168條。

圖1 原始數據長度分布圖

圖2 質控后序列長度分布圖

2.2 OTU多樣性分析

通過繪制OTU數目變化與聚類similarity值之間的關系圖，從中選擇最佳的similarity值進行OTU分析和分類學分析，分析使用的similarity值為97%的序列相似性(圖3)[12]。樣品的基因組 DNA 經 PCR 擴增 16S rDNA V4-V5 區后進行高通量測序，經統計分析獲得序列的 alpha 多樣性，結果見表 1。Shannon值越大，說明群落多樣性越高，而Simpson值越大說明群落多樣性越低[13]，由表1 可知樣本群落多樣性較高。此外，樣本的測序覆蓋率在1.00以上，表明樣品中序列未被測到的概率較低。

圖3 OTU數目與聚類相似度值關系圖

采用對測序序列進行隨機抽樣的方法，以抽到的序列數與它們所能代表OTU的數目構建曲線，即稀釋曲線(Rarefaction Curve)。它可以用來比較測序數據量不同的樣本中物種的豐富度，也可以用來說明樣本的測序數據量是否合理[14]。當曲線趨向平坦時，說明測序數據量合理，更多的數據量只會產生少量新的OTU，反之則表明繼續測序還可能產生較多新的OTU。樣本的Shannon-Wiener曲線(圖4)， 發現曲線已到達平臺期， 說明更大的測序量不會引起物種多樣性的顯著增長，基于現有數據量的分析結果準確可靠。

表1 年糕腐敗樣品菌群的alpha多樣性比較

注：橫坐標為樣本中隨機抽取序列數，縱坐標為所得相應的Alpha指數圖4 Alpha指數稀疏曲線圖

2.3 細菌種類多樣性和豐度分析

利用blastn將OTU序列與對應數據庫進行比對，篩選出OTU序列的最佳比對結果，并對比對結果進行過濾，默認滿足相似度>90%且coverage>90%的序列被用來后續分類，不滿足條件的序列則被歸為unclassified。根據分類學分析結果，統計樣本在各個分類層級水平上的群落組成如表2所示。以OTU的物種分類結果為依據， 分別在門、綱和屬3個分類級別上對樣本中的細菌種類和相對豐度進行統計分析。樣本的細菌種類覆蓋5個門，其中以變形桿菌門(Proteobacteria)為主， 占總數的69.23%，厚壁菌門(Firmicutes)28.26%次之，其中變形桿菌門genus水平上樹狀圖如圖5所示，厚壁菌門genus水平上樹狀圖如圖6所示。在綱的分類水平上， 包含8個綱，以變形菌綱(Betaproteobacteria，37.21%)和芽孢桿菌綱(Bacilli，28.25%)為主。在目的分類水平上，包含13個目，以伯克氏菌目(Burkholderiales，37.21%)和芽孢桿菌目(Bacillales，24.38%)為主，在科的分類水平上，包含了25個科，以叢毛單胞菌科(Comamonadaceae，36.80%)和芽孢桿菌科(Bacillaceae，24.38%)為主。在屬的分類水平上， 包含29個屬，以玫瑰色半光合菌屬(Roseateles，36.75%)和芽胞桿菌屬(Bacillus，24.38%)為主。豐度大于1%的還有寡養單胞菌屬(Stenotrophomonas，9.36%)、藤黃色桿菌屬(Luteibacter，1.45%)、乳球屬(Lactococcus，2.83%)。而劉青梅等[15]通過傳統培養方法獲得散裝年糕中疑似腐敗菌，并采用分子生物學方法對細菌的進行擴增結合數據庫比對等方法鑒定出優勢腐敗菌為芽孢桿菌和短芽孢桿菌兩大屬。胡慶松等[16]采用VITEK-32自動化微生物分析儀鑒定系統、 API微生物鑒定系統和手工法鑒定引起散裝年糕腐敗的典型菌株為波茨坦短芽胞桿菌、巨大芽孢桿菌、乙酰短桿菌、希氏短桿菌。俞科偉等[2]通過對發生腐敗變質的真空包裝年糕中微生物進行分離、純化和反證，經過形態觀察、 生理生化、微生物鑒定系統等方法,初步鑒定 2株細菌分屬于枯草芽孢桿菌、不動桿菌屬。黃麗金等[17]采用分子生物學方法對真空年糕的腐敗菌進行分離和鑒定，確定其優勢菌為枯草芽孢桿菌。陳挺[18]對真空年糕制品表面和內部的微生物進行分離純化,并通過形態特征的觀察和一些生理生化實驗對其進行鑒定得出主要優勢菌株為醋酸鈣不動桿菌、乙酰短桿菌、枯草芽孢桿菌、蠟狀芽孢桿菌、巨大芽孢桿菌和短芽孢桿菌。由此可以看出利用傳統純培養方法等和高通量測序得到的年糕腐敗菌屬具有基本一致的結果，但又不完全相同。此外，在屬的水平上，盡管不同方法的分析數據有差異，造成分析數據差異的原因可能在于年糕的不同包裝方式、不同培養時間和不同的腐敗時期菌群結構鑒定也會造成優勢菌屬及豐度的不同，但不同方法所得出豐度較高的腐敗菌屬中均含有芽胞桿菌屬。這可能是由于芽孢桿菌屬的微生物能夠形成芽孢的特性使得它們能夠抵抗各種極端環境如高溫、極酸、極鹽，而且對各種殺菌劑具有抵抗力，因此，它們在各種自然環境中都能被分離到[19]，并且引起年糕脹袋的微生物也很有可能是芽胞桿菌屬(Bacillus)[20]。所以，如何抑制芽胞桿菌屬對于真空年糕的貯運保鮮具有重要意義后續研究應重點關注。除此之外，新鮮年糕中另一種優勢菌屬玫瑰色半光合菌屬(Roseateles)是伯克氏菌目(Burkholderiales)分支的一個屬，而伯克氏菌目是一類能夠利用自身固氮的微生物，在水稻與環境的互作過程中可能發揮著重要的作用，是水稻的內生有益菌[21]。故由稻米制作出的新鮮年糕可能攜帶該菌，造成該菌在年糕微生物豐度中所占比例較大。

圖5 樣本genus水平上變形桿菌門樹狀圖組成

圖6 樣本genus水平上厚壁菌門樹狀圖組成

表2 genus水平上樣本主要rank reads數目

3 結論

本研究通過高通量測序技術分析了新鮮年糕經真空包裝后置于37 ℃下腐敗后的總菌群結構和優勢菌屬。結果表明，年糕內生細菌主要分布于玫瑰色半光合菌屬(Roseateles)和芽胞桿菌屬(Bacillus)、寡養單胞菌屬(Stenotrophomonas)，這3個屬是年糕的優勢菌屬。