戴云山物种多样性与系统发育多样性海拔梯度分布格局及驱动因子

李梦佳,何中声,江 蓝,谷新光,晋梦然,陈 博,刘金福

1 福建农林大学林学院,福州 350002 2 福建农林大学海峡自然保护区研究中心,福州 350002 3 生态与资源统计福建省高校重点实验室,福州 350002

生物多样性海拔梯度分布格局是多样性研究的重要问题[1]。海拔作为一个综合性生态因子,会引起光、温、水与土壤因子等的综合变化,进而影响到植物群落的物种、系统发育与功能多样性的空间分布[2]。海拔梯度的环境异质性产生了多种物种多样性格局,主要有单调递减和中间高度膨胀格局[3]。由此衍生出多种探讨海拔梯度格局的内在机制假说,如水热动态假说[4]、中域效应假说[5]和环境能量假说等[6]。与物种多样性相比,系统发育多样性可以从进化历史积累的特征来阐述生物多样性,能更好地解释群落构建、多样性维持与生态系统功能[7-8]。由于性状演化的保守性,系统发育多样性往往和功能多样性密切相关,加之其避免了功能多样性测定时性状选择的主观性,尤其在关键性状缺失条件下,其代表性优于功能多样性[9]。

虽然系统发育多样性与物种多样性高度相关[10],而目前系统发育多样性及结构的海拔梯度格局并无统一规律。在生态位保守前提下,影响群落系统发育结构的主要因素为生境过滤和种间竞争。如果群落中生境过滤占主导地位,将筛选出适应能力相似、亲缘关系相近物种,群落呈现系统发育聚集；相反,竞争排斥作用会使生态位相似的物种无法共存,群落内物种亲缘关系较远,表现为系统发散[11]。海拔梯度增加了群落环境异质性,而环境异质性是否会影响群落在海拔梯度上的系统发育结构,以及导致系统发育结构产生差异的驱动因子有哪些？解决这些问题将进一步丰富多样性海拔梯度格局理论。

戴云山国家级自然保护区位于福建省中部,属于南亚热带和中亚热带过渡带[12]。戴云山植被覆盖率93.4%,且具有明显的垂直带谱,植被分层结构明显,是研究山地局部尺度系统发育多样性的良好平台。目前戴云山海拔梯度研究主要集中在土壤空间变异、功能性状对比、群落特征以及物种多样性等方面,而海拔梯度上系统发育多样性驱动因子研究尚未报道[12- 14]。为此,以戴云山南坡不同海拔梯度植物群落为研究对象,探讨物种多样性、系统发育多样性和系统发育结构的海拔变化及其影响因素。旨在解决：(1)物种多样性与系统发育多样性的海拔分布格局；(2)物种多样性与系统发育多样性间相互关系；(3)环境如何驱动物种与系统发育多样性海拔分布？通过研究戴云山系统发育与物种多样性海拔梯度格局及其驱动因子,对探讨环境因子在群落构建与多样性维持中的作用具有重要意义。

1 研究区概况

图1 戴云山样地设置示意图Fig.1 The sample plots of Daiyun MountainDYS-900至DYS-1600：戴云山海拔 900 m至1600 m

戴云山国家级自然保护区位于福建省泉州市德化县内(118°05′22″—110°20′15″E,25°38′07″—25°43′40″N),最高海拔1856 m,最低海拔650 m,是闽中生物多样性中心,具有南亚热带与中亚热带过渡特征[12]。保护区属亚热带海洋性季风气候,年平均气温15.6—19.5 ℃,最冷月均温在6.5—10.5 ℃(1月),最热月均温23—27.5 ℃(7月)。年降水量1700—2000 mm,年平均相对湿度≥80%,年平均日照时数1875.4 h,年均雾日220 d。土壤主要为赤红壤、山地红壤、山地黄红壤和泥炭沼泽土。垂直海拔梯度大,随海拔升高依次出现垂直带谱为：常绿阔叶林(650—1000 m)、针阔混交林(1000—1200 m)、温性针叶林(1200—1650 m)、苔藓矮曲林(1650—1765 m)与草甸(1750 m以上)[15]。

2 研究方法

2.1 植物调查

为减少人为干扰影响,于2017—2018年自海拔900 m起,每间隔100 m设置一个海拔梯度,共计八个梯度(图1)。根据CTFS(热带林业科学中心)标准设置植物群落固定样地,每个海拔设置一个20 m×30 m样地,每个样地划分3个10 m×20 m样方,样地选择林相整齐且具有代表性的森林群落。采用GPS记录每个样方的经纬度、海拔、坡度、坡向等地形因子[16]。调查样方内胸径≥ 1 cm的树木种类与多度。记录种名,测量其胸径、树高和冠幅等。物种定名参考《中国在线植物志》(http://frps.iplant.cn/)。

表1 每个样地内重要值前五的树种

2.2 环境因子调查

土壤因子选取10个土壤指标,即土壤全氮(TN)、全磷(TP)、全碳(TC)、有效磷(AP)、全钾(TK)、碳氮比(C/N)、碳磷比(C/P)、氮磷比(N/P)、土壤pH。每个样方设置三个土壤取样点,分别位于样方中心、东北角和西南角。在每个取样点采集0—20 cm土层土壤样品,并使用环刀收集土壤样品,带回实验室,自然风干,测定其理化性质。各指标测定方法参考《土壤农化分析》[17]。土壤全氮、全碳含量测定采用碳氮分析仪测定;全磷、全钾使用等离子电感耦合装置测定;速效磷采用碳酸氢钠浸提-钼锑抗比色法;酸碱度采用电极法测定。空气与土壤温度监测使用美国MAXIM iButton(Maxim Intergrated, iButton, San Jose, USA)DS1923-F5型(空气温湿度)和DS1922L-F50型(土壤温度)记录器,在每个海拔放置1个空气与土壤温度记录器,每日0点开始每隔4 h记录一次土壤与空气温度,以平均值作为各海拔空气与土壤温度最终数据。生长季土壤温度(GST)和空气温度(GAT)分别为3月—8月土壤和空气温度的平均值。休眠季土壤温度(DST)和空气温度(DAT)分别为9月—次年2月土壤和空气温度的平均值。

2.3 物种多样性计算

物种多样性指数选择用基于物种重要值(如表1),来表示的Margalef、Simpson、Shannon-Wiener和Pielou指数。Margalef指数仅考虑群落的物种数量和总个体数,将一定大小的样本中的物种数量定义为多样性指数[18]。

式中,N代表总个体数量,S为总物种数量。

Simpson指数是假设在对无限个群落随机取样,样本种两个不同种的个体相遇的概率为多样性的测度[19]。公式如下：

式中,N为总个体数,Ni为第i个种的个体数量。

Pielou指数表示与物种丰富度无关的均匀度指数[19],公式如下:

式中,S为物种数,N为全部种的个体总数。Ni是物种i的个体数。

Shannon-Wiener指数是假设在无限大的群落中对个体随机取样,而且样本包含了群落中所有的物种,个体出现的机会即为多样性指数[20]。公式如下:

式中,Pi表示第i个种的多度比例。

2.4 系统发育多样性与系统发育结构计算

利用在线软件phylomatic (http://phylodiversity.net/phylomatic/) 在线生成基于被子植物分类系统APGIII且具有分化时间的谱系树结构,再利用R中picante软件包计算系统发育多样性和系统发育结构。系统发育多样性采用PD(Phylogenetic diversity)指数[21],该指数用于量化群落系统发育总枝长[22]。系统发育结构指数选用净谱系亲缘关系指数(Net relatedness index,NRI)和净最近种间亲缘关系指数(Net nearest taxa index, NTI)。NRI指数侧重于物种之间相似度,NTI侧重于相似物种间的影响[23]。

假设群落调查的物种随机分布构成了一个物种库,根据平均发育系统(MPD)计算平均系统发育距离[24]。在物种中提取999种提取物构成随机零模型种物种的MTD分布。NRI和NTI计算公式如下：

式中，MPDs和MNTDs分别表示平均谱系距离观测值和最近相邻谱系距离平均观测值,MPDmds和MNTDmds表示通过软件随机999次模拟出平均谱系距离和最近相邻谱系距离的平均值,SD为标准差。

通过计算NTI和NRI分析各个群落的系统发育结构。若NRI>0,且NTI>0表明系统发育结构聚集,代表此群落更趋向于亲缘相近的物种组成；若NRI<0,且NTI<0,代表此群落物种间物种群落更趋于分散；若NRI=0,NTI=0,表明群落系统发育结构随机。

2.5 数据分析

2.5.1海拔梯度分布格局

采用单因素方差(one-way ANOVA)分析各海拔梯度间物种多样性Margalef指数、Simpson指数、Shannon-Wiener指数和Pielou指数的差异性,以及系统发育多样性(PD)、系统发育结构(NRI与NTI)的差异性。

2.5.2相关性分析

检验数据是否符合正态分布,符合正态分布相关性分析采用Pearson相关性分析物种多样性与系统发育多样性、系统发育结构之间的关系,采用一般线性模型对物种多样性与系统发育多样性进行拟合,根据赤池信息量准则(AIC)选择最优模型。

2.5.3驱动力分析

采用Pearson相关性分析地形因子有海拔(ELE)；气候因子有生长季土温(GST)、休眠季土温(DST)、生长季气温(GAT)和休眠季气温(DAT)；土壤因子包括pH、土壤含水量(SWC)、土壤全氮(TN)、全磷(TP)、全碳(TC)、有效磷(AP)、钾(TK)、碳氮比(C/N)、碳磷比(C/P)和氮磷比(N/P)等影响物种多样性与系统发育多样性的环境因子。对显著性因子采用广义线性模型(Generalized linear model)中基于正态分布的Gaussian模型来探讨影响系统发育指数与物种多样性指数的驱动因子。根据AIC准则采用逐步回归对模型进行优化,同时通过方差膨胀因子VIF对回归模型的共线性进行诊断。当VIF>10时,剔除VIF值较大的变量,重新对模型进行拟合和优化,得到物种多样性、系统发育多样性与系统发育结构对环境因子的最优响应模型。

数据分析主要采用R 3.5.3软件,使用 mass 包进行正态分布检验、单因素方差分析,corrplot 包进行相关性分析和 broom 包进行广义线性模型分析[24]。

3 结果与分析

3.1 物种多样性与系统发育指数的海拔梯度分布格局

8个海拔梯度森林群落的物种多样性、系统发育多样性与系统发育结构结果见表2。森林群落系统发育多样性指数范围在1270.979—3091.344。物种数越多,系统发育多样性PD值越大,Margalef、Pielou、Shannon-Wiener和Simpson指数也相对较大。PD指数随着海拔升高呈先升高再下降的单峰格局。Margalef指数沿海拔梯度的趋势与PD指数一致。而Pielou、Shannon-Wiener和Simpson指数在海拔梯度上没有明显规律。系统发育结构呈现为高海拔聚集,低海拔发散的格局。

3.2 相关性与回归分析

物种多样性Shannon-Wiener指数与生物多样性指数均具有显著相关性(图2)。Shannon-Wiener指数与物种多样性指数呈极显著正相关关系(P<0.01),与系统发育多样性指数PD之间具有显著相关性(P<0.05),与系统发育结构NRI(P<0.01)和NTI(P<0.05)也呈显著相关性。系统发育多样性指数PD仅与物种多样性指数中的Margalef指数及Shannon-Wiener指数有显著相关性。Margalef指数与PD指数之间相关性呈极显著,相关系数高达0.898。系统发育结构指数NRI、NTI指数与Simpson指数、Shannon-Wiener指数、Pielou指数显著相关(P<0.05),其中Simpson指数与NRI、NTI指数线性模型拟合度较高(表2)。

对具有显著关系的指数进行回归分析(表2),PD指数与物种丰富度Margalef指数间线性关系拟合最佳(R2=0.8072,AIC=328.0105)。系统发育多样性PD随着物种丰富度Margalef指数升高呈单调递增模式。其中Shannon-Wiener指数与系统发育多样性、系统发育结构都存在线性相关关系,Simpson指数与系统发育结构线性关系拟合较好。

表2 森林群落物种多样性与系统发育多样性的海拔梯度变化

同一列中不同的字母表示物种多样性指数、系统发育指数在海拔间的显著差异性(P<0.05)

图2 森林群落物种多样性指数与系统发育指数Pearson相关性Fig.2 Pearson correlation between species diversity and phylogenetic diversity indices for forest communities*表示P<0.05, **P<0.01, ***P<0.001. 蓝色表示正相关,红色表示负相关;颜色越深,圆形越大,表示相关性越强;PD: 系统发育多样性指数Phylogenetic diversity index NRI:净谱系亲缘关系指数Net relatedness indexNTI: 净最近种间亲缘关系指数Net nearest taxa indexS.W: Shannon-Wiener指数;

3.3 驱动因子分析

土壤因子中全磷含量是影响物种丰富度Margalef指数的重要因素(表3、表4),随着土壤全磷含量的增加,Margalef指数呈现下降趋势。Simpson指数、Shannon-Wiener指数、Pielou指数与土壤含水量呈显著正相关(P<0.05)。Simpson指数、Shannon-Wiener指数和Pielou指数随着土壤含水量增加而增加。其中土壤含水量随着海拔上升呈下降趋势,Simpson指数、Shannon-Wiener指数和Pielou指数也随之下降。PD指数与Margalef指数一致,与土壤全磷含量高度相关。随着海拔升高NRI与NTI指数也逐渐升高,海拔作为最主要驱动因子,改变了温度与土壤理化性质进而影响NRI与NTI。海拔是影响系统发育结构变化的主要环境因子,而土壤因子是影响物种多样性与系统发育多样性的主要因素。

表3 森林群落物种多样性指数与系统发育指数的回归

表4 森林群落多样性指数与环境因子相关系数

4 讨论

4.1 群落多样性和系统发育多样性的海拔梯度变化

戴云山物种多样性与系统发育多样性均随海拔梯度分布呈单峰格局,均在1200 m处出现峰值,表明戴云山南坡群落物种丰富度格局与系统发育多样性符合中间高度膨胀格局,即物种丰富度与系统发育多样性最大值出现在中海拔区域[6]。植被类型和植被演化历史是影响物种丰富度格局重要因素[25]。从植被演化角度,福建戴云山自晚全新世来已出现常绿阔叶林,地质条件使戴云山具有鲜明山地特色,并保有丰富久远的植物区系。有研究表明,植物群落多样性与植物系统发育的关系实质上是资源生产力、系统发育和群落物种多样性三者的关系[26]。从物种垂直地带性分布格局中,资源生产力高的亚热带常绿阔叶林[14]。高的资源生产力环境有利于古老类群的生存[26]。资源生产力从高到低的环境梯度上,随植物群落多样性降低的同时,植物类群由古老类群向年轻进化的类群演化[14, 15, 27, 28]。戴云山由于未受第四纪冰川直接影响,戴云山目前有许多古老蕨类植物分布于低海拔处,长梗薹草(Carexglossostigma)与九仙山薹草(C.jiuxianshanensis)常分布于海拔1200 m以下的林下阴湿地与溪边,九仙莓(Rubusyanyunii)多分布于中海拔处[15]。1200 m的中海拔地区拥有古老的植物类群及较高的资源生产力,且人为干扰较少,表现出其群落物种多样性与系统发育多样性呈现出中间高度膨胀格局。

表5 广义线性模型分析各指数与环境因子的关系

群落系统发育结构是生态环境相互作用,环境过滤和历史进化相结合的结果[29]。戴云山系统发育结构随着海拔梯度上升由系统发散逐渐演变为聚集。在低海拔地区(900—1200 m)水分与温度更适宜植物生长,群落主要通过生态位分化来驱动,那些占据相似生态位的物种之间竞争激烈,物种间生态位逐渐分化,共存物种之间亲缘关系相对较远,从而导致系统发育呈现发散结构[30]。高海拔地区较为恶劣的生境条件使环境过滤占据主要地位,环境筛选出戴云山高海拔地区以黄山松为优势的群落,表现为相对聚集。

在海拔1000 m与1400 m的植物群落中,NRI指数离散,而NTI指数聚集,两个指数不一致可能与植被稳定程度有关[31]。群落中存在一些易受胁迫影响的物种,当这些物种丧失后,群落极其脆弱,稳定程度下降,致使系统发育结构没有明显的聚集和发散趋势。从物种重要值角度看,均匀度和优势度可视为戴云山群落稳定程度的判定标准之一。优势度与重要值较低的群落中,如海拔1000 m与1400 m处,优势乔木所占优势度最大,在群落中最为重要,当优势种丧失后,群落结构与群落生境(微环境)发生改变,优势种的去除可导致受优势种压制的多个物种得到释放,而中间竞争与生境过滤的作用还未得到体现,这可能是导致群落物种系统发育结构无明显趋势的一个重要原因。戴云山1000 m与1400 m处NRI指数与NTI指数不一致与植被群落稳定健康程度相关,表明仅从优势度与重要值角度推测植被群落的健康稳定程度不够全面,有待增加取样数量进一步研究。

4.2 物种多样性指数与系统发育指数相关性分析

戴云山森林群落系统发育多样性与物种丰富度之间存在极显著正相关关系[32]。其中Shannon-Wiener指数与物种多样性指数、系统发育多样性指数与系统发育结构指数高度相关,在一定程度上可替代上述指数对森林群落结构进行定量描述。系统发育多样性受物种丰富度(Margalef)、多样性(Shannon-Wiener)的影响,而物种优势度(Simpson)、均匀度(Pielou)对系统发育多样性影响不显著。已有多数研究证明了两者之间的相关性,但在进化树不平衡或独特物种狭域分布的情况下,二者相关性会下降[33, 34]。如南非好望角发现PD值与物种丰富度所指示的物种多样性区域不吻合,与进化树不平衡相关[10]。在戴云山南坡群落中,未出现进化树不平衡或独特物种分布区域,群落系统发育多样性与物种多样性之间呈现极显著相关性。同时,系统发育结构NRI指数、NTI指数与Simpson、Shannon-Wiener、Pielou指数显著相关,说明了物种均匀度与优势度对系统发育结构具有一定影响,物种优势度指数与均匀度指数不仅能描述物种多样性,在某一方面也能描述群落物种结构。Simpson、Shannon-Wiener、Pielou指数与系统发育结构的显著关系,进一步说明了物种多样性与系统发育多样性高度相关[35]。因此,在制定保护政策时,有必要结合物种多样性和系统发育多样性进行综合考虑。

4.3 物种多样性与系统发育指数驱动因子分析

物种多样性指数和系统发育指数与环境因子均存在相关性。土壤磷元素作为红壤区植物生长关键元素,是影响戴云山植物物种多样性和系统发育多样性海拔分布的主要环境因子。戴云山土壤全磷和有效磷含量处于较低水平[36]。Margalef指数、PD指数与土壤全磷含量存在显著负相关关系。土壤磷含量与群落中丛枝菌根(Arbuscular mycorrhiza,AM)的比例有关[37],丛枝菌根比例的变化可能导致不同海拔梯度土壤磷含量的变化,从而影响了群落系统发育多样性。戴云山南坡菌根群落结构从900 m的75%外生菌根(Ecto-mycorrhiza, ECM)+21.5%丛枝菌根转变为1100—1300 m的2%外生菌根+63%丛枝菌根+17.6%杜鹃花类菌根(Ericoid mycorrhiza ERM)最后转变为海拔1400—1600 m处37%杜鹃花类菌根+35.2%外生菌根。因此,海拔梯度上树种变化导致菌根变化,使得群落系统发育多样性与土壤全磷含量相关。系统发育结构与有效磷含量显著相关。随着海拔梯度升高,有效磷含量显著下降,戴云山高海拔有效磷不足的生境条件相当于环境筛,生境过滤作用筛选出适应磷不足的物种,使得种间亲缘相近物种增加,物种间亲缘关系接近,这也是导致高海拔地区系统发育结构聚集的重要原因。系统发育多样性与土壤磷含量显著相关,与程毅康等的研究结果一致[38],在这两个地区中土壤磷元素含量的缺少对于物种丰富度和系统发育多样性的限制大于土壤含水量与降水的限制。此外,内蒙古草原坡位对高寒灌丛群落谱系多样性和谱系结构的研究也发现有效磷是决定系统发育结构的主导因子[39]。

系统发育多样性与水热条件关系也十分密切[40]。随着海拔升高,土壤温度逐渐降低,系统发育多样性也随之下降。土壤温度主要影响物种根系生长、分布以及土壤微生物结构,进而影响到系统发育多样性分布格局。NTI指数与NRI指数与温度也显著相关,随海拔升高,温度是影响植物分布最重要因素之一。 Mastrogianni 等研究表明土壤裸露导致土温调节能力下降可能使群落在系统发育上呈现聚集状态[31]。戴云山高海拔地区空气温度较低,土层较薄且温度调节能力差,温度因子作为主要环境筛,使得戴云山森林群落系统发育结构随海拔上升呈聚集趋势。戴云山脉作为海峡西岸的第一道山脉,降雨量丰富[41]。与优势物种密切相关的Shannon-Wiener、Simpson、Pielou指数与土壤含水量显著正相关。土壤含水量受降水影响,对群落中优势物种产生一定作用,从而改变了植物多样性与密度[42]。不同生境中气候和土壤条件的异质性导致了群落建群机制差异[43]。在具有明显海拔梯度分化的生境中,海拔变化所引起的土壤、水分、温度,特别是温度的变化是导致群落受到竞争排斥或是生境过滤支配的原因[43]。在戴云山南坡八个海拔梯度中呈现出从低海拔到高海拔,系统发育结构由发散到聚集的结构。戴云山森林群落在低海拔受竞争排斥的主导,高海拔受生境过滤的影响,在中间海拔呈现出竞争排斥与生境过滤相互平衡所导致的随机模式。

5 结论

戴云山南坡植物群落物种多样性与系统发育多样性在900—1600 m海拔梯度格局下呈现中间高度膨胀格局,其系统发育多样性与物种多样性之间高度相关,海拔与土壤磷含量是物种丰富度、系统发育多样性和系统发育结构沿海拔梯度变化的重要驱动因子之一。土壤含水量是影响群落物种优势度、均匀度的最显著因子。海拔是影响系统发育结构变化的主要环境因子,而土壤因子是影响物种多样性与系统发育多样性的主要因素。