郫县豆瓣后发酵过程中微生物群落与理化指标相关性研究

刘 平，王雪梅，向 琴，孙文佳，林俊帆，车振明

（西华大学食品与生物工程学院，四川 成都 610039）

传统郫县豆瓣后发酵采用天然晒露法，开放式的发酵环境造就了发酵体系中微生物群落的复杂性和多样性。郫县豆瓣风味物质的形成与微生物群落的代谢活动密切相关[1]，微生物或通过代谢合成产生风味或风味前体物质，或通过分泌酶类将大分子物质降解成小分子风味物质[2]。而发酵体系的水分、盐分、总酸、氨基酸态氮和还原糖等为微生物提供了不可或缺的生长代谢环境，在发酵过程中的作用不可小觑。目前，已有研究大多是运用各种分析技术探究郫县豆瓣中微生物群落的动态变化[3-5]，推断风味物质的形成机理[6]，探讨风味物质与理化指标的关联性[7]，或探究细菌群落与风味物质的相关性[8]，但在郫县豆瓣后发酵过程中微生物群落演替与发酵体系理化指标间关联性方面的研究未见报道。

本文采用高通量测序技术对郫县豆瓣后发酵过程中微生物的变化规律进行监测和多样性分析，并在明晰群落组成和核心微生物后，将其与动态监测的理化指标进行PLSR与Pearson相关性分析，通过探究后发酵过程中微生物群落演替及其与理化指标之间的关联，进一步加深对郫县豆瓣发酵机制的认识，以期对郫县豆瓣的质量控制与调控提供理论支撑。

1 材料与方法

1.1 材料与试剂

郫县豆瓣样品取自四川省郫县豆瓣股份有限公司，取样时间为：自然发酵0天和1、2、3、4、5、6、9、12个月，并依次编号T-0D、T-1M、T-2M、T-3M、T-4M、T-5M、T-6M、T-9M和T-12M。甲醇为色谱纯，其他试剂均为分析纯，购于成都市科龙化工试剂厂。

1.2 仪器与设备

PHS-320显数式pH计，成都世纪方舟科技有限公司；K9860全自动凯氏定氮仪，济南海能仪器股份有限公司；WFJ7200分光光度计，上海尤尼柯仪器有限公司；SB-5 200 DTN超声清洗机，宁波新芝生物科技股份有限公司；VG3S025涡旋混匀器，德国IKA公司；TD-5 M台式低速离心机，四川蜀科仪器有限公司。

1.3 实验方法

1.3.1 郫县豆瓣理化指标的测定

参考已有报道[7]，总酸含量采用酸碱滴定法测定，氨基酸态氮含量采用甲醛滴定法测定，还原糖含量采用直接滴定法测定，水分含量采用直接干燥法测定，色价采用丙酮超声法结合紫外分光光度计测定，pH值用pH计直接测定，盐分含量采用硝酸银滴定法[9]测定，可溶性氮含量采用凯氏定氮法[10]测定。

1.3.2 郫县豆瓣微生物群落的测定

郫县豆瓣后发酵过程中微生物群落采用高通量测序技术测定，测序及数据处理方法见文献[11]。

1.4 数据处理

采用SPSS 23.0软件，对数据进行单因素方差分析（ANOVA）、邓肯多重比较检验（Duncan）和皮尔森相关性分析（Pearson）；采用The Unscrambler软件对实验数据进行PLSR相关性分析；采用Excel 2016对实验数据进行处理和图形绘制。

2 结果与分析

2.1 郫县豆瓣后发酵过程中的理化指标变化规律

郫县豆瓣后发酵过程中的各项理化指标检测结果如图1所示。由图1（a）可知，整个后发酵阶段，豆瓣含水量呈整体轻微下降趋势，但均维持在50%水平以上，这与日晒夜露和复水勾兑的传统发酵工艺密切相关。发酵期间，豆瓣中的水不断挥发散失，需人工补水以维持湿润状态，从而帮助微生物生长繁殖和形成豆瓣独特风味[12]。豆瓣中的盐主要在甜瓣子发酵和椒醅盐渍阶段加入，大量的食盐抑制了腐败微生物的生长[13]，高盐环境下，部分耐受力和适应能力较差的菌群逐渐消亡，其中既包括有益微生物，又包括腐败微生物。在发酵体系中，盐分无明显变化，均控制在产品标准15～22 g/100 g之间[14]。色价呈持续下降趋势，相比于发酵初期，发酵12 月的豆瓣其色价降幅高达42%，且发酵前2 月是主要下降阶段，此后下降速率放缓。郫县豆瓣的色泽主要由辣椒中的红色素、发酵过程中的美拉德反应和酶促褐变综合作用形成[15]，而辣椒红色素受光照和氧气的影响，极易发生分解产生褪色，因此色价在前2 月迅速下降可能与辣椒红色素的分解有关。

图1 郫县豆瓣后发酵过程中的理化指标变化规律

从图1（b）可以看出，豆瓣的pH值呈缓慢下降趋势，而总酸的趋势则刚好相反。在发酵的前6月，总酸明显增加，可能是因为酵母菌和乳酸菌等微生物生长代谢产生了有机酸[16]，导致总酸含量增加。而6—12 月期间，总酸上下波动，无明显增加，可能是因为微生物代谢生成的有机酸，与自我挥发的有机酸或和醇类物质发生酯化反应消耗的有机酸达到了平衡状态所致[15]。

由图1（c）可知，在整个后发酵阶段，氨基酸态氮的含量维持在0.26～0.32 g/100 g之间，基本保持不变，可能是因为积累氨基酸态氮的主要阶段是前期甜瓣子发酵阶段。豆瓣中可溶性氮的含量维持在0.598～0.645 g/100 g，整体变化幅度较小。还原糖含量在后发酵4月时达到最高值。发酵中米曲霉代谢所产生的淀粉酶不断地将淀粉分解成还原糖，伴随着发酵的进行，一部分还原糖作为碳源物质为体系中的微生物所利用，另一部分还原糖与氨基酸类物质发生美拉德反应，从而产生消耗。在发酵前期（0—4 月），还原糖的生成速率高于消耗速率，而越往后发酵（4—12 月），还原糖的形成阻力不断加大，含量急剧下降。原因可能在于：乙醇、有机酸等次级代谢产物的积累，使得产淀粉酶的微生物受到抑制[16]，导致还原糖生成量减少；而愈发庞大的微生物菌系在维持生长时，消耗的还原糖也越来越多；此外，还原糖与其他物质之间相互反应生成风味物质，也会产生消耗[17]。

2.2 郫县豆瓣后发酵过程中的微生物群落分析

2.2.1 郫县豆瓣样品测序数据质量评估

根据FLASH[18]双端序列拼接所得原始序列发现，郫县豆瓣后发酵过程中细菌、真菌序列的平均长度分别为414、394 bp，且各时期样本中平均原始序列数均大于30 000条，测序深度足够。使用Trimmomatic软件[19]对原始序列进行质控[3]及去除嵌合体后，通过97%相似水平的OTU聚类分析，共获得6 280个细菌，2 799个真菌，各样本中OTU数分别为233、104个。

2.2.2 郫县豆瓣后发酵过程中的微生物群落α-多样性分析

在微生物群落生态学的研究中，α-多样性用来表征样本的物种丰富程度[20]，常见的有Chao1、Simpson和Shannon指数，它们可反映微生物群落的丰度、均匀性及多样性等。郫县豆瓣后发酵过程中的微生物群落α-多样性如图2所示。为了比较组间差异，采用检duncan多重比较，显著水平为0.05。

从图2（a）可以看出，在郫县豆瓣后发酵阶段，发酵体系中细菌群落的相对丰度和多样性呈现出动态变化趋势。Chao1指数随发酵时间的延长呈现不规律波动，在发酵至6 月时其值最高。Simpson指数和Shannon指数发酵前3 月和发酵4—6 月两阶段基本保持不变，且后一阶段数值略有降低，发酵9 月时分别出现最低值0.61和2.16。由此可知，细菌多样性除在发酵9 月时出现低谷外，整个发酵初期基本维持稳定，发酵9 月之后开始回升。

由图2（b）可知，与细菌群落相比，真菌群落的多样性指数均整体呈升高趋势，Chao1、Simpson和Shannon在发酵开始阶段均显著低于细菌，至发酵12 月时，上升幅度均超过2倍。由此可知，郫县豆瓣后发酵阶段，体系中真菌群落的丰富度和多样性具有明显的规律性，发酵前2 月略微降低，此后匀速升高。

细菌和真菌群落出现不同变化趋势，其原因可能在于微生物对环境的适应能力不同。众所周知，当发酵成熟的甜瓣子和椒醅按照一定比例混合均匀后，才开始郫县豆瓣后发酵[16]。甜瓣子发酵时，需加入食盐水，此时体系盐浓度约为16%，此阶段米曲霉和环境中落入的微生物不断适应高盐发酵体系，逐渐形成稳定的耐盐性高的真菌群落。对椒醅进行盐渍时，体系盐浓度约为18%，此时也完全依赖高盐体系对微生物进行抑制，辣椒携带的微生物与环境中落入的微生物相互竞争，最终形成以乳酸菌为主的生态群落[1]。混合发酵后，椒醅中的稳定性菌群在混合发酵初期并未出现明显下降，同时，体系中酵母菌代谢产生了大量乙醇及其衍生物[12]，对体系中的细菌表现出一定抑制，使得细菌趋势始终处于波动状态。对于真菌而言，豆瓣发酵体系的盐分高于甜瓣子体系且低于椒醅体系，从甜瓣子发酵到豆瓣发酵，高盐环境使得部分真菌的生长暂时受到抑制，因此，出现发酵前2 月丰度及多样性下降的现象。此时，下降速率逐渐增大，其原因可能在于：发酵初期真菌细胞壁的几丁质成分对真菌具有一定保护作用，使得真菌群落下降缓慢，随着发酵的进行，某些微生物代谢产生的酸类物质逐渐累积，发酵体系酸性增强，部分适应性强的真菌继续生长，而部分耐酸性较弱的真菌出现消亡，当消亡速率高于生长速率时，整个体系中真菌群落多样性出现降低[2]。经过较长时间的适应，耐盐耐酸真菌迅速繁殖，使得酱醅中真菌多样性持续升高。

图2 后发酵阶段郫县豆瓣样品α-多样性指数

2.2.3 郫县豆瓣后发酵过程中的微生物群落β-多样性分析

为获取郫县豆瓣后发酵过程中的微生物群落结构差异，采用β-多样性进行表征[20-21]，基于物种丰度数据（OTU表）计算样本间的UniFrac距离[22]。考虑到物种丰度，本文采用加权法计算各样本间距离，并以Bray-Curtis距离为基础，使用主坐标分析（PCoA）[21]，观察样本间的差异程度以及差异变化规律，样品的群落组成越相似，则它们在图中的距离越接近。

由图3可知，后发酵时间与发酵体系中微生物群落结构的分布具有较强的相关性，细菌群落变化存在3个明显的区分阶段，即后发酵初期（发酵0—3 月）、后发酵中期（发酵4—6 月）、后发酵后期（发酵9—12 月）。在这3个阶段内部，细菌的丰度及多样性十分接近，群落组成相似度较高，这一发现与前人的研究报告较为一致[23]。同样，真菌群落变化也被划分为3个阶段，即0—1、2—6 和9—12月。造成上述变化趋势的原因可能在于：后发酵初期，甜瓣子与椒醅混合后，甜瓣子作为能量供体，被椒醅中的细菌充分利用，急速增殖；发酵2 月后，混合体系中酵母菌等真菌大量生长，代谢产物中的乙醇具有抑菌作用，使得酱醅中的细菌多样性略有降低；经过为期3 月发酵，体系中乙醇由于消耗逐渐趋于稳定，同时某些适应能力强且对乙醇耐受力强的细菌，甚至以乙醇为碳源，不断增殖[2]。此外，豆瓣发酵体系的盐分略高于甜瓣子，甜瓣子带入的米曲霉等真菌在高盐体系中存在适应过程，生长受限，造成真菌群落组成发生明显变化。适应力强的高耐盐耐酸真菌在为期2—6 月发酵过程中持续累积，加上发酵环境中不断掉落的真菌，使得真菌群落多样性逐渐回升，最后趋于稳定。经过长达9 月的后熟发酵，酱醅中的能量物质逐渐消耗殆尽，乙醇浓度再度升高，细菌及真菌的生长再次受到抑制。

图3 基于加权Unifrac的PCoA分析

2.3 郫县豆瓣后发酵过程中的微生物群落与理化指标相关性分析

鉴于风味是影响人们选择郫县豆瓣的主要因素之一，故在确定豆瓣的核心微生物菌群时，需考虑其对风味物质的影响效果。因此，以Lefse分析具有显著差异的物种（LDA值≥5）、或在发酵体系中始终处于优势地位（即平均相对丰度＞0.1%）、具有较大变异系数（CV＞50%）、与至少2类风味组分强烈相关（|ρ|＞0.8且P＜0.05）等为判定条件，对豆瓣的细菌群落与真菌群落进行鉴定，得出满足条件的核心微生物有Kosakonia、Kazachstania、Debaryomyces、Lactobacillus、Myroides、Stenotrophomonas、Ochrobactrum、Wohlfahrtiimonas和Lactococcus属，以上9种微生物菌群与风味物质高度相关。为探究上述核心微生物菌群与理化指标的关系，进一步加深对郫县豆瓣发酵机制的认识，将其丰度值（见表1）与理化指标进行相关性分析。

表1 核心微生物菌群

经偏最小二乘回归法（PLSR）进行分析，前二维主成分相关载荷图如图4所示。以微生物相对丰度值作为自变量X，可解释的X变量方差贡献率为PC1=54%、PC2=28%。以风味指标作为因变量Y，可解释的Y变量方差贡献率为PC1=45%、PC2=12%。对X、Y变量解释良好，表明微生物菌群与理化指标具有一定的相关性。

图4 核心微生物菌群与理化指标的相关性载荷图

可以观察到，在坐标轴的最右侧Kazachstania、Kosakonia、Lactobacillus、Wohlfahrtiimonas属与色价相关性较好。发酵体系中氨基酸与还原糖的含量影响着美拉德反应的进行，进而影响豆瓣色泽。结合研究数据发现[11]，Kazachstania属的真菌和Kosakonia属的细菌均与天冬酰胺、谷氨酰胺、脯氨酸等氨基酸呈正相关，Wohlfahrtiimonas属的细菌与天冬酰胺、谷氨酰胺呈正相关。这说明该类微生物生长时，并未大量消耗氨基酸，故对美拉德反应负面影响较小，使得体系色价较高，也有可能是该类微生物能够分泌蛋白酶，使得体系中的蛋白质不断分解成氨基酸，为美拉德反应提供原料。而细菌属中的Lactobacillus属是整个后发酵阶段的优势菌[8]，也是多种发酵制品的主要菌群[24-26]，由此可知，Lactobacillus属在生长时也未大量消耗氨基酸。

Wohlfahrtiimonas、Myroides、Debaryomyces、Kosakonia、Lactobacillus、Lactococcus属与水分呈正相关性，说明水分充足的发酵环境更适合此类微生物的生长，比如后发酵初期豆瓣发酵体系含水量较大，氧含量较低，Lactobacillus属的厌氧或兼性厌氧菌快速生长繁殖。从表1也可看出，Lactobacillus属在发酵前3 月的相对丰度的确较高。该6种微生物菌群与总酸呈负相关性，说明环境中酸类物质的增加会抑制它们的生长。结合之前的研究数据发现，Wohlfahrtiimonas属的细菌与酒石酸呈负相关性，Myroides属的细菌和Debaryomyces属的真菌均与苹果酸呈强烈的负相关性，Kosakonia属的细菌与乳酸和琥珀酸存在强烈的负相关性，Lactococcus属的乳球菌和Lactobacillus属的乳杆菌都属于产酸的乳酸菌，但产酸过多会影响它们的生长[27]，故同样表现为负相关性。

由于豆瓣发酵体系中的盐分无明显变化，因此微生物的演替与盐分含量相关性不大。pH虽未在50%～100%方差解释椭圆之间，但与总酸含量仍呈负相关关系，总酸含量增加，pH减小。氨基酸态氮、可溶性氮、还原糖含量也未在50%～100%方差解释椭圆之间，与核心微生物的相关性不明显。其原因可能在于：一方面，蛋白酶与淀粉酶不断水解蛋白质与淀粉生成氨基酸态氮、可溶性氮和还原糖；另一方面，氨基酸态氮、可溶性氮和还原糖作为碳源和氮源不断被微生物消耗。故在整个后发酵过程中，该3项指标看似与微生物演替毫无关系，实际上却是因为受限于分析方法，不能对它们之间的关系做更明确的讨论。

因此，对核心微生物菌群与理化指标做进一步分析，由Pearson相关系数分析得到定量系数表[8]，如表2所示。可以看出，Lactobacillus、Kosakonia、Myroides属与水分显著相关，Wohlfahrtiimonas、Lactococcus属与水分呈极显著相关性；Lactobacillus、Wohlfahrtiimonas属与色价显著相关，Kazachstania、Kosakonia属与色价呈极显著相关性；Lactococcus属与总酸显著相关，Myroides、Wohlfahrtiimonas属与总酸呈极显著相关性；9种核心微生物菌群均与pH呈正相关，与还原糖、氨基酸态氮、可溶性氮含量呈负相关，该结论与PLSR分析基本吻合。

表2 核心微生物菌群与理化指标的Pearson系数表

综上所述，Lactobacillus、Kosakonia、Myroides、Wohlfahrtiimonas、Lactococcus属与水分含量呈显著正相关，Kazachstania、Kosakonia、Lactobacillus、Wohlfahrtiimonas属与色价呈显著正相关，Lactococcus、Myroides、Wohlfahrtiimonas属与总酸含量呈显著负相关，9种核心微生物菌群均与pH呈正相关，与还原糖、氨基酸态氮及可溶性氮含量呈负相关，与盐分含量相关性不大。

3 结论

对郫县豆瓣后发酵过程中的理化指标和微生物群落进行动态监测，发现体系的水分、色价、pH值呈下降趋势，总酸、可溶性氮则相反；盐分、氨基酸态氮含量基本不变，还原糖则先增后减；细菌的相对丰度先升高再降低再升高最后趋于平稳，其多样性除在9 月出现低谷外，整个发酵阶段基本维持不变；真菌丰富度和多样性则持续升高，升幅均超过2倍；后发酵时间与微生物群落结构的分布具有较强的相关性。

将核心微生物菌群与理化指标进行PLSR分析与Pearson相关系数分析，发现核心微生物菌群与水分、色价、pH值呈正相关，与总酸、还原糖、氨基酸态氮、可溶性氮呈负相关，与盐分相关性不大。微生物群落演替与理化指标相关性的探讨，可作为郫县豆瓣后发酵过程中质量调控的参考。