浙江地區中溫大曲貯存過程細菌群落多樣性及功能預測分析

陳樂樂，王乙伊，汪怡寧，虞 敏

（1.寧波職業技術學院 化學工程學院，浙江 寧波 315800；2.寧波大學 食品與藥學學院，浙江 寧波 315800）

大曲是白酒生產過程中重要的糖化發酵劑，是在開放環境下以未經高壓滅菌的原料自然接種生產而成的，對白酒產品品質有重要影響[1-2]。生產過程中的微生物主要來自未滅菌的原料和生產環境[3]。大曲的制備需經過谷物潤濕、粉碎混合、固定成型、堆積發酵、貯存成熟等多個步驟[4]。通常，大曲的貯存過程所需時間較長，需存放在通風良好的開放環境下經歷長達數個月之久[5]。目前，大曲成熟時間仍取決于個體生產者的經驗，導致大曲產品品質參差不齊且穩定性差。因此，了解大曲在貯存期間的微生物變化規律極為重要，已成為行業關注重點。

近年來，高通量測序技術已被廣泛用于研究大曲的微生物群落組成[6]。陳曉茹等[7]通過高通量測序技術發現，魏斯氏菌屬（Weissella）、刺糖多孢菌屬（Saccharopolyspora）和乳酸桿菌屬（Lactobacillus）是大曲貯存過程中的優勢細菌，但細菌群落結構變化顯著。GUAN T等[1]分析了濃香型大曲貯存6個月期間的細菌群落變化規律，魏斯氏菌屬（Weissella）、乳酸桿菌屬（Lactobacillus）、片球菌屬（Pediococcus）、芽孢桿菌屬（Bacillus）和醋酸桿菌屬（Acetobacter）是大曲貯存期間的優勢細菌，且不同貯存時間的細菌群落結構差異較大。由此可知，貯存過程是大曲最復雜、最關鍵的環節，大量微生物在自然環境變化下不斷演替而趨于平衡，為后續釀酒的多邊發酵提供重要的微生物基礎[7]。然而，當前研究多集中于大曲貯存過程中的微生物群落結構組成，忽視了大曲貯存過程的功能演替，尤其是利用生物信息學對大曲貯存過程細菌群落的功能預測研究還鮮見報道。

當前關于大曲貯存過程細菌群落多樣性及功能的研究還相對較少，導致大曲貯存過程的微生物生態學機制研究仍存在不足，尤其是在南方沿海地區生產的大曲。大部分研究者都是以我國比較著名的白酒生產所用大曲為研究對象，但是不同地域由于氣候環境的不同，環境中微生物群落差異較大，而大曲通常是在相對干燥、通風良好，溫濕度恒定的曲房進行貯藏，不可避免的會受到當地環境的影響。目前，針對沿海地區生產的大曲還尚未報道。本研究以浙江地區生產的中溫大曲為研究對象，通過高通量測序技術，分析中溫大曲貯存過程中細菌群落結構及多樣性的演替規律，并基于多元統計分析揭示影響細菌群落結構組成變化的理化因素，同時，結合生物信息學對細菌功能進行預測，從微生物生態學視角全面深入理解大曲貯存過程的重要性，旨在為生產品質穩定的成熟中溫大曲提供理論和數據參考。

1 材料與方法

1.1 材料與試劑

樣品：取自浙江省某白酒大曲生產車間，成曲常溫貯存90 d，每隔30 d取樣。

鄰苯二甲酸氫鉀、酚酞、氫氧化鈉（均為分析純）：國藥集團化學試劑有限公司；DNeasyPowerSoil Pro Kit：德國QIAGEN公司；BMS1000-2核酸純化磁珠：無錫百邁格生物科技有限公司；Taq脫氧核糖核酸（deoxyribo nucleic acid，DNA）聚合酶：上海生工生物工程股份有限公司；細菌16S rRNA V3-V4區擴增引物（343F：5′-TACGGRAGGCAGCAG-3′和806R：5′-GGACTACVVGGGTATCTAATC-3′）：上海歐易生物醫學科技有限公司；Nextera XT DNA Sample Preparation Kit：美國Illumina公司。

1.2 儀器與設備

AB-50型電子分析天平、FE28型pH計：瑞士Mettler-Toledo公司；GL-21M型高速冷凍離心機：長沙平凡儀器廠；VORTEX-3型旋渦混勻儀：德國IKA公司；NanoDrop 2000型紫外微量分光光度計：美國Thermo公司；DYCP-31DN型瓊脂糖水平電泳儀電泳槽：北京六一儀器廠；SEDI G型聚合酶鏈式反應（polymerase chain reaction，PCR）儀、KETA G型全自動凝膠成像系統：美國Wealtec公司；Illumina PE250測序儀：美國Illumina公司。

1.3 試驗方法

1.3.1 樣品收集

為消除曲房空間差異和大曲個體差異，分別從曲房中間的曲堆上、中和下層取樣，并粉碎混勻，將粉碎混勻的大曲分為兩份，分別置于-20 ℃和-80 ℃保存，一份用于理化檢測，一份用于DNA提取，每份樣品設置三個平行樣本。

1.3.2 理化特性分析

參考QB/T 4257—2011《釀酒大曲通用分析方法》[8]測定樣品的pH、酸度和含水量。

1.3.3 DNA提取及高通量測序

1.3.4 生物信息學分析

原始序列使用Qiime2 pipeline[11]進行處理，具體如下：利用Demux插件對原始序列進行解復用，使用Cutadapt插件去除測序引物，并過濾低質量堿基，隨后采用DADA2插件合并序列并去除嵌合體，通過預訓練SILVA 138數據庫[12]對序列進行注釋，得到非單擴增序列變異體（amplicon sequence variants，ASV），將葉綠體、線粒體和不能注釋到門水平的ASV刪除。以最少樣本序列數為標準進行抽平處理，分析不同樣品的多樣性差異。利用R軟件實現基于ASVs的柱狀圖、樣本聚類分析、主坐標分析和冗余分析。

根據樣品中各ASV代表序列及其豐度表，使用Phylogenetic Investigation of Communities by Reconstruction of Unobserved States 2（PICRUSt2）算法的標準集成基因組數據庫推斷16S rRNA基因的功能組成[13]。該方法是將樣品的16S rRNA基因序列定位到京都基因與基因組百科全書（Kyoto Encyclopedia of Genes and Genomes，KEGG）數據庫進行注釋，并通過PICRUSt2算法預測群落的代謝途徑豐度。

1.3.5 數據處理及統計分析

所有實驗重復進行三次，數據結果以“平均值±標準偏差”表示。用SPSS 19.0軟件進行方差分析，P＜0.05表示差異顯著。采用GraphPad Prism 8.0和Excel 2017繪制圖表。

2 結果與分析

2.1 貯存期間大曲理化特征變化

大曲在自然環境變化脅迫下會發生復雜的理化反應，故從理化角度揭示大曲貯存過程中的各理化特征的差異。大曲貯存過程中pH、酸度和含水量的變化見圖1。

圖1 大曲貯存期間pH（A）、酸度（B）和含水量（C）的動態變化Fig.1 Dynamic changes of pH (A),acidity (B) and moisture (C)during Daqu storage

由圖1可知，在90 d貯存過程中，大曲樣品的pH和含水量隨貯存時間的不斷延長而逐漸降低，而酸度則逐漸增加。成曲的pH和含水量分別為6.69和11.97%。貯存90 d后，大曲樣品的pH和含水量分別下降至5.93和9.98%，顯著低于成曲樣品的pH和含水量（P＜0.05）。此外，大曲的酸度在貯存90 d后可達0.96 mmol/10 g，比貯存前的成曲樣品提高了95.92%，二者的酸度差異顯著（P＜0.05）。這一結果與GUAN T等[1]的研究結果相一致，其研究發現濃香型白酒大曲貯存6個月后會引起含水量和pH的降低以及酸度的升高。造成這一現象的原因可能是由于每塊大曲磚塊之間都留有一定的縫隙，而曲房是一個相對干燥、通風良好、溫濕度恒定的空間[14]。pH的降低及酸度的增加與乳酸菌、醋酸菌等細菌的相對豐度持續較高有關[1]。故貯存會造成大曲的pH、酸度和含水量發生顯著變化。

2.2 貯存期間大曲中細菌群落多樣性變化

使用Observed指數、Shannon指數和Simpson指數評估大曲貯存過程中細菌群落的α-多樣性變化，結果見表1。

表1 大曲貯存期間細菌群落序列數量和α-多樣性分析結果Table 1 Bacterial community sequence numbers and alpha-diversity analysis results during Daqu storage

由表1可知，12個大曲樣品的Illumina高通量測序共得到663 332條有效序列，平均42 188～65 079條/大曲樣品。各大曲樣品的測序覆蓋率均較高（＞98%），說明這些有效測序序列可真實反映大曲中細菌群落的物種組成和多樣性。成曲細菌群落的Observed指數、Shannon指數和Simpson指數分別為348、3.69和0.77。隨著貯存時間的不斷延長，大曲樣品細菌群落的Observed指數和Shannon指數呈先下降后逐漸增加的變化趨勢，而Simpson指數則呈逐漸增加的變化趨勢。貯存90 d后，大曲細菌群落的Observed指數、Shannon指數和Simpson指數分別為411、4.22和0.90，與成曲細菌群落的α-多樣性結果相比差異顯著（P＜0.05）。這與陳曉茹等[7]研究結果相一致。大曲樣品細菌群落α-多樣性的變化與大曲貯存環境有關，通常pH和酸度會顯著影響微生物群落的α-多樣性，表現為一些不耐酸的微生物不能適應該環境而消亡[15]。由于曲房具有良好的通風性，環境中耐酸的微生物可附著在大曲表面，利用高粱、小麥中的糖類物質進行生長繁殖，從而提高貯存后期大曲樣品的α-多樣性。

基于Bary-Curtis算法計算各大曲樣品的距離，隨后采用主坐標分析（principal coordinate analysis，PCoA）比較不同貯存時間大曲樣品中細菌群落的差異，結果見圖2。

圖2 大曲貯存期間細菌群落結構的主坐標分析結果Fig.2 Principal coordinate analysis results of bacterial community structure during Daqu storage

由圖2可知，主坐標第一軸方差貢獻率為49.83%，而主坐標第二軸方差貢獻率為31.11%，由此可知，物種差異的累計方差貢獻率為80.94%，可利用前兩個主成分來解釋該模型的變量。不同貯存時間的大曲樣品在主坐標前兩軸的距離較遠，說明各貯存時間的大曲樣品中細菌群落結構差異顯著，尤其是成曲樣品和其他貯存時間的大曲樣品差異最大。這一結果與陳曉茹等[7]研究結果相一致，表明大曲樣品中細菌群落結構會隨著貯存時間延長而不斷變化。

2.3 貯存期間大曲中細菌群落結構組成演替

不同貯存時間大曲樣品中細菌群落門水平的演替規律見圖3。由圖3可知，經物種分類注釋后，12個大曲樣品共檢測到來自11個門水平的物種，其中厚壁菌門（Firmicutes）、變形菌門（Proteobacteria）和放線菌門（Actinobacteria）是大曲樣品貯存期間相對豐度最高的三個門。厚壁菌門在大曲貯存期間占絕對優勢地位，其相對豐度范圍在86.83%～97.61%，表明貯存期間的大曲樣品具有相對穩定的門水平細菌群落結構。這一結果與唐賢華等[16]研究結果相似，說明來自厚壁菌門的微生物普遍具有耐酸、分解高粱和小麥中淀粉和蛋白等特性。

圖3 大曲貯存期間細菌門的相對豐度動態變化Fig.3 Dynamic changes of relative abundance of bacterial phyla during Daqu storage

為進一步更為準確的描述大曲樣品在不同貯藏時間的細菌群落結構，從屬水平揭示其組成差異。不同貯存時間大曲樣品中細菌群落屬水平的演替規律見圖4。由圖4可知，所有大曲樣品共檢測到79個屬，其中乳酸桿菌屬（Lactobacillus）、魏斯式菌屬（Weissella）、片球菌屬（Pediococcus）、明串珠菌屬（Leuconostoc）、假單胞菌屬（Pseudomonas）、芽孢桿菌屬（Bacillus）、醋酸桿菌屬（Acetobacter）、克羅彭斯特菌屬（Kroppenstedtia）、不動桿菌屬（Acinetobacter）、鏈霉菌屬（Streptomyces）和乳球菌屬（Lactococcus）等是大曲貯存期間的優勢屬。此外，在屬水平的聚類分析也發現同一貯存時間的大曲樣品被聚在一起，且整個貯存期間的大曲樣品可主要分為兩簇，成曲樣品簇和其他貯存時間大曲樣品簇，這一結果與圖2的結果相一致，說明成曲和其他貯存時間的大曲樣品在屬水平組成上差異較大。成曲樣品中乳酸桿菌屬和魏斯式菌屬的相對豐度最大，分別為47.56%和42.77%。隨著貯存時間的不斷延長，乳酸桿菌屬的相對豐度逐漸升高，而魏斯式菌屬的相對豐度則逐漸降低。貯存90 d后，大曲樣品中乳酸桿菌屬和魏斯式菌屬的相對豐度分別為63.44%和25.16%。乳酸桿菌屬是一類兼性厭氧的微生物，是谷物類產品發酵所必需的微生物之一，廣泛存在于大曲及白酒釀造中，其可通過自身代謝系統將谷物中淀粉轉化為乳酸[17-18]。這可能是促使大曲樣品在貯存期間酸度增加的主要因素。魏斯式菌屬也被確定為大曲中最豐富的細菌類群之一，其能夠降解纖維二糖并在酸性發酵環境中生長代謝[17]。片球菌屬和明串珠菌屬的相對豐度在貯存30 d前逐漸升高，分別從2.71%和2.22%升高至5.57%和3.49%，之后隨著貯存時間的增加逐漸降低，貯存90 d后，二者的相對豐度分別僅為1.89%和0.89%。這一結果與GUAN T等[1]研究結果相一致，表明片球菌屬和明串珠菌屬受其他微生物之間的競爭抑制，不能在貯存期間大量生長繁殖，從而達不到相對較高的相對豐度。據報道，從葡萄酒中分離出的片球菌具有酯化能力[19]。明串珠菌屬通常附著在谷物原料表面，可將原料中糖類物質轉化為甘露糖和葡聚糖，為其他乳酸菌的生長繁殖提供物質基礎[20]。芽孢桿菌屬和不動桿菌屬在貯存期間的相對豐度范圍分別為0.64%～2.07%和0～6.28%，這些少量的芽孢桿菌屬具有產生蛋白水解酶和淀粉酶等特征，可提高大曲中有機酸和吡嗪類化合物的含量[21-23]。而不動桿菌屬在大曲中的作用目前尚不清楚。此外，醋酸桿菌屬和克羅彭斯特菌屬的相對豐度則隨貯存時間的不斷延長而逐漸升高，二者分別從成曲樣品的0.11%和0.04%升高至貯存90 d后大曲樣品的2.04%和2.35%。醋酸桿菌屬是一類好氧的微生物，可將谷物原料中糖類和酒精氧化成醋酸，從而提高大曲樣品的酸度。克羅彭斯特菌屬廣泛存在各種不同香型的大曲中，可作為大曲成熟的生物標志微生物[24-25]。

圖4 大曲貯存期間細菌屬的相對豐度動態變化Fig.4 Dynamic changes of relative abundance of bacterial genera during Daqu storage

2.4 貯存期間大曲理化因子與細菌群落的冗余分析

冗余分析（redundancy analysis，RDA）是一種回歸分析結合主成分分析的排序方法，可用于揭示響應變量和解釋變量的線性關系。不同貯存時間大曲樣品理化因子與優勢細菌屬的冗余分析見圖5。

圖5 大曲貯存期間理化因子與優勢細菌屬的冗余分析Fig.5 Redundancy analysis of physicochemical factors and dominant bacterial genera during Daqu storage

由圖5可知，貯存期間大曲理化因子對優勢細群屬的總解釋度為60.44%，其中冗余分析第一軸的解釋度為45.76%，可較好的區分不同貯存時間的大曲樣品。魏斯式菌屬相對豐度與含水量和pH呈正相關，與酸度呈負相關，而酸度則與乳酸桿菌屬、醋酸桿菌屬和克羅彭斯特菌屬相對豐度呈正相關。此外，通過條件限制分析闡明每個理化因子對大曲貯存期間細菌群落結構的影響程度，發現pH、含水量和酸度對大曲細菌群落主要優勢屬演替的貢獻度相差無幾，說明多個理化因子共同調控貯存期間大曲樣品細菌群落結構的演替變化。綜上可知，貯存期間大曲樣品細菌群落組成和理化因子變化顯著，這些理化因子共同影響了大曲細菌群落多樣性及結構穩定，從而促進了大曲的成熟。

2.5 貯存期間大曲細菌群落功能預測分析

利用PICRUSt2算法預測細菌的代謝功能，是宏基因組分析的一種代替方法，具有成本低，精度較高等特點。不同貯存時間大曲樣品細菌群落的功能預測分析見圖6。

由圖6可知，貯存期間大曲樣品共涉及6類一級功能層，包括代謝（L1）、遺傳信息處理（L2）、環境信息處理（L3）、細胞過程（L4）、有機系統（L5）和人類疾病（L6），其中代謝、遺傳信息處理和環境信息處理是主要功能，貯存期間的相對豐度分別為60.11%～61.82%、20.21%～21.05%和14.22%～14.84%。此外，這6類一級功能層的相對豐度在不同貯存時間的大曲樣品中均無顯著差異（P＞0.05）。進一步對二級功能層進行分析，發現貯存期間大曲樣品共涉及27個二級功能層，其中碳水化合物代謝（24.79%～26.85%），翻譯（10.52%～11.02%）和膜轉運（9.86%～10.21%）是最主要的三個二級功能層。通過聚類分析發現27個二級功能層可被分為兩簇，其中運輸和代謝、消化系統、碳水化合物代謝、氨基酸代謝、能力代謝等11個二級功能層的相對豐度在貯存90 d時相對高于其他貯存時間。而折疊、分選和降解、衰老、脂類代謝、萜類和聚酮類化合物代謝、翻譯等9個二級功能層在貯存30 d時富集；膜轉運、信號轉導、核苷酸代謝、復制和修復及細胞運動5個二級功能層則在成曲樣品中富集。

圖6 大曲貯存期間細菌群落的PICRUSt2功能預測分析Fig.6 PICRUSt2 functional predictive analysis of bacterial communities during Daqu storage

3 結論

本研究利用16S rRNA高通量測序技術探究浙江地區大曲貯存過程中理化特征和細菌群落結構的演替變化，并結合PICRUSt2預測貯存期間細菌群落功能差異。研究發現大曲的pH和含水量隨著貯存時間的不斷延長逐漸降低，而酸度呈相反的變化趨勢。貯存90 d后，大曲樣品的細菌群落結構和多樣性差異顯著（P＜0.05）。門水平上，貯存期間的絕對優勢菌門為厚壁菌門（Firmicutes），其相對豐度＞85%。屬水平上，乳酸桿菌屬（Lactobacillus）和魏斯式菌屬（Weissella）是大曲貯存期間的絕對優勢屬，其相對豐度范圍分別為47.56%～63.44%和25.16%～42.77%。冗余分析結果表明魏斯式菌屬相對豐度與含水量和pH呈正相關，而酸度與乳酸桿菌屬、醋酸桿菌屬和克羅彭斯特菌屬相對豐度呈正相關。此外，PICRUSt2預測結果表明大曲貯存期間的基因功能主要與代謝、遺傳信息處理和環境信息處理等功能相關。本研究旨在為生產品質穩定的成熟大曲提供理論和數據參考。