利用高通量测序分析云南豆豉中细菌群落多样性

李晓然， 李 洁， 刘晓峰， 冯 阳， 柳陈坚

（昆明理工大学 生命科学与技术学院，云南 昆明 650500）

发酵食品是指利用微生物加工制造而成的一类具有独特风味与营养价值的食品，诸如酸奶、奶酪、泡菜、豆豉等。其中豆豉是以大豆或黄豆为主要原料，经由原料处理、制曲及后发酵三阶段酿制而成，它不仅风味独特，而且还富含多种生理活性物质，因此豆豉具有开胃增食、消积化滞、驱风散寒及预防心脑血管疾病等机能。豆豉品种繁多，根据发酵微生物的不同可以分为细菌型、曲霉型、毛霉型、根霉型和脉孢菌型等五大类[1]。由于云南省独特的地理位置、气候条件和民族文化传统，生物资源繁多，细菌型豆豉中具有独特和宝贵的菌种资源，由于缺乏对这些细菌全面的了解，这些宝贵的资源仍未得到很好的开发。

在传统发酵食品中，通常采用富集纯分离培养与鉴定方法对其中微生物群落结构进行研究，然而该方法不但操作繁琐、费时费力，尤为重要的是在现有技术条件下，它仅能分离少数参与发酵的微生物，而大多数微生物得不到相应的分离培养与鉴定，通过纯培养方法仅能得到传统发酵样品中微生物群落的极少部分信息，而不能真正反映传统发酵食品中微生物群落的真实信息，从而给充分发掘与利用发酵食品中的微生物资源带来一定困难。李祥研究了参与细菌型豆豉发酵的细菌群落，发现主要有豆豉芽孢杆菌 （Bacillus douchi）、枯草芽孢杆菌（B.subtilis）、 乳酸杆菌 （Lactobacillus） 和微球菌（Micrococcus）[2]。基于 16S rRNA 基因克隆文库或变性梯度凝胶电泳 （DGGE）能够规避培养方法的弊端，使得对于样品中的未培养细菌的分析成为可能。Chen等采用了PCR-DGGE的方法分析了豆豉中微生物群落结构，发现大多数豆豉中都检测到了枯草芽孢杆菌、解淀粉芽孢杆菌（B.amyloliquefaciens）、假单胞菌（Pseudomonas）、酿酒酵母 （Saccharomyces cerevisiae）和粉状毕赤酵母（Pichia farinose），同时也检测到了潜在的病原微生物如葡萄球菌（Staphylococcus）、泛菌（Pantoea）和肠杆菌（Enterobacter）等[3]，暗示了在传统发酵食品中潜在致病菌存在的广泛性，因此使用分子生物学方法对其进行检测是极为必要的。

随着测序技术发展，被称为“第二代测序（nextgeneration sequencing[NGS]）”技术在过去几年中不断涌现与发展。第二代测序技术的共同特点是：高通量、低成本、可以同时对多个单独的DNA分子进行测序分析。基于Pyrosequencing[4-5]测序法一次能够测定100万个长度达到400 bp的核苷酸序列，其费用只相当于使用Sanger法测定几千个序列的费用；由于该方法不需要克隆，可解决构建克隆文库时所呈现的偏向性问题。2008年，Hamady等[6]展示一种同时可以平行分析多个样品的焦磷酸测序法。在对不同样品进行PCR扩增时，在引物5′末端添加由8个碱基组成的标记（bar-code），应用这种方法，在同一次测序反应中可以同时分析1 544个样品。这一方法在多种类型样品中的微生物多样性及群落结构研究中得到极为广泛应用。然而，该方法在传统发酵食品中的应用还极为有限，仅有少数发酵食品使用该方法进行相应研究[7-8]。

作者通过高通量测序和序列分析，以期获得云南细菌型豆豉中群落结构和多样性，为进一步利用和开发豆豉中为细菌资源提供研究基础，并为食品质量和安全提供一定信息。

1 材料与方法

1.1 材料与试剂

1.1.1 样品采集 试验用豆豉：购于云南省玉溪市易门县菜市场，为当地百姓手工制作。购买后置于冰盒中运送到实验室，-80℃保存。

1.1.2 主要缓冲液 DNA提取裂解缓冲液：0.1 mol/L Tris/HCl，0.1 mol/L EDTA，0.75 mol/L 蔗糖。

1.2 方法与步骤

1.2.1 DNA提取 DNA提取参考Schmidt等[9]的方法并做了部分改变，步骤如下：取约0.2 g样品置于无菌的2 mL离心管中，加入450 μL裂解缓冲液（0.1 mol/L Tris/HCl，0.1 mol/L EDTA，0.75 mol/L 蔗糖），加入 10 μL、50 mg/mL 溶菌酶（lysozyme）和 10 μL、20 mg/mL 溶壁酶（lyticase），避免剧烈振荡，充分混合均匀后，37℃水浴30 min。加入25 μL 20%SDS和5 μL 20 mg/mL蛋白酶K，避免剧烈振荡，充分混合均匀后，55℃水浴2 h以上，加入0.1 g玻璃珠（直径 212～300 μm，美国 SIGMA），于 labnet230 V EU涡旋仪（美国labnet）以最高转速振荡5 min，将管子放置于研钵中，加入液氮覆盖，待液体全部冻住后放入55℃水浴锅中融化，重复3次。加入80 μL 5 mol/L NaCl，小心混匀后，加入60 μL 10%CTAB/NaCl，小心混匀，65℃水浴20 min。加入700 μL 苯酚：氯仿：异戊醇（体积比 25∶24∶1），颠倒混匀后于4 000 r/min离心20 min，将上层水相转移至新的2 mL离心管中，加入等体积的氯仿∶异戊醇（24∶1，体积比），颠倒混匀后于4 000 r/min离心20 min。将上层水相转移至新的2 mL离心管中，加入0.6体积的异丙醇，-20℃沉淀过夜后于12 000 r/min离心15 min，倒掉液体，用500 μL预冷的70%乙醇洗沉淀，去掉液体后，将离心管打开盖子，置于60℃烘干箱干燥，待液体全部挥发后取出，加入120 μL TE（10 mmol/L Tris-HCl（pH 8.0），1 mmol/L EDTA（pH8.0））缓冲液溶解DNA。提取好的DNA存放于-20℃。

1.2.2 细菌16S rRNA基因扩增和焦磷酸测序 使用细菌的16S rRNA基因通用引物对提取的DNA进行扩增。扩增使用通用引物338F（5’-ACH YCT ACG GGA GGC HGC-3’） 和 907R （5’-CCG TCA ATT CMT TTG AGT TT-3’），并在引物338F的5’末端添加8个碱基的标记和焦磷酸测序的AdapterA，在907R的5’末端添加AdapterB。为了减少PCR产物中非特异性扩增的异源双链DNA的含量，将PCR产物稀释10倍作为模版后按照原来的条件再做5个循环[10]。PCR扩增采用Premix Ex TaqTM（大连TaKaRa），具体体系及程序如下：

1）PCR 反应体系：ddH2O 18 μL；Forward 引物（5 μmol/L）2 μL；Reverse 引 物 （5 μmol/L）2 μL；Premix Ex TaqTM25 μL；模板（15 ng）3 μL；总体积50 μL。

2）PCR扩增程序：95℃预变性 5 min；94℃变性1 min；56℃退火1 min；共30个循环。72℃延伸1 min 30 s；72 ℃延伸 7 min。

将PCR产物在1 g/dL的琼脂糖凝胶中电泳，电泳后将凝胶浸泡于含有20 μg/mL EB的TAE缓冲液中，10 min后在紫外灯（216 nm）照射下呈现橙色条带，目标条带割下后使用QIAquick Gel Extraction Kit（德国Qiangen）对PCR产物进行回收。使用nanodrop ND-1000分光光度计（美国Thermo）对完成纯化的PCR扩增产物定量。使用GS-FLX（瑞士Roche）测序系统对DNA进行测序。

1.2.3 测序结果分类鉴定 对于焦磷酸测序的测序结果分析质量文件，去掉测序质量不好的序列，并对序列长度进行筛选，删掉长度小于300 bp的序列，根据barcode将序列确定为本样品的序列，并去除barcode和引物序列，得到有效的序列文件。使用Mallard v1.02[11]来检测PCR产物序列中的嵌合体。使用mothur v1.25.1[12]的classify.seqs命令，采用Silva的核糖体小亚基（SSU）rRNA序列数据库V102[13]的分类信息，得到序列的分类信息，分别采用80%的阈值（threshold）。

1.2.4 多样性指数和系统发育分析 对所有序列使用mothur v1.25.1[12]进行比对，以0.03为cutoff值划分可操作分类单元 （operational taxonomic unit，OTU）。计算cutoff为0.03时的ACE和Chao1值，并绘制rarefaction曲线。从细菌16S rRNA基因序列中每个OTU选取代表序列，将其在NCBI的GenBank进行Blast比对，将结果中具有确定分类信息的序列用ClustalX 1.83[14]进行比对，将比对的结果用Mega v4.0[15]中的Neighbour-joining方法构建系统发育树。

2 结果与分析

2.1 序列数目和多样性指数

经过对焦磷酸测序所得序列进行质量和长度删选，得到356条高质量序列。对这些序列以0.03为cutoff进行分析，共得到62个OTU，其中，singleton42个，doubleton7个。 ACE值为 387，Chao1值为170，rarefaction曲线见图1。可见OTU数目随着序列数上升的趋势还远未达到平台期。

图1 根据rarefaction方法估计豆豉样品的多样性Fig.1 Estimated OTU numbers，according to the rarefaction method，depending on OTUs identified with a 3%cut-off

2.2 细菌群落多样性分析

根据Silva的SSUrRNA序列数据库对序列进行分类，在356条序列中，有91%（324条）的序列属于厚壁菌门（Firmicutes），9%（32 条）的序列属于变形菌门（Proteobacteria）。在厚壁菌门中，占总序列72%（258条）的序列属于乳杆菌科（Lactobacillaceae），10%（36 条） 的序列属于芽孢杆菌科（Bacillaceae），4%（15 条）的序列属于肉杆菌科（Carnobacteriaceae）；在变形菌门中，所有的序列都属于肠杆菌科 （Enterobacteriaceae），2%的序列不能确定分类信息（图2（a））。在属这一分类水平上，总序列的72%（256条）的序列都属于乳酸杆菌属（Lactobacillus），8%（28条） 的序列属于芽孢杆菌属（Bacillus），而肠杆菌科的序列则大多属于肠道菌簇（Enteric_Bacteria_cluster），3%的序列不能确定分类信息（图 2b）。

图2 豆豉中细菌在不同分类水平的分布Fig.2 Major families and genera in fermented soybean

2.3 系统发育分析

以0.03为cutoff，对356条序列进行序列比对和距离矩阵分析，确定了62个OTU，将其中的singleton和doubleton去除后，还有13个OTU，将其与GenBank数据库中的序列进行比对，并下载相近序列共同构建系统发育树，见图3。丰度最高的OTU是YM_QT，占了总序列数的64%，经过比对，它与乳酸杆菌科的L.acidifarinae和L.zymae都有大于97%的相似度。丰度第二高的是OTU YM_DT，有11条序列，在系统发育树上并未与已知分类信息的细菌呈现较高的相似度，但是可以确定分在芽孢杆菌属中。有10条序列属于OTU YM_P8，也属于芽孢杆菌属，与B.thermoamylovorans有98%的相似度。另有两个OTU，YM_FM和YM_DP，均含有9条序列，虽然在系统发育树上可以归入乳酸杆菌属，但是不能与已知分类信息的乳酸杆菌聚类在一起。在变形菌门中，有两个OTU，其中YM_W2中有8条序列，与Erwinia persicina有99%以上的相似度，而YM_JS只有3条序列，与Shigella flexneri较为相似。虽然在属的水平只有3%的序列不能确定分类信息，但是对于丰度较高的OTU，还是有约占总序列数10%的序列不能确定到种的水平。

3 结语

近年来，食品发酵中的微生物变化受到关注并得到了研究，科学家们通过rRNA基因序列分析的方法研究了如泡菜[16]、发酵的芥末酱[17]和发酵乳品[18]等世界各地不同的发酵食品的微生物群落多样性，作者对云南特色的豆豉中细菌群落多样性进行了研究，使用高通量焦磷酸测序技术，得到了356条序列，较为充分地展示了豆豉中细菌群落结构。通过多样性指数和rarefaction曲线的结果可以看出，这一样品细菌群落的多样性是非常高的，因此，为了分析豆豉中的稀有群落（rare biosphere），使用高通量测序分析群落结构是必须的。在本研究中，约10%的序列不能确定具体的分类信息，有2%的序列甚至不能确定到科，这说明豆豉虽然为富营养环境，含有较多的可培养细菌，但是仍有一部分细菌为目前不能培养的细菌类型。

在他人的研究中，传统豆豉的生产有很大的地域性差异，大多豆豉为真菌（霉菌）型豆豉[1]，而在本研究中，虽然对易门豆豉采用溶壁酶共同提取了真菌的DNA，在对真菌转录间隔区（ITS）序列扩增时则未得到相应条带，因此，易门豆豉应该是以细菌起主要作用的细菌型豆豉。

在本研究分析的豆豉样品中，大多数序列都集中在乳酸杆菌中，作为主要的乳酸菌种类，乳酸杆菌属的细菌在多种发酵食品中检测到并发挥重要作用。作为乳酸杆菌中丰度最高的OTU，YM_QT与L.acidifarinae和L.zymae都显示了很高的相似度，这两株菌都是从小麦酵母中分离到的[19]，也证明这两株菌对于植物性底物的发酵具有较好的作用，是豆豉的生产过程中最为主要的细菌。丰度第二高的OTU为YM_P8，可认为是B.thermoamylovorans，这是一株中度嗜热和淀粉的细菌[20]，在豆豉制作过程中也起着重要作用。这一结果与以往的研究认为细菌型豆豉以芽孢杆菌为主[1-2]并不相符，主要原因为：一方面，由于以前的研究大多是采用培养的方法，这使得芽孢杆菌更容易在培养基中培养出来，因此在一定程度上高估了芽孢杆菌的含量；另一方面，在基于PCR的方法的研究中，也有以芽孢杆菌为主的结果[3]，这应该与豆豉的类型有关。在Chen的另一项研究中也发现主要的细菌为乳酸杆菌[21]。

图3 豆豉中细菌的系统发育树Fig.3 Phylogenetic tree of the 16S rRNA gene sequences and their phylogenetic relatives

除了乳酸杆菌和芽孢杆菌这类在发酵中起作用的细菌，在这个豆豉样品中，也发现了潜在的致病菌，YM_JS虽然只有3条序列，却与导致人类腹泻的病原菌S.flexneri显示出高度的相似性。S.flexneri为公认的食源致病菌，在多种类型的食品中发现其存在，并可能导致大规模腹泻的爆发[22-25]。另一个OTU YM_W2则与E.persicina高度一致，而E.persicina则是植物的病原菌[26]，会导致豌豆苗萎蔫[27]，暗示了用来制作豆豉的大豆可能有过细菌性疾病。在市场随机购买的豆豉样品中发现的这两个属于变形菌门的OTU，一个是人类的致病菌，一个是食品原料的致病菌，无不显示着食品安全问题的紧迫性，对于致病菌快速准确的检测是非常必要的[28]。

[1]胡会萍，李秀娟，黄贤刚.传统豆豉微生物学研究综述 [J].中国调味品，2012，37（6）：4-13.HU Hui-ping，LI Xiu-juan.HUANG Xian-gang.Research of traditional douchi microbiology[J].China Condiment，2012，37（6）：4-13.（in Chinese）

[2]李祥.细菌型豆豉生产的研究 [J].中国调味品，1999（10）：14-17.LI Xiang.Study on production of bacterial douchi[J].China Condiment，1999（10）：14-17.（in Chinese）

[3]CHEN T，JIANG S，XIONG S，et al.Application of denaturing gradient gel electrophoresis to microbial diversity analysis in Chinese douchi[J].J Sci Food Agric，2012，92（10）：2171-2176.

[4]EDWARDS R A，RODRIGUEZ-BRITO B，WEGLEY L，et al.Using pyrosequencing to shed light on deep mine microbial ecology[J].BMC Genomics，2006，7：57.

[5]NYR N P，PETTERSSON B，UHL N M.Solid phase DNA minisequencing by an enzymatic luminometric inorganic pyrophosphate detection assay[J].Anal Biochem，1993，208（1）：171-175.

[6]HAMADY M，WALKER J J，HARRIS J K，et al.Error-correcting barcoded primers for pyrosequencing hundreds of samples in multiplex[J].Nat Methods，2008，5（3）：235-237.

[7]LI X R，MA E B，YAN L Z，et al.Bacterial and fungal diversity in the traditional Chinese liquor fermentation process[J].Int J Mol Sci，2011，146（1）：31-37.

[8]KIM Y-S，KIM M-C，KWON S-W，et al.Analyses of bacterial communities in meju，a Korean traditional fermented soybean bricks，by cultivation-based and pyrosequencing methods[J].J Microbiol，2011，49（3）：340-348.

[9]SCHMIDT T M，DELONG E F，PACE N R.Analysis of a marine picoplankton community by 16s rRNA gene cloning and sequencing[J].J Bacteriol，1991，173（14）：4371-4378.

[10]THOMPSON J R，MARCELINO L A，POLZ M F.Heteroduplexes in mixed-template amplifications：formation，consequence and elimination by'reconditioning PCR'[J].Nucleic Acids Res，2002，30（9）：2083-2088.

[11]ASHELFORD K E，CHUZHANOVA N A，FRY J C，et al.New Screening Software Shows that Most Recent Large 16s rRNA Gene Clone Libraries Contain Chimeras[J].Appl Environ Microbiol，2006，72（9）：5734-5741.

[12]SCHLOSS P D，WESTCOTT S L，RYABIN T，et al.Introducing mothur：open-source，platform-independent，communitysupported software for describing and comparing microbial communities[J].Appl Environ Microbiol，2009，75（23）：7537-7541.

[13]PRUESSE E，QUAST C，KNITTEL K，et al.Silva：a comprehensive online resource for quality checked and aligned ribosomal RNA sequence data compatible with arb[J].Nucleic Acids Res，2007，35（21）：7188-7196.

[14]THOMPSON J，GIBSON T，PLEWNIAK F，et al.The clustal_X windows interface：flexible strategies for multiple sequence alignment aided by quality analysis tools[J].Nucleic Acids Res，1997，25：4876-4882.

[15]TAMURA K，DUDLEY J，NEI M，et al.Mega4：molecular evolutionary genetics analysis（Mega） software version 4.0[J].Mol Biol Evol，2007，24（8）：1596-1599.

[16]CHANG H W，KIM K H，NAM Y D，et al.Analysis of yeast and archaeal population dynamics in kimchi using denaturing gradient gel electrophoresis[J].Int J Food Microbiol，2008，126（1-2）：159-166.

[17]CHAO S H，WU R J，WATANABE K，et al.Diversity of lactic acid bacteria in suan-tsai and fu-tsai，traditional fermented mustard products of Taiwan[J].Int J Food Microbiol，2009，135（3）：203-210.

[18]剧柠，张家超，孙志宏，等.制作云南乳扇用酸乳清中乳杆菌的多样性分析[J].食品与生物技术学报，2010，29（5）：123-129.JU Ning，ZHANG Jia-chao，SUN Zhi-hong，et al.The polymorphism of Lactobacilli from acid whey for dairy fan in Yunnan area[J].Journal of Food Science and Biotechnology，2010，29（5）：123-129.（in Chinese）

[19]VANCANNEYT M，NEYSENS P，WACHTER M D，et al.Lactobacillus acidifarinae sp.nov.and Lactobacillus zymae sp.nov.，from wheat sourdoughs[J].Int J Syst Evol Microbiol，2005，55（Pt2）：615-620.

[20]COMBET-BLANC Y，OLLIVIER B，STREICHER C，et al.Bacillus thermoamylovorans sp.nov.，a moderately thermophilic and amylolytic bacterium[J].Int J Syst Evol Microbiol，1995，45（1）：9-16.

[21]CHEN T，XIONG S，JIANG S，et al.Molecular identification of microbial community in Chinese douchi during post-fermentation process[J].Food Sci Biotechnol，2011，20（6）：1633-1638.

[22]LAU H K，CLOTILDE L M，LIN A P，et al.Comparison of ims platforms for detecting and recovering Escherichia coli O157 and Shigella flexneri in foods[J].Jala，2013，18（2）：178-183.

[23]MAHBUB M M，AHSAN C R，YASMIN M，et al.Analysis of different prognostic Indicators for malnutrition and Shigella flexneri infection among the children in bangladesh[J].Indian J Microbiol，2012，52（3）：400-405.

[24]MAHMOUD B S M.Effects of X-ray radiation on Escherichia coli O157：H7，Listeria monocytogenes，Salmonella enterica and Shigella flexneri inoculated on shredded iceberg lettuce[J].Food Microbiol，2010，27（1）：109-114.

[25]李长军，肖桂芝.志贺菌的菌型分布与耐药性研究[J].航空航天医学杂志，2011，22（8）：1008-1009.LI Chang-jun，XIAO Gui-zhi.Serotype distribution and drug resistance of Shigella[J].Journal of Aerospace Medicine，2011，22（8）：1008-1009.（in Chinese）

[26]STARR M P，CHATTERJEE A K.The genus erwinia：enterobacteria pathogenic to plants and animals[J].Annu Rev Microbiol，1972，26：389-426.

[27]崔汝强，胡学难，廖金铃，等.进境澳大利亚豌豆细菌性萎蔫病病原菌的分离鉴定[J].植物病理学报，2009，39（3）：238-242.CUI Ru-qiang，HU Xue-nan，LIAO Jin-ling，et al.Isolation and identification of bacterial pathogen causing wilt on pea seed imported from Australia[J].Acta Phytopathologica Sinica，2009，39（3）：238-242.（in Chinese）

[28]孙秀兰，张芳，张银志.多种致病菌同步检测方法的研究进展[J].食品与生物技术学报，2012，31（5）.SUN Xiu-lan，ZHANG Fang，ZHANG Yin-zhi.Research advance study of the synchronous detection of multiplex foodborne pathogens[J].Journal of Food Science and Biotechnology，2012，31（5）.（in Chinese）