多菌种参与的固态发酵被认为是一种古老、经济且对人体有益的食品发酵技术^[1]。由于多种微生物共同作用并具有复杂的代谢活力，这种类型的发酵可以产生独特的食品风味和口感特征^[2]。通过寻找核心微生物群有利于人工构建合成微生物群落(synthetic microbial communities)并获得易操作且可重复的发酵系统^[3-4]。目前从系统水平上寻找核心微生物并构建模拟的微生物群落来研究群落的性能和稳定性是常用的方法^[5]。有研究者发现，不同国家的奶酪中具有稳定且丰度较高的微生物属，如Staphylococcus、Brevibacterium、Brachybacterium、Penicillium和Candida等。通过筛选相应的物种并人工构建奶酪皮，证明了奶酪发酵过程中微生物群落聚集的动态变化具有一定的可操控性和可重复性^[6]。

寻找复杂发酵系统中的核心微生物群是构建模拟微生物群落的关键问题之一。目前寻找核心微生物群的方法可以归纳为以下几种：(1)功能相关的微生物(functional microbiota)^[7-8]，如通过微生物群落和代谢物的相关性分析在食醋发酵中找到了由7个属组成的功能微生物群^[8]；(2)丰度界定的优势微生物(dominant microbiota)^{[6, 9]}，如Hu把平均相对丰度＞1%且至少出现在一个样品以上的属定义为窖泥中的核心微生物属，共鉴定出17个核心属^[9]；(3)在样本中分布频率较大的普遍微生物群(ubiquitous microbiota)^[10-11]，如Rui等将至少分布在90%的家用沼气池中的OTU定义为核心OTU^[10]；(4)共现性网络分析出共现的微生物群(co-occurring microbiota)^[12-13]，如Barberan等根据微生物的同现率(incidence of co-occurrence)找出了土壤群落中的2个关键门为Verrucomicrobia和Acidobacteria^[14]。

本研究以中国传统的白酒发酵为研究载体，试图尝试多种方法以更全面地揭示该发酵系统的核心微生物群(core microbiota)并探究环境因素对其的影响。通过高通量扩增子测序深入挖掘发酵过程中的微生物群落结构，GC-MS测定样品中的挥发性化合物含量，并使用多元统计学分析解析发酵过程中微生物群落和代谢轮廓的变化规律，找出了白酒发酵过程中的核心微生物群，并量化了环境因素对核心微生物群的影响，为白酒生产的定向调控提供一定的理论基础。

1 材料和方法 1.1 样品的采集

样品从江苏北部某酒厂采集，分别从2个窖池发酵的0、5、8、10、15、20、30、40 d取样，共计16份。取样时间为2015年4月。样品采集，一式2份，立即密封。一份–20 ℃保存用来测定理化指标和物质含量；一份–80 ℃保存用来提取DNA。

1.2 理化指标、乙醇、乳酸和挥发性化合物的测定

使用温度计插入窖池1 m深处测定酒醅温度。酒醅中的水分采用恒温干燥法测定。酸度采用酸碱中和滴定法。还原糖采用3, 5-二硝基水杨酸(DNS)法测定^[15]。乙醇采用液相色谱仪(Agilent 1200)测定，色谱柱为AminexHPX-87H (300 mm×7.8 mm，Bio-Rad)，示差折光检测器(RI)为Schambeck SFD GmbH，测定条件参考Wu等^[16]的方法。乳酸采用反相高效液相色谱仪(Waters)测定，色谱柱为Atlantis T3 C₁₈柱(250 mm×4.6 mm)，流动相为0.5%的NaH₂PO₄ (pH 2.7)，测定条件：进样量10 μL，流速0.25 mL/min，检测器波长210 nm，柱温40 ℃，保留时间14 min。挥发性化合物的测定采用气相色谱-质谱联仪(GC-MS；Agilent 6890N气相色谱仪和Agilent 5975质谱检测器)，色谱柱为DB-Wax (30 m×0.25 mm×0.25 μm；J & W Scientific)，内标为薄荷醇(106.25 mg/L)，检测条件参考Kong等^[17]的报道。

1.3 DNA的提取、PCR扩增和Miseq测序

称取7 g样品，根据E.Z.N.A.^® soil DNA Kit (Omega Bio-Tek, Norcross，GA，USA)的操作说明提取各样品中的DNA。对于细菌，针对16S的V3-V4区进行扩增，引物为338F (5′-ACTCCTACG GGAGGCAGCAG-3′)和806R (5′-GACTACHVGG GTWTCTAAT-3′)。对于真菌，扩增区域为ITS2区，引物为ITS3 (5′-GCATCGATGAAGAACGCAGC-3′)和ITS4 (5′-TCCTCCGCTTATTGATATGC-3′)。这些引物针对不同的样品都包含有一个8核苷酸barcode的序列。25 μL的反应体系包含2.5 μL的10×Pyrobest buffer，2 μL的dNTPs (2.5 mmol/L)，正反向引物各1 μL (10 μmol/L)，0.4 μL的Pyrobest DNA polymerase (Takara)和15 ng的模板DNA。细菌的PCR扩增程序：95 ℃预热5 min；95 ℃变性30 s，55 ℃退火30 s，72 ℃延伸40 s，共25个循环，最后72 ℃保持10 min；真菌的PCR扩增程序：95 ℃预热2 min；95 ℃变性30 s，61 ℃退火30 s，72 ℃延伸45 s，共32个循环，最后72 ℃保持10 min。PCR产物根据之前的报道进行纯化^[18]。纯化后的产物经NanodropND-1000 UV-Vis定量，然后将PCR产物等分子量混合，并根据Low Sample Protocol进行文库的制备，最后在Illumina Miseq (PE300)平台上机测序。测序后所有的原始序列保存在DDBJ数据库，接受号为DRA005468和DRA005469。

1.4 序列处理

原始序列使用QⅡME pipeline进行处理^[19]：(1)序列归类，去除标签、引物和接头序列；(2)去除长度小于110 bp，模糊碱基大于2和序列平均质量小于30分的序列；(3) UCHIME去除嵌合体。最后使用UCLUST把高质量的序列根据97%的序列相似度聚成不同的操作性分类单元(OTU)，并计算相应的Chao1和Shannon指数。

1.5 数据分析

1.5.1 冗余分析(RDA): 为寻找出环境因子和微生物群落或风味轮廓之间的关系，在R语言中使用Vegan包分析环境因子对群落结构的影响，并进行蒙特卡洛置换检验。

1.5.2 Heatmap绘制和聚类分析: 为找出风味化合物随发酵时间的动态变化规律和样本之间的相似度，在R语言中使用ComplexHeatmap包进行数据标准化、样品和变量的聚类以及热图的绘制。

1.5.3 网络绘制: 计算微生物与代谢物之间的Spearman相关系数(ρ)，找出p>0.5且P < 0.05的连接作为可视化对象。对于微生物之间的共现性(co-occurrence)分析，首先计算微生物之间的Spearman相关系数，然后选取|ρ|>0.8且P < 0.05的作为互作对象。使用Cytoscape 3.4.0 (http://www.cytoscape.org)对微生物和代谢物以及微生物与微生物之间的相互关系进行可视化，以表征微生物对代谢物的贡献和微生物与微生物之间的共现性关系。

2 结果和分析 2.1 发酵过程中微生物群落结构组成及其动态变化

通过Illumina Miseq测序来确定发酵过程中的微生物群落组成。对于细菌，质控后共得到450273条高质量的序列，平均每个样本具有28142±5001条；并得到3824个OTUs，每个样品共得到239±95个OTUs。对于真菌，质控后共得到450273条高质量的序列，平均每个样本具有89898±39131条；并得到1799个OTUs，每个样品共得到112±15个OTUs。每个样品的覆盖率均在99%以上，表明足够的测序深度和可靠的数据质量^[20]。具体的测序信息如Chao1和Shannon指数等见表 1。

通过统计发现，在真菌群落中，0–8 d主要由3个真菌目主导，包括Saccharomycetales、Eurotiales和Mucorales；8 d之后，则由Saccharomycetales和Trichosporonales替代(图 1-A)。对于细菌，在最初的10 d内，Lactobacillales和Bacillales占优；之后Lactobacillales迅速增长并占据绝对优势。总体来讲，Saccharomycetales (51.12%–95.21%)和Lactobacillales (39.39%–94.65%)在整个发酵过程中具有很高的优势性。在属水平共得到69个真菌属和215个细菌属。在这些属中，Pichia (11.27%–49.51%)、Dipodascus (4.03%–41.38%)、Candida (1.48%–25.88%)、Saccharomyces (1.19%–13.48%)和Aspergillus (1.19%–13.84%)在整个发酵过程中具有一定的优势。Bacillus、Enterococcus和Lactococcus在0–10 d占优，10 d后Lactobacillus迅速从6.03%增长到92.80% (图 1-B)。

2.2 发酵过程中代谢物的变化

在发酵过程中，共鉴定出61种代谢物，包括5种醇类、10种酸类、25种酯类、17种芳香族化合物和4种杂环类化合物(图 2)。基于这些化合物相对含量的聚类分析表明，样本被分成了2个发酵阶段(阶段Ⅰ：0–20 d；阶段Ⅱ：20–40 d)。代谢物则可以分成3组：组Ⅰ包含21种化合物(9种酯类，6种酸类，5种芳香类和1种杂环类)；组Ⅱ包含37种物质(4种酸类，11种芳香类，16种酯类，5种醇类，1种杂环类)；组Ⅲ则包含有3种物质，分别是2, 4-二叔丁基苯酚、糠醛和四甲基吡嗪。组Ⅰ的化合物在前期(10–15 d)和后期(20–40 d)均有产生；组Ⅱ中的化合物几乎在后期才产生；组Ⅲ在发酵初始(0 d和5 d)和发酵后期(20 d之后)的含量相对较高。这表明20 d之后的发酵期是白酒发酵过程中的风味物质高产期。

2.3 揭示酿造核心微生物群

我们选取丰度在前20的细菌属和前20的真菌属，计算这40个属与61种代谢物之间的Spearman相关系数(p)，选取p>0.5且P < 0.05作为有效的网络连接并绘图如图 3。可以发现，Lactobacillus和Saccharomyces是连接数最大的2个属(36和22)，相关的化合物主要是酯类、酸类、醇类和芳香族化合物，尤其是这2个属都贡献重要物质乙酸(VC10)和乙醇(E)，且Lactobacillus与乳酸(VC61)也有显著的正相关。Candida与辛酸乙酯(VC44)，乙酸异戊酯(VC3)、2-甲基-1-丁醇(VC1)、苯乙醇(VC54)、异戊醇(VC2)和乙醇(E)都有一定的正相关性。还原糖(S)和糠醛(VC24)则与多个属具有正向的联系。其次是Dekkera、Pseudeurotium、Trichosporon与一些酯类、芳香类和杂环类化合物具有显著的相关性。以上结果表明在白酒酒醅发酵过程中，比较重要的代谢产物贡献者是Lactobacillus、Saccharomyces和Candida。具体的网络图信息见表 1。

为探究微生物之间的相互作用关系，通过共现性网络分析(co-occurrence network analysis)，共得到38个有效的连接节点(nodes)和267个边(edges)。网络的聚类系数为0.622，较强的连接节点(每个节点有超过12条边的连接)大多分布在Lactobacillales、Enterobacteriales和Bacillales，包括Lactobacillus、Bacillus、Enterococcus、Lactococcus、Pantoea、Rhizopus、Carnobacterium、Acinetobacter、Leuconostoc、Oceanobacillus、Paenibacillus、Saccharomycopsis、Staphylococcus、Streptococcus、Alkaliphilus、Thermoascus、Clostridium、Wickerhamomyces、Penicillium、Corynebacterium和Trichomonascus，各菌属的共现率(incidence of co-occurrence)在2%以上(图 4)，尤其是Lactobacillus与其他15个属具有显著的负相关性。此外，Lactobacillus、Rhizopus、Bacillus、Lactococcus、Saccharomycopsis、Alkaliphilus、Thermoascus、Enterococcus和Pantoea的平均相对丰度在1%以上，且共现率均在3.5%以上。Pichia和Dipodascus的平均相对丰度尤其较高，分别为28.91%和25.95% (表 2)。Alkaliphilus、Enterococcus和Pantoea被认为是来自于动物肠道、植物根部或土壤且对食品发酵具有污染性的菌属^[21-23]。

根据风味贡献程度，相互作用关系强度和相对丰度(图 3，图 4，表 2)，得到发酵过程中的核心微生物群共包含10个属：Lactobacillus、Saccharomyces、Candida、Rhizopus、Saccharomycopsis、Pichia、Dipodascus、Bacillus、Thermoascus和Lactococcus。为验证所鉴定的核心微生物群是否推动了风味的演变，通过提取二者的第一主成分并进行相关性分析后发现，核心微生物群与风味轮廓的演变具有极显著的相关性(图 5，r²=0.627，P < 0.001)，说明这10个属所组成的核心微生物群是白酒发酵过程中风味演替的主要推动者。

2.4 环境因素对核心微生物群的影响

为揭示由10个属组成的核心微生物群和环境因子之间的关系，使用冗余分析进行定量化分析并进行蒙特卡洛置换检验。首先将环境因子分为两类：物理因素包括水分、酸度和温度，化学因素包括还原糖、乙醇、乳酸和乙酸(表 2)。冗余分析表明7个因素能够解释绝大部分的核心微生物群的动态变化(74.38%)，说明环境因素的确是驱动群落演变的主要因素(图 6)。蒙特卡洛置换检验的结果表明，对核心微生物群的演替具有重要影响的因素共有5个：水分、酸度、温度、还原糖和乙醇，因为这5个因素与该微生物群的分布具有显著的相关性(r²>0.5，P < 0.01) (表 3)。化学因素对群落变化的解释性(21.46%)要高于物理因素(3.37%) (图 6)，表明核心微生物群的演变主要是由微生物产生的代谢物引起的化学因素所驱动的。此外，物理因素和化学因素的互作对其演变也具有非常重要的影响(49.55%)。如表 4所示，还原糖几乎与除Lactobacillus和Saccharomyces外的所有菌属具有显著的正相关；其他6个环境因素则与Lactobacillus、Candida和Saccharomyces具有一定的正相关性，而与Bacillus、Lactococcus、Rhizopus、Thermoascus和Saccharomycopsis具有一定的负相关性。

3 讨论

近年来，有关发酵食品的研究越来越多地集中在微生物群落结构的挖掘上。据报道，许多分子手段尤其是二代测序技术已经广泛应用到食品发酵的研究当中，为我们了解食品发酵的本质提供了更广阔的视角^[24]。由于大曲和新粮的加入以及酒醅暴露在空气中，在白酒的固态发酵前期，在酒醅入窖池时，会有多种微生物混入其中，从而导致前期的微生物群落结构比较离散。随着发酵的进行，由于温度、酸度和乙醇的升高以及可利用还原糖的下降，窖池内的生态环境逐渐变成一个相对极端的环境，导致后期发酵微生物群落尤其是细菌群落的相对集中。

不同的发酵食品由不同的核心微生物组成，如Acetobactor和Lactobacillus对于醋的发酵^[25]，Lactobacillus和Streptococcus对于奶酪的发酵^[26]。核心微生物群对于保证发酵食品的风味多样性和产品质量的批次稳定性具有重要的作用^[27]。本研究从白酒发酵过程中共揭示出10个核心微生物属。Saccharomyces是白酒酿造过程中的优势酵母属之一^[16]，且能够产生多种风味酯并贡献乙醇和其他有机酸^[28]。Lactobacillus是该过程的优势细菌菌属，如Lactobacillus rhamnosu、Lactobacillus fermentum、Lactobacillus acetotolerans、Lactobacillus kefiri等乳酸菌在不同的发酵食品中能够贡献乙醇、乙酸、乳酸和其他的重要风味物质^[29-31]，由于这些物质(乙醇、乙酸和乳酸等)能够较快地改变发酵过程中的生态环境，使得其他不耐受微生物的生长受到抑制，这也是Lacobacillus与很多其他菌属具有负相关，并且使得化学因素成为驱动菌群变化的重要原因。有报道指出，Lactococcus具有一定的氨肽酶活性和蛋白质水解力，并且有一定的风味贡献功能^[32]。Bacillus、Candida、Pichia、Saccharomycopsis以及Rhizopus是大曲成熟过程中的优势菌属之一^[33]，前三者在白酒酿造中具有潜在的风味贡献功能^[34-35]，如 Bacillus能够产生多种酚类和含氮类杂环化合物。Saccharomycopsis和Rhizopus则具有较高的淀粉酶和糖苷酶活性^[36]，Thermoascus能够产生热稳定的水解酶如纤维素酶和木聚糖酶^[37]，这些菌属在白酒酿造过程中能够降解谷物中的淀粉或纤维素，是白酒发酵不可或缺的重要菌属。而目前关于Dipodascus在发酵食品中的报道较少，但在本研究中由于所检测到的比例较高，所以也将其列为核心微生物属之一。以上10个属所组成的核心微生物群与风味的演变具有显著的相关性，说明这10个属组成的核心微生物群是白酒发酵过程中的关键微生物群。

本研究从微生物组和代谢组的角度，结合多元统计学分析更加全面地揭示了白酒发酵过程中的核心微生物群，并量化了各环境因素对核心微生物群的影响，为下一步白酒酿造的定向调控提供参考。