秦岭山地天然次生林群落MRT数量分类、CCA排序及多样性垂直格局

物种多样性是表征群落学的重要指标, 在反映植物群落生境差异、群落结构组成和稳定性等方面具有重要的意义。因此, 进行植物物种多样性研究可以更好地评价群落组成、结构及其发展变化, 同时测度群落及其环境的保护状态, 对控制和减少濒危物种的丧失具有重要意义, 是群落生态学乃至整个生态学研究的核心之一[1-2]。生物多样性沿环境梯度的变化格局一直是生物多样性研究的一个重要议题, 海拔梯度的变化综合了温度、湿度和光照等多种环境因子, 在一定程度上能够显著地影响植物种群的空间分布格局[3]。坡位、坡度和坡向也是影响植物生长的重要地形因子[4-5]。大量研究表明, 植物物种多样性随着海拔梯度的变化呈现出比较复杂的规律[6-11], 结论不尽相同。

植被格局的数量分析是植被生态学的主要方法和核心内容之一, 包括分类、排序、格局分析、植被与环境因子关系分析等主要方面[12-14], 成为近几十年来植被生态研究中较为活跃的领域之一。植被分类可采用人为分类方法, 也可采用自然分类方法。自然分类的重要原则是尽量减少分类过程中的人为性, 数量分类方法比较客观, 是自然分类的重要手段[15], 可以深刻地揭示植物种、群落与生境之间的关系, 是研究群落分类与排序的重要途径。数量分析方法的目的是揭示植被空间分布格局, 并提出相应的生态学解释。目前, 常用的植被群落数量分类方法有聚类分析(cluster analysis)、双向指示种分析(two-way indicator species analysis, TWINSPAN)和多元回归树(multiple regression trees, MRT)[16]。其中多元回归树是一种较新的群落数量分类方法, 比其他的分析方法更加客观[15]。排序是分析群落之间的连续分布关系, 将样方或植物种排列在一定的空间, 使得排序轴反映一定的生态梯度, 从而解释植物物种、植物群落的分布与环境因子之间的生态联系[16]。目前生态学上常用的排序分析方法有除趋势对应分析(detrended correspondence analysis, DCA)、典范对应分析(canonical correspondence analysis, CCA)和除趋势典范对应分析(detrended canonical correspondence analysis, DCCA)等。CCA在研究过程中将植物物种数据结合多个环境因子数据, 包含大量信息, 每一步的计算结果都与环境因子进行回归, 揭示环境因子对物种组成和群落分布格局的影响, 从而对群落和环境因子间的关系的做出更加直观和容易的生态解释[17]。实际上, 分类和排序都是以不同的方式、从不同的角度来研究植被和物种与环境的关系, 每一类方法主要反映问题的某一方面。因此, 只有不同类型的方法结合起来才能对植被格局及其成因获得较为完整的认识。

秦岭作为长江和黄河两大水系的分水岭, 为我国南北方的地理分界线, 也是我国暖温带和亚热带气候的分界线[18], 属于我国气候变化敏感区域, 地貌复杂, 物种多样性丰富, 其山地气候变化与植被覆盖变化等研究较为丰富。由于受到人类活动的影响, 秦岭植被景观类型发生着结构性的改变[19]。因此, 本文以秦岭北坡中段的周至国家级自然保护区天然次生林群落为研究对象, 探讨其群落数量分类、排序以及物种多样性沿海拔梯度的垂直格局, 以期为该地区的生物多样性保护和天然次生林的可持续发展提供理论依据。

调查区域处于周至国家级自然保护区, 该保护区位于秦岭北坡中段的周至县境内, 地理位置在107°39′—108°19′E, 33°41′—33°57′N, 海拔1000—2900 m。该区属暖温带大陆性季风气候区, 气候特征是夏短热, 冬长冷, 秋多雨。年平均气温6.4—8.4℃, 极端最高温29.7℃, 极端最低温-10℃, 无霜期150 d左右, 年均降水量约1000 mm。土壤自下而上依次为山地棕壤、山地暗棕壤、山地草甸土。区内植被垂直分布明显, 自下而上依次分布有暖温带落叶阔叶林、温性针阔混交林、温性针叶林、寒温性针叶林、亚高山灌丛和草甸, 森林覆盖率高达98.3%[20]。多样的地貌类型和气候孕育了保护区内丰富的植物种类和多样的植被类型, 该区域植被在秦岭区域有一定的典型性和代表性。

2016年7月在周至国家级自然保护区内通过全面踏查, 选择相对成熟的天然次生林作为调查样地。从海拔1000 m到2000 m沿海拔梯度取样, 海拔每升高200 m设一条样带, 共计5条样带, 编号分别为E1、E2、E3、E4、E5, 样带基本信息如表 1所示。在每条样带上沿山脊两侧选择具有代表性的群落分别设置5个面积为20 m×25 m的样地[21], 同一坡面样方按照坡位分别分布在上、中、下坡位, 共获得25个样地, 同一样地内保证坡形、坡向、土壤母质等因子一致。在每个样地的中心及四角分别设立2 m×2 m的灌木样方和1 m×1 m的草本样方, 共设置乔、灌、草样方275个。分别记录每个样地内高度高于150 cm、胸径大于5 cm的乔木层的树种名、胸径、树高、冠幅及生长状况等；记录灌木层和草本层的物种名称、高度、冠径、盖度和株数等[22]。用手持GPS和地质罗盘测量各个样地的经纬度、海拔、坡度、坡向, 用土钎和钢卷尺测定土壤厚度和枯枝落叶层厚度, 同时记录样地的生境条件。

分别计算每个样地中各物种的重要值, 以及群落各层次的物种丰富度(S)、Shannon-Wiener多样性指数(H′)和Pielou均匀度指数(J)。

重要值计算公式[12]如下:

乔木重要值=(相对多度+相对优势度+相对高度)/3

灌木层及草本层物种重要值=(相对盖度+相对高度)/2

物种多样性的相关指标计算公式[11]如下：

S=样方内物种总数；

式中, Pi为物种i的相对重要值。

环境数据中海拔高度、坡度、凋落层厚度、土壤厚度以实际观测记录值表示, 坡向、干扰程度采用数字等级表示。坡向以数字等级表示, 以北为0°, 将坡向按照45°的夹角, 顺时针方向旋转分为8个坡向等级, 以数字表示各等级：1表示北坡(0°—22.5°, 337.5°—360°), 2表示东北坡(22.5°—67.5°), 3表示西北坡(292.5°—337.5°), 4表示东坡(67.5°—112.5°), 5表示西坡(247.5°—292.5°), 6表示东南坡(112.5°—157.5°), 7表示西南坡(202.5°—247.5°), 8表示南坡(157.5°—202.5°), 数字越大, 表示坡向越向阳, 得到的光照越多[23]。干扰包括人类活动和动物活动, 干扰情况分为4个等级：无干扰为1级, 轻微干扰为2级, 中度干扰为3级, 强度干扰为4级[24]。共计6种环境因子。

天然次生林群落的物种多样性指数均采用Excel 2007进行统计计算；采用SPSS 21.0软件, 对次生林群落各指标进行单因素方差分析(one-way ANOVA)。

多元回归树(MRT)是一种较新的研究物种分布与环境因子之间关系的数量分类方法。它将环境因子梯度作为分类节点, 利用递归划分法, 将样地划分为尽可能同质的类别, 减少样地中群落的非相似性, 基于预测误差的测度(CVRE)值最小代表预测能力最佳的原则对回归树进行裁剪, 获得预测能力最好的回归树。同时采用交叉验证来确定分类结果, 能很好地反映异质环境下群落的结构特征[25]。本研究选用25个天然次生林群落样地中6个生境因子为自变量, 195种植物在样地内的重要值为因变量进行群落分类。使用R语言中的mvpart和MVPART-wrap程序包运行MRT并结合IndVal函数来寻找每个群落的指示植物种。

典范对应分析(CCA)是一种非线性多元直接梯度分析方法, 它把分析与多元回归结合起来。将25×195的物种重要值矩阵与25×6的环境因子矩阵结合, 用CCA排序分析群落分布格局与环境因子的关系；使用前向选择法和蒙特卡罗显著性检验, 研究影响群落生境分布格局的主要环境因子, 同时根据环境变量的方差膨胀因子(VIF)检验变量的多重共线性程度。CCA排序使用生物统计学软件R软件中的Vegan软件包进行统计分析。

根据275个样方的调查结果统计, 调查区域共有种子植物195种, 隶属61科128属。蔷薇科(Rosaceae)、菊科(Compositae)是该地区最具优势的两个科。其中蔷薇科21种, 占总种数的10.71%；菊科12种, 占总种数的6.12%。单种属有91个, 占总属数的75.83%。研究结果表明, 调查区域天然次生林群落的物种丰富度、Shannon-Wiener多样性指数和Pielou均匀度指数沿海拔梯度的变化趋势基本一致, 除灌木层的Shannon-Wiener物种多样性指数和Pielou均匀度指数外其余均为先上升后下降的趋势, 呈单峰模型, 在1600—1800 m海拔范围内达到最大值, 而灌木层的Shannon-Wiener物种多样性指数随着海拔的升高而升高, 在1800—2000 m范围内达到峰值；Pielou均匀度指数沿海拔梯度变化趋势较平缓。图 1表明, 物种丰富度表现为乔木层>草本层>灌木层, 乔、灌、草3个层次的物种丰富度在5个海拔梯度间差异显著。由图 2可以看出, 多样性指数表现为乔木层最大, 灌木层和草本层没有明显差别。其中乔木层中各样带间没有明显差异；灌木层样带E1分别与E4、E5差异显著；草本层中, 样带E2分别与其他4个样带差异显著, 样带E1和E4差异显著。由图 3可以看出, 均匀度指数表现为灌木层>草本层>乔木层, 并且群落3个层次物种均匀度指数沿海拔梯度的变化波动较小, 趋势平稳且无显著差异。

基于CVRE值最小代表预测能力最佳的原则, MRT以坡度24为节点(图 4), 经交叉验证并依据植物群落分类和命名原则, 将25个样地划分为2个群落。

群落Ⅰ冬瓜杨(Populus purdomii)+陇东海棠(Malus kansuensis)+蛇莓(Duchesnea indica)群落。该群落包括6、7、9、12、14、15、17、20、23、24共10个样地, 乔木层优势种为冬瓜杨, 伴生种为鹅耳枥(Carpinus turczaninowii)、油松(Pinus tabuliformis)、白桦(Betula platyphylla)、毛樱桃(Cerasus tomentosa)等；灌木层优势种为陇东海棠, 主要伴生种为胡颓子(Elaeagnus pungens)、绿叶胡枝子(Lespedeza buergeri)等；草本层优势种为蛇莓, 伴生种为荩草(Arthraxon hispidus)、华蟹甲(Sinacalia tangutica)、野棉花(Anemone vitifolia)、日本羊茅(Festuca japonica)等。此类型群落大部分分布在山谷、平地或背阴的缓坡, 生境地势较平坦, 水分比较充足而热量较低, 土壤为棕壤和棕壤性石渣土。

群落Ⅱ锐齿槲栎(Quercus aliena var. acuteserrata)+黄栌(Cotinus coggygria)+茜草(Rubia cordifolia)群落。该群落包括1、2、3、4、5、8、10、11、13、16、18、19、21、22、25共计15个样地, 乔木层优势种为锐齿槲栎, 主要伴生种为山杨(Populus davidiana)、华山松(Pinus armandii)、刺榛(Corylus ferox)、红桦(Betula albosinensis)等；灌木层优势种为黄栌, 伴生种为胡枝子(Lespedeza bicolor)、华北绣线菊(Spiraea fritschiana)、荚蒾(Viburnum dilatatum)等；草本层优势种为茜草, 主要伴生种为鹅观草(Elymus kamoji)、牛尾蒿(Artemisia dubia)、求米草(Oplismenus undulatifolius)、日本羊茅(Festuca japonica)等。该类型群落分布较广, 占据山体地势陡峭的半阳坡或阴坡, 土壤水分中等偏湿润, 土壤为褐土和棕壤性石渣土。

通过使用前向选择法和蒙特卡罗显著性检验, 从6个环境因子出筛选出包括海拔(Altitude)、坡向(Aspect)、凋落层厚度(Litter layer thickness)和干扰情况(Disturbance)4个环境代理变量, 同时计算环境变量的方差膨胀因子(VIF), 其VIF值均小于10, 说明无明显的共线性问题。

通过对调查区域的25个样地进行典范排序(CCA)得到二维排序图。结果如图 5所示, 箭头表示环境因子, 箭头和排序轴的夹角代表该环境因子与排序轴的相关性大小, 箭头的长短表示该环境因子与群落分布的相关程度的大小。4个环境因子的解释量为20.73%, 前4个排序轴的特征值分别为0.5494、0.4796、0.3999、0.3078, 表明CCA分析在一定程度上能解释植物群落的分布格局与环境因子的关系。结合图 5和表 2可以看出, 与第一排序轴正相关最高的是凋落层厚度, 呈极显著相关, 相关系数为0.99491, 其次是海拔, 也呈极显著相关, 相关系数为0.94465；与第一排序轴负相关的是干扰情况和坡向, 相关系数为0.80342和0.03173, 说明第一轴排序主要反映了植物群落分布格局在凋落层厚度、海拔和干扰情况梯度上的变化, 即沿排序轴第一轴从左到右凋落层厚度呈逐渐增加趋势、海拔越来越高, 人类或动物活动逐渐减弱。与第二排序轴呈现正相关的是坡向和凋落层厚度, 相关系数分别为0.99950和0.10080；与第二排序轴呈现负相关的是干扰情况和海拔, 相关系数为0.59541和0.32808, 说明第二轴排序主要反映的是坡向、干扰情况和海拔梯度上的变化, 具体情况为沿排序第二轴从下往上, 坡向由阴坡转为阳坡、人类或动物活动逐渐减弱, 海拔越来越低。综合前两轴, 海拔、坡向、凋落层厚度和干扰情况对研究区域次生林群落的分布有显著的影响。从箭头长短来看, 箭头最长的是海拔, 说明海拔高低可能是影响该地区植物群落分布最重要的因素。

结合群落分类和排序结果可以看出两个群落的分布规律, 群落物种组成是沿着环境因子梯度连续渐变的, 两个群落有个别样方在排序轴中重叠。群落Ⅰ分布在生境地势较为平坦的背阴缓坡, 位于排序图的下方中部。群落Ⅱ位于排序图的上方偏左, 该群落位于凋落层较厚、受人类或动物活动干扰较小的阳坡或半阳坡。

在研究区域调查的275个样方中, 种子植物包括61科128属195种, 以蔷薇科、菊科为主。植物区系中单种科、单种属所占比例较大, 可能是因为研究区域海拔较低, 人为干扰情况严重造成了该地区植物总体上呈现出少数科属多样性突出、植物群落组成比较分散的特点。有些研究者认为低海拔地区受人为干扰较为严重, 所以其物种多样性较低。人为干扰所造成的生境破碎对物种多样性会造成很大的影响, 导致生长力较弱的物种种群生长受阻、发育过慢, 逐渐衰落甚至消亡[26]。

结果表明, 调查区域天然次生林群落不同层次的调查区域天然次生林群落的物种丰富度、Shannon-Wiener多样性指数和Pielou均匀度指数沿海拔梯度的变化趋势基本一致, 除灌木层的Shannon-Wiener物种多样性指数和Pielou均匀度指数外其余均为先上升后下降的趋势, 呈单峰模型, 在1600—1800 m范围内达到峰值, 这可能是因为中海拔地区温、湿状况、水热条件较高海拔地区适宜, 并且相对于低海拔地区干扰较少。分布在海拔范围为1600—1800 m的群落光照相对充足, 土壤和气候都较为湿润, 水热组合最佳, 资源的可利用性可能是最高的, 因此形成了最高的物种丰富度。这与Liu等[27]的研究结果基本一致。灌木层的Pielou均匀度指数波动较小, 趋势平稳且无显著差异, 可能是因为海拔变化后, 群落内小生境依然呈现均质性, 群落各层次物种分布聚集度变化较小。群落不同层次的物种多样性对垂直海拔梯度的响应也不一样。乔木层的物种丰富度和Shannon-Wiener多样性指数明显大于灌木层和草本层, 可能是因为乔木层物种丰富, 林间郁闭度大, 林下植被光照不足, 限制了林下植物物种的发育。Pielou均匀度指数呈现为灌木层>草本层>乔木层, 说明乔木层植物分布不均匀, 且优势种现象明显, 这与相关研究结果一致[28]。

采用MRT分类方法将研究区域25个样地分成两个群落类型, 与样方在CCA二维排序图上的分布格局基本相同, 25个样地仅包含两个群落的原因可能是本研究所选取的样方分布在一个连续的局部区域, 生境因子较为接近。群落分类不仅要求物种组成和结构相对一致, 占据的生境也是很重要的特征, 并且有研究表明不同植物具有不同的生境偏好[29]。群落Ⅱ分布在地势陡峭的山坡, 光照相对充足, 出现了黑水亚麻(Linum amurense)、阴行草(Siphonostegia chinensis)等喜干燥气候及土壤的植物种。而群落Ⅰ分布在山谷, 平地或缓坡地带, 水分条件较好, 多生长蛇莓等喜湿的物种。有研究表明[30], 坡度较缓地带, 土壤湿度较大, 多生长喜湿植物；而坡度较大的地带, 在重力作用和地表径流的作用下水分不断向下渗透, 导致该区域的土壤含水量低于坡度小的生境, 生长对水分响应较小的植物, 这与本研究的结果也一致。

CCA二维排序图表明海拔梯度、坡向、凋落层厚度和干扰情况4个环境因子对该区域天然次生林群落的分布有显著的作用, 主要反应了海拔的垂直变异和坡向、凋落层厚度等方面的水平变异, 影响其分布格局, 形成不同的植被类型。CCA第一轴排序主要与海拔和凋落层厚度有关, 植物对海拔梯度变化的响应, 本质上反映了其对水热条件需求的差异, 而凋落层越厚, 表明土壤养分越丰富。第二排序轴主要与坡向有关, 坡向对物种组成也有较显著的影响, 坡向其实反映的是光照条件的变化, 直接影响了太阳辐射和水肥的空间再分配, 代表着光照、水分、养分等环境因素的常年、综合生态梯度变化[31-32]。同时, CCA排序结果表明, 所调查的环境因子只能解释部分次生林群落的群落分布, 而未能解释部分超过50%, 暗示了还存在很多未调查的因子可能也是显著影响该地区次生林群落分布格局的重要因子。

数量分类和排序方法在群落学分析中结合使用, 往往效果更好。这是因为植物群落既有间断性又有连续性, 分类侧重于揭示群落的间断性, 增加排序的有效性；而排序侧重于揭示群落的连续性, 能够更好地理解分类的结果[33]。分类和排序后的环境解释可以深刻地揭示环境和植物群落间的数量关系。MRT和CCA的结果均表明本研究所调查的环境因子对群落的数量分类和排序的贡献都不大, 所以, 在今后的研究中, 需测量更多的环境因子, 挖掘更多的数据, 从而给出更深入全面关于群落分布与环境因子的关系的解释。

物种多样性与海拔梯度的关系错综复杂, 研究天然次生林群落组成和物种多样性随海拔梯度垂直变化格局对揭示生物多样性环境梯度的变化规律具有重要意义, 该结果有助于更好地揭示研究区域次生林群落空间分布格局特征及各物种的生态适应性, 可为深入了解海拔梯度的垂直格局与生物多样性之间的关系累积数据；MRT数量分类和CCA的排序研究能够反映群落分布与环境因子之间的关系, 为该地区植被恢复以及保护工作提供相应的理论依据。