青柠（Citrus aurantiifolia）隶属于芸香科柑橘属，在我国云南和广西等地均有栽培，年产量可达3.5～4.0 t，是世界柑橘类水果中的四大主要栽培种之一[1]。青柠中富含丰富的柠檬酸和维生素C，具有助消化、预防感冒以及其他抗疾病的功效[2]。酸柠檬是以青柠为原料，并添加食盐等辅料，经过晒干、混合和密封发酵等过程加工而成的[3]。酸柠檬，既保留了青柠原本的清香与大部分营养，又赋予了成品咸爽的口感，能够生津健胃，消食解乏[4]。广西地处亚热带，具有四周多山和高原、气候温暖、降雨量大等特征，果蔬品类较多，因此独具特色的“酸嘢”进入广西人民的餐桌上，当地各家都有腌制柠檬的习惯，并常将其用作烹饪肉禽鱼类等的辅料，其中尤以柠檬鸭最为出名，这种制作方式不仅可以降低原料的腥膻味，还能保留酸柠檬酸鲜口感[5-6]。酸柠檬加工工艺简单，整个制作过程一般都在常温常压下进行，成品品质易受多种因素影响，尤其是参与发酵的微生物。目前，对于酸柠檬的研究主要集中在酸柠檬食品的烹饪[7]和生产工艺的研究[8]等方面，而对于酸柠檬中微生物多样性的研究尚少。
随着测序技术的迅猛发展，Illumina MiSeq高通量测序技术以其快速高效的特点，成为微生物群落多样性解析的重要工具，目前已在发酵食品领域得到广泛应用[9]。王成等[10]通过高通量测序技术对Back-slopping发酵酸粥研究发现，乳酸杆菌属（Lactobacillus）和醋酸杆菌属（Acetobacteria）为优势细菌属；安飞宇等[11]通过高通量测序技术对自然发酵豆酱研究发现，乳杆菌属和四联球菌属（Tetragenococcus）为优势细菌属，青霉菌属（Penicillium）和异常威克汉姆酵母属（Wickerhamomyces anomalus）为优势真菌属；贾晶晶等[12]通过高通量测序技术发现，腌渍白菜液中，变形菌门（Proteobacteria）和普罗威登斯菌属（Providencia）分别为优势细菌门和优势细菌属。高通量测序技术的发展及在发酵食品领域的应用为解析酸柠檬中微生物群落结构组成提供了重要的技术支撑。
本研究从广西壮族自治区南宁市农贸市场共采集酸柠檬样品7份，采用Illumina MiSeq高通量测序技术对其微生物群落多样性进行了解析，同时采用PICRUSt软件对酸柠檬样品中微生物的基因功能进行预测，以期为后续南宁地区酸柠檬产品品质的提升提供理论依据和数据支撑。
1 材料与方法
1.1 材料与试剂
1.1.1 样品
从广西自治区南宁市秀厢农贸市场和淡村农贸市场共采集酸柠檬样品7份，编号为NM1～NM7，放置于含有冰袋的采样箱中运回实验室后进行分装，-20 ℃保存备用。
1.1.2 主要试剂
脱氧核糖核苷三磷酸（deoxy-ribonucleosidetriphosphate，dNTPs）Mix、FastPfu Fly DNA Polymerase（酶活500 U/g）、5×Trans StartTMFastPfu Buffer：北京全式金生物技术有限公司；琼脂糖（生化试剂）：索来宝生物技术有限公司；克隆载体pMD18-T：大连宝生物工程有限公司；QIAGEN DNeasy mericon Food Kit DNA基因组提取试剂盒：德国QIAGEN公司；聚合酶链式反应（polymerase chain reaction，PCR）清洁试剂盒：Axygen生物技术（杭州）有限公司。
1.2 仪器与设备
SW-CJ-2D型双人单面净化工作台：苏州安泰空气技术有限公司；VeritiFAST梯度PCR仪：美国ABI公司；UVPCDS8000凝胶成像分析系统：美国ProteinSimple公司；Illumina MiSeq高通量测序平台：美国Illumina公司；R920机架式服务器：美国Dell公司。
1.3 实验方法
1.3.1 酸柠檬样品中微生物宏基因组DNA提取
称取2.0 g酸柠檬样品，采用QIAGEN DNeasy mericon Food Kit DNA对样品中微生物宏基因组DNA进行提取，将琼脂糖凝胶电泳检测合格的宏基因组DNA放置于-20 ℃保存，备用。
1.3.2 细菌16S rRNA的PCR扩增
使用加入核苷酸标签（barcode）的成对引物338F（5'-ACTCCTACGGGAGGCAGCAG-3'）和806R（5'-GGACTACHVGGGTWTCTAAT-3'），参照文献[13]中的PCR扩增体系和条件对酸柠檬样品中细菌的16S rRNA V3-V4区基因序列进行PCR扩增。PCR扩增体系：5×PCR缓冲液4 μL，dNTPs Mix（2.5 mmol/L）2 μL，正反向引物（5 μmol/L）各0.8 μL，DNA聚合酶（5 U/μL）0.4 μL，DNA模板10 ng，补充无菌双蒸水（ddH2O）至20 μL。PCR扩增条件：95 ℃预变性3 min；95 ℃变性30 s，55 ℃退火30 s，72 ℃延伸45 s，循环30次；72 ℃再延伸10 min。
1.3.3 高通量测序及生物信息学分析
将经琼脂糖凝胶电泳检测合格的PCR扩增产物，通过Illumina MiSeq PE300高通量测序平台进行测序。将返回序列拼接后进行质控。参照沈馨等[14]的方法，通过QIIME（v1.90）平台对质控合格的序列进行生物学信息分析，将使用PyNAST软件进行对齐后的序列，按照100%和97%相似度进行两步UCLUST法构建操作分类单元（operational taxonomic units，OTU）矩阵；通过ChimeraSlayer软件对含有嵌合体的OTU进行选取与去除，其后选取代表性OTU序列在核糖体核糖核酸（ribonucleic acid，RNA）数据库RDP、SILVA和Greengenes数据库中进行同源性比对，确定其不同分类学地位；计算超1指数和香农指数从而对样品中微生物的丰度和多样性进行评价。
1.3.4 PICRUSt功能预测
将通过Greengenes数据库比对后的高质量序列进行UCLUST划分与注释，通过PICRUSt软件对酸柠檬样品中微生物的基因功能进行预测[15]，并依照直系同源蛋白数据库（clusters of orthologous groups of proteins，COG）进行功能注释[16]。
1.3.5 多元统计学分析
使用R软件进行气泡图与瀑布图的绘制；使用Origin 2017软件进行柱状图的绘制；使用MEGA 5.0软件的邻接（neighbor joining，NJ）法构建系统发育树。
2 结果与分析
2.1 酸柠檬细菌序列丰度和多样性分析
对本研究采集的7份酸柠檬样品进行Illumina MiSeq高通量测序，其16S rRNA V3-V4区基因序列测序结果、各分类地位数量统计及α-多样性分析结果见表1。
表1 基于测序结果酸柠檬样品细菌不同分类学地位数量统计及α-多样性分析结果
Table 1 Results of quantitative statistics of different taxonomic status and α-diversity analysis of bacteria in acidic lemon samples based on sequencing results
注：计算细菌超1指数和香农指数时，样品测序量均为33 010条序列。
由表1可知，7个酸柠檬样品共产生312 150条高质量16S rRNA序列，平均每个样品产生44 593条。按照97%和100%的相似度对所得序列进行划分且去除嵌合体后共得到8 958条代表性序列和5 223个细菌OTU。由表1亦可知，当测序量达到33010时，样品NM3的超1指数最高，为1116.79，样品NM5的香农指数最高，为4.13。有研究报道，超1指数越大，样本菌群丰富度越高；香农指数越小，群落多样性越高[17]。由此可见，样品NM3的物种丰富度最高，样品NM5的群落多样性最高。
2.2 基于微生物分类学地位酸柠檬样品细菌菌群多样性分析
本研究进一步从质控合格的序列中挑选代表性序列与数据库进行比对，共注释到32个门、87个纲、122个目、198个科和356个属，仅有0.076%和5.63%的合格序列无法注释到门和属水平。将7份酸柠檬样品中平均相对含量＞1.0%的门和属定义为优势细菌门和优势细菌属，将相对含量＜1.0%的细菌门和属定义为others，细菌门种类及平均相对含量气泡图见图1。
图1 基于门水平酸柠檬样品中细菌菌群结构分析结果
Fig.1 Analysis results of bacterial community structure in acidic lemon samples based on phylum level
由图1可知，7份酸柠檬样品在门的分类水平上，平均相对含量＞1%的细菌门有4个，分别为变形菌门（Proteobacteria）（80.61%）、硬壁菌门（Firmicutes）（9.24%）、放线菌门（Actinobacteria）（4.25%）和拟杆菌门（Bacteroidetes）（4.13%）。由图1亦可知，变形菌门在所有样品中均存在且相对含量均＞70%，为酸柠檬中的绝对优势菌门，不同酸柠檬样品优势细菌门的组成与含量无明显差异。葛东颖等[18]同样采用MiSeq高通量测序技术对与酸柠檬制作工艺相似的来凤地区盐渍藠头盐水中细菌多样性进行了解析，结果表明，优势细菌门为硬壁菌门（Firmicutes）、变形菌门（Proteobacteria）和放线菌门（Actinobacteria），与本研究结果相似。
本研究进一步对酸柠檬样品中优势细菌属的种类及平均相对含量进行统计，结果见图2。
图2 基于属水平酸柠檬样品中细菌菌群结构分析结果
Fig.2 Analysis results of bacterial community structure in acidic lemon samples based on genus level
由图2可知，7份酸柠檬样品在属分类水平上，平均相对含量＞1%的细菌属有5个，分别为劳尔氏菌属（Ralstonia）（66.13%）、伯克霍尔德菌属（Burkholderia）（8.18%）、诺卡氏菌属（Nocardia）（2.87%）、假单胞菌属（Pelomonas）（1.27%）和别样杆菌属（Alistipes）（1.01%）。值得一提的是劳尔氏菌属在样品NM1和NM2中含量较高，分别为71.14%和72.22%，远高于其他细菌属，为绝对优势细菌属。然而隶属于上述5个属的部分菌种为条件致病菌，如劳尔氏菌属中大部分菌种都具有将硝酸盐还原为亚硝酸盐的能力，从而引起人体中毒乃至死亡，有的则是植物病原菌，从而导致烟草青枯病[19-20]；伯克霍尔德菌属对多种抗菌药物具有耐药性，其中洋葱伯克霍尔德菌（Burkholderia cepacia）是引起囊性纤维化及慢性肉芽肿患者感染的最重要条件致病菌[21]。由此可知，自制的酸柠檬样品存在一定安全隐患。究其原因，可能与原料特性或农户制作酸柠檬开放的发酵条件有关，因此在后续对酸柠檬样品的研究中，应该在改善加工条件的基础上筛选出更优良的菌株，从而提高产品的安全性；或是采用纯培养技术收集菌株，发掘其潜在价值，如大多数伯克霍尔德氏菌属物种已证明对多种作物的不同植物病原体具有良好的生物防治作用；同时该属的部分菌株显示出巨大的生物技术潜力，可以作为新型抗生素和生物活性次生代谢物的来源[22]。
2.3 OTU统计分析
本研究进一步统计了样品中平均相对含量＞1%的核心OTU构成及平均相对含量，结果见图3。
图3 细菌核心OTUs平均相对含量分析结果
Fig.3 Analysis results of average relative content of bacterial core OTUs
由图3可知，7份酸柠檬样品中核心细菌OTU有6个，仅占总OTU数的0.11%，其中有2个隶属于仅鉴定到科的瘤胃菌科（Ruminococcaceae），4个分别隶属于瘤胃球菌属（Ruminococcus）、诺卡氏菌属、Reyranella和伯克霍尔德菌属。基于6个核心OUT序列及其模式菌株进行系统发育树的构建，结果见图4。
图4 酸柠檬样品中核心OTUs的系统发育树
Fig.4 Phylogenetic tree of core OTUs in acidic lemon samples
由图4可知，6个核心OTU共分为2大分支，其中OTU2529和OTU4497与纤维素分解瘤胃梭菌（Ruminiclostridium cel lulolyticum）聚为一类，OTU7750与扭链瘤胃球菌（Ruminococcustorques）聚为一类，OTU6193和产腔诺卡氏菌（Nocardia coeliaca）聚为一类，OTU1885与胰岛素伯克霍尔德氏菌（Burkholderia insulsa）聚为一类，OTU2043与马西里雷拉菌（Reyranella massiliensis）聚为一类。此外，6个核心OTU累计平均相对含量达到了89.24%，尤其是OTU4497（隶属于Ruminococcaceae）在7份酸柠檬样品中含量分别为76.06%、78.60%、72.80%、77.41%、69.82%、77.75%、72.93和75.05%，由此可见，OTU4497为绝对优势核心OTU。
由图4亦可知，OTU4497和OTU2529与R.cellulolyticum的相似度≥99%，由此可见，这2个核心OTU可能具有与R.cellulolyticum相似的功能。许多研究学者发现，R.cellulolyticum一方面与其他厌氧纤维素分解的微生物一样，能产生胞外多酶复合物，其中包含Cel9V和Cel9E等家族9糖苷水解酶，从而使纤维素分解，释放纤维二糖和纤维糊精，进一步降解后，为细胞提供能量[23-24]；另一方面因其可以分解纤维素的特点，会对未来生物物质转化成为生物燃料提供菌株支持[25]。由此可知，酸柠檬中亦可能存在大量有利于生物燃料合成的菌株，从而对其工业化的发展起到积极作用。
2.4 PICRUSt功能预测
本研究进一步使用PICRUSt软件对酸柠檬样品中细菌菌群的基因功能进行预测，在此基础上，根据蛋白质直系同源簇数据库对样品预测到的功能类别进行注释，结果见图5。
图5 酸柠檬样品中细菌功能类别分析结果
Fig.5 Analysis results of functional categories of bacteria in acidic lemon samples
7份酸柠檬样品中细菌序列注释到的COG有4 173个，分别隶属于23个功能类别。由图5可知，有关氨基酸转运与代谢、能量生产和转换功能的序列在酸柠檬细菌序列中占比分别为7.19%和9.04%，而细胞周期控制、细胞分裂、染色体分割、细胞运动、RNA的加工与修饰、染色质结构与动力学和细胞骨架的潜在表达则相对较低。由图5亦可知，7份酸柠檬样品共分为2大类，样品NM1、NM2、NM3、NM5、NM6和NM7中细菌功能之间存在一定的相似性，而样品NM4与其他样品之间细菌功能存在较大的差异。高秀芝等[26]研究发现，在豆酱发酵过程中，蛋白质被分解成氨基酸，从而产生独特的风味与优良的品质。由此推测，氨基酸转运与代谢功能在酸柠檬样品中的高表达，会对其品质与风味产生一定的积极作用。
3 结论
本研究利用Illumina MiSeq高通量测序技术对7份广西南宁地区酸柠檬样品的细菌多样性进行分析。结果表明，酸柠檬样品的细菌优势菌群组间差异不显著，且分布较为均匀。优势细菌门为Proteobacteria、Firmicutes、Actinobacteria和Bacteroidetes；优势细菌属为Ralstonia、Burkholderia、Nocardia、Pelomonas和Alistipes。基因功能及表型预测结果表明，酸柠檬中具有氨基酸转运与代谢、能量生产和转换功能的细菌序列占比较高。此外，潜在致病菌的检出指示，应该改善发酵条件或筛选出更优良的菌株控制发酵进程以提高产品的安全性。
[1]李素清，秦文，杨雨薇，等.响应面法优化青柠汁的果胶酶澄清工艺[J].饮料工业，2015，18（2）：39-43.
[2]沈月.柑橘家族，营养各有千秋[J].农家致富，2018（4）：59-59.
[3]刘志全.咸酸组合让柠檬菜出彩[J].四川烹饪，2007（8）：28.
[4]SILVA A F R,MONTEIRO M,RESENDE D,et al.Inclusion complex of resveratrol with γ-cyclodextrin as a functional ingredient for lemon juices[J].Foods,2020,10(1):16.
[5]李云腾.九个柠檬品种在广西南宁试种初报[J].广西农学报，2016，31（6）：21-24.
[6]广西壮族自治区卫生健康委员会.DBS 45/056—2018 食品安全地方标准柠檬鸭[S].南宁：广西出版社，2018.
[7]刘春荣，吴雪珍，郑雪良，等.柑桔资源制作美食现状与产业化开发建议[J].中国果业信息，2016，33（9）：13-15.
[8]张玉莲.广西靖西“酸嘢”腌制工艺考察与文化流变研究[D].南宁：广西民族大学，2012.
[9]SHAO J L,LI H,ZHANG X R,et al.Identification of serum exosomal microRNA expression profiling in menopausal females with osteoporosis by high-throughput sequencing[J].Curr Med Sci,2021,40(6):1161-1169.
[10]王成，肖晗，郭浩楠，等.自然发酵酸粥中细菌多样性分析及优势乳酸菌的分离与鉴定[J].食品科技，2020，45（10）：46-52.
[11]安飞宇，姜静，武俊瑞，等.自然发酵豆酱的滋味特性与微生物多样性分析[J].中国食品学报，2020，20（7）：207-215.
[12]贾晶晶，易洋，韩睿媛，等.高通量测序技术分析腌渍白菜液中原核微生物群落结构[J].中国调味品，2020，45（12）：87-90.
[13]王玉荣，折米娜，刘康玲，等.基于MiSeq 高通量测序技术内蒙古地区酸粥细菌多样性研究[J].食品工业科技，2018，39（19）：124-129.
[14]沈馨，王艳，代凯文，等.基于MiSeq 高通量测序技术的辣椒酱核心细菌类群研究[J].食品研究与开发，2018，39（10）：151-157.
[15]CHEN Y,TIAN W,SHAO Y,et al.Miscanthus cultivation shapes rhizosphere microbial community structure and function as assessed by Illumina MiSeq sequencing combined with PICRUSt and FUNGUIld analyses[J].Arch Microbiol,2020,202(5):1157-1171.
[16]HUERTA-CEPAS J,SZKLARCZYK D,HELLER D,et al.eggNOG 5.0:a hierarchical,functionally and phylogenetically annotated orthology resource based on 5090 organisms and 2502 viruses[J].Nucleic Acids Res,2018,47(D1):309-314.
[17]郭倩倩，卢彪.基于高通量测序不同生产工艺晴隆酸菜细菌多样性分析[J].中国酿造，2019，38（8）：73-76.
[18]葛东颖，仇港，赵楠，等.盐渍藠头盐水细菌多样性及基因功能预测分析[J].中国酿造，2020，39（10）：157-161.
[19]柴枝楠.亚硝酸盐中毒急救治疗和预防[J].中老年保健，2020，323（9）：36-36.
[20]向立刚，周浩，汪汉成，等.健康与感染青枯病烟株根际土壤与茎秆细菌群落结构与多样性[J].微生物学报，2019，59（10）：1984-1999.
[21]CUNHA M V,PINTO-DE-OLIVEIRA A,MEIRINHOS-SOARES L,et al.Exceptionally high representation of Burkholderia cepacia among B.cepacia complex isolates recovered from the major portuguese cystic fibrosis center[J].J Clin Microbiol,2007,45(5):1628-1633.
[22]ELSHAFIE H S,CAMELE I.An overview of metabolic activity,beneficial and pathogenic aspects of Burkholderia spp.[J].Metabolites,2021,11(5):321.
[23]VITA N,RAVACHOL J,FRANCHE N,et al.Restoration of cellulase activity in the inactive cellulosomal protein Cel9V from Ruminiclostridium cellulolyticum[J].FEBS Lett,2018,592(2):190-198.
[24]ROMAIN B,MINH-UYEN D T,HENRI-PIERRE F,et al.Catalytic subunit exchanges in the cellulosomes produced by Ruminiclostridium cellulolyticum suggest unexpected dynamics and adaptability of their enzymatic composition[J].FEBS J,2020,287(12):2544-2559.
[25]LIU N,FOSSES A,KAMPIK C,et al. In vitro and in vivo exploration of the cellobiose and cellodextrin phosphorylases panel in Ruminiclostridium cellulolyticum:implication for cellulose catabolism[J].Biotechnol Biofuels,2019,12(1):208.
[26]高秀芝，王小芬，李献梅，等.传统发酵豆酱发酵过程中养分动态及细菌多样性[J].微生物学通报，2008，35（5）：748-753.