广陈皮为芸香科植物橘（Citrus reticulata Blanco）的栽培变种茶枝柑（Citrus reticulate Chachi）的干燥成熟果皮[1]，有理气健脾、燥湿化痰之功。在我国有着悠久的历史，具有很高的药用价值，同时作为传统的香料和调味佳品，享有盛誉[2]。传统理论认为陈皮“陈久者良”，即随着陈化时间的增加，陈皮的质量会提高。目前，关于陈皮“陈久者良”的研究主要集中在不同陈化年限陈皮化学成分、药理作用的强弱比较方面[3-4]，而对其科学内涵的现代研究相对薄弱。
本课题组前期报道陈皮表面的曲霉属（Aspergillus）是引起黄酮类和挥发油含量与类型变化的重要因素[5]。陈聪聪[6]研究发现陈皮表面芽孢杆菌属（Bacillus）、乳酸链球菌属（Lactococcus）、假单胞菌属（Pseudomonas）等优势细菌可加速其挥发油的生物转化。这些研究证实了陈皮表面微生物会引起陈皮质量的变化。但未对陈皮陈化过程中微生物种群结构变化与气候条件的关联进行研究。文献报道，不同地区特定的气候环境影响微生物群落结构的变化[7-8]。课题组前期也通过对四川地区自然陈化1年川陈皮表面微生物与环境因子相关性研究，初步证实陈皮表面微生物群落变化与环境因子密切相关[9]。
高通量测序技术具有测序深度高、利于鉴定低丰度群落物种以及费用低的特点[10-11]。不仅能获得样品中微生物群落组成，还能将其相对含量进行数字化，已在白酒[12-13]、
醋[14]、腐乳[15]、发酵饮料[16]等的微生物多样性分析得到广泛应用，为环境样品微生物群落结构多样性分析提供了一种可靠有效的方法[17]。因此，本实验采用高通量测序技术首次对温带季风性气候（北京）、亚热带季风湿润气候（四川）、亚热带季风气候（广东）、亚热带高原季风气候（云南）自然陈化1年广陈皮表面细菌16S rDNA V4区、真菌内转录间隔区（internal transcribed spacer，ITS）1区基因进行测序，并通过α和β多样性分析对细菌、真菌群落结构组成进行分析，比较不同地区广陈皮中微生物群落结构的差异，为探究适合陈皮陈化的气候条件，同时为深入研究陈皮陈化机制，阐释陈皮“陈久者良”科学内涵提供新思路。
1 材料与方法
1.1 材料与试剂
1.1.1 材料
为保证样品的代表性和均一性，2020年1月初在广东省江门市新会区三江镇新马单村采集当年产的橘皮，样品均为芸香科植物橘的栽培变种茶枝柑（Citrus Reticulate Chachi）的干燥成熟果皮。以麻布袋将样品包装后分别将其放置于北京市丰台区西四环中路100号解放军302医院、四川省成都市温江区柳台大道1166号成都中医药大学、广东省江门市蓬江区堤东路208号首层、云南省昆明市呈贡区雨花路1076号云南中医药大学的阴凉干燥通风的仓库内，自然陈化1年后取样，样本编号分别为BJ、SC、GD、YN，每个样本取5个生物学重复。
1.1.2 试剂
Gene JET胶回收试剂盒：美国Thermo Scientific公司；Phusion
High-Fidelity 聚合酶链式反应（polymerase chain reaction，PCR）Master Mix with GC Buffer、高效高保真酶：美国New England Biolabs公司；Illumina TruSeq 脱氧核糖核酸（deoxyribonucleic acid，DNA）PCR-Free Library Preparation Kit：美国Illumina公司。
1.2 仪器与设备
SW-CJ-2FD型洁净工作台：苏州安泰空气技术有限公司；DZX-50KBS型立式高压灭菌锅：上海申安医疗器械厂；Bio-rad T100梯度PCR仪：美国Bio-Rad公司；NovaSeq 6000测序系统：美国Illumina公司。
1.3 方法
1.3.1 样品总DNA提取和PCR扩增
参照文献[18]提取广陈皮样品的总DNA，以其为模板，采用带有Barcode序列的引物对515F（5'-GTGCCAGCMGCCGCGGTAA-3'）/806R（5'-GGACTACHVGGGTWTCTAAT-3'）和ITS1-1F-F（5'-CTTGGTCATTTAGAGGAAGTAA-3'）/ITS1-1F-R（5'-GCTGCGTTCTTCATCGATGC-3'）分别对细菌16S rDNA V4区、真菌ITS1-1F区基因进行PCR扩增[19]。PCR扩增体系：Phusion Master Mix（2×）15 μL，正反向引物（2 μmol/L）各1.5 μL，gDNA（1 mg/L）10 μL，双蒸水（ddH2O）2 μL。细菌PCR扩增条件：98 ℃预变性1 min；98 ℃变性10 s，50 ℃退火30 s，72 ℃延伸30 s，30个循环；72 ℃再延伸5 min。真菌PCR扩增条件：98 ℃预变性30 s；98 ℃变性10 s，56 ℃退火30 s，72 ℃延伸15 s，35个循环；72 ℃再延伸7 min。利用2%琼脂糖凝胶电泳检测PCR扩增产物，然后根据Thermo Scientific公司GeneJET胶回收试剂盒说明书切胶回收PCR扩增产物。
1.3.2 文库构建和上机测序
使用Illumina公司的TruSeqDNAPCR-FreeLibraryPreparation Kit进行文库构建，构建好的文库经Qubit定量和文库检测合格后，采用NovaSeq 6000上机测序。
1.3.3 高通量测序数据处理与分析
去除Barcode序列和引物序列，使用FLASH（V1.2.7）软件进行拼接，得到原始序列数据。原始序列数据进行严格的过滤处理[20]后参照QIIME（V1.9.1）[21]的序列质量控制流程，与物种注释数据库比对剔除嵌合体，得到有效序列。以97%的一致性进行操作分类单元（operationaltaxonomic units，OTUs）聚类和物种分类分析，在聚类过程中去除嵌合体，得到OTUs的代表序列。根据OTUs聚类结果，细菌用Mothur方法与SSU rRNA数据库[22]对OTUs代表序列进行物种注释分析（设定阈值为0.8～1.0）；真菌通过QIIME软件中的BLAST方法与Unit数据库进行物种注释分析，得到对应的物种信息和基于物种的丰度分布情况。同时，使用QIIME软件计算Alpha多样性（Chao1、Shannon、Simpson、ACE、Coverage指数），采用R软件（V2.15.3）分别进行Alpha多样性指数组间差异分析、Venn图与非度量多维尺度（nonmetric multidimensional scaling，NMDS）图绘制，以得到样品内物种丰富度和均匀度信息、不同样品或分组间的共有和特有OTUs信息等。
2 结果与分析
2.1 测序数据统计分析
不同地区自然陈化1年的广陈皮中细菌、真菌的有效序列及OTUs数分别见表1和表2。由表1可知，不同地区广陈皮样品中，细菌共获得1 279 520条有效序列，平均碱基长度为252～253 bp，共得到可分类的OTUs数目约40 267个。由表2可知，真菌共获得1 284 394条有效序列，平均碱基长度为226～233 bp。共得到可分类的OTUs数目约6 163个。
表1 不同地区广陈皮中细菌的有效序列及OTUs数目
Table 1 Effective sequences and OUTs numbers of bacteria in Chachi samples from different regions
表2 不同地区广陈皮中真菌的有效序列及OTU数目
Table 2 Effective sequences and OUT numbers of fungi in Chachi samples from different regions
2.2 α-多样性分析
Coverage指数代表测序对物种的覆盖度；Shannon指数和Simpson指数用于衡量微生物群落的多样性，Chao1指数和ACE指数用于估算样品中微生物的丰富度[23]。不同地区自然陈化1年广陈皮表面细菌及真菌的α-多样性分析结果分别见表3和表4。
表3 不同地区广陈皮中细菌群落α-多样性分析结果
Table 3 Analysis results of α-diversity of bacteria community of Chachi samples from different regions
表4 不同地区广陈皮中真菌群落α-多样性分析结果
Table 4 Analysis results of α-diversity of fungi community of Chachi samples from different regions
由表3和表4可知，所有样本的Coverage指数均达到0.99，表明测序结果能完整反映陈皮表面微生物群落组成信息。就细菌而言，不同地区自然陈化1年广陈皮样品Shannon指数、Simpson指数无显著差异（P＞0.05）；云南地区样品的Chao1指数及ACE指数高于其他地区样品，分别为2 088.45、2 157，但不存在显著差异性（P＞0.05）；就真菌而言，不同地区自然陈化1年广陈皮样品Shannon指数、Simpson指数无显著差异（P＞0.05）；广东地区样品的Chao1指数及ACE指数高于其他地区样品，分别为326.59、338.77，但不存在显著差异性（P＞0.05）。说明在不同地区自然陈化1年的广陈皮表面微生物的多样性及丰富度无显著性差异，这可能与陈化时间较短有关。
2.3 微生物群落组成分析
基于门水平，不同地区广陈皮样品中细菌及真菌群落结构见图1。
图1 基于门水平不同地区广陈皮样品中细菌（a）及真菌（b）群落结构
Fig.1 Community structure of bacteria (a) and fungi (b) in Chachi samples from different regions based on phylum level
Others为10个门之外的其他所有门的相对丰度之和。
由图1a可知，在门水平上，不同地区自然陈化1年的所有样本中细菌主要包括变形菌门（Proteobacteria）、蓝藻门（Cyanobacteria）、厚壁菌门（Firmicutes）、放线菌门（Actinobacteria）、Acidobacteriota等。不同地区自然陈化1年广陈皮样品中细菌群落结构组成基本相同，但相对丰度差异较大。其中Proteobacteria、Cyanobacteria在4个地区广陈皮样品中相对丰度较高，Proteobacteria在北京、四川、广东和云南地区的广陈皮样品中的相对丰度分别为47.62%、49.48%、48.69%、41.94%，Cyanobacteria的相对丰度分别为37.99%、31.44%、36.44%和43.45%。
由图1b可知，子囊菌门（Ascomycota）和担子菌门（Basidiomycota）在不同地区自然陈化1年广陈皮表面真菌中占主导地位。不同地区自然陈化1年广陈皮样品中真菌群落结构组成基本相同，但相对丰度差异较大。Ascomycota在云南地区广陈皮样品中相对丰度最高，为83.23%，在四川地区相对丰度最低，为54%。Basidiomycota在北京、四川、广东和云南地区广陈皮样品中的相对丰度分别为18.83%、37.34%、22.86%和14.12%。此外，罗兹菌门（Rozellomycota）仅在四川地区存在，相对丰度为0.96%。
基于属水平，不同地区自然陈化1年广陈皮样品中细菌及真菌群落结构见图2。
图2 基于属水平不同地区广陈皮样品中细菌（a）及真菌（b）群落结构
Fig.2 Community structure of bacteria (a) and fungi (b) in Chachi samples from different regions based on genus level
Others为30个属之外的其他所有属的相对丰度之和。
由图2a可知，在属水平上，不同地区自然陈化1年广陈皮样本中细菌主要包括unidentified Chloroplast、unidentified Mitochondria、甲基杆菌属（Methylobacterium）、鞘氨醇单胞菌属（Sphingomonas）、中慢生根瘤菌属（Mesorhizobium）、芽孢杆菌属（Bacillus）等。其中unidentified Chloroplast和unidentified Mitochondria在4个地区广陈皮样品中相对丰度较高，在云南地区两者丰度之和最大，为62.27%，广东地区两者丰度之和最小，为42.18%。说明不同地区广陈皮中细菌多样性有待进一步挖掘和发现。Methylobacterium在四川地区广陈皮样品中相对丰度最大，为15.59%，在广东地区广陈皮样品中相对丰度最小，为4.99%，在北京和云南地区广陈皮样品中相对丰度分别8.97%和6.81%。Sphingomonas、Mesorhizobium、Bacillus在四川地区广陈皮样品中相对丰度最高，分别为2.86%、2.04%和1.22%。值得注意的是，研究报道Bacillus、Lactococcus可加速陈皮中挥发油的生物转化[6]。
由图2b可知，不同地区自然陈化1年广陈皮样本中真菌属主要包括平脐疣孢属（Zasmidium）、枝孢霉属（Cladospo rium）、短梗霉属（Aureobasidium）、酵母菌属（Symmetrospora）、unidentified Cystobasidiomycetes sp.、镶刀菌属（Fusarium）、红酵母属（Rhodotorula）、担孢酵母属（Erythrobasidium）、炭疽菌属（Colletotrichum）、赤霉菌属（Gibberella）、拟棘壳孢属（Pyrenochaetopsis）和曲霉属（Aspergillus）等。其中Zasmidium在云南地区相对丰度最大，为45.77%；Cladosporium、Aureobasidium和Colletotrichum在北京地区相对丰度最大，分别为16.80%、27.27%、4.18%；Symmetrospora、Aspergillus在四川地区相对丰度最大，分别为13.00%、1.16%。unidentified Cystobasidiomycetes sp.在广东地区相对丰度最大，为7.99%。Rhodotorula仅在北京、四川地区中存在，且在四川地区相对丰度最大，为8.77%。外瓶霉属（Exophiala）仅在四川地区存在，其相对丰度为1.32%。值得注意的是，课题组前期报道陈皮表面黑曲霉（Aspergillus niger）是引起黄酮类和挥发油含量与类型变化的重要因素[5]。但不同地区自然陈化1年广陈皮样品表面Aspergillus检出较少，这可能与陈皮陈化时间较短相关。
2.4 不同地区广陈皮样品微生物群落比较
为了说明不同地区自然陈化1年广陈皮样本物种组成的相似性情况[24]，对不同地区广陈皮中细菌及真菌的OTU进行Venn图绘制，结果见图3。
图3 基于OTU水平不同地区广陈皮中细菌（a）及真菌（b）Venn图
Fig.3 Venn diagrams of bacteria (a) and fungi (b) of Chachi samples from different regions based on OTU level
由图3可知，在OTUs组成上，四川、北京、云南、广东地区自然陈化1年的广陈皮表面细菌群落分别得到585、748、794和697个特有OTU，共有OTU数为1 826个；真菌群落分别得到282、154、71和159个特有OTU，共有的OUT数为272个，且不同地区自然陈化1年广陈皮表面，真菌种类明显少于细菌种类。
非度量多维标度（nonmetric multidimensional scaling，NMDS）法是一种将多维空间的研究对象（样本）简化到低维空间进行定位、分析和归类，同时又保留对象间原始关系的数据分析方法[25]。对不同样本间的差异程度，通过点与点间的距离进行体现，距离越近，则样本组成越相似。基于属水平，不同地区广陈皮中细菌及真菌的NMDS分析结果见图4。
图4 基于属水平不同地区广陈皮中细菌（a）及真菌（b）的NMDS结果
Fig.4 NMDS results of bacteria (a) and fungi (b) of Chachi samples from different regions based on genus level
由图4可知，不同地区所有广陈皮样品距离较近，说明不同地区自然陈化1年广陈皮样品中细菌及真菌群落组成相似。
为进一步比较不同地区自然陈化1年广陈皮样品在属水平上群落组成的相似性和差异性，根据所有样品在属水平的物种注释及丰度信息，选取相对丰度排名前35的属进行相似性聚类，得到聚类分析热图，结果见图5。一个样品在某个分类上的Z值为样品在该分类上的相对丰度和所有样品在该分类的平均相对丰度的差与所有样品在该分类上的标准差的比值，将聚类分析热图中Z＞1的物种定义为优势菌属。
图5 基于属水平细菌（a）及真菌（b）群落物种相对丰度聚类热图
Fig.5 Cluster heatmap of relative abundance of bacteria (a) and fungi (b) community species based on genus level
由图5a可知，北京地区自然陈化1年广陈皮表面优势细菌属为Methylocella、苔藓杆菌属（Bryobacter）和短小杆菌属（Curtobacterium）；四川地区自然陈化1年广陈皮表面优势细菌属为Methylobacterium、Sphingomonas、Bacillus、类芽孢杆菌属（Paenibacillus）、副球菌属（Paracoccus）、膨胀芽胞杆菌属（Tumebacillus）、溶菌杆菌（Lysinibacillus）、Tundrisphaera、乳球菌属（Lactococcus）、嗜酸菌属（Acidiphilium）；广东地区自然陈化1年广陈皮表面优势细菌属为Rubellimicrobium、乳杆菌属（Lactobacillus）、放线孢菌属（Actinomycetospora）；云南地区自然陈化1年广陈皮表面优势细菌属为unidentified Chloroplast、罗尔斯顿菌属（Ralstonia）、Glutamicibacter、马赛菌属（Massilia）、罗河杆菌属（Rhodanobacter）。
由图5b可知，北京地区自然陈化1年广陈皮表面优势真菌属为Aureobasidium、Colletotrichum、Pyrenochaetopsis、Roussoella、Sarcopodium、Wallemia、掷孢酵母属（Sporobolomyces）；四川地区自然陈化1年广陈皮表面优势真菌属为Symmetrospora、Rhodotorula、Gibberella、unidentified Pucciniomycetes sp.、红冬孢酵母属（Rhodosporidiobolus）、Aspergillus、Exophiala、Pseudoramichloridium、Naganishia；广东地区自然陈化1年广陈皮表面优势真菌属为unidentified Cystobasidiomycetes sp.、Pseudocercospora、毛壳菌属（Chaetomium）、耐干霉菌属（Xeromyces）、Cystobasidium、汉纳酵母属（Hannaella）、Ectophoma、隐球酵母属（Cryptococcus）、链格孢属（Alternaria）、弯孢属（Curvularia）；云南地区自然陈化1年广陈皮表面优势真菌属为Zasmidium、拟暗球腔菌属（Phaeosphaeriopsis）、小粉孢属（Microidium）和多臂菌属（Trichomerium）。文献报道，陈皮表面Methylobacterium、Sphingomonas、Lactococcus、Bacillus等细菌可加速挥发油的生物转化[6，26]，Aspergillus是引起黄酮类和挥发油含量与类型变化的重要因素。Methylobacterium、Sphingomonas、Lactococcus、Bacillus、Aspergillus为四川地区的优势菌属，说明四川地区的气候条件比较适合广陈皮的自然陈化。
3 结论
本研究通过高通量测序技术分析了不同地区自然陈化1年广陈皮样品表面微生物群落多样性变化及优势菌群分布。分析结果表明，不同地区自然陈化1年广陈皮样品丰富度和多样性不存在显著差异（P＞0.05），微生物群落结构相似，但相对丰度存在差异，优势菌属不同。Methylocella、Bryobacter和Curtobacterium为北京地区优势细菌，Methylobacterium、Sphingomonas、Bacillus等为四川地区优势细菌，Rubellimicrobium、Lactobacillus、Actinomycetospora为广东地区优势细菌，unidentified Chloroplast、Ralstonia、Glutamicibacter、Massilia、Rhodanobacter为云南地区优势细菌。Aureobasidium、Colletotrichum、Pyrenochaetopsis等为北京地区优势真菌，Symmetrospora、Rhodotorula、Gibberella等为四川地区优势真菌，Pseudocercospora、Chaetomium、Xeromyces等为广东地区优势真菌，Zasmidium、Phaeosphaeriopsis、Microidium和Trichomerium为云南地区优势真菌。四川地区优势细菌属Methylobacterium、Sphingomonas、Lactococcus、Bacillus可加速挥发油的生物转化，优势真菌属Aspergillus是引起黄酮类和挥发油含量与类型变化重要因素。说明四川地区的气候条件比较适合广陈皮的自然陈化。
[1]国家药典委员会.中华人民共和国药典[M].北京：中国医药科技出版社，2020：199-200.
[2]傅曼琴，肖更生，吴继军，等.广陈皮促消化功能物质基础的研究[J].中国食品学报，2018，18（1）：56-64.
[3]王智磊，张鑫，刘素娟，等.陈皮“陈久者良”历史沿革和研究现状[J].中华中医药学刊，2017，35（10）：2580-2584.
[4]杨放晴，何丽英，杨丹，等.不同陈化时间广陈皮中黄酮类成分的UPLCQ-Orbitrap HRMS 分析[J].食品与发酵工业，2021，47（15）：267-275.
[5]王福，张鑫，卢俊宇，等.陈皮“陈久者良”之黄酮类成分增加原因探究[J].中国中药杂志，2015，40（24）：4890-4896.
[6]陈聪聪.广陈皮陈化过程中微生物群落多样性解析及代谢物成分变化分析[D].广州：华南理工大学，2017.
[7]程伟，李娜，张杰，等.江淮地区自然环境对柔和型白酒品质的影响[J].酿酒科技，2019（6）：59-64.
[8]王美溪，刘珂艺，邢亚娟.气候变化背景下土壤微生物与植物物种多样性关联分析[J].中国农学通报，2018，34（20）：111-117.
[9]刘素娟.陈皮陈化过程“环境因子-真菌-药效物质基础”相关性研究[D].成都：成都中医药大学，2018.
[10] DEGNAN P H,OCHMAN H.Illumina-based analysis of microbial community diversity[J].ISME J,2012,6(1):183-194.
[11] CAPORASO J G,LAUBER C L,WALTERS W A,et al.Ultra-highthroughput microbial community analysis on the Illumina HiSeq and MiSeq platforms[J].ISME J,2012,6(8):1621-1624.
[12]孙利林，李立郎，胡萍，等.酱香型白酒第四轮次酒酿造过程中细菌多样性分析[J].中国酿造，2020，39（5）：35-39.
[13]杨磊，程铁辕，吴树坤，等.沉香型白酒酒醅理化特性及细菌群落结构研究[J].中国酿造，2020，39（2）：109-114.
[14]王俊奇，黄卫红，李双彤，等.永春老醋不同生产阶段细菌和真菌多样性动态变化特征分析[J].食品与发酵工业，2021，47（2）：38-44.
[15]WAN H F,LIU T,SU C W,et al.Evaluation of bacterial and fungal communities during the fermentation of Baixi sufu,a traditional spicy fermented bean curd[J].J Sci Food Agr,2020,100(4):1448-1457.
[16]PEREZ-CATALUNA A,ELIZAQUIVEL P,CARRASCO P,et al.Diversity and dynamics of lactic acid bacteria in Atole agrio,a traditional maize-based fermented beverage from South-Eastern Mexico,analysed by high throughput sequencing and culturing[J].Anton Leeuw,2017,111(3):385-399.
[17]TRINGE S G,MERING C V,KOBAYASHI A,et al.Comparative metagenomics of microbial communities[J].Science,2005,308(5721):554-557.
[18]张鑫，刘素娟，王智磊，等.橘皮表面真菌群落结构多样性分析[J].微生物学通报，2017，44（5）：1089-1098.
[19]WANG X Y,ZHU X G,BI Y L,et al.Dynamics of microbial community and changes of metabolites during production of type I sourdough steamed bread made by retarded sponge-dough method[J].Food Chem,2020,330:127316.
[20]BOKULICH N A,SUBRAMANIAN S,FAITH J J,et al.Quality-filtering vastly improves diversity estimates from Illumina amplicon sequencing[J].Nat Methods,2013,10(1):57-59.
[21]CAPORASO J G,KUCZYNSKI J,STOMBAUGH J,et al.QIIME allows analysis of high-throughput community sequencing data[J].Nat Methods,2010,7(5):335-336.
[22] EDGAR R C.UPARSE:highly accurate OTU sequences from microbial amplicon reads[J].Nat Methods,2013,10(10):996-998.
[23]武亚婷，杜木英，何欢欢，等.基于高通量测序技术分析新疆不同地区自然发酵辣椒酱微生物群落多样性[J].食品与发酵工业，2019，45（21）：221-228.
[24] SHADE A,HANDELSMAN J.Beyond the Venn diagram:the hunt for a core microbiome[J].Environ Microbiol,2012,14(1):4-12.
[25]李靖宇，张琇.腾格里沙漠不同生物土壤结皮微生物多样性分析[J].生态科学，2017，36（3）：36-42.
[26]何静.广陈皮细菌群落和植物化学物独特性及其相关性分析[D].广州：华南理工大学，2019.