基于高通量测序技术奶豆腐核心细菌类群的研究

雷炎，周亚澳，单春会，蔡文超，宋凯，郭壮

（1.湖北文理学院食品科学技术学院，湖北襄阳441053；2.石河子大学食品学院，新疆石河子832003）

0 引 言

奶豆腐作为新疆地区的传统发酵食品之一，因其口味和营养价值而广受欢迎[1]。越来越多的研究人员也针对奶豆腐展开了各种研究工作[2]。结果表明地域[3]、气候[4]和制作工艺[5]等均会影响奶豆腐的品质，而奶豆腐中的微生物在其发酵过程中起着至关重要的作用，并直接决定了奶豆腐的质量[6]，因而开展奶豆腐中微生物菌落结构的解析具有重要的意义。而以Illumina MiSeq 为代表的第二代测序技术的迅速发展也为解析环境中的微生物结构[7-8]，特别是为解析传统发酵食品中的微生物群落结构提供了有效的手段[9]。

本研究采用Illumina MiSeq 高通量测序技术对采集自新疆维吾尔自治区塔城地区的15 份奶豆腐样品中微生物构成和物种多样性进行了解析，以期于为奶豆腐的产业化发展、优良菌株的筛选以及传统奶豆腐的安全提供一定的参考依据。

1 实 验

1.1 试剂及仪器设备

1.1.1 材料

奶豆腐，采集自新疆维吾尔自治区塔城地区的乌苏市、额敏县和裕民县牧民家中。

5× TransStartTMFastPfu Buffer，FastPfu Fly DNA Polymerase，dNTPs Mix，DNA 1000 试剂盒，QIAGEN DNeasy mericon Food Kit DNA基因组提取试剂盒。

1.1.2 设备

R930 机架式服务器，美国DELL 公司；TGL-16M冷冻离心机，德国Eppendorf 公司；UVPCDS8000 凝胶成像分析系统，美国Protein Simple 公司；ND-2000C微量紫外分光光度计，美国TherimoFisher 公司；Illumina MiSeq 高通量测序平台，美国Illumina 公司；vetiri 梯度基因扩增仪，美国AB 公司。

1.2 样品的采集

使用在酒精灯上灼烧过的手术刀从牧民家自制的奶豆腐上切下一块样品，将奶豆腐的表面部分去除，修正成长×宽×高约为10 cm×10 cm×5 cm 的立方体，装入无菌袋中置于含有冰袋的采样箱中运回实验室，在超净台中使用灭菌的手术剪将其剪碎后装入灭菌的50 mL样品管中备用。

1.3 样品微生物宏基因组DNA 提取

准确称取2.0 g 奶豆腐样品，按DNA 基因组提取试剂盒说明书所约束的步骤进行宏基因组DNA 的提取，通过分光光度法和琼脂糖凝胶电泳法对提取DNA 的质量和浓度进行检测。

1.4 细菌16S rRNA V3-V4区PCR 扩增及测序

本研究对奶豆腐中细菌微生物的16S rRNA V3-V4区进行扩增，其正反向引物分别为338F（5'-ACTCCTACGGGAGGCAGCAG-3'）和806R（5'-GGACTACHVGGGTWTCTAAT-3'）。同时将成对的核苷酸标签（barcode）加入引物中用以区分归并每个样品的序列。参照王玉荣的方法进行PCR 扩增[10]，并将检验合格的PCR 产物送至上海美吉生物医药科技有限公司使用Illumina MiSeq 高通量测序平台进行测序。

1.5 生物信息学分析

本研究首先对下机后的双端序列依据其重叠关系进行拼接，并依据barcode 信息将所有序列归并到不同的样品中，同时参照Illumina MiSeq 标准数据质控流程进行质控[11]。将所有样品的序列合并后，参照Hou 的方法[12]，使用QIIME（v1.70）平台对所有有效序列进行微生物物种注释和多样性分析。其大致流程为100%和97%相似性序列的划分、构建分类操作单元矩阵（operational taxonomic units，OTU）、序列的同源性比对和多样性分析和基因功能预测。

1.6 核酸登录号

本研究中所有序列数据已提交至MG-RAST 数据库，登录号为mgp90776。

1.7 多元统计学分析

采用Pearson 相关性分析法对奶豆腐样品中优势细菌属（样品的平均相对含量大于1.00%）的相关系数和P 值进行计算，并用热图（Heatmap）对其进行可视化。使用Origin 8.5 软件（OriginLab Corp，MA，USA）进行绘图；使用R 软件（v3.3.2）进行相关性和显著性检验，加载waterfall包绘制瀑布图。

2 结果与分析

2.1 序列丰富度和多样性分析

采用MiSeq 高通量测序技术对奶豆腐中细菌微生物的构成和多样性进行了解析，并对奶豆腐的测序结果和分类学地位数量信息进行了统计，其结果如表1所示。

全部奶豆腐样品经高通量测序共产生了470210条高质量的16S rRNA 序列，每个样品平均可得到31347 条序列。按97%的相似性对序列进行划分后去除嵌合体OTU 共得到8 525 个OTU，进一步删除单序列OTU 得到4 209 个OTU 用于后续分析。由表1可知，A8 奶豆腐样品的OTU 数和超1 指数最大，A13奶豆腐样品的辛普森指数最高，且不同奶豆腐样品之间的OTU 数、超1 指数和辛普森指数均存在较大差异。超1 和辛普森指数常用来评估样品中的微生物丰度和多样性，因而A8 样品中微生物的丰度最高，而A13样品中微生物的多样性最高。

表1 样品16S rRNA测序情况及各分类地位数量

2.2 门属水平奶豆腐中细菌组成分析

本研究对OTU 的代表性序列进行了同源性比对，15 个奶豆腐样品中共鉴定出了32 个细菌门，540个细菌属。优势细菌门（平均相对含量大于1.00%的门）和优势细菌属（平均相对含量大于1.00%的属）在各个样品中的相对含量如图1所示。

由图1 可知，32 个细菌门中有4 个优势菌门，分别为Proteobacteria（变形菌门）、Frimicutes（硬壁菌门）、Actinobacteria（放线菌门）和Bacteroidetes（拟杆菌门），其平均相对含量分别为47.98%，43.84%，3.24%和1.90%。共鉴定出10 个优势属，分别是Lactobacillus（乳杆菌属，19.66%）、Acetobacter（醋酸杆菌属，16.69%）、Lactococcus（乳球菌属，13.90%）、Enterobacter（肠杆菌属，11.42%）、Stenotrophomonas（寡养单胞菌属，3.37%）、Staphylococcus（葡萄球菌属，3.33%）、Acinetobacter（不动杆菌属，2.79%）、Moraxella（莫拉菌属，1.78%）、Paenibacillus（类芽孢杆菌属，1.66%）和Streptococcus（链球菌属，1.14%）。进一步对奶豆腐中优势菌属之间的相互关系进行了探究，其相关性热图如图2所示。

图1 奶豆腐中优势细菌门和属相对含量的比较分析

图2 奶豆腐样品中优势细菌属之间的秩相关热图

由图2 可知，尽管大多数优势菌属之间存在一定的相关性，但其相关性并不强，经差异性检验可知，仅不动杆菌属与葡萄球菌属和链球菌属之间存在显著正相关性（P＜0.05,r＞0.5），其他优势菌属之间均不存在显著相关性，由此也可说明优势菌属之间的相互作用强度较弱。本研究进一步将在每个样品中均存在，且平均相对含量大于1.0%的属定义为优势核心属，其在每个样品中的相对含量如图3所示。

图3 奶豆腐中优势核心属相对含量的比较分析

由图3 可知，奶豆腐样品中共有7 个优势核心属，分别为乳杆菌属、醋酸杆菌属、乳球菌属、肠杆菌属、葡萄球菌属、不动杆菌属和链球菌属，其累计平均相对含量为68.94%。由图3 亦可知，A12 样品中核心优势菌属占比最多为98.13%，而A9 样品中核心优势菌属占比最少仅为12.49%。由此可知，尽管优势核心属存在于所有样品中，但其在不同样品中的分布却存在着较大的差异。

2.3 OTU 水平奶豆腐中细菌组成分析

本研究进一步对奶豆腐样品中OTU 出现的频率进行了分析，结果如图4所示。

图4 奶豆腐样品中OTU出现次数的统计分析

表1 结果显示不同样品中OTU 的数量存在较大差异。由图4 可知，仅出现1 次的OTU（只在一个样品中出现），即样品中的特有OTU 占总数的67.23%，其所包含的序列数为12 353 条，占总序列数的2.64%。出现15 次的OTU（在所有样品中均出现），即奶豆腐的共有OTU仅为总数的0.09%，但却包含250 078条序列，占全部序列数的53.50%。由此说明，虽然15个奶豆腐样品共有大量的核心菌群，但部分样品具有较为独特的细菌类群。在所用样品中均存在且平均相对含量大于0.10%的核心OTU 如图5所示。

由图5 可知，在所有样品中均存在且平均相对含量大于0.10%的OTU 有8 个，所包含序列数占总序列数的53.50%，其中OTU8262、OTU1070 和OTU2416均隶属于乳酸杆菌属，而OTU549，OTU1164，OTU1698，OTU7521 和OTU6656 分别隶属于Macrococcus（巨球菌属）、乳球菌属、肠杆菌属、不动杆菌属和醋酸杆菌属。由此可见，奶豆腐中的核心细菌类群主要以乳酸杆菌、醋酸杆菌属和不动杆菌属为主。

乳酸菌作为发酵乳制品中最常见的细菌之一，能发酵碳水化合物产生大量的乳酸，对乳制品的发酵进程有着重要的影响，已有研究表明乳酸菌能产生多种生物活性物质，且这些活性物质在发酵完成后仍然保留于乳制品中[13]。DUGAT-BONY E 研究了12 种法国奶酪发现其均以乳酸乳球菌（Lactococcus lactis）为主[9]，DOMINGOS-LOPES M F P 也发现乳酸菌为手工奶酪细菌菌群的重要组成部分[14]，上述研究结果均与本文中的结论相一致。

乳杆菌属和乳球菌属是奶豆腐中最常见的益生菌，与人体具有共生关系保护宿主免受病原体的潜在入侵[15]，它在维持人类健康方面的应用广泛，可以帮助治疗腹泻、阴道感染和湿疹等疾病[16]。醋酸杆菌属的特点是能够将乳酸盐和乙酸盐氧化成二氧化碳和水，从而赋予奶豆腐样品特殊的口感[17]，这也赋予了新疆奶豆腐样品特有的功能与风味。值得注意的是，奶豆腐样品中存在着一定量的肠杆菌属，可能导致有机体的多种疾病，如败血症和腹泻等[18]。调查发现，牧民制作奶豆腐的环境较为开放和复杂，制作过程中奶豆腐与空气和模具等接触频繁，这也可能是导致奶豆腐中存在条件致病菌的原因之一。

图5 奶豆腐样品中平均相对含量大于0.1%的OTU

2.4 PICRUSt基因预测

牛奶中富含蛋白质和碳水化合物等营养物质，其中的微生物具有代谢和利用这些营养物质的基因。细菌作为奶豆腐发酵进程的主要参与者，通过利用乳源中的营养物质进行新陈代谢，从而赋予奶豆腐特有的风味。本研究使用PICRUSt 对奶豆腐样品中细菌微生物的基因功能进行了预测，并参考蛋白质直系同源簇数据库（COG）进行功能注释[19]。研究共注释到4 792 个COGs，所有COGs 分别隶属于20 个功能大类。奶豆腐样品中细菌功能类别如图6所示。

由图6 可知，奶豆腐中细菌的功能基因主要与氨基酸、碳水化合物和无机盐离子运输和代谢相关，与能量的产生和转化相关的基因进行了高表达，而与细胞运动和RNA 加工修饰相关的基因表达量却相对较低。乳糖是牛奶中的主要碳水化合物，也是细菌能量的重要来源。奶豆腐中存在的大量乳酸菌通过其相关基因的高表达将乳源中的乳糖转化为乳酸，并产生一些风味物质。相关报道显示，尽管不同的乳酸菌都能代谢乳糖转化为乳酸，但此过程中生成的一些风味物质却并不相同，这也使得奶豆腐的品质存在较大的差异。通过对奶豆腐中微生物多样性进行解析，探讨了乳酸菌在奶豆腐发酵过程中的重要作用，对奶豆腐的工业化生产和菌株筛选具有一定的参考意义。

图6 奶豆腐中细菌功能类别

3 结 论

结果表明，塔城地区奶豆腐中优势细菌门分别为变形菌门、硬壁菌门、放线菌门和拟杆菌门，优势细菌属主要为乳杆菌属、醋酸杆菌属、乳球菌属、克雷伯氏菌属、寡养单胞菌属、葡萄球菌属、不动杆菌属、莫拉菌属、类芽孢杆菌属和链球菌属。塔城地区奶豆腐中主要以乳酸菌为主，通过生成乳酸和其他代谢产物，从而赋予奶豆腐特殊的风味。