华北落叶松林土壤细菌群落结构与表型的环境异质性响应*

赵鹏宇 白 雪 燕平梅 赵晓东 武晓英 柴宝峰

(1.太原师范学院生物系 晋中 030619 2.山西大学黄土高原研究所 晋中 030006)

同一基因型受不同环境因子影响而产生不同的表型，这一现象称为表型可塑性(West-Eberhardetal.，1989)。表型可塑性是生物个体表型变异和生态适应的基础，其可使物种在资源充足条件下达到时空上的“表现最大化”，并且在资源受限、存在胁迫的条件下仍然维持必要的生理功能，以确保对环境异质性的适应度(West-Eberhardetal.，1989；王文晓等，2020)。这种机制使得物种具有更宽的生态幅和更高的耐受性，可以适应更广阔的地理范围和更加多样化的生境条件(王文晓等，2020；赵孟良等，2020)。生物体表型的变化可以体现在行为上、生理上、形态上、生长上或生活史上等，也可以跨越世代。目前，关于生物体表型可塑性的研究常见于大型动植物领域的报道(王小艺等，2015；Ackermannetal.，2015)，如昆虫的翅多型现象(王文晓等，2020)；黑壳白点的螺生活在低盐分地区，而白壳黑纹的螺生活在高盐分地区等(刘彩云等，2020；魏莹等，2020；王文晓等，2020)。因此，开展表型可塑性研究对理解种群的适应性十分重要(朱璟等，2019)；也对保护生态环境、维持物种共存有很重要的生态意义(王小艺等，2015)。

表型可塑性在生物界普遍存在，但在微生物领域的研究却相对较少。微生物具有个体微小、代际周期短、生长极其迅速、种类繁多、存在横向基因转移等独特的生理特性。有研究指出，金黄色葡萄球菌(Staphylococcusaureus)能够依据环境变化在生理表型方面主动做出对环境信号的多样性响应(朱璟等，2019)。即使生存在同一生境的微生物，也可能产生表型的差异(Ackermannetal.，2015)。因此，以微生物群落为研究对象，可以为表型可塑性研究提供更多的理论线索。土壤理化性质不但能够驱动微生物群落的多样性、种类和结构发生巨大变化(赵鹏宇等，2018；Zhaoetal.，2019b；2019c)，而且会对微生物群落功能产生深远影响(赵鹏宇等，2018；Zhaoetal.，2019c；王瑞宏等，2018；段鹏等，2020)。但是，关于土壤理化性质对微生物群落结构和表型变化的影响研究还很少报道。因此，研究土壤理化因子对微生物群落结构和表型的影响，能为理解微生物的环境适应机制提供一定的理论支持。

基于此，为揭示土壤细菌群落的结构与表型的关系及其对土壤异质性的响应机制，本研究做出以下3个假设：1)五台山、芦芽山和关帝山的环境因子存在显著差异；2)细菌群落的结构和表型对环境异质性的响应程度各不相同；3)细菌群落的结构和表型在空间分布格局上具有一致性。对所采取的土壤样品通过高通量测序、BugBase分析、数据统计分析等，旨在探明不同研究站点中土壤细菌群落的结构、表型变化，及二者对土壤异质性的具体响应机制，从而有助于理解微生物适应环境的机制。

1 研究区概况

研究区位于山西省境内的华北落叶松集中分布区，属北温带大陆性气候，气候寒冷干燥，多大风，四季分明，冬季漫长，无霜期短，昼夜温差大。植被类型为寒温性针叶混交林、针阔混交林和温性阔叶混交林。依建群种大致可分为5种类型：华北落叶松(Larixprincipis-rupprechtii)天然次生林纯林、云杉(Piceaasperata)林、油松(Pinustabulaeformis)林、辽东栎(Quercusliaotungensis)林、山杨(Populusdavidiana)、红桦(Betulaalbosinensis)阔叶林。研究区土壤以褐土、山地棕壤、亚高山草甸土为主。

2 研究方法

2.1 样方设置

研究样方选择华北落叶松纯林内，目的是减少植被类型等未知因子的干扰。具体样方数量：五台山叶斗峰8个、关帝山庞泉沟自然保护区10个、芦芽山5个，各站点的地理信息见表1。

表1 研究区样地位置信息①Tab.1 The location of sampling plots

为选择较强的环境异质性，将样方建立在特定海拔梯度上，这样可在较短距离内有较大的气候梯度，从而为细菌群落的结构和表型响应环境异质性的研究提供理想的试验条件。

野外采样时，在每个1 m×1 m的样方中，采用S型采样法，随机选取5个土芯(深度15 cm)，然后混合为一个样本。混合土样被密封在自封袋中并冷藏，立即送往实验室。土壤样本用筛子滤过植物根系和碎石块，一部分土样储存在-80 ℃冰箱中，用于后续的高通量测序；另一部分样本放置在阴凉的实验台上，自然风干，用于测定土壤理化性质。

2.2 土壤环境因子分析

使用仪器HANNA HI3221 pH仪(加拿大)测定土壤pH值(1∶2.5的土/水质量比的土壤悬浮液)(赵鹏宇等，2018)；使用vario MACRO cube元素分析仪(德国)测定土壤总氮、总碳和总硫含量(Zhaoetal.，2019a)；使用间断式元素分析仪(德国；CleverChem 380)测定土壤硝态氮、铵态氮和亚硝态氮含量。采用重铬酸钾-外加热法测定土壤有机碳含量(赵鹏宇等，2018)。

2.3 DNA提取和PCR扩增及测序

使用E.Z.N.A.®Soil DNA Kit试剂盒(OMEGA,USA)提取土壤DNA。每个样本采用25 μL的扩增体系：PCR扩增体系主要包括：0.5 μL的PlantiumTaq(5 U·μL-1)、5.0 μL的10 × EasyTaq Buffer(Mg2+Plus)、4.0 μL的dNTPs、0.5 μL的引物(50 μmol·L-1)，加入5.0 μL的DNA模板和无菌水至50 μL的反应体系。PCR扩增的反应条件：95 ℃预变性3.0 min；95 ℃变性30 s，30个循环；62 ℃退火30 s；72 ℃延伸30 s；72 ℃最后延伸5 min。使用酶标仪Infinite200 PRO(TECAN,Switzerland)测定DNA提取物的纯度和浓度。只有测定值的A260/A280吸光度位于1.8～2.0之间时，才能用于下一步试验分析。从每个土样中提取3个DNA样品，充分混合并送至上海派森诺生物技术工程有限公司进行高通量测序。通过Illumina Miseq测序平台，对16S rDNA高变区V3-V4区进行测序。

2.4 生物信息学分析

使用QIIME软件(v1.8.0,http:∥qiime.org/)对测序所得数据进行分析，具体分析步骤如下。首先对测得的原始数据进行拼接、质量控制，然后去除序列中的嵌合体，最终得到高质量序列。使用抽平的方法，避免在分析过程中由于测序深度不同而造成的分析结果偏差。使用 MOTHUR 软件，采用97%的阈值，对高质量序列进行物种注释分析。

BugBase分析能基于OTU表和Mapping文件对菌群进行表型预测。它可预测的表型类型包括革兰氏阳性、革兰氏阴性、生物膜形成、致病性、移动元件含量、氧需求(包括有氧、厌氧、兼性厌氧)及氧化胁迫耐受等7类(West-Eberhardetal.，1989；赵孟良等，2020；魏莹等，2020)。

2.5 统计分析

使用单因素方差分析(ANOVA)评估不同研究站点的土壤理化性质、细菌群落组成和菌群表型的差异。利用非度量多维尺度分析(NMDS)对细菌群落结构和表型进行排序分析。在NMDS分析中使用Stress值检验NMDS分析结果的优劣。通常认为stress<0.2时可用，其图形有一定的解释意义。

为了评估环境因子的多重共线性，计算每个环境变量(数据标准化后)的方差膨胀因子(variance inflation factor，VIF)：使用R语言分析软件Vegan程序包中的vif.cca函数，并将VIF >10的环境因子剔除，最终筛选出pH值、土壤有机碳、总碳和总氮含量作为环境因子，用于进一步的数据统计分析。然后通过冗余分析(RDA)研究细菌群落的组成结构、表型与环境因素之间的关联关系。以上数据统计分析都是使用R语言的Vegan、ggplot2等程序包下的多个函数完成(Ginestetetal.，2011)。

使用普鲁克(Procrustes)分析评估不同站点的细菌物种丰度组成和表型可塑性的空间分布格局的潜在一致性(Zhaoetal.，2019；Peres-Netoetal.，2001)。线段一端实心形状代表高通量测序数据的排序结果，另一端代表表型数据的排序结果，连线代表2排序构型的残差，可评价二者间的变异情况；M2是普鲁克分析中2排序结果的拟合优度统计量，用于评估2排序结果的相关性，大小范围为0～1，越接近于0说明2排序结果的相关性越大；P-value由Monte Carlo模拟生成，用于检验M2的显著性。细菌群落的结构、表型和环境因子之间的相关分析使用R语言ggcor程序包分析并美化。

3 结果与分析

3.1 土壤化学性质的站点差异

土壤化学性质的单因素方差分析(图1)表明：土壤pH值在庞泉沟最高(6.38)，且和芦芽山有显著差异(P<0.05)；最高的铵态氮含量在芦芽山(36.91 mg·kg-1)，最低的亚硝态氮含量在五台山(0.16 mg·kg-1)，并且2个地点均呈显著差异(P<0.05)；硝态氮含量(6.45 mg·kg-1)、土壤有机碳含量(70.29 mg·kg-1)和总氮含量(0.57 mg·kg-1)均在五台山最高，在芦芽山最低，且庞泉沟和五台山之间均有显著差异(P<0.05)；总碳和总硫含量在所有站点间均差异不显著。这表明3个研究站点环境异质性比较高。

图1 不同研究站点的土壤理化性质Fig.1 Soil physicochemical properties at different sampling sites不同字母表示土壤理化因子在不同的研究站点的差异达到0.05的显著水平。WT：五台山叶斗峰；PQG：关帝山庞泉沟自然保护区；LY：芦芽山。下同。Different letters mean the significant difference (P>0.05)for soil physicochemical properties at different sites.WT：Yedoufeng of Wutai Mountain；PQG：Pangquangou Nature Reserve of Guandi Mountain；LY：Luya Mountain.The same below.

3.2 细菌群落结构和表型分析

通过高通量测序，总共测得1 062 241个高质量序列，共比对出4 258个OTUs。测序的覆盖度指数为95.19%～99.75%，表明所比对的序列涵盖了绝大部分土壤样品中的细菌序列。

相对丰度大于0.01%的细菌门总共比对出15个(图2)。在所有研究站点中，变形菌门(Proteobacteria)(平均相对丰度=30.59%)的相对丰度最高；其次是酸杆菌门(Acidobacteria，19.63%)、放线菌门(Actinobacteria，16.51%)和绿弯菌门(Chloroflexi，13.22%)。其中，变形菌门的相对丰度在庞泉沟最高(34.39%)，放线菌门的相对丰度在芦芽山最高(26.29%)。酸杆菌门(28.68%)和绿弯菌门(16.09%)的相对丰度在五台山最高。另外，相对丰度大于0.01%的细菌纲有20个。

图2 优势细菌门在各研究站点的分布Fig.2 Distribution of dominant bacterial phylum in each site

Bugbase菌群表型预测分析表明，不同站点的菌群表型有很大差异。单因素方差分析(图3)表明，需氧表型在芦芽山与五台山间有显著差异(P<0.05)，但在芦芽山和庞泉沟、五台山与庞泉沟之间无显著差异(P>0.05)。厌氧表型在芦芽山与五台山、庞泉沟与五台山之间有显著差异(P<0.05)，而在芦芽山与庞泉沟之间无显著差异(P>0.05)。移动元件含量在芦芽山和五台山、庞泉沟和五台山之间均有显著差异(P<0.05)。生物膜形成在五台山与庞泉沟间有显著差异(P<0.05)。另外，革兰氏阴性、革兰氏阳性、兼性厌氧、胁迫耐受和致病潜力等表型在不同站点间无显著差异(P>0.05)。

图3 不同研究站点的细菌表型的单因素方差分析Fig.3 ANOVA of bacterial phenotype at different sampling sites不同字母表示细菌表型在不同的研究站点间的差异达到0.05的显著水平。Different letters indicate that the difference in bacterial phenotype between different sampling sites reaches a significant level of 0.05.

非度量多维尺度分析(NMDS)表明，在二维空间上，土壤细菌群落的结构(图4a)和表型(图4b)同一站点内的点离得相对较近，而不同站点间的点离得相对较远，表明群落相异性在同一站点较小，不同站点间较大。进一步通过ANOSIM分析表明，微生物群落的结构和表型在不同站点间差异显著(P<0.01)。

图4 不同研究站点细菌群落的非度量多维尺度分析Fig.4 NMDS analysis(Non-metric multidimensional scaling；NMDS)of bacterial community(a)and phenotype(b)at sampling sites不同颜色或形状的点代表不同分组的样本，2样本点越接近，表明2样本物种组成越相似。横纵坐标表示相对距离，无实际意义。Stress：检验NMDS分析结果的优劣。通常认为Stress<0.2时可用，其图形有一定的解释意义。Points of different colors or shapes represent samples in different groups.The closer the two sample points are,the more similar the species composition of the two samples.The abscissa and ordinate indicate the relative distance and have no practical meaning.Stress:Test the pros and cons of the NMDS analysis results.It is generally considered that it is available when stress<0.2,and its graphics have a certain explanatory significance.

3.3 环境因子与土壤细菌群落的结构和表型的关联分析

细菌群落结构与土壤理化性质的RDA分析结果(图5a)表明，土壤的pH值和有机碳、总碳和总氮含量对细菌群落结构的影响显著(蒙特卡洛置换检验，P<0.01)。RDA的2轴总共解释了88.16%的总变异。其中，变形菌门和蓝藻门(Cyanobacteria)等受pH值影响较多，而土壤有机碳，总氮和总碳含量主要影响浮霉菌门(Planctomycetes)。

细菌群落表型与土壤理化性质的RDA分析结果(图5b)表明，土壤的pH值和有机碳、总碳和总氮含量对细菌群落表型也有显著影响(蒙特卡洛置换检验，P<0.01)。RDA的2轴总共解释了95.43%的总变异。其中，胁迫耐受、需养、革兰氏阳性和兼性厌氧等表型受总碳和总氮含量的影响较大，而革兰氏阴性和治病潜力等表型则受到pH值影响较大，需氧、厌氧和移动元件含量表型受土壤有机碳含量影响较大。

图5 冗余分析揭示土壤理化性质和土壤细菌群落之间的关系Fig.5 Redundancy analysis(RDA)reveals the relationship between soil physicochemical properties and soil bacterial communities图中横纵坐标表示土壤环境因子对细菌的群落结构(a)或表型(b)的解释度。图中箭头表示环境因子，箭头连线的长度表示该环境因子与细菌的群落结构(a)或表型(b)相关程度的大小，连线越长，相关性越大，反之相关性越小。箭头的连线与排序轴以及箭头连线之间的夹角表示相关性，夹角越小相关性越高。锐角表示正相关关系，钝角表示负相关关系。The horizontal and vertical coordinates in the figure indicate the degree of explanation of the soil environmental factors to the bacterial community structure(a)or phenotype(b).The arrow in the figure indicates the environmental factor,and the length of the arrow line indicates the degree of correlation between the environmental factor and the bacterial community structure(a)or phenotype(b).The longer the line,the greater the correlation,and vice versa..The angle between the line of the arrow and the sort axis and the line of the arrow indicates the correlation.The smaller the angle,the higher the correlation.An acute angle indicates a positive correlation,and an obtuse angle indicates a negative correlation.

3.4 土壤细菌群落和表型对环境因子的一致性分析

基于细菌群落的结构和表型的主成分分析(PCA)所获得的2组PCA坐标，进行普鲁克分析(图6)，表明不同站点的细菌物种丰度组成和表型间存在很强的关联(M2=0.217 8，且P<0.001)，具有类似的生物地理分布格局，从而进一步推断细菌群落的结构和表型的空间分布格局具有潜在的一致性，即土壤中细菌群落结构与表型对环境响应具有一致性。

图6 细菌群落结构和表型的普鲁克分析Fig.6 Procrustes analysis of bacterial species abundance composition and phenotypic plasticity

通过对表型可塑性与环境因子和相关优势菌群的分析(图7)表明，土壤总氮含量和酸杆菌门、绿弯菌门极显著正相关，与需氧表型显著负相关(P<0.05)，这从另一方面佐证了土壤理化性质与土壤细菌群落的结构和表型具有很强的关联，且菌群结构和表型之间表现出比它们与土壤理化性质之间更明显的相关。

图7 表型可塑性与环境因子和优势菌群的相关分析Fig.7 Correlation analysis of phenotypic plasticity with environmental factors and dominant bacterial phylum图中的*表示显著水平；在相关分析中，只有呈现显著的相关关系，才会以线条的形式显示出来。* in the figure indicates the level of significance;in correlation analysis,only significant correlations will be displayed in the form of lines.

4 讨论

本研究中3个站点的土壤pH值和总氮、有机碳、硝态氮的含量都比较高，表明土壤肥力较好，与前人研究结果类似(乔沙沙等，2017)，可能因华北落叶松林有大量凋落物、土壤深厚、肥沃、湿润且排水良好(Zhaoetal.，2019a)。本研究中所有站点的变形菌门的相对丰度最高，其次是酸杆菌门、放线菌门和绿弯菌门等，这可能因土壤中变形菌门和拟杆菌门的菌群大多会富集在弱酸性、土壤湿度大、资源有效性高的土壤中，并可以将易分解的碳作为碳源(Zhaoetal.，2019；Pangetal.，2011；李毳等，2018)。

本研究中，与氧需求相关的表型(需氧、厌氧)、生物膜形成和与移动元件含量相关的表型在不同站点间差异显著，这可有以下几个原因：

1)与氧需求相关的表型(需氧、厌氧、兼性厌氧)变化可能与土壤盐度有关(Bronicketal.，2005)。有研究指出，厌氧表型的相对丰度与盐度呈显著正相关，而兼性厌氧和生物膜形成表型的相对丰度与盐度呈显著负相关(Zhangetal.，2019)。一方面，土壤中较高的盐分含量可能导致土壤颗粒分散(Roger-Estradeetal.，2009)和结构变化(Wongetal.，2009)，从而迫使微生物产生表型的差异。另一方面，本研究区中3个站点的土壤pH值和总碳、总氮、有机碳、氨态氮、硝态氮和亚硝态氮含量差异显著。在山地环境下可能形成很多不同的微气候、微生境。在这种情况下，土壤微生物响应外界胁迫的时候，需消耗更多能量，导致更快的新陈代谢速率，更快的电子交换速率。因此，与氧需求相关的表型和生物膜形成表型、移动元件表型的含量得以大大加强。

2)也可能和研究区的海拔梯度有关。有研究指出，一些嗜冷细菌、真菌和革兰氏阴性细菌主要出现在高海拔地区(Margesinetal.，2009)。本研究中，五台山站点的样方的海拔跨度最大，由此推断五台山样地有较大环境异质性梯度。土壤理化性质分析表明，土壤的pH值和总氮、有机碳，硝态氮等含量在五台山最高，由此推断土壤最肥。而群落表型分析表明，需氧、厌氧、兼性厌氧和生物膜形成等表型在五台山的丰度也是最高的。也就是说，环境异质性梯度越大的地方，土壤微生物群落的表型变化越明显。由此推断，微生物能通过调整其表型变化，加强其环境适应性。这可能与2方面有关：一方面微生物能够调整其在生长、存活、繁殖等方面资源分配的权衡性，最直接的外在表现是微生物群落表型发生变化；另一方面，生物有机体能通过表型的适应性改变，产生选择性的性状或功能变化或反应，以应对环境的改变。

3)土壤表层的微生物群落在对抗极端环境过程中可能会形成特有的生存适应能力和群落组成差异。除此之外，林相和地形的差异，可能也会增大土壤表层的含氧量差异，从而迫使微生物群落调整其在生长、存活、繁殖等方面资源分配的权衡性，从而导致菌群与氧需求相关的表型变化。。

本研究表明，土壤细菌的物种丰度组成和表型具有潜在的环境适应一致性(M2=0.217 8，且P<0.001)。有学者通过普鲁克分析指出，在活性污泥中的细菌群落的物种丰度组成和抗生素抗性基因(antibiotic resistance genes,ARGs)丰度组成具有环境适应的潜在一致性(Zhaoetal.，2019)。这可能是由于细菌群落的优势菌群和表型之间表现出更明显的相关，这可能是两者空间分布格局类似的主要原因。由此推断微生物菌群结构和表型在应对外界环境变化时表现为协同变化。

本研究的不足之处是，Bugbase表型分析仅仅是基于OTU数据的预测分析，并不能真实反映出菌群全部表型的变化。但是，本研究结果仍可为土壤微生物的表型可塑性研究提供一些理论参考。另外，在以后研究中，对土壤的性状的异质性分析应加入土壤粒径、孔隙度等物理性质的分析。

5 结论

在山西省华北落叶松林自然集中分布区山地的3个站点，土壤理化性质差异显著，并对土壤细菌群落的结构和表型产生了很大影响，尤其是与氧需求相关的表型(需氧、厌氧)与移动元件含量。土壤细菌群落的优势菌群的物种丰度和表型具有潜在的环境适应一致性，表现出对环境因子的协同变化和相似的生物地理分布格局。这些结果可为研究土壤微生物的表型可塑性、理解物种共存机制等提供更多理论依据。