基于高通量测序技术的不同产地腐乳微生物多样性分析

付瑞敏，刘春雷，徐 良，张 红，夏铁骑，陈五岭

（1.河南财政金融学院体育与健康管理学院，河南郑州 450046；2.陕西师范大学食品学院，陕西西安 710069）

腐乳，又称为“东方奶酪”，是我国一种传统的发酵豆制品，它起源于我国大约两千年前的魏朝，因其独特的营养和功能成分、芳香和风味特性而被国人广泛用作调味品和开胃菜[1]。经过微生物发酵后的腐乳以其低胆固醇、高蛋白和高钙等营养特点而被人们认为是一种健康食品。根据参与发酵的微生物类型不同，可将腐乳分为细菌发酵、霉菌发酵和自然发酵等三种类型。其中，以放线菌、毛霉和根霉为发酵剂的霉菌发酵法最为普遍[2]。

腐乳中的微生物可代谢合成产生多种风味物质，其群落构成在调节腐乳产品风味方面发挥着重要作用，被认为是决定腐乳最终质量的关键因素[3-5]。腐乳发酵前期，葡萄球菌属、明串珠菌属和毛霉属占据了优势，其酶系分泌活跃，使得产品中谷氨酸、缬氨酸、苯丙氨酸和亮氨酸等含量有较大幅度增长，苯丙氨酸和亮氨酸可进一步形成醇类化合物，这些醇类化合物可和霉菌分泌的脂肪酶降解豆腐中粗脂肪所产生的脂肪酸以及细菌代谢产生的乙酸、乳酸等发生酯化反应，进而产生乳酸乙酯、乙酸乙酯等风味物质[6]。芽孢杆菌属和假单胞菌属也可分解蛋白和脂肪，在腐乳后发酵阶段产生醇类、酯类和吡嗪类等风味物质；这些微生物通过合成代谢产生各种风味物质。不同的微生物组成受其产地的环境群落及气候特点等多种因素的影响，这导致不同产地的腐乳中微生物群落存在复杂多样的差异[7]。因此，阐明腐乳发酵过程中微生物的组成及其与食品性质的关系，并在此基础上对腐乳发酵进行合理的菌群设计是进一步提高腐乳发酵质量的必要条件。当前，有许多研究使用纯培养技术来揭示腐乳中的优势菌群，但这种技术受限于纯培养，对于含量较低、没有菌株分离、培养和鉴定过程的情况无法对腐乳中的整体微生物概况进行全面分析[8]。

近年来，随着测序技术的不断发展，高通量测序（high-throughput sequencing，HTS）作为新一代测序技术以其价格低廉、可读量大等优点被广泛用于探索各种豆豉、醋和米酒等生态系统中微生物的多样性，徐琼等[9]采用高通量测序技术分析不同地区红腐乳的细菌多样性，然而关于腐乳细菌及真菌的微生物群落同其产地和化学成分之间关系的研究还鲜见报道。

基于此，本研究以我国东部、南部和北部等不同地区所产的16种腐乳为研究对象，采用高通量测序（high-throughput sequencing，HTS）技术对腐乳中的细菌16S rRNA基因和真菌内部转录间隔区2（ITS2）基因进行测序，并针对腐乳中微生物多样性、关键微生物种属与腐乳化学性质、区域分布和调味处理的相关性展开研究，以期为提高传统发酵调味品的品质提供理论支持。

1 材料与方法

1.1 材料与仪器

16种腐乳样品 通过京东、淘宝等线上平台购自国内华东、华南和华北地区，具体编号、品牌、产地及调味方式如表1所示；实验中所用试剂 均购自郑州优吉生物科技有限公司，所有试剂纯度均为分析纯；磁珠法通用型基因组DNA提取试剂盒（DP705） 天根生化科技有限公司。

表1 腐乳16种样品的产地信息Table 1 Origin information of 16 sufu samples

欧莱博9830凯氏定氮仪 济南欧莱博科技有限公司；TGL-16M台式高速冷冻离心机 上海卢湘仪离心机有限公司；Sub-Cell GT电泳仪 美国伯乐；NanoPhotometer N60分光光度计 德国implen；Gene Amp 9700PCR仪 成都东胜科创；GelSMART凝胶成像仪 北京海鹏翔；Illumina Miseq高通量测序平台 上海派森诺生物科技有限公司。

1.2 实验方法

1.2.1 DNA提取和高通量测序

1.2.1.1 DNA提取 从16种不同的腐乳产品中各取1.0 g腐乳样品，使用基因组DNA提取试剂盒进行基因组DNA的提取和纯化。参考余巨全等[10]的方法，对所得到的基因组DNA进行纯度和浓度检测。吸取1.0 μL的DNA将其加至分光光度计中以检测所得基因组DNA中A260与A280的比值，以检测所得DNA的纯度；吸取2.5 μL的DNA将其加至1%的琼脂糖凝胶中，并将所得结果置于凝胶成像系统下进行成像以检测所得DNA的浓度。

1.2.1.2 PCR扩增 针对细菌16S rDNA基因的V3和V4区域合成带barcode的引物338F 5'-ACTCC TACGGGAGGCAGCAG-3'、806R 5'-GGACTACHV GGGTWTCTAAT-3'，针对真菌rRNA ITS2的保守序列合成带barcode的引物5'-TCCGTAGGTGGA ACCTGCGG-3'、5'-TCCTCCGCTTATTGATATGC-3'。以所提取的DNA为模板，使用Gene Amp 9700PCR仪进行PCR扩增，PCR扩增体系：基因组模板DNA 2 μL，5×PCR Master Mix 5 μL，正反引物各0.5 μL，用双蒸水补足至25 μL。细菌PCR扩增条件：94 ℃预处理3 min；30循环（94 ℃变性25 s，56 ℃退火40 s，72 ℃延伸40 s），72 ℃延伸5 min。真菌PCR扩增条件：94 ℃预处理3 min；32循环（94 ℃变性30 s，56 ℃退火35 s，72 ℃延伸40 s），72 ℃延伸5 min。用0.01 g/mL琼脂糖凝胶电泳进行检测扩增出的PCR产物并使用PCR产物纯化试剂盒纯化回收PCR产物，而后使用Illumina Miseq进行高通量测序。

1.2.2 腐乳样品的化学成分测定

1.2.2.1 样品处理 取腐乳样品20 g并将其分别置于研钵中，用杵磨成糊状，然后用80 mL去离子水在60 ℃下溶解。将混合物轻轻搅拌并煮沸，然后转入200 mL的容量瓶中，用水稀释至总体积200 mL。将制备好的腐乳悬浮液过滤，取滤液进行化学成分的测定。

1.2.2.2 氨基酸态氮含量的测定 以蒸馏水为空白对照，吸取10 mL的样品滤液，将其加入盛有90 mL蒸馏水的三角瓶（250 mL）中，将其充分混匀后，再加入10 mL的甲醛溶液（浓度为0.36 g/mL），用浓度为0.05 mol/L的氢氧化钠溶液进行滴定，直至pH为9.2，通过计算氢氧化钠溶液的消耗量来计算腐乳样品中的氨基氮含量[11]。

1.2.2.3 水溶性蛋白质含量的测定 以蒸馏水为空白对照，吸取10 mL的样品滤液，将其加入盛有50 mL蒸馏水的三角瓶（150 mL）中，充分混匀后100 ℃水浴加热10 min，待温度降至室温后将其转入容量瓶（100 mL）中，参照GB 5009.5-2016中的第一法进行水溶性蛋白质的测定，其中蛋白质的换算系数为5.71[12]。

1.3 数据处理

用测序数据质控软件Trimmomatic对原始FASTQ文件进行质量控制。使用UPARSE7.1软件（http://drive5.com/uparse/）对操作分类单元（operational taxonomic unit，OTUs）进行预聚类，其序列一致性截断率为97%。真菌ITS和细菌16S rRNA的所有操作分类单元均采用核糖体数据库程序（RDP）分类器（http://rdp.cme.msu.edu/）对Silva（SSU128）数据库进行分类，置信阈值为70%。

使用美吉生物云平台的（www.majorbio.com）的免费在线平台对所获序列进行生信统计分析。使用QIIME（Quantitative Insights Into Microbial Ecology）的α多样性模块对样品中Shannon、Simpson、Coverage和Chao1指数进行alpha多样性统计与分析，其中Shannon指数和Simpson指数代表样品中的菌落多样性，coverage代表样品文库的覆盖率、Chao1指数代表菌落丰度。使用Mann-Whitney检验对样本间各操作分类单元（OTU）相对丰度进行统计差异分析，P＜0.05表示差异显著。使用Majorbio Cloud platform（www.majorbio.com）的免费在线平台进行主坐标分析（Principal coordinate analysis，PCoA）和层次聚类分析（Hierarchical Cluster Analysis，HCA）。最后使用R语言进行微生物群落组成与腐乳化学成分之间的典范对应分析（Canonical correspondence analysis，CCA）和冗余分析（Redundancy analysis，RDA）。

2 结果与分析

2.1 腐乳样品的文库构建及测序数据分析

使用高通量测序技术，对来自不同地区的16份腐乳样品进行微生物多样性分析。共获得1032656个真菌群落的读数（OTU范围从22到239）以及689359个细菌群落的读数（OTU范围从68到289），样品中真菌ITS和细菌16S rRNA序列的覆盖度均大于99.90%，具有较好的微生物群落代表性。对所得样本序列进行随机取样，将抽到的样本序列与其代表的OTU数目构建稀释曲线，并以此来对测序数据量的合理性予以说明。所得结果如图1所示，由图1可知，16个腐乳样本的稀释曲线均趋于平缓，说明测序数据量合理，文库构建合理，测序深度可较好地反映微生物群落的信息[13]。

图1 16种腐乳样品的稀释性曲线Fig.1 Dilution curves in the 16 sofu samples

2.2 腐乳样品的真菌多样性分析

α多样性分析可反映测序样本中微生物群落的物种丰度和多样性，本研究采用ACE、Chao1、Shannon指数、Simpson指数对对腐乳样本中的真菌α-多样性进行评价，所得结果如表2及图2～图3所示。

表2 16种腐乳样品中真菌的α多样性分析Table 2 α diversity analysis of fungi in 16 sofu samples

其中，表2中的ACE指数和Chao1指数代表的是样本中微生物菌落的丰度，Shannon指数和Simpson指数代表的是样本中微生物菌落的多样性，但具体代表的值有所差别。Shannon指数值越大，代表微生物群落的多样性越高；Simpson指数越大，却代表微生物群落多样性越低。Coverage指的是所测样本文库的覆盖率，其数值的高低反映了所得到的测序结果是否完全代表样本的微生物群落[9]。从表2中可看到，所有样本的覆盖率均大于99%，说明文库的覆盖率基本涵盖样本中所有序列，本次测序结果可以代表样本中微生物的真实情况。

在真菌多样性方面，16个样本中的BJ18以最高的ACE（246.280）、Chao1（242.667）、Shannon指数（4.231）和最低的Simpson指数（0.046）而展现出最高的真菌丰度和多样性。BJ19和SH8分别以较低的ACE（27.057和46.035）和Chao（25.667和45.333）以及较高的Simpson指数（0.557和0.706）展现出较低的真菌丰度和多样性。在腐乳样品中，数量最多的5个真菌属分别为土赤壳属（Ilyonectriasp.）、刚毛藻属（Cadophorasp.）、曲霉属（Aspergillussp.）、指趾藤属（Dactylonectriasp.）和瓶霉菌属（Phialophorasp.）（图2），平均占真菌群落组成的56.08%。此外，不同的腐乳样品在真菌属成分上存在显著差异，这可能与腐乳产地、发酵原料和发酵工艺的不同有关[14]。

高通量测序（high-throughput sequencing，HTS）腐乳16个样品所得到的真菌属相对丰度百分比如图2所示。16种腐乳样品中的微生物群落主要由土赤壳属（Ilyonectriasp.）、曲霉属（Aspergillussp.）、酵母属（Saccharomycopsissp.）、放射毛霉属（Actinomucorsp.）、假丝酵母属（Candidasp.）、根霉属（Rhizopussp.）等胞外水解酶产生菌组成。此外，还存在刚毛藻属（Cadophorasp.）、枝孢属（Cladosporiumsp.）、束毛藻属（Trichomeriumsp.）、瓶霉属（Phialophorasp.）和镰刀菌属（Fusariumsp.）等大量的大豆共生体和病原菌。在产生水解酶的属中，放射毛霉属（Actinomucorsp.）在除了SH8和BJ19的其余14个腐乳样本中均有分布，而根霉（Rhizopussp.）分布比较有限，仅在ZJ1（2.11%）中被检测出；综合对比各样本中的真菌多样性，推测放射毛霉属（Actinomucorsp.）比其他丝状属更适合作为腐乳生产发酵剂中的水解酶生产菌[15]。此外，在16个腐乳样品中均检测到了曲霉属（Aspergillussp.），平均占所有属的8.65%，在样本SH8中更是以85.32%的比例占绝对优势。有报道称曲霉在腐乳生产过程中的后发酵中非常重要并占主导地位[16]，这与本研究的结果一致。

图2 16种腐乳样品中真菌属的相对丰度百分比（%）Fig.2 Relative abundance percentages (%) of the fungal genera in 16 sufu samples

利用主坐标分析（Principal coordinate analysis，PCoA）和层次聚类分析（Hierarchical Cluster Analysis，HCA）对16种腐乳中真菌组成的异同进行了分析，所得结果如图3所示。根据PCoA图，主坐标成分PC1和PC2分别占真菌群落方差的36.98%和18.86%。中国东部地区的真菌（ZJ2、ZJ4、ZJ5、ZJ6、SH9、AH11和JX12）可以归为一组，说明腐乳的真菌群落可能受地理因素的影响，因为来自同一地区的样品在进行腐乳发酵制备时通常具有相似的温度和水分条件[17]。除此之外，课题组研究发现，中国南方（GD14、GL15、LGM16）和中国北方（JL17、BJ18和BJ20）各样本间的距离均比较近，表明这些样本之间的菌群差异性较小，可归为一类，这表明除了地理因素之外，还有其它因素可影响腐乳真菌的多样性。这些腐乳产品可能使用类似的接种发酵剂、大豆原料或预发酵工艺[18]，这还需要进一步确认。另一方面，花香调味的ZJ1、ZJ4、ZJ5、SH9、BJ19以及酒香调味的ZJ2、ZJ6、SH8、AH11、GL15、JL17、BJ20和香辛调味的JX12、GD14、GZ16、BJ18这3组不同调味香料的样品在PCoA图中呈现不规则的分布，这表明调味香料不是影响腐乳中真菌群落组成的主要因素。

图3 16种腐乳样品中真菌组成的主坐标分析（A）和层次聚类分析（B）Fig.3 Principal coordinate analysis (A) and hierarchical cluster analysis (B) of the fungal genus compositions among 16 sufu samples

2.3 腐乳样品中的细菌多样性分析

16种腐乳样品的细菌α多样性检测结果如表3所示。在细菌多样性方面，16个样本中的ZJ1以最高的ACE（230.465）、Chao（186.895）、Shannon指数（1.117）和最低的Simpson指数（0.297）而展现出最高的细菌丰度和多样性。BJ19和GZ16分别以较低的ACE（39.697和86.113）和Chao（39.500和85.857）以及较高的Simpson指数（0.523和0.565）展现出较低的细菌丰度和多样性。对比真菌的α多样性结果，BJ19的细菌α多样性与真菌的α多样性结果一致，说明该产品在腐乳的发酵过程中应该是选用了特殊的发酵或加工方法[19]。此外，ZJ2和JX12的Shannon值最高，分别为2.187和2.162，Simpson值最低，分别为0.321和0.153，说明其菌群分布较其他样本更为均匀。而AH11的Shannon值最低（0.844），Simpson值最高（0.677），表明该样品中以少数物种为主[20]。

表3 16种腐乳样品中细菌的α多样性分析Table 3 α diversity analysis of bacteria in 16 sufu samples

16个腐乳样品中细菌属的相对丰度百分比如图4所示，16份腐乳样品中平均菌群含量排名前6位的是乳球菌属（Lactococcussp.）（33.83%）、芽孢杆菌（Bacillussp.）（11.79%）、四联球菌（Tetragenococussp.）（7.72%）、肠杆菌（Enterobactersp.）（6.36%）、链球菌（Streptococcuasp.）（6.10%）和明串珠菌（Leuconostocsp.）（5.72%），占菌群组成的71.52%，这与alpha多样性分析结果一致，即细菌群落倾向于以少数物种为主[21]。

从图4中还可以看到，腐乳样品中的细菌群落由可水解蛋白的芽孢杆菌属（Bacillussp.）、乳杆菌属（Lactobacillussp.）、魏斯氏菌属（Weissellasp.）、梭菌属（Clostridiumsp.）、肠球菌属（Enterococcussp.）、棒状杆菌属（Corynebacteriumsp.）等增味菌和不动杆菌属（Acinetobactersp.）、葡萄球菌属（Staphylococcussp.）等潜在的致病性或腐败菌株组成。所有腐乳产品中，乳球菌属（Lactococcussp.）最占优势，在腐乳样本ZJ1、GL15、SH9、JL17、ZJ2、ZJ5和GZ16中分布范围高达42.28%～76.69%。链球菌（Streptococcussp.）在SH8和BJ19中分布范围高达59.04%和16.36%。明串珠菌属（Leuconostocsp.）在样品SH9、GZ16和GD14中分别为35.21%、21.47%和13.67%。魏斯氏菌属（Weissellasp.）在BJ19、JX12和BJ20中的分布占12.68%、10.14%和10.00%。样品JL17中以乳酸菌（Lactobacillussp.）含量较高，为32.66%。这些不同样品中的微生物组成支持了乳酸菌在腐乳发酵中的重要性和复杂机制。此外，从图4中看到，芽孢杆菌菌株广泛存在于样品AH11、JX12和GD14中，但在许多其他腐乳样品中分布较少。

图4 16种腐乳样品细菌属的相对丰度百分比（%）Fig.4 Relative abundance percentages (%) of the bacterial genera in 16 sufu samples

利用主坐标分析（Principal coordinate analysis，PCoA）和层次聚类分析（Hierarchical Cluster Analysis，HCA）对16种腐乳中细菌组成的异同进行了分析，所得结果如图5所示。与真菌群落相似，不同腐乳样品的细菌群落差异显著。16种腐乳样品大体可分为乳酸杆菌类进化支（ZJ2、ZJ5、SH9、ZJ1、GL17、BJ20）和芽孢杆菌进化支（ZJ4、ZJ6、AH11、JX12、GD14）。这与之前对黄酒发酵的研究一致，8种发酵剂中有4种菌群以乳酸菌为主，而其余4个以芽孢杆菌为主。芽孢杆菌能提供多种水解酶，加速发酵物质中蛋白质、脂类和碳水化合物的降解，但芽孢杆菌对腐乳发酵的起源、影响及意义有待进一步研究。与真菌群落不同的是，样品中的细菌群落没有按区域进行聚类，说明腐乳样品中的细菌群落具有较强的变异性，可能受到多种因素的影响。与真菌群落相似，16种不同产品的菌群在发酵后也不能很好地按照调味方式进行聚类。因此可得出以下结论：腐乳中的细菌群落比真菌群落变化更大，其群落多样性在发酵后可能受到气候变化、浸泡和盐渍等高渗处理等多种因素的影响。

图5 16种腐乳样品中细菌属组成的主坐标分析（A）和层次聚类分析（B）Fig.5 Principal Coordinate Analysis (A) and Hierarchical Cluster Analysis (B) of the bacterial genus compositions among 16 sufu samples

2.4 腐乳样品中的化学成分分析

对选自华东、华南、华北等不同产地的16份腐乳样品中可溶性蛋白和氨基氮含量等化学成分进行检测，所得结果如表4所示。从表中可以看出，不同产地对于腐乳中的化学成分具有一定的影响，比如，华东产区中，产自浙江的腐乳样品（ZJ1、ZJ2、ZJ4、ZJ5、ZJ6）和产自上海的腐乳样品（SH8和SH9）中的可溶性蛋白含量和氨基氮含量彼此间比较接近，说明这些腐乳样品中所使用的发酵剂菌种具有较为相近的产蛋白酶能力。同时，产自安徽和江西的腐乳样品（AH11）和（JX12）可溶性蛋白含量和氨基氮含量明显高于其它华东地区的腐乳样品，结合前边的细菌丰度检验结果，推测是由于这两种样品中含量较为丰富的芽孢杆菌作用所致，为验证该推测，课题组又检测了也具有高丰度的芽孢杆菌分布的产自华南地区GD14样品中的可溶性蛋白和氨基氮含量，发现该腐乳样品也具有较高的可溶性蛋白和氨基氮含量，这说明芽孢杆菌在腐乳发酵的蛋白质分解中扮演极为重要的角色。

表4 16种腐乳样品的化学成分测定Table 4 Chemical determination of 16 different sufu samples

3 讨论

研究表明，酱油、醋、腐乳等发酵调味品中的微生物组成是决定其最终质量的关键因素[22]。本研究基于高通量测序技术，将腐乳发酵过程中微生物群落、化学成分、产地及风味相关性进行了探讨分析，以期揭示腐乳发酵过程中微生物多样性的影响因素及潜在后果。

研究中课题组发现华东地区所产的腐乳样本中的真菌群落可大致分为一类，这说明了地理因素在腐乳发酵过程中的重要作用。然而，腐乳样品中的细菌群落就不能很好地按照地域进行聚类，这表明细菌群落在腐乳发酵过程中的变异性更大，除去地理条件之外，还可能受到加工技术及发酵环境因素的影响，这和前人的报道相一致[23]。

微生物多样性分析有助于揭示腐乳发酵的关键菌群和功能菌群，为今后腐乳发酵菌群的合理设计提供依据。结果表明，放射毛霉属（Actinomucorsp.）作为产水解酶的真菌，在腐乳发酵过程中起着极为重要的作用，通过分析16种腐乳样品的真菌多样性，发现放射毛霉属（Actinomucorsp.）作为优势菌以较高的丰度百分比广泛分布于华南、华东和华北等多个产地的腐乳样品样本（GD14、ZJ1、AH11、BJ20）中；此外，课题组也在16种腐乳样品中检测到大量的可产蛋白酶、脂肪酶和淀粉酶等水解酶的曲霉属（Aspergillussp.），这些曲霉因其可广泛增加可溶性蛋白含量和氨基氮水平等特点被广泛用于豆制品的食品发酵中，腐乳样品的真菌多样性分析结果说明毛霉和曲霉对于腐乳的发酵生产具有极为重要的作用。陈卓等发现放射毛霉属是红腐乳发酵前期的主要菌种之一，红曲霉在发酵后占主导地位，该观点和本研究所得结果一致[24]。

本研究中根据细菌属的丰度，可将16种腐乳样品的细菌群落分为两大进化支：乳酸菌进化支和芽孢杆菌进化支。这两个分支的存在引起了课题组对于芽孢杆菌在豆制品发酵中的潜在作用的关注。芽孢杆菌是对食品发酵至关重要的菌种[25]，它可以产蛋白酶等水解酶，因而在腐乳的发酵和成熟中发挥重要作用；本研究通过对比AH11、JX12和GD14这三个产地的细菌多样性结果和化学成分的关联，发现腐乳样品中可溶性蛋白和氨基氮的含量同芽孢杆菌的含量呈现正相关，这为后期针对不同用户需求设计生产腐乳发酵剂提供了理论依据。此外，腐乳样品中的细菌多样性分析结果显示，包括乳球菌属（lactococcussp.）和链球菌（Streptococcuasp.）在内的乳酸菌广泛地存在于16种腐乳样品中，研究显示，这些菌群腐乳pH及风味的维持及生物防腐等方面发挥着重要影响[26]。

HTS技术可有效地揭示腐霉样品中细菌和真菌生态系统的整体状况。采用高通量测序，不仅检测出毛霉、曲霉和芽孢杆菌以及乳酸菌等优势菌群在各样品中的分布丰度，也检测出短梗霉（Aureobasidiumsp.）、聚孢霉（Clonostachyssp.）和金黄杆菌属（Chrysobacteriumsp.）等腐乳中的一些稀有属在AH11、ZJ2、SH9、ZJ5、ZJ6和ZJ4等样本中具有丰富的含量。研究表明，这些稀有属均来自于大豆材料[27]，且其存在会影响腐乳发酵微生物群落的组成[28-29]，具体影响机制还需要对腐乳发酵过程中微生物多样性及其特性的动态变化做进一步的研究。