基于16S rDNA扩增子测序分析灵芝连作覆土细菌群落的变化

袁 源，李 琳，黄海辰，刘国辉，谢福泉，傅俊生，吴小平

（1福建农林大学生命科学学院，福州350002；2福建农林大学菌物研究中心，福州350002；3仙芝科技（福建）股份有限公司，福建南平353400；4福建省食用菌技术推广总站，福州350003）

0 引言

灵芝是一种天然的药用真菌，目前已报道400多个分类单元的灵芝科，国内普遍栽培的灵芝通常是指‘赤芝’，崔宝凯等[1]认为‘赤芝’的拉丁学名为Ganoderma lingzhi Sheng H.Wu,Y.Cao&Y.C.Dai，而国内外误用Ganoderma lucidum(Curtis)P.Karst.作为‘赤芝’的拉丁学名已经100余年。其富含多糖，三萜和灵芝酸等活性成分，具有显著的抗肿瘤、抗氧化和降血脂等活性，具有良好的药用价值和保健功能，深受众多消费者的喜爱[2]。为满足消费者对灵芝日益增长的市场需求，人工栽培目前已成为灵芝供应的主要来源。到目前为止，中国人工栽培灵芝的历史已超过了60年，灵芝段木覆土栽培是一种重要的灵芝的栽培模式，段木为灵芝的生长提供丰富的营养成分，覆土能稳定料中的湿度，维持良好的通气性，该模式生产上具有显著的低污染，易控制，产量高，并且子实体营养丰富和药效成分全面等特点，但连作2～3年后灵芝覆土依然存在常见的连作障碍问题，往往容易引发灵芝子实体畸形，病害严重以及产量急剧下降等现象，将会导致大面积的灵芝基地不得不停用，制约了灵芝产业的良性发展，也不利于土地资源的循环利用[3-5]。

研究发现，连作障碍的主要原因包括微生物群落的改变、化感作用和自毒作用等[6-7]。同时，多项报道显示灵芝连作障碍与覆土中微生物群落的变化具有较强的相关性，研究微生物群落变化有助于阐明灵芝连作障碍的机制[4-5]。微生物是土壤生态系统中关键的构成要素，也是一种重要的土壤质量评价指标[5]。一方面，在灵芝段木覆土栽培过程中，段木可以更新，大棚内的温度和湿度可进行人工调节，但是灵芝的自毒作用、化感作用以及土壤理化性质的改变等均会导致覆土内微生物群落发生变化，这种菌群的改变反映了土壤质量的变化；另一方面，当土壤的质量发生变化，灵芝覆土中微生物的群落的丰度和数量又会相应的改变，从而对覆土中微生物群落进行研究就显得尤为重要[8]。当前对灵芝连作障碍的相关研究报道主要集中在覆土中的真菌群落，如Kang等[9]基于形态学和系统发育树分析，发现木霉属与灵芝连作过程中产生的黄腐病有关；马红梅等[10]通过平板对峙培养法，发现子囊菌门的木霉属、青霉属和链孢霉属等真菌是抑制灵芝生长主要的病原菌，而对灵芝连作覆土中细菌群落系统变化的研究却鲜有报道。因此，文章对灵芝覆土中的细菌群落进行深入的探究。

近年来，基因测序在有关土壤微生物方面的研究中发挥着十分重要的作用，能对土壤中不同类群的微生物进行全面系统的分析[11]。16S rDNA扩增子测序技术是研究土壤微生物群落多样性变化的基本手段之一，常用于反映样本中细菌类群的种间差异[12]。该研究基于16S rDNA扩增子技术，对不同栽培年限灵芝覆土中细菌群落进行系统的分析，以期为灵芝连作障碍机制的研究奠定基础。

1 材料和方法

1.1 材料与试剂

试验菌株为‘赤芝119’（‘沪农灵芝1号’），由上海市农业科学院食用菌研究所进行新品种认定，审定编号为沪农品认食用菌(2009)第003号。土壤DNA提取试剂盒（OMEGA公司）。

1.2 仪器与设备

22R冷冻离心机（Microfuge）；DK-8D型三温三控恒温水浴锅（上海博迅实业有限公司）；WH-861型旋涡震荡仪（华利达公司）；PCR仪（杭州朗基科学仪器有限公司）；EPS300型电泳仪（上海天能科技有限公司）等。

1.3 试验基地概况

试验地位于仙芝科技（福建）股份有限公司有机灵芝博览园（福建省南平市浦城县，28°01’N，118°32’E），地处福建省最北端，属于亚热带季风气候，四季分明，全年平均气温17.4℃[13]。基地一共设有41个标准大棚，试验土壤均源自基地同一区域内品质优良的沙壤土，pH呈弱酸型(5～6)，适于灵芝生长出菇，现已成功培植赤芝等14个优良品种。

1.4 试验方法

采用灵芝段木覆土栽培模式对‘赤芝119’进行栽培，将同一批生长状态一致的‘赤芝119’脱袋后，2018年3月中旬待其菌丝长满段木，及时移入大棚开始栽培，分别埋入通气较好的栽培1年（30号大棚）、连作1年（38号大棚）、2年（6号大棚）和3年（17号大棚）的覆土中，使覆土能覆盖住段木顶部2～3 cm，进一步进行出芝管理。在此过程中，4个大棚内的光照、温度、湿度、通风以及操作过程等均控制一致，并且每年采集一次灵芝孢子粉和采收一次灵芝。

根据浦城县有机灵芝博览园的灵芝生产安排，8月至9月为灵芝采收阶段，9月底灵芝采收结束。为了更好的反映灵芝连作覆土中细菌群落的变化，灵芝栽培1年、连作1年、2年和3年的覆土以及基地内未栽培过灵芝的邻近沙壤土，于2018年9月25日灵芝采收期间取自有机灵芝博览园。试验设计如下：在4个实验组中，每个分组面积为150 m²的大棚分别均匀栽培100棒灵芝。以邻近沙壤土为对照组，将栽培1年、连作1年、2年、3年的灵芝覆土采用S形取样法，按照半径20～25 cm、深度0～5 cm的范围使用铁铲进行取样，每个分组取5个点的混合土样，分别记为GL0、GL1、GL2、GL3和GL4组。每个分组混合土样大约200 g，采集后分别装入自封袋并且编号，置于冰盒内带回实验室，于-20℃冰箱保存。

1.5 土壤总DNA的提取

准确称量0.2 g土壤样品，加入含500 mg玻璃珠的15 mL离心管中，参照OMEGA土壤DNA提取试剂盒(D5625-02)说明书，对灵芝土样的DNA进行提取，通过1%琼脂糖凝胶电泳和核酸微量测定仪对所获得的DNA进行质检和定量，每个分组设定3个生物学重复，一共15个土壤样本（FS_1、FS_2、FS_3、ONE_S_1、ONE_S_2、ONE_S_3、TWO_S_1、TWO_S_2、TWO_S_3、THREE_S_1、THREE_S_2、THREE_S_3、FOUR_S_1、FOUR_S_2和FOUR_S_3）。

1.6 PCR扩增及测序

参考赵力的方法，采用特异性引物338F(5'-ACTCCTACGGGAGGCAGCAG-3')和 806R(5'-GGACTACHVGGGTWTCTAAT-3')对16S rDNA可变保守区(V3+V4)进行PCR扩增[14]，PCR反应体系和反应条件如实验室已有报道所示[15]，对PCR扩增产物进行检测和回收。在Qbit荧光定量系统上采用PicoGreen双链DNA荧光定量测定试剂盒对纯化后的PCR产物文库进行定量。将合格的文库（浓度＞2 nmol/L）梯度稀释后，按一定要求进行比例混合和变性处理，变性为单链后可进行上机测序；于MiSeq测序仪上使用MiSeq Reagent Kit试剂盒进行2×300 bp的双端测序。

1.7 生物信息学分析流程

将灵芝样品送至杭州联川生物技术股份有限公司(LC-Bio)在Illumina MiSeq平台上进行测序。依据双端序列之间的重叠关系，采用FLASH(v1.2.8)软件将序列拼接成长的tag，去除建库引入的barcode和引物序列，采用Vsearch(v2.3.4)过滤嵌合体。预处理之后的clean data使用Vsearch将序列相似度大于97%的有效序列定为一个OUT(operational taxonomic unit)，选出位于几何中心的代表序列。多样性通过QIIME(v1.8.0)分析。通过blast进行序列比对，将OTU代表序列与RDP(ribosomal database project)以及UNITE数据库对每个代表性序列进行物种注释。其余图片均是通用R包(V3.2.5)实现。并将序列数据和基因注释信息上传到NCBI网站数据库。

1.8 数据统计学分析

采用SPSS17.0统计软件对测序数据进行分析处理。组间比较采用t检验，以P＜0.05表示差异显著，P＜0.01表示差异极显著。

2 结果与分析

2.1 OTU聚类分析

对原始数据进行质控、双端拼接、去除低质量以及嵌合体过滤后，对高质量数据进行统计，一共得到有效序列477556条。其中，400～500 bp序列约占总序列的97.05%，小于200 bp，200～300 bp和300～40 bp的序列分别占2.20%，0.64%和0.11%。将上述15个土壤样本的477556条有效序列按照碱基相似度大于97%聚类成为不同的OTUs，得到每个cluster的代表序列则为一个OTU[16]。如图1所示，Venn图可直观表示5个分组之间共有和特有的OTU序列和丰度，不同颜色的椭圆表示不同的分组[15]。5个分组经聚类后一共获得9129个OTU，其中，GL0、GL1、GL2、GL3和GL4组OTU的数目分别为4689、5129、6506、5186和4281个。连作1年的覆土中OTU数目达到最大值，即6506个；连作3年的灵芝覆土中OTU数目最少，为4281个。表明随着灵芝连作年限的增加，其覆土中的细菌类群数量呈先增后减的趋势。同时，5个分组共有的OTU数目为851个，特有的OTU数目从GL0至GL4组依次为214、101、121、510和311个，灵芝连作2年和3年覆土细菌中特有的OTU数目大幅增加，表明特有细菌类群增加。

图1 覆土细菌群落OTU聚类维恩图

2.2 覆土细菌物种的Alpha多样性

Alpha多样性分析常用于衡量群落中物种丰富度、多样性以及测序深度，通常包括chao1、observed species、Shannon、Simpson和 goods_coverage5个指数。chao1指数和observed species指数反映样品中物种的丰富度；Shannon和Simpson指数分别反映样品中物种的多样性以及不同种类的物种丰度的均匀性；goods_coverage指数是指微生物覆盖率，用于反映测序的深度，其数值越接近1，则样本中没有被测出的新物种概率越低。如图2A和2B所示，通过对5个分组的chao1指数和observed species指数进行分析，发现灵芝连作1年覆土中细菌物种数目极显著增加，随后开始下降。如图2C和2D的Shannon和Simpson指数所示，灵芝连作前2年覆土中细菌群落多样性水平极显著增加，随后开始减少，且物种丰度的均匀性较好。如图2E和2F所示，5个分组中goods_coverage指数均接近1，说明测序深度已足够，基本覆盖土壤样品中的所有物种。由此表明，灵芝覆土中细菌的丰富度和多样性随着连作年限的推移均呈先增后减的趋势，并且测序深度足够。

图2 覆土细菌群落的Alpha多样性指数的分析

2.3 覆土细菌物种的Beta多样性

Beta多样性在时空尺度上衡量物种组成的变化，以反映不同环境中群落间的物种组成差异[17-18]。如图3A所示，采用非加权平均距离法(UPGMA)对15个样本进行聚类，聚类树展现了样品间的相似度。其中，5个分组一共聚为两个大分枝，GL0和GL1组灵芝覆土距离最近，表明邻近沙壤土与栽培1年灵芝覆土之间细菌群落组成结构相似度最高，而GL2，GL3和GL4组聚为一个分枝并与GL0组的距离逐渐加大，表明连作逐渐降低了灵芝覆土中细菌物种组成的相似度。主坐标分析(PCoA)是一种常用的数据降维方法，可反映各个样本群落组成的相似性或差异性[19]。样本间距离越小，细菌群落组成的相似度越高；反之，细菌群落组成的差异性越大。如图3B所示，与UPGMA聚类分析结果类似，随着灵芝连作年限的增加，覆土中细菌组成的差异逐渐加大。此外，采用箱线图对组内和组间Unifrac距离的差异进行比较如图4A和4B所示，组间差异距离大于组内差异距离，表明5个分组多组比较的组间具有更显著的差异。GL0_GL4组、GL1_GL4组和GL0_GL3组组间差异较大；GL0_GL1组组间差异较小，且GL2_GL4组存在异常值，其可能是因为两组不同样本之间比较的次数较多，若有个别数据与平均值差别较大则容易出现异常值，同样表明连作逐步加大了邻近沙壤土与连作覆土间的细菌组成差异。

图3 覆土细菌群落的UPGMA聚类和PCoA图

图4 覆土细菌群落多组比较箱线图

2.4 覆土细菌物种组成分析

2.4.1 门水平组成的分析 从15个灵芝土样中一共获得1050种细菌，分别隶属于28门、89纲、169目、287科、624属。如图5A所示，变形菌门(Proteobacteria)、酸杆菌门(Acidobacteria)、放线菌门(Actinobacteria)和绿弯菌门(Chloroflex)的相对丰度约占细菌总序列的80%，表明5个分组中主要的优势菌门为变形菌门、酸杆菌门、放线菌门和绿弯菌门。如图5B的热图所示，随着连作年限的增加，变形菌门、放线菌门和拟杆菌门(Bacteroidetes)的相对丰度显著增 加 ；绿 弯 菌 门 (Chloroflex)、Bacteria_unclassified、candidate_division_WPS-2、糖化念珠菌(Candidatus_Saccharibacteria)和螺旋体门(Spirochaetes)的相对丰度显著减少；而芽单胞菌门(Gemmatimonadetes)、厚壁菌门(Firmicutes)、硝化螺旋菌门(Nitrospirae)和疣微菌门(Verrucomicrobia)的相对丰度呈现先增加后减少的趋势。其中，厚壁菌门和硝化螺旋菌门的相对丰度在GL2组开始下降，表明连作会影响灵芝覆土中优势菌群的分布情况。

图5 细菌群落门水平的相对丰富度

2.4.2 属水平组成的分析 在属水平上一共获得624个属，其相对丰度排名前20的优势菌属如图6A所示。从图中可以看出，在属水平上检出优势细菌的相对丰度不断减少，其他细菌相对丰度逐年增加。如图6B所示，在灵芝栽培0、1、2、3和4年覆土中，变形菌门的溶杆菌属(Lysobacter0.01%，0.11%，1.05%，1.14%，3.58%)与酸杆菌门的Gp6和Gp16的相对丰度随连作年限的增加而不断上升；变形菌门鞘氨醇单胞菌属(Sphingomonas 5.33%，1.12%，3.09%，1.33%，2.71%)、Anaeromyxobacter（2.23% ，2.36% ，1.33% ，0.71% ，0.15%）、慢生根瘤菌属(Bradyrhizobium 1.63%，1.54%，0.88%，0.75%，0.53%)和绿弯菌门脱氯菌属(Dehalococcoides 8.89%，3.70%，2.02%，0.55%，0.27%)、Ktedonobacter（5.48%，3.66%，2.95%，0.09%，0.14%）的相对丰度不断减少；变形菌门地杆菌属(Geobacter)、假单胞菌属(Pseudomonas)、红游动菌属(Rhodoplanes)，酸 杆 菌 门 Gp1、Gp7、Edaphobacter、Gp3，芽单胞菌门芽单胞菌属(Gemmatimonas)以及硝化螺旋菌门硝化螺旋菌属(Nitrospira)相对丰度随灵芝连作年限的增加呈现先增加后减少的趋势。其中，地杆菌属(Geobacter)、假单胞菌属和Gp3的相对丰度在GL2组开始下降，表明长期连作会引起灵芝覆土中大量优势细菌群落相对丰度的减少。

图6 细菌群落属水平的相对丰富度

3 讨论

随着国内健康产业的快速发展，灵芝产品因其较高的营养保健价值备受消费者喜爱，但人工栽培灵芝深受连作障碍现象的制约。以往研究发现，土壤微生物群落与连作障碍关系密切，常常表现为由“细菌型”转向“真菌型”[8]。为了探索灵芝连作障碍的形成机制和防治措施，有人发现子囊菌门的木霉属和青霉属等真菌是灵芝主要的病原菌，并分析灵芝覆土微生态系统，提出对灵芝进行轮作、仿野生栽培、有益菌富集以及有机质使用的栽培方式[8,10]。吴晓明等[4]发现液氨熏蒸能有效降低灵芝覆土中的绿色木霉菌，不仅显著提高了灵芝的生长速度，而且产量大约增长20.95%。在前期研究中，实验室通过对灵芝连作覆土中的真菌群落多样性进行分析，发现真菌群落的相对丰度逐年增加，并且在GL4组中青霉属的相对丰度相对于GL0组约增加了57.92%[15]。李艳春等[20]发现在离茶树大约25～30 cm的范围内间作灵芝，能有效提高土壤中鞘氨醇单胞菌属和伯克氏菌属等有益细菌群落的相对丰度，显著改善土壤的理化性质。然而，本研究发现灵芝连作覆土中鞘氨醇单胞菌属等细菌群落相对丰度随着连作年限的增加而递减。

细菌能促进土地资源有效的循环利用，而连作障碍形成了特殊的土壤环境，会使细菌群落发生重要变化，因此解析灵芝连作覆土中细菌群落的变化情况对防治灵芝连作障碍具有重要指导意义[21]。在灵芝覆土中，变形菌门是5个分组中相对丰度最高的一类优势菌群，约占细菌总量的34.77%～40.45%，并且其相对丰度随着连作年限的增加而显著增加。类似的，陈政在生姜连作土壤中也发现变形菌门的相对丰度位居第一，约占细菌类群总量的49.46%[22]。以往研究表明，变形菌门鞘氨醇单胞菌属能促进芳香族化合物的降解，对于保护环境具有潜在的应用价值[23]。Hwang等[24]通过分析变形菌门Anaeromyxobacter Fw109-5的基因组序列，发现该菌株对地下污染物具有原位生物修复的潜力。而变形菌门慢生根瘤菌属能较好的与降香黄檀等植物共生，具有高效的固氮能力[25]。研究结果同样发现，变形菌门的鞘氨醇单胞菌属、Anaeromyxobacter和慢生根瘤菌属等作为灵芝生长过程中的优势菌属，随连作的年限增加，其相对丰度均呈下降的趋势。然而，变形菌门的溶杆菌属是一种新型的生物防治药物，在植物中具有防治病害和促进生长的功能[26]。其相对丰度在灵芝覆土中却随连作年限的增加而急剧增加，其在GL4组中大约提高了452倍，表明连作可能引发灵芝覆土中多数变形菌门优势细菌属的相对丰度下降。同时，在以往研究中，马红梅等[5]通过平板分离法发现放线菌相对丰度随着连作年限的增加而呈现先增后减的趋势，本研究发现放线菌相对丰度在栽培0、1、2、3和4年灵芝覆土中(7.86%，5.29%，6.72%，11.34%，23.90%)先减后增，其可能原因有待进一步探讨。

酸杆菌门是灵芝覆土中仅次于变形菌门的第二优势菌门，在土壤的循环利用和生态组建方面发挥着十分重要的作用[27]。结果发现酸杆菌门的Gp1、Gp7、Gp3和Edaphobacter的相对丰度呈现出先增后减的现象，表明长期连作会显著降低灵芝覆土中酸杆菌门优势菌属的数量，不利于土壤资源的循环利用。在厌氧条件下，绿弯菌门的脱氯菌属可对高氯代有机污染物进行脱氯并降解成为低氯代的产物，能有效的对污染的土壤进行修复和治理[28]。研究结果发现，灵芝覆土中脱氯菌属的相对丰度逐年递减，在GL4组中大约减少了96.96%。此外，芽单胞菌属和硝化螺旋属也是净化环境常见的有益菌，田美等通过高通量测序技术发现芽单胞菌属在污水净化系统中发挥生物除磷作用[29]。Daims等[30]研究表明硝化螺旋属可将亚硝酸盐氧化成硝酸盐，广泛运用于水体净化以及生态修复领域。然而，本研究中发现，覆土中芽单胞菌属和硝化螺旋属的细菌类群数量由于长期连作后也呈现大幅降低的趋势。

此外，马红梅等[10]基于对峙平板培养法发现灵芝覆土中分离的分支杆菌属、梭杆菌属和产碱杆菌属等细菌对灵芝菌丝体具有较强的化感作用，能显著抑制灵芝菌丝体的生长。本研究进一步基于扩增子测序技术对灵芝覆土中的细菌群落的多样性和结构组成进行系统的探究，其相比对传统的平板分离法具有覆盖面广、通量高、准确和灵敏等特点，能全面系统反映覆土中细菌群落的结构[12]。实验结果发现连作会加大灵芝覆土与邻近沙壤土中细菌群落的物种组成差异，并且细菌群落的丰富度和多样性呈现先增加后减少的趋势。在属水平上，灵芝覆土中的细菌总的相对丰度随着连作年限的增加逐渐减少。其中，鞘氨醇单胞菌属、Anaeromyxobacter、慢生根瘤菌属和脱氯菌属等土壤有益细菌的相对丰度逐年递减，表明灵芝连作障碍与土壤中细菌群落的减少存在明显的关联性。本研究目前尚未对灵芝连作覆土中有益细菌群落的作用机制进行探讨，有待进一步对其进行分离、回接和验证。

4 结论

灵芝连作会逐年加大覆土中细菌群落的物种组成差异，其可能与灵芝栽培0、1、2、3和4年覆土中鞘氨醇单胞菌属(Sphingomonas 5.33%，1.12%，3.09%，1.33%，2.71%)、Anaeromyxobacter(2.23% ，2.36% ，1.33% ，0.71%，0.15%)、慢生根瘤菌属(Bradyrhizobium 1.63%，1.54%，0.88%，0.75%，0.53%)和脱氯菌属(Dehalococcoides 8.89%，3.70%，2.02%，0.55%，0.27%)、Ktedonobacter(5.48%，3.66%，2.95%，0.09%，0.14%)等多种有益细菌相对丰度的减少有关，以期为灵芝连作障碍生态防控奠定理论基础。