高寒草地植物群落关键种对不同放牧家畜组合放牧的响应

刘玉祯,孙彩彩,刘文亭,杨晓霞,冯 斌,时 光,张 雪,李彩弟,杨增增,高 婕,张小芳,俞 旸,张春平,董全民

青海大学畜牧兽医科学院, 青海省高寒草地适应性管理重点实验室 西宁 810016

作为陆地最大的生态系统之一[1],草地不仅是畜牧业可持续发展的基础,同时也发挥着涵养水源、维持碳平衡和保护生物多样性等重要生态服务功能[2]。青藏高原高寒草地分布广泛,由于其特殊的生态序列位置,使之成为我国重要生态安全屏障的同时,也使其对外界变化响应更为敏感[3—4]近几十年来,受气候变化与人为过度干扰双重因素的影响,高寒草地生态系统结构与功能变化明显[5—8]。放牧作为青藏高原高寒草地最主要的利用方式之一,也是该系统所受干扰最主要的来源[9]。相关研究表明,我国约90%的草地退化过程由过度放牧主导[10]。草地在退化过程中会逐渐破碎化并形成退化斑块[11],在退化斑块的介导下植被群落特征与土壤理化特征会在局部尺度下呈现出高度的空间变异性[12],并进一步导致草地生产力降低、生物多样性下降甚至生态系统功能的丧失[13]。因此,如何合理利用高寒草地、实现其生态系统的可持续发展、保障其生态系统安全是当前亟待解决的问题。

群落稳定性作为衡量草地生态系统功能与健康程度的重要指标,同时受诸多生物因素与非生物因素的共同调控[14]。伴随空间尺度大小的改变,调节群落稳定性的关键指标亦随之变化,较大空间尺度下,温度、降水等气候因子是决定群落稳定性的关键因素[15]；较小空间尺度下,生物因素、土壤可利用资源等则会成为决定群落稳定性的关键因素[16]。大型草食动物作为小空间尺度下放牧草地食物网的主要消费者,其长期选择性采食会在一定程度上改变相邻植被群落内物种间的竞争关系,进而导致群落内的物种发生变化[17]。与此同时,大型草食动物的践踏、排泄物返还等生命活动不但会改变草地土壤物理性质,也会对一些土壤限制性资源如氮等起到调控作用,进而影响群落的稳定性[18]。

1969年,生态学家Paine[19]首次提出了关键种的概念。随着对关键种研究的逐渐深入,生态学家们发现,关键种在群落结构、功能与稳定性中发挥着重要的作用,且群落中的优势种并不一定是关键种[20]。关键种的确定通常有2种方法[21],即移除法与对比法。物种移除法破坏性较强,且需要一个先决条件,即需假设生态系统中有无限多物种,而真实的生态系统难以满足这一先决条件,另外,不同试验区域的空间异质性也会在一定程度上降低移除法的可靠性[22]；对比法实质上是利用空间替代时间的方法,在一定程度上克服了移除法的局限性,但是由于历史观测过程中的主观性太强且不同观测地点本底差异的存在,该方法的准确性依旧有待商榷[23]。网络分析作为评价复杂群落稳定性的有力方法[24],已经被广泛应用于微生物群落结构的评估以及关键种的识别中(在网络中具有最高的中心度)[25]。Zhang等[26]应用网络分析系统分析了不同恢复年限高寒草地植物群落结构稳定性及关键种演替状况,明确了高寒草地最佳围封年限(7年),并通过进一步分析发现基于物种重要值构建的贡献网络最具代表性。目前,网络分析在微生物群落结构特征研究方面应用已十分广泛,但相比之下,将其应用于高寒草地植物群落结构特征研究的报道还较为少见。鉴于此,本研究选取青海省海北藏族自治州海晏县西海镇的高寒草地为研究对象,设置中等放牧强度下不同放牧家畜组合放牧(牦牛单牧、藏羊单牧、牦牛藏羊1∶2混牧、牦牛藏羊1∶4混牧、牦牛藏羊1∶6混牧)以及围封试验,系统分析其植物群落特征、关键种演替特征与土壤理化性质,着重探讨以下2个科学问题：(1)高寒草地植物群落关键种对不同放牧家畜组合放牧的响应；(2)明确驱动不同放牧家畜组合放牧下关键种演替的关键理化因子,以期为高寒草地的管理与可持续发展提供科学依据。

1 材料与方法

1.1 研究区域概括

研究区位于青海省海北藏族自治州海晏县西海镇(36°92′ N,100°93′ E),平均海拔3100 m以上。气候类型隶属于大陆高原寒温带季风气候,试验期间年均降水量为424.82 mm,多为地形雨且集中于5—9月。全年无绝对无霜期,仅以冷热两季区分,年均温为-1.69℃(图1)。最冷月平均气温为-24.8℃,最热月平均气温为12.5℃。草地类型为草原化草甸,土壤类型为高山草甸土[27]。主要植物有天山针茅(Stipatianschanica)、矮生嵩草(Kobresiahumilis)、星毛委陵菜(Potentillaacaulis)和赖草(Leymussecalinus)等。

1.2 研究方法

1.2.1样地设置

于2014年6月选取基况一致的冬春季放牧高寒草地建立放牧样地,样地小区设置采取完全随机区组设计,分别设置了牦牛单独放牧、藏羊单独放牧、牦牛和藏羊不同比例混合放牧以及围封对照,每个处理均设有3个重复,小区面积及放牧家畜比例详见表1。自放牧样地建立之日起,于每年高寒草地生长季(6—10月)进行放牧,放牧强度水平控制在中等放牧强度,以牧草利用率为50%左右界定,每月放牧时间约为10 d,其余时间将家畜转移至临近草地进行放牧[28]。

表1 放牧试验设计

1.2.2植物样品采集方法

植被调查工作于2020年8月15日进行,每个小区随机调查3个50 cm×50 cm的样方,分别记录每个样方中出现的物种种类、盖度(针刺法测定)、高度(每种植物测量5次求平均值作为原始数据)、地上生物量(齐地面分种刈割后在烘箱中70℃下烘干至恒重,随后用精度为0.001 g的分析天平称量),每个处理合计调查9个样方,所有处理共计调查54个样方。

1.2.3土壤样品采集方法

土壤样品采集工作于样方调查结束后进行,在刈割后的样方中,首先使用体积为100 cm3的环刀采集0—15 cm土层土壤样品,用于土壤容重(Bulk density,BD)与土壤含水量(Soil moisture,SM)的测定；随后使用直径为7 cm的土钻采集0—15 cm土层土壤样品,每个样方中取3钻混合为一个样本。将采集好的样品置于阴凉通风处阴干后,剔除多余杂质并碾碎过不同孔径筛(0.15、0.25、1.40 mm)用于以下化学性质的测定：土壤全氮(Total nitrogen,TN)；土壤全磷(Total phosphorus,TP)；土壤全碳(Total organic carbon,TOC)；土壤速效氮(Available nitrogen,AN)；土壤速效磷(Available phosphorus,AP)；土壤速效钾(Available potassium,AK)；土壤pH。所有指标测定均参照鲍士旦[29]《土壤农化分析》中方法进行。

1.3 数据分析

首先,计算每个物种的重要值(Important Value,IV)大小,IV采用如下公式计算[30]：

IV=(RHi+RCi+RABi)/3

式中,RH为物种i在群落中的相对高度；RC为物种i在群落中的相对盖度；RAB为物种i在群落中的相对生物量。每个样方中所有物种的重要值和为1,因此,具有较高重要值的物种会在群落中占据一定优势。与此同时,还将每个样方的植被按照其经济型划分为禾本科(Gramineae)、莎草科(Cyperaceae)、豆科(Leguminosae)以及杂类草(Forbs),并分别计算各处理下各功能群的重要值。为进一步探讨植物群落变化特征,基于物种重要值进行非度量多维尺度分析(Non-metric multidimensional scaling,NMDS)分析,NMDS分析在Vegan包[31]中进行。

随后,对各处理下的植物群落进行网络分析。进行网络分析时,将每个处理下出现的所有物种纳入,用于构建该处理下网络图的节点；基于物种重要值计算两两物种间的spearman相关性,用于构建该处理下网络图的边[32]。为了使构建的网络图清晰明了,在进行边的构建时剔除显著性水平大于0.05且相关性系数绝对值小于0.8的边。基于上述的节点和边,使用Gephi 0.9.2绘制各处理下的无向随机网络图。随机网络中,特征向量中心度大小可以衡量节点在网络中的重要性[33],因此,本研究中具有最高中心度的节点(物种)即为该处理下植被群落中的关键种。平均连通度大小能够衡量群落中物种的互作强度,平均连通度越高说明植被群落复杂度越高、稳定性越强[34]。平均聚类系数大小能够衡量群落组织水平,平均聚类系数越大说明群落组织水平越高[35]。

最后,计算每个处理下关键种重要值与其土壤理化性质间的Spearman相关性,通过adespatial包[36]的forward.sel()函数执行前向选择以降低共线性,使用lm()函数拟合环境变量与关键种重要值的多元线性回归以获取最优模型,再使用relaimpo包[37]的calc.relimp()函数进行方差分解,求得各环境变量对关键种重要值变化的相对重要性,随后使用ggplot2包[38]进行绘图。对各处理间盖度、物种丰富度、地上生物量、功能群重要值、关键种重要值、土壤理化性质进行单因素方差分析。上述所有统计分析、绘图均在R 4.0.3(R Development Core Team,2020)中进行。

2 结果

2.1 不同放牧家畜组合下植物群落特征变化

不同放牧处理下,NG的禾本科重要值显著高于其他的放牧处理(P<0.05)；MG 1∶6的莎草科重要值显著高于NG(P<0.05),但与其他放牧处理下的莎草科重要值差异不显著；豆科与杂类草重要值在各放牧处理下差异均不显著(表2)。同一放牧处理下,YG、SG的豆科重要值显著低于其他科重要值(P<0.05)；MG1∶2、MG1∶4、MG1∶6的各功能群重要值大小排序为杂类草、莎草科>禾本科>豆科；NG的各功能群重要值大小排序为禾本科>莎草科、杂类草>豆科(表2)。由此可以看出,放牧会显著降低禾本科的重要值,增加莎草科的重要值。

不同放牧处理下,NG的盖度以及高度显著高于其它放牧处理(P<0.05)；物种丰富度以及多样性指数在各处理间差异均不显著；地上生物量则沿着YG、SG、MG1∶2、MG1∶4、MG1∶6的方向呈现出逐渐递减的趋势(表3)。

表2 不同放牧家畜组合下功能群重要值大小

表3 不同放牧家畜组合下高度、盖度、物种丰富度、地上生物量变化

图2 不同放牧家畜组合下植被群落非度量多维尺度(NMDS)分析Fig.2 A non-metric multidimensional scaling (NMDS) ordination of the plant community among different herbivore assemblage YG：牦牛单牧 Yak grazing alone；SG：藏羊单牧 Tibetan sheep grazing alone；MG1∶2：牦牛、藏羊1∶2混牧 Yak and Tibetan sheep 1∶2 mixed grazing；MG1∶4：牦牛、藏羊1∶4混牧 Yak and Tibetan sheep 1∶4 mixed grazing；MG1∶6：牦牛、藏羊1∶6混牧Yak and Tibetan sheep 1∶6 mixed grazing；NG：围封 No grazing

NMDS结果分析表明,YG、SG、MG1∶4、MG1∶6的植物群落结构特征较为相似,而MG1∶2的植被群落结构则与NG的较为相似。相似的群落结构在一定程度上意味着,MG1∶2与NG这两个处理下的植被群落间互作强度、群落稳定性以及组织水平均较为一致(图2)。

2.2 不同放牧家畜组合下植物群落网络结构变化

通过网络分析发现(图3),在不同放牧条件下,NG的植被群落网络共有28个节点91条边,表现出较高水平的群落间关联性；MG1∶2的植被群落网络共有19个节点65条边,群落间关联性次之；YG、SG、MG1∶6、MG1∶4的植被群落网络结构较为相似,群落间关联性也较为一致,分别拥有19个节点28条边、18个节点26条边、17个节点32条边、17个节点27条边。

不同放牧家畜组合下群落间平均连通度与平均聚类系数如图4所示,MG1∶2拥有最高的平均连通度,表明其群落复杂度最高、稳定性最强；NG的平均连通度次之,但依旧明显高于其它放牧处理下的平均连通度；YG、SG、MG1∶4、MG1∶6的平均连通度大小差异不大。不同放牧处理下的平均聚类系数大小排序如下,MG1∶2>MG1∶6>MG1∶4>NG>SG>YG,说明混牧能够明显提高植被群落的组织水平与群落稳定性,而围封相较于牛、羊单牧也能较好的提高植被群落的组织水平与群落稳定性。

图3 不同放牧家畜组合下植被群落网络分析Fig.3 The plant community networks under different herbivore assemblageStipa tianschanica：天山针茅；Stipa sareptana：紫花针茅；Leymus secalinus：赖草；Poa annua：早熟禾；Koeleria macrantha：洽草；Agropyron cristatum：扁穗冰草；Elymus nutans：垂穗披碱草；Kobresia humilis：矮生嵩草；Carex aridula：干生薹草；Oxytropis subfalcata：白花棘豆；Astragalus membranaceus：斜茎黄芪；Medicago archiducis-nicolai：青海苜蓿；Thermopsis lanceolata：披针叶黄华；Potentilla acaulis：星毛委陵菜；Potentilla bifurca：二裂委陵菜；Potentilla multicaulis：多茎委陵菜；Dracocephalum heterophyllum Benth：异叶青兰；Artemisia scoparia：猪毛蒿；Artemisia frigida：冷蒿；Taraxacum mongolicum：蒲公英；Sibbaldia adpressa：伏毛山莓草；Ajuga lupulina：白苞筋骨草；Plantago asiatica：车前；Euphorbia fischeriana：狼毒；Gentiana squarrosa：鳞叶龙胆；Heteropappus hispidus：阿尔泰狗娃花；Thalictrum aquilegifolium var. sibiricum：唐松草；Gentianopsis paludosa：湿生扁蕾；Geranium sibiricum：鼠掌老鹳草；Achnatherum splendens：芨芨草；Gentiana dahurica：达乌里秦艽；Achnatherum inebrians：醉马草；Rhodiola bupleuroides：柴胡；Leontopodium leontopodioides：香芸火绒草

图4 不同放牧家畜组合下植被群落平均连通度与平均聚类系数的变化Fig.4 Changes of the average connectivity and average clustering coefficient of the plant community among different herbivore assemblage

2.3 不同放牧家畜组合下植物群落关键种变化

通过网络分析中各处理下物种中心度大小可得(图3),YG与MG1∶2下的关键种为矮生嵩草；SG下的关键种为天山针茅；MG1∶4下的关键种为星毛委陵菜；MG1∶6与NG下的关键种为赖草(图5)。

图5 不同放牧家畜组合下的关键种变化Fig.5 The shifting of keystone species among different herbivore assemblage

不同放牧处理下天山针茅与星毛委陵菜的重要值差异均不显著；NG的矮生嵩草重要值显著低于MG1∶6的重要值(P<0.05),但与其它放牧处理下的重要值差异均不显著；MG1∶2与NG的赖草重要值显著高于除SG外的其它放牧处理(P<0.05),YG、SG、MG1∶6间赖草重要值差异不显著(表4)。

2.4 不同放牧家畜组合下植物群落关键种变化的土壤理化驱动因子

各放牧处理下的土壤全氮含量与土壤pH值大小差异均不显著；SG的土壤全磷含量显著高于NG和MG1∶4下的土壤全磷含量(P<0.05)；MG1∶4、NG的土壤全碳含量显著高于MG1∶2和YG(P<0.05)；SG的土壤速效氮含量显著高于其他放牧处理下的土壤速效氮含量(P<0.05),YG、MG1∶4和MG1∶6的土壤速效氮含量显著高于MG1∶2下的土壤速效氮含量,而NG的土壤速效氮含量则显著低于其他放牧处理下的土壤速效氮含量(P<0.05)；SG的土壤速效磷含量显著高于其它放牧处理,而MG1∶4和NG的土壤速效磷含量则显著低于其它放牧处理(P<0.05)；NG的土壤速效钾含量显著高于其它放牧处理(P<0.05)；MG1∶6和YG的土壤含水量显著高于其它放牧处理(P<0.05)；除MG1∶4外,SG的土壤容重大小显著高于其它放牧处理,而NG的土壤容重大小则显著低于其它放牧处理(P<0.05)(表5)。

表4 不同放牧家畜组合下关键种重要值大小

表5 不同放牧家畜组合下土壤理化性质变化

通过相关性分析、多元线性回归以及方差分解发现(图6),土壤速效氮含量、土壤含水量以及土壤容重是驱动不同放牧家畜组合下关键种重要值变化的关键环境变量。SG下土壤速效氮与土壤容重对关键种重要值变化的贡献尤为明显,且其关键种重要值大小与土壤速效氮含量呈正相关,与土壤容重大小呈负相关；而NG下土壤速效氮含量、土壤含水量以及土壤容重对关键种重要值变化的贡献均较低,且其关键种重要值大小与土壤速效氮含量呈负相关,与土壤含水量以及土壤容重近乎无相关关系。

图6 不同放牧家畜组合下土壤因子对关键种重要值变化的贡献Fig.6 Contributions of soil physic-chemical properties to the keystone species important value among different herbivore assemblageTN：全氮 Total nitrogen；TP：全磷 Total phosphorus；TOC：全碳 Total organic carbon；AN：速效氮 Available nitrogen；AP：速效磷 Available phosphorus；AK：速效钾 Available potassium；pH：氢离子浓度指数 Hydrogen ion concentration；SM：土壤含水量 Soil moisture；BD：土壤容重 Bulk density

3 讨论

3.1 植物群落特征变化与放牧家畜组合的联系

诸多研究表明,与围栏封育相比,中等强度放牧对物种丰富度的影响并不显著[39—41]。本研究中,连续7年的不同放牧方式下,各放牧处理物种丰富度与禁牧处理相比较差异也不显著,说明在中等放牧强度下,不同放牧家畜组合对物种丰富度也不会造成显著影响。究其原因,不同放牧家畜组合下,植被群落长期受选择性采食的影响,其生活型与繁殖方式也会出现分异与权衡。Tilman[42]提出的资源比率/异质性假说认为,相较于匀质环境的单一限制性资源,异质性环境中的多种限制性资源比率会使生态位分化更加丰富。相较于禁牧单一的资源限制,放牧时家畜施加于草地的胁迫会提高资源比率/异质性,虽然会发生物种的更替,但是其总丰富度将保持在一个较为稳定的状态。此外,土壤种子库的存在也能够使群落组成变化应对外界干扰响应时起到有效的缓冲作用[43]。物种丰富度的变化同时受到多种因素的共同调控,致使其在不同放牧家畜组合下并未呈现出显著差异。

不同家畜的长期选择性采食会改变物种间的竞争关系,进而改变植物群落组成[9]。本研究中,NMDS与网络分析均表明牦牛藏羊1∶2混合放牧下的植物群落结构与围封处理下的植物群落结构有较高的相似性。Gao等[44]的研究表明,对于较为湿润的高寒草地而言,中等放牧强度能够有效增加植物群落间物种的互作强度,这与本研究结果一致。物种间互作强度的增大能够明显提高群落的稳定性[45],而较为稳定的丰富度与相似的群落组织水平表明植被群落间的竞争状况也较为相似,因此,单从群落稳定性方面而言,牦牛藏羊1∶2混合放牧是青藏高原高寒草地较为理想的利用方式。

3.2 驱动关键种变化的土壤理化因子与放牧家畜组合的联系

土壤理化性质作为调控高寒草地生态系统植被生长发育的关键生态因子,亦是维系其群落稳定性的重要因素之一[3]。放牧过程中,家畜主要通过采食、践踏以及排泄物返还3个途径共同作用于草地生态系统。其中,选择性采食会改变群落物种间的竞争关系；践踏则能作用于草地多组分且效果持久,对土壤容重的作用尤为明显[46]；排泄物返还会导致草地局部含水量与氮素的空间变异[47]。本研究中,驱动不同放牧家畜组合下关键种变化的主要理化因子分别为土壤速效氮、含水量以及容重,这与放牧过程中家畜对草地的直接作用关系密切。相较于单独放牧,混合放牧能够显著降低牦牛和藏羊的采食口数与活动量[9],而藏羊的活动量本就大于牦牛,因此,藏羊频繁的践踏活动导致藏羊单牧处理下的容重显著高于其它处理。另外,藏羊单牧的土壤速效氮含量之所以最高是因为藏羊排泄物氮含量本就较高,加之藏羊活动频繁,排泄物分布更加分散且更易分解,致使藏羊单牧处理下的土壤速效氮含量增加。柴锦隆等[46]研究表明,家畜的践踏作用能够显著促进氮的矿化作用,这也是为什么速效氮对禁牧处理下关键种变化的贡献最低的原因。Ganjurjav等[48]研究表明,土壤水分是青藏高原中部高寒草地植被生长的限制性因素。本研究中,土壤水分在所有放牧处理下均是驱动关键种变化的物理因子。围封处理下,土壤水分对植被关键种变化的贡献最低,这是由于长期禁牧导致土壤容重降低,表层土壤变得疏松,从而导致蒸发强度的增加[49]。同时,围封会导致禾本科植被增多,而禾本科植被相较于其它科植被而言水分利用量及利用效率更高[50],进而导致围封处理下的土壤含水量低于其它放牧处理。

3.3 青藏高原高寒草地利用方式的综合考量

一直以来,青藏高原高寒草地的退化均被归因于过度放牧[51],也正因为如此,围栏封育作为较为经济有效的恢复措施被广泛应用于青藏高原退化高寒草地的修复中[52]。然而,许多研究结果表明,青藏高原高寒草地的恢复水平与围栏封育时间并非是线性关系[53]。与此同时,也有研究表明,在关键期(如春季或秋季)进行临时禁牧或轮牧能够更加有效的促进植被生产力的恢复[54]。因此,合理的利用方式能够保障且有利于维持高寒草地生态系统的健康。本研究中,不同放牧处理下,各处理的关键种重要值并非最高,这也说明优势种不一定是关键种。而矮生嵩草通常被认为是青藏高原高寒草地的关键种[55],因此,从关键种演替角度出发,结合各处理下群落的稳定性与组织水平,牦牛藏羊1∶2混合放牧依旧是青藏高原高寒草地较为有效的利用方式。此外,随着社会与经济的发展,虽然游牧这种青藏高原高寒草地传统且可持续放牧方式的生产效率已经大打折扣,但依旧有其借鉴意义[56]。鉴于此,我们建议利用高寒草地时,将关键时期休牧或轮牧与最佳混牧比例相结合。

4 结论

(1)中等放牧强度时,不同比例的牦牛、藏羊混牧组合均会显著降低植物的高度、盖度以及禾本科的重要值,显著增加莎草科的重要值,但对物种多样性并不会产生显著影响。

(2)中等放牧强度时,从植物群落结构稳定性角度出发,NG与MG 1∶2下的植物群落组织水平与稳定性均最佳。因此,牦牛藏羊1∶2混牧是青藏高原高寒草地中等放牧强度下较为理想的利用方式。

(3)中等放牧强度时,随着牦牛、藏羊混牧比例的不同,植物群落关键种会发生明显演替。其中,YG与MG 1∶2下的关键种为矮生嵩草；SG下的关键种为天山针茅；MG 1∶4下的关键种为星毛委陵菜；MG 1∶6与NG下的关键种为赖草。土壤速效氮、含水量、容重是驱动关键种演替的重要土壤理化因子。