EN中文

地理科学 >

2023 , Vol. 43 >Issue 3: 398 - 410

DOI: https://doi.org/10.13249/j.cnki.sgs.2023.03.003

青藏高原土地利用变化及其生态环境效应专栏

基于地理加权随机森林的青藏地区放牧强度时空格局模拟

  • 李兰晖 , 1, 2 ,
  • 黄聪聪 1 ,
  • 张镱锂 , 2 ,
  • 刘林山 2 ,
  • 王兆锋 2 ,
  • 张海燕 2 ,
  • 丁明军 3 ,
  • 张华敏 3
展开
  • 1.厦门理工学院计算机与信息工程学院，福建 厦门 361024
  • 2.中国科学院地理科学与资源研究所陆表格局与模拟重点实验室，北京 100101
  • 3.江西师范大学地理与环境学院，江西 南昌 330022
张镱锂。E-mail:

李兰晖（1989—），男，江西赣州人，副教授，主要从事人类活动对地表覆被变化的影响。E-mail:

收稿日期: 2022-06-20

  修回日期: 2022-08-20

  网络出版日期: 2023-03-20

基金资助

中国科学院战略性先导科技专项(XDA20040201)

第二次青藏高原综合科学考察研究(2019QZKK0603)

国家自然科学基金项目(42101099)

版权

版权所有，未经授权，不得转载、摘编本刊文章，不得使用本刊的版式设计。
收起

Mapping the multi-temporal grazing intensity on the Qinghai-Tibet Plateau using geographically weighted random forest

  • Li Lanhui , 1, 2 ,
  • Huang Congcong 1 ,
  • Zhang Yili , 2 ,
  • Liu Linshan 2 ,
  • Wang Zhaofeng 2 ,
  • Zhang Haiyan 2 ,
  • Ding Mingjun 3 ,
  • Zhang Huamin 3
Expand
  • 1. School of Computer and Information Engineering, Xiamen University of Technology, Xiamen 361024, Fujian, China
  • 2. Key Laboratory of Land Surface Pattern and Simulation, Institute of Geographic Sciences and Natural Resources Research, Chinese Academy of Sciences, Beijing 100101, China
  • 3. School of Geography and Environment, Jiangxi Normal University, Nanchang 330022, Jiangxi, China

Received date: 2022-06-20

  Revised date: 2022-08-20

  Online published: 2023-03-20

Supported by

Strategic Priority Research Program of Chinese Academy of Sciences(XDA20040201)

the Second Tibetan Plateau Scientific Expedition and Research(2019QZKK0603)

National Natural Science Foundation of China(42101099)

Copyright

Copyright reserved © 2023.
Fold

摘要

精确刻画放牧强度时空格局是深刻理解青藏地区高寒草地服务功能动态及其驱动机制的重要基础，对辅助制定区域生态安全和可持续发展战略具有意义。基于牲畜年末存栏量、牧区人口密度、生长季NDVI、年降水量、年平均气温和居民点分布等数据，采用地理加权随机森林模型，模拟了2000年、2010年和2020年3期青藏地区放牧强度的空间格局，并定量评价了环境因子对放牧强度空间分布解释性的区域差异。结果表明：① 地理加权随机森林模拟结果可精细地刻画青藏地区放牧强度的空间特征，与经典随机森林模型相比，判定系数更高，平均绝对值误差和均方根误差更低。② 青藏地区放牧强度呈现东南高、西北低的基本特征，其中，西北部地区放牧强度低于25羊单位/km2的区域约占青藏地区面积的1/2。③ 与2000年和2010年相比，2020年青藏地区放牧强度呈现总体下降、局部抬升的态势；其中，较2010年，2020年放牧强度下降超过1羊单位/km2的区域占牧区面积的61.69%。④ 牧区人口密度是解释放牧强度空间异质性最主要的因素，其相对重要性呈现西高东低的特征，而降水量和生长季NDVI的相对重要性则呈现西北高、东南低的特征。研究结论可为青藏地区草地可持续管理和生态安全屏障建设提供科学参考。

关键词： 放牧强度; 地理加权随机森林; 牧区人口密度; 青藏地区

本文引用格式

李兰晖 , 黄聪聪 , 张镱锂 , 刘林山 , 王兆锋 , 张海燕 , 丁明军 , 张华敏 . 基于地理加权随机森林的青藏地区放牧强度时空格局模拟[J]. 地理科学, 2023 , 43(3) : 398 -410 . DOI: 10.13249/j.cnki.sgs.2023.03.003

Abstract

Accurately quantifying the spatiotemporal pattern of grazing intensity on the Qinghai-Tibet Plateau (QTP) is crucial to improving our understanding of the driving mechanism of alpine grassland change, and is of great significance for maintaining regional ecological security and promoting sustainable development policies. Based on the data of livestock inventory at the end of the year and environmental covariates (e.g. population density in pastoral areas, growing season NDVI, annual precipitation, annual mean temperature, and settlements), the gridded grazing intensities on the QTP in 2000, 2010, and 2020 were simulated by geographically weighted random forest (GRF), and the regional differences in the interpretability of environmental variables were then analyzed. The results showed that grazing intensity maps predicted by the GRF could mirror the spatial distribution of grazing intensity on the QTP, compared with the classic random forest model, the R2 was higher, and both the mean absolute value error (MAE) and root mean square error (RMSE) were lower. The grazing intensity was generally higher in the southeast and lower in the northwest of the QTP. The areas with grazing intensity of less than 25 sheep units/km2 in the northwest part accounted for about half of the QTP. Compared with 2000 and 2010, the grazing intensity of the QTP in 2020 showed a trend of overall decrease but local increase. For example, compared with 2010, the areas of grazing intensity decreased by higher than 1 sheep unit/km2 in 2020 accounting for 61.69% of the pastoral area. The population density in pastoral areas was the most important factor explaining the spatial heterogeneity of grazing intensity, and its relative importance was higher in the western QTP and lower in the eastern QTP. On the contrary, the relative importance of both precipitation and growing season NDVI was higher in the northwest part and lower in the southeast part of QTP. Our results provide scientific references for sustainable grassland management and ecological safety barrier construction on the QTP.

Key words： grazing intensity; geographically weighted random forest; population density in pastoral areas; the Qinghai-Tibetan Plateau

草地约占陆地表面的1/3，维持着全球数亿人的生计[1-2]。放牧活动是人类利用草地最直接方式[3]，不仅关系到畜产品的生产和供应，也影响着全球碳循环和生物多样性[4]。放牧强度是反映放牧活动对草地扰动程度的重要指标[5]，单位面积上的牲畜数量能有效地表征草地的放牧强度[6]，往往被作为衡量放牧强度的代用指标[4]。随着人口的快速增加，畜产品需求不断上升，草地的放牧强度随之增大，全球大量草地因过度放牧导致了不同程度的退化[7-8]，进而影响到草地生态服务功能的稳定供给[2,6]。因此，掌握详细的放牧强度变化信息对于草地资源的合理利用、生态系统影响评估以及畜牧业可持续管理具有重要意义。
传统上，牲畜数量按行政范围（单元）统计，难以详细地反映放牧强度的空间信息，也难以与网格尺度上的生态环境空间数据集成[6,9]。鉴于此，联合国粮农组织（FAO）推出了世界网格化牲畜数量（Gridded Livestock of the World，GLW）项目[10]，并在2014年和2018年分别采用多层回归和随机森林（Random forest，RF）模型生产了空间分辨率约为1 km×1 km（GLW2）和10 km×10 km（GLW3）的全球牲畜网格化数据产品[11-12]。近些年，放牧强度空间化研究也逐渐受到国内学者的关注。有学者基于植被生产力或人口密度单一因子数据对全国或局部地区的放牧强度进行了空间化研究[13-14]；也有学者基于人口密度、地表覆盖类型、降水量等多因子数据，采用支持向量机和RF等模型对中国西北部牲畜数量进行了空间化，并指出了各环境因子的相对重要性[9]。在当前的研究中，环境因子由单变量到多变量，模型由线性回归模型到非线性的机器学习模型，模拟精度逐渐提高。然而，这些研究使用的模型均为全局校准的，未顾及到放牧强度与环境因子之间的空间异质性，即未考虑各环境因子的相对重要性在空间上的差异。为处理变量的空间异质性问题，Georganos等[15]在RF模型的基础上，借鉴了地理加权回归的思路，提出了局部校准的地理加权随机森林（Geographical random forests，GRF）模型。在已有的GRF模型应用研究中，如人口数量和小麦产量空间分布模拟，其性能均比RF模型更佳[15-16]；并且GRF还能有效地展示环境因子解释性的区域差异[17-19]。这为进一步精确模拟放牧强度的时空格局以及评估其环境因子的重要性提供了新的契机。
青藏地区包括青海省和西藏自治区，作为青藏高原的主体，是国家生态安全屏障的关键地区[20]，有近2/3的面积为草地，是全球著名的放牧生态系统之一[21]。放牧活动是青藏地区最主要的人类活动类型之一[22]。21世纪以来，放牧强度分布及其变化已成为青藏地区生态安全评估和草地可持续研究重要内容[23-24]。已有研究将青藏地区县域牲畜数量统计数据直接用于分析近几十年来放牧强度的变化特征[22,25]，或将牲畜统计资料简单地平均到行政单元范围内后，作为量化人类活动强度的重要因子[26]以及用于分析放牧活动对生态环境的影响[27]，然而，这些研究难以刻画青藏地区县域内放牧强度的空间差异。也有学者通过修正FAO的GLW数据后开展青藏及其周边地区的放牧强度研究[28-29]，然而，在中国范围内，由于该网格数据是基于中国2001年的牲畜量统计值研制而成，需通过插值等外推方式获得其他年份数据[28,30]。
综上所述，虽然青藏地区放牧强度数据具有广泛的用途，但缺乏详细的时空格局信息。因此，有必要增强青藏地区放牧强度空间模拟研究，并进一步掌握放牧强度空间分布的影响因素。本研究基于2000年、2010年和2020年青藏地区牲畜年末存栏量统计数据及牧区人口、年降水量和年平均气温等12个环境因子数据，对比RF和GRF模型在放牧强度方面的模拟精度，并生成3期高精度的放牧强度公里网格数据。在此基础上，分析青藏地区放牧强度的空间格局，并探讨环境因子的相对重要性及其区域分异特征。研究结论可为青藏地区草地可持续管理和生态安全屏障建设提供科学参考。

1 数据来源与方法

1.1 数据来源和预处理

本研究使用的数据包括2000年、2010年和2020年的县级牲畜年末存栏量统计数据、人口密度等栅格数据以及保护区等矢量数据[11-12,31-35]（表1）。为了实现放牧强度空间化辅助数据空间分辨率的一致性，使用双线性插值法将气温和降水等环境因子数据重新采样到1 000 m空间分辨率，并转为Albers_Conic_Equal_Area投影。
表1 数据来源汇总

Table 1 Summary of the datasets in this study

数据名称 时间粒度（年份） 数据类型/分辨率 来源/参考文献
　　注：区县、地级市和自治州牲畜年末存栏统计数据来源于《西藏统计年鉴》和《青海统计年鉴》（ https://data.cnki.net/Yearbook/）及当地政府的统计数据；资源环境科学与数据中心网址 https://www.resdc.cn/；“–”为无此项。
县级牲畜年末存栏量 2000，2010，2020 Excel 统计数据
乡镇牲畜年末存栏量（部分） 2010 Excel 统计数据
地表覆盖类型 2000，2010，2020 1 000 m Liu等[31]、Zhang等[32]
夜间灯光指数 2000，2010，2020 15 arcsec （～500 m） Chen等[33]
人口密度 2000，2010，2020 1 000 m Li等[34]
到居民点距离 2000，2012，2016 1 000 m Li等[34]
居民点密度 2000，2012，2016 1 000 m Li等[34]
到道路的距离 2000，2012，2016 1 000 m Li等[34]
河流 – 矢量 资源环境科学与数据中心
保护区 – 矢量 资源环境科学与数据中心
行政区划 – 矢量 资源环境科学与数据中心
DEM – 1 000 m 资源环境科学与数据中心
生长季NDVI均值 2000—2016 1 000 m 资源环境科学与数据中心
年平均气温均值 2000—2016 0.1° He等[35]
年降水量均值 2000—2016 0.1° He等[35]
GLW2牲畜密度 2001 ~1 000 m Robinson等[12]
GLW3牲畜密度 2001 ~10 000 m Gilbert等[11]
区县、地级市和自治州牲畜年末存栏统计数据来源于《西藏统计年鉴》和《青海统计年鉴》（ https://data.cnki.net/Yearbook/）及当地政府的统计数据。由于青海省各县牲畜年末资料已收集的年份截至2017年，为保持时间系列的一致性，采用2020年《青海统计年鉴》的地级市和自治州牲畜年末存栏量对各县数据进行校准。由于乡镇牲畜年末存栏数据收集难度较大，收集的乡镇数据涉及2010年西藏阿里地区5县（日土县、普兰县、札达县、革吉县和噶尔县）、那曲市2县（安多县和聂荣县）以及青海的泽库县，这些乡镇大多分布于青藏地区的西北部，乡镇牲畜年末存栏数据仅用于模拟结果的辅助验证。计算放牧强度时，由于统计数据未涉及牲畜重量和年龄，在参照NY/T 3647-2020草食家畜羊单位换算标准和相关研究[36]的基础上，将各类牲畜折算为羊单位进行简化换算，即1绵羊/山羊=1羊单位，1大牲畜（牛和马等）=5羊单位。然后，将羊单位数据匹配至对应的县域/乡镇行政区划矢量数据中。需要指出的是，青藏地区的牧区存在季节性放牧和常年放牧等多种放牧方式，为保持可比性，以一个完整年为单位对放牧强度进行统计。
数据预处理主要包括两个步骤：一是提取草地面积比重、计算牧区人口密度以及将矢量数据转为栅格数据等；二是在这基础上，采用县级行政区划矢量数据对所有环境因子进行区域平均值统计。预处理结束后，将进行模型训练与预测以及精度评估。放牧强度空间化模拟流程如图1所示。
图1 青藏地区放牧强度空间化模拟流程图

Fig. 1 Flowchart illustrating the algorithm used for mapping grazing intensity on the Qinghai-Tibet Plateau

由于牲畜数量分布与农村人口密度分布高度相关[14]，将牧区人口密度作为环境因子有助于提升放牧强度空间模拟的精度[37]，而获取高精度的牧区人口密度数据仍面临挑战[38]。首先，需要获得精确可靠的人口密度栅格数据。2000年和2010年青藏地区人口密度栅格数据来自文献[34]，该数据基于公开发布的2000年和2010年青藏地区乡镇（街道）的人口普查数据和地表覆被类型等因子数据，采用随机森林模型研发而成；基于研究区各县网站公布的2020年第七次人口普查公报中的乡镇人口数据（其中，丁青县等20个县的乡镇人口数据缺失）和同期土地覆被类型等因子数据，采用文献[34]相同的流程研制2020年人口密度栅格数据。其次，基于研究区野外调查和牧户走访发现，牲畜分布的区域一般需满足2个条件：① 有人口分布的草地区域，即无人区的动物则视为野生动物；② 非城区，即城区为非牧区，城区人口一般不是牧民，即不从事放牧活动。基于上述两个条件，以夜间灯光指数等于0表征非城区，以人口密度大于0表征有人口分布区。其次，基于30 m空间分辨率的土地覆被类型数据提取草地分布数据，再使用ArcGIS的聚合工具生成1 km×1 km像元中草地的面积比重数据，以草地面积比重大于0表征有草地分布的区域。最后，综合这3个指标提取牧区人口密度，进一步提取牧区人口分布的草地区域（即放牧活动区）。本研究将县域放牧强度限定在放牧活动区开展空间模拟。
基于道路、河流和居民点矢量数据，采用ArcGIS空间分析中的欧氏距离分析工具，分别生成到道路距离、到河流距离和到居民点距离的栅格数据，再采用核密度工具生成居民点密度栅格数据[34]。考虑到保护区内核心区、缓冲区和试验区保护要求的差异，分别赋值3、2和1，非保护区则赋值为0，并转为栅格数据。基于DEM数据获得坡度数据，并基于2000—2016年气候数据以及NDVI数据，获得年平均气温、年降水量和生长季NDVI多年平均值。主要环境因子的空间分布如图2所示。
图2 2010年青藏地区主要环境因子的空间分布特征

Fig. 2 Spatial distribution of nine environmental factors for grazing intensity on the Qinghai-Tibet Plateau in 2010

1.2 随机森林模型

经典的随机森林（RF）是由Breiman[39]提出的一种基于分类回归树的非参数机器学习模型，能从数据中学习复杂的非线性关系和交互特征。该模型可用于样本数量相对较少的高维变量，对多重共线性不敏感，也不易发生欠拟合和过拟合现象。在模型中，构建分类回归树的样本量约占训练数据集的2/3，剩下的1/3称为袋外（Out-of-Bag，OOB）数据。这些袋外数据则用于内部交叉验证建立无偏估计，以测评模型的性能，也用于评估变量相对重要性的得分[40]。该分值越大，其重要性越高，对因变量的解释性也越强。

1.3 地理加权随机森林模型

地理加权随机森林（GRF）是经典RF在地理学领域的扩展[15]。GRF的主要思想类似于地理加权回归（GWR），是局部校准的。GWR是通过允许局部估计将普通线性回归进行改造升级而成，已被广泛应用于空间数据分析[41]；Georganos等[15]则将类似的原理应用于非线性的RF，改进升级为GRF。具体而言，GRF模型是通过引入数据观测点的空间信息生成空间权重矩阵，然后和RF集成到局部回归分析框架中[19]。与全局的RF相比，GRF可提供更卓越的模拟性能，也可作为处理因变量和自变量空间异质性的有效工具。在本研究中，对于每个县域行政单元i，仅用单元i相邻单元的放牧强度统计值来运行一个RF，即构建出一个子模型，该子模型提供其预测性能以及各环境因子重要性得分（%IncMSE）等信息。GRF的简化公式如下：
$ {y}_{i}=a({lat}_{i},{lon}_{i}){x}_{i}+e,i=1, 2, 3\cdots n $
式中， $ {y}_{i} $ 是第i个县域的放牧强度， $ a({lat}_{i},{lon}_{i}){x}_{i} $ 是基于位置i上校准的RF预测值， $ ({lat}_{i},{lon}_{i}) $ 是县域行政单元i的质心坐标，e是误差项，n为县级行政单元的数量。
本文借鉴GWR原理[41]，基于交叉验证法（CV）自动确定模型的邻近点数量。在预测时，将局部子模型中提取的局部异质性信号（低偏差）与使用更多数据的全局模型（低方差）融合，可增强模型的预测效果[15]。鉴于研究区东南部和西北部县域面积的较大差异，本文采用权重系数融合全局和局部估计的方法[15]，经多次调试，权重系数设置为0.5时获得较理想的模拟效果。更多GRF信息可参考R语言中“SpatialML”包[42]。
最后，将公里网格的环境因子数据输入模型，生成公里网格的放牧强度权重模拟值，再利用该权重将县级放牧强度统计值映射至每个像元，获得青藏地区公里网格放牧强度预测值，公式如下：
$ {I}_{p}=\dfrac{{I}_{c}\times {W}_{p}} {{W}_{c}} $
式中，Ip为像元放牧强度预测值，Ic为县域放牧强度统计值，二者分别是公里网格的放牧强度预测值和县域单元的放牧强度统计值；Wp为像元权重模拟值，Wc为县域权重模拟值，二者分别是GRF输出的像元权重模拟及其在县域单元的总和。

1.4 模拟精度评估

为对比GRF和RF模型的拟合优度和整体模拟性能，拟采用以下3个稳健的误差测量指标：判定系数（R2）、平均绝对误差（MAE）、以及均方根误差（RMSE）：
$ {R}^{2}=1-\dfrac{{\displaystyle\sum }_{i}{({y}_{i}-\hat {{y}_{i}})}^{2}}{{\displaystyle\sum }_{i}{({y}_{i}-\bar {y})}^{2}} $
$ MAE=\dfrac{{\displaystyle\sum }_{i=1}^{n}\left|{y}_{i}-\hat {{y}_{i}}\right|}{n} $
$ RMSE=\sqrt{{\displaystyle\sum} _{i=1}^{n}\dfrac{{({y}_{i}-\hat {{y}_{i}})}^{2}}{n}} $
式中， $ {y}_{i} $ 是行政单元i的放牧强度统计值， $ \hat {{y}_{i}} $ 是行政单元i各像元的权重模拟值之和， $ \bar {y} $ 是各行政单元的放牧强度平均值。

2 研究结果

2.1 模拟精度评估

基于县域统计值的GRF模型训练时的评估结果显示，2000年、2010年和2020年青藏地区放牧强度与环境因子数据拟合结果的判定系数（R2）均高达0.97。同时，基于袋外（OOB）数据交叉验证建立的无偏估计得出，2000年、2010年和2020年放牧强度与环境因子数据的判定系数（OOBR2）分别达到0.83、0.83和0.82，这表明当前模型中环境因子对放牧强度的解释度可达82%。
基于像元尺度的预测结果显示，GRF预测的放牧强度权重模拟值与县域放牧强度统计值具有较好的一致性，R2达到0.75（图3）。其中，2000年、2010年和2020年分别为0.78、0.75和0.75，均高于RF模拟结果，并且回归系数也更倾向于1。进一步统计显示，GRF预测的放牧强度权重模拟值平均绝对值误差（MAE）和均方根误差（RMSE）均低于RF预测结果（图3）。这些结果表明，GRF模型能有效提升放牧强度空间模拟的精度。因此，在此基础上，采用公式（2）将县级放牧强度统计值映射至每个像元，获得2000年、2010年和2020年青藏地区公里尺度的放牧强度预测值，并用于后续的分析。
图3 2000年、2010年和2020年青藏地区县域放牧强度权重模拟值与统计值对比

Fig. 3 Comparisons of predicted grazing intensity weights and statistical data at the county level on the Qinghai-Tibet Plateau in 2000, 2010 and 2020

为进一步验证预测结果的精度，将2010年放牧强度预测值与安多县和日土县等地区乡镇统计值进行对比，结果显示，乡镇预测值与统计值的回归系数达到0.84，R2达到0.79（图4），说明当前公里网格的放牧强度数据具有较高的精度。
图4 2010年乡镇放牧强度预测值与统计值对比

Fig. 4 Comparison of predicted grazing intensity and statistical data at the town level in 2010

2.2 青藏地区放牧强度空间分布和变化特征分析

基于GRF模型的结果精细地刻画了2000年、2010年和2020年青藏地区放牧强度东南高、西北低的分布特点（图5）。放牧强度低于25羊单位/km2的地区约占全区面积的3/5，主要分布于青藏地区的西北部以及藏东南森林区，其中，研究区西北部的面积约占全区面积的1/2，并且大致形成一条以25羊单位/km2为界的放牧强度东北–西南分界线。该分界线主要位于250~350 mm等降水量分布区，该线西北和东南两侧的人口和居民点分布、降水量、NDVI和草地类型（高寒草甸和高寒草原）均呈现明显的差异。在该线的西北部，近1/5的地区无牲畜或极少牲畜分布，主要分布在羌塘–可可西里–阿尔金山无人区。在该线东南侧，超过100羊单位/km2的地区约占全区的1/4，集中分布于青海省东部和西南部地区、西藏的一江两河地区（雅鲁藏布江中游及其支流年楚河、拉萨河中部流域地区）以及那曲、昌都等地区。
图5 2000—2020年青藏地区放牧强度的空间分布格局

Fig. 5 Spatial pattern of grazing intensity on the Qinghai-Tibet Plateau from 2000 to 2020

基于3个年份放牧强度空间分布数据，分析了2000—2020年青藏地区放牧强度的变化特征（图6）。2000—2010年和2010—2020年放牧强度上升超过1羊单位/km2的区域面积分别占放牧活动区面积的46.24%和31.31%，前者明显高于后者；而下降超过1羊单位/km2的区域面积分别占放牧活动区面积的30.86%和61.69%，后者约为前者的2倍（图6 a-b）。具体而言，2000—2010年和2010—2020年上升超过10羊单位/km2的区域分别占放牧活动区面积的21.60%和12.94%，前者主要分布在放牧强度超过100羊单位/km2的地区，如青海湖以南地区和玉树南部地区、西藏的一江两河地区以及那曲、昌都等地区；而后者则分布在青海湖周边地区和西藏的那曲部分地区（如尼玛县和申扎县）。2000—2010年和2010—2020年下降超过10羊单位/km2的区域分别占放牧活动区面积的8.95%和24.02%，后者约为前者的3倍，前者主要分布在青海的三江源南部地区以及西藏零星地区，后者不仅分布在上述地区，还分布在西藏的南部和东部大片地区。此外，2010—2020年下降5~10羊单位/km2的区域占比12.74%，约为2000—2010年的2倍。由上述分析可知，放牧强度由2000—2010年增强态势转为2010—2020年的减弱态势，这可能与2011年6月以来政府有效实施草原生态保护补助奖励机制等政策有关[43]。
图6 2000—2020年青藏地区放牧强度的空间变化特征

Fig. 6 Spatial pattern of the change in grazing intensity on the Qinghai-Tibet Plateau from 2000 to 2020

相比于2000年，2020年放牧强度增幅和降幅大于1羊单位/km2的区域面积分别占放牧活动区面积的39.69%和52.20%（图6c）。除青海省西南部和西藏一江两河地区等局部区域与2000—2010年相近外，2000—2020年青藏地区放牧强度上升和下降的分布区域大致与2010—2020年相近。上述结果表明，近20 a来青藏地区放牧强度大致呈现总体下降、局部上升的态势。

2.3 青藏地区放牧强度空间分布影响因子分析

基于县域统计值，采用RF和GRF模型均可得到环境因子的重要性得分（%IncMSE）。从图7可看出，环境因子重要性得分（RF）或得分的平均值（GRF）从高到低的排序中，牧区人口密度均位于第一位，并且降水量、生长季NDVI、到道路的距离、居民点密度和到居民点距离等5个环境因子的重要性得分均位居前六。RF和GRF模型评估的环境因子重要性得分排序也存在差异，如在GRF中，到居民点距离和居民点密度的重要性得分由RF模型中的第六和第五位分别转为第二和第三位，而降水量、生长季NDVI则由第二和第三位分别转为第四和第六位。可见，在全区尺度，研究结果验证了青藏地区由牧民管理的牲畜分布呈现“逐水草而居”的特征，在更细的尺度，也能反映出放牧强度的空间分布与居民点和道路等社会因素的密切程度。
图7 RF模型中青藏地区放牧强度各环境因子的重要性和GRF模型中各环境因子的局部重要性均值

Fig. 7 Variable importance of grazing intensity on the Qinghai-Tibet Plateau in RF model and averaged local importance of per factor in GRF model

在空间上，牧区人口密度对放牧强度空间分布解释性最高的县域占比近2/3，主要分布在青藏地区的中东部（图8）。到居民点距离作为解释性最高变量的县域占15.65%，主要分布在西藏西部地区，尤其是阿里地区；居民点密度也达到13.91%；二者累计占比约为1/3，这可能与该地区人口稀少、居民点分散有关。综合来看，牧区人口密度和居民点分布对青藏地区放牧强度分布起关键作用。
图8 GRF模型中重要性得分排第一位的环境因子空间分布特征

Fig. 8 Spatial pattern of importance score of key environmental factors in GRF model

将重要性得分前六位的环境因子空间可视化（图9），结果显示牧区人口密度、到居民点距离和居民点密度等3个变量的重要性得分大致呈现西高东低的特征，其中，在青藏地区中西部地区，这3个变量的重要性得分普遍高于15。降水量、生长季NDVI的重要性得分总体呈现西北高东南低的特征，其中东南地区的重要性得分大多低于10，这与降水量和NDVI东南偏高和西北偏低的特征相反（图2）一定程度上反映了青藏高原西北地区降水少，制约植被生长，进而成为影响放牧活动的重要因子。
图9 GRF模型中主要环境因子重要性得分（%IncMSE）的空间分布特征

Fig. 9 Spatial pattern of the importance scores of main environmental factors （%IncMSE） in GRF model

3 讨论

3.1 模拟精度对比

本文基于县域统计的均值数据，采用GRF自带的袋外（OOB）数据无偏估计得出放牧强度与环境因子数据的判定系数（OOBR2）超过0.82，高于Li等[9]基于RF和深度神经网络模型对中国西北地区牲畜数量空间模拟时的判定系数（R2），也接近甚至高于Reed等[44]基于RF模型开展人口空间化研究中的解释度。同时，将GRF模拟的放牧强度权重值与县域统计值拟合时的判定系数（R2）达到0.75，均优于RF的预测结果（图3），也接近甚至高于Li等[9]模拟中国西北部牲畜数量空间分布的判定系数（R2）。在以往人口空间模拟研究中，若人口数量的预测权重值与县域统计值之间的判定系数（R2）达到或高于0.70左右，能较好地刻画人口数量的空间特征[45-46]。此外，经公式（2）计算后得出的2010年放牧强度公里网格数据与乡镇统计值的回归系数达到0.84（图4）。上述结果说明GRF能较好地模拟青藏地区放牧强度的空间特征。
为进一步验证本研究的结果，与GLW2和GLW3放牧强度（各类牲畜数量转为羊单位）数据进行对比[11-12]（图3、图5、图10和图11）。在青藏地区，GLW2和GLW3放牧强度预测值与县域统计值拟合的判定系数（R2）分别为0.06和0.12，明显低于本研究的结果（图10和图3a），这与本研究采用更为详尽的环境因子数据以及采用更优化的模型有关。GLW产品主要基于气候、海拔、NDVI和人口密度作为自变量，而本研究除这些因子外，引入了居民点和道路等与放牧密切相关的因子[5]。与GLW产品的全球或区域尺度模型相比，本文不仅研究尺度更小，而且采用了GRF模型，更充分地考虑了区域差异性（图9），进而有助于提升预测数据的精度。在空间上，3套数据存在明显差异。虽然三者均显示放牧强度的高值区集中分布于青海东部和西藏一江两河等地区，但在这些地区，GLW2和GLW3放牧强度预测值要高于本研究结果（图5a和图11）。由于在模型训练过程时是采用行政单元的统计值，而预测则采用像元值，难以顾及局部区域的牲畜承载总量。本文先将牛羊等牲畜类型转换为放牧强度（羊单位），再进行训练和预测，而GLW2和GLW3分别直接模拟牛、羊等牲畜数量，这可能导致放牧强度在高值区更高，反之亦然。此外，GLW2数据显示青藏地区中西部地区的放牧强度多为0；而GLW3数据则显示几乎所有地区都大于100羊单位/100 km2；本研究利用人口密度和草地分布等信息提取放牧活动区，避免了雅鲁藏布江流域上游等大量牧区模拟出无放牧强度，也避免了在藏北无人区模拟出放牧强度的错误现象。综上所述，与GLW数据相比，本研究结果能更精确地刻画青藏地区放牧强度的空间分布特征。
图10 2001年GLW2和GLW3放牧强度县域预测值与统计值的对比（数据引自Robinson等[12]和Gilbert等[11]）

Fig. 10 Comparisons of predicted grazing intensity and statistical data at the county level and their spatial distribution in GLW2 and GLW3 in 2001 (dataset from Robinson et al [12] and Gilbert et al[11])

图11 2001年GLW2和GLW3放牧强度的空间分布特征（数据引自Robinson等[12]和Gilbert等[11]）

Fig. 11 Spatial distribution of predicted grazing intensity in GLW2 and GLW3 in 2001 (dataset from Robinson et al [12] and Gilbert et al[11])

3.2 牧区人口密度是牲畜数量空间分布关键因子

基于RF和GRF模型获得的环境因子重要性得分均表明，牧区人口密度是最重要的环境因子（图7和图8），这与Li等[9]得出人口密度是牦牛空间分布的主要因素的结论基本一致。与Li等[9]的研究相比，本研究基于夜间灯光数据和草地分布比重来剔除城区人口，提取牧区人口密度，这有利于减少放牧强度在人口稠密的城区被高估而在人口稀少的牧区被低估的现象[37]。此外，由于牲畜年末数据多为县域统计值，而目前公开发布的人口密度数据往往也基于县域统计值模拟生成，如果直接使用该人口密度数据，可能影响空间模拟效果。本文采用乡镇级人口统计值模拟生成更高精度人口密度栅格数据[34]，有助于提高放牧强度空间模拟精度。
除牧区人口密度这一因子外，到居民点距离、居民点密度以及到道路的距离等与人类活动有关的环境因子位居重要性得分的前六位，也是部分区域放牧强度空间分布的关键环境因子（图7和图8），因此，合理地引入与放牧强度紧密相关的人类活动因子，有利于提高放牧强度空间模拟精度。此外，放牧强度空间分布影响因子的相对重要性在空间上存在明显的差异（图9），表明这些因子对放牧强度空间分布的影响是非线性的，在局部区域使用的线性赋值方法拓展到大尺度放牧强度空间化时需谨慎[5]。

4 结论与展望

4.1 结论

基于2000年、2010年和2020年青藏地区牲畜年末存栏量统计数据和12个环境因子数据，GRF模拟的结果能精确地刻画青藏地区放牧强度的时空特征及其环境因子相对重要性的区域差异。青藏地区放牧强度的基本特点是“东南多、西北少”，其中，低于25羊单位/km2的地区集中分布于青藏地区西北部，约占全区面积的1/2。过去20 a放牧强度呈现总体下降、局部上升的态势，且由2000—2010年的放牧强度增强态势转为2010—2020年的减弱态势，在一定程度上反映了草原生态保护补助奖励机制等政策实施的成效。牧区人口密度是放牧强度空间分布解释性最高的环境因子，与降水量、生长季NDVI、到道路的距离、居民点密度和到居民点距离位居重要性得分的前六位，且这些环境因子的相对重要性呈现明显区域差异。

4.2 局限性和展望

在青藏高原普遍存在冬季牧场和夏季牧场转场的现象，并且自2004年以来，青藏高原逐步实施退牧还草等一系列生态工程，然而，受牧场、禁牧区和育肥场等数据资料的限制，本研究尚未考虑草畜平衡面积[5]、牧场转场以及育肥场[14]的情况，这可能导致当前放牧强度模拟在局部地区与真实场景存在偏差。在青藏地区，县级行政区划的面积相差较大，尤其是西北部地区，其中尼玛县（2012年行政区划调整为尼玛县和双湖县）和改则县的面积超过13万km2，相当于福建省等东部省级行政区的面积。随着输入单元的面积增加，空间模拟精度也随之下降[11,47]，因此，青藏地区西北部等地区的模拟精度低于总体水平。本文选择的乡镇样本大多位于西北部，这意味着青藏地区放牧强度空间分布的总体精度可能高于图4的结果。此外，西藏地区的尼玛县和申扎县放牧强度持续增强，与周边县域之间出现明显界限，也说明西北部地区放牧强度的模拟精度仍需进一步提高。
虽然本文顾及空间异质性的GRF模型的预测效果高于RF模型，但仍以截面数据作为输入变量，仅研究了青藏地区2000年、2010年和2020年3期的放牧强度，难以全面反映研究区放牧强度的连续变化。最近研究还指出深度学习有助于提高放牧强度空间模拟的精度[9]。随着深度学习算法的发展，顾及空间异质性的时空加权深度神经网络算法已运用于环境要素的空间模拟[48]，有望提高放牧强度空间模拟的精度，并用于研制长时间连续系列的栅格数据，以进一步分析放牧强度的时空变化特征及其影响因素。

参考文献

原文顺序 | 文献年度倒序 | 文中引用次数倒序
[1]
Scurlock J M O, Hall D O. The global carbon sink: A grassland perspective[J]. Global Change Biology, 2002, 4(2): 229-233.

[2]
Bardgett R D, Bullock J M, Lavorel S et al. Combatting global grassland degradation[J]. Nature Reviews: Earth & Environment, 2021, 2(10): 720-735.

[3]
文雯, 周丁扬, 李广泳, 等. 基于牧民放牧行为的青海湖流域放牧路线研究——以刚察县为例[J]. 干旱区资源与环境, 2016, 30(7): 95-100.

DOI

Wen Wen, Zhou Dingyang, Li Guangyong et al. The research on grazing routes in Qinghai Lake watershed based on the grazing behavior of herdsmen—A case of Gangcha County. Journal of Arid Land Resources and Environment, 2016, 30(7): 95-100.

DOI

[4]
Hankerson B R, Schierhorn F, Prishchepov A V et al. Modeling the spatial distribution of grazing intensity in Kazakhstan[J]. PLOS ONE, 2019, 14(1): e210051

[5]
胡晓阳, 王兆锋, 张镱锂, 等. 青藏高原放牧强度空间化方法与应用[J]. 地理学报, 2022, 77(3): 547-558.

DOI

Hu Xiaoyang, Wang Zhaofeng, Zhang Yili et al. Spatialization method of grazing intensity and its application in Tibetan Plateau. Acta Geographica Sinica, 2022, 77(3): 547-558.

DOI

[6]
Dara A, Baumann M, Freitag M et al. Annual landsat time series reveal post-Soviet changes in grazing pressure[J]. Remote Sensing of Environment, 2020, 239: 111667.

DOI

[7]
Hilker T, Natsagdorj E, Waring R H et al. Satellite observed widespread decline in Mongolian grasslands largely due to overgrazing[J]. Global Change Biology, 2014, 20(2): 418-428.

DOI

[8]
李文龙, 薛中正, 郭述茂, 等. 基于3S技术的玛曲县草地植被覆盖度变化及其驱动力[J]. 兰州大学学报(自然科学版), 2010, 46(1): 85-90,95.

Li Wenlong, Xue Zhongzheng, Guo Sumao et al. Vegetation coverage changes and analysis of the driving forces in Maqu county based on 3S technology. Journal of Lanzhou University (Natural Science), 2010, 46(1): 85-90,95.

[9]
Li X, Hou J, Huang C. High-resolution gridded livestock projection for Western China based on machine learning[J]. Remote Sensing, 2021, 13(24): 5038.

DOI

[10]
Robinson T P, Franceschini G, Wint W. The food and agriculture organization's gridded livestock of the world[J]. Veterinaria Italiana, 2007, 43(3): 745

[11]
Gilbert M, Nicolas G, Cinardi G et al. Global distribution data for cattle, buffaloes, horses, sheep, goats, pigs, chickens and ducks in 2010[J]. Scientific Data, 2018, 5: 180227.

DOI

[12]
Robinson T P, Wint G R W, Conchedda G et al. Mapping the global distribution of livestock[J]. PLOS ONE, 2014, 9(5): e96084.

DOI

[13]
Sun Y, Liu S, Liu Y et al. Grazing intensity and human activity intensity data sets on the Qinghai-Tibetan Plateau during 1990–2015[J]. Geoscience Data Journal, 2022, 9(1): 140-153.

DOI

[14]
Zhu Z, Zhang X, Dong H et al. Integrated livestock sector nitrogen pollution abatement measures could generate net benefits for human and ecosystem health in China[J]. Nature Food, 2022, 3(2): 161-168.

DOI

[15]
Georganos S, Grippa T, Niang Gadiaga A et al. Geographical random forests: A spatial extension of the random forest algorithm to address spatial heterogeneity in remote sensing and population modelling[J]. Geocarto International, 2021, 36(2): 121-136.

DOI

[16]
Khan S N, Li D, Maimaitijiang M. A geographically weighted random forest approach to predict corn yield in the US Corn Belt[J]. Remote Sensing, 2022, 14(12): 2843.

DOI

[17]
Alampi Sottini V, Barbierato E, Capecchi I et al. Assessing the perception of urban visual quality: An approach integrating big data and geostatistical techniques[J]. Aestimum, 2022, 79: 75-102.

DOI

[18]
Grekousis G, Feng Z, Marakakis I et al. Ranking the importance of demographic, socioeconomic, and underlying health factors on US COVID-19 deaths: A geographical random forest approach[J]. Health & Place, 2022, 74: 102744

[19]
Luo Y, Yan J, Mcclure S. Distribution of the environmental and socioeconomic risk factors on COVID-19 death rate across continental USA: A spatial nonlinear analysis[J]. Environmental Science and Pollution Research, 2021, 28(6): 6587-6599.

DOI

[20]
孙鸿烈, 郑度, 姚檀栋, 等. 青藏高原国家生态安全屏障保护与建设[J]. 地理学报, 2012, 67(1): 3-12.

DOI

Sun Honglie, Zheng Du, Yao Tandong et al. Protection and construction of the national ecological security shelter zone on Tibetan Plateau. Acta Geographica Sinica, 2012, 67(1): 3-12.

DOI

[21]
贺有龙, 周华坤, 赵新全, 等. 青藏高原高寒草地的退化及其恢复[J]. 草业与畜牧, 2008(11): 1-9.

He Youlong, Zhou Huakun, Zhao Xinquan et al. Alpine grassland degradation and its restoration on Qinghai-Tibet Plateau. Prataculture & Animal Husbandry, 2008(11): 1-9.

[22]
吴雪, 刘峰贵, 刘林山, 等. 青藏高原牲畜数量变化及其空间特征[J]. 生态科学, 2021, 40(6): 38-47.

Wu Xue, Liu Fenggui, Liu Linshan et al. Changes and spatial characteristics of livestock population in Tibetan Plateau. Ecological Science, 2021, 40(6): 38-47.

[23]
Chen B, Zhang X, Tao J et al. The impact of climate change and anthropogenic activities on alpine grassland over the Qinghai-Tibet Plateau[J]. Agricultural and Forest Meteorology, 2014, 189–190: 11-18.

[24]
Zhang L, Fan J, Zhou D et al. Ecological protection and restoration program reduced grazing pressure in the Three-River Headwaters Region, China[J]. Rangeland Ecology & Management, 2017, 70(5): 540-548.

[25]
王小青. 2000-2018年青藏高原放牧活动时空变化及草畜平衡分析[D]. 西安: 长安大学, 2021: 72.

Wang Xiaoqing. Analysis on the temporal and spatial changes of grazing activities and the balance between forage and livestock in the Qinghai-Tibet Plateau from 2000 to 2018. Xi'an: Chang'an University, 2021: 72.

[26]
Li S, Zhang Y, Wang Z et al. Mapping human influence intensity in the Tibetan Plateau for conservation of ecological service functions[J]. Ecosystem Services, 2018, 30: 276-286.

DOI

[27]
Wei Y, Lu H, Wang J et al. Dual influence of climate change and anthropogenic activities on the spatiotemporal vegetation dynamics over the Qinghai-Tibetan Plateau from 1981 to 2015[J]. Earth's Future, 2022, 10(5): e2021E-e2566E.

[28]
张文强, 罗格平, 郑宏伟, 等. 基于随机森林模型的内陆干旱区植被指数变化与驱动力分析: 以北天山北坡中段为例[J]. 植物生态学报, 2021, 44(11): 1-14.

Zhang Wenqiang, Luo Geping, Zheng Hongwei et al. Analysis of vegetation index changes and driving forces in inland arid areas based on random forest model: A case study of the middle part of northern slope of the north Tianshan Mountains. Chinese Journal of Plant Ecology, 2021, 44(11): 1-14.

[29]
Han Q, Luo G, Li C et al. Simulated grazing effects on carbon emission in Central Asia[J]. Agricultural and Forest Meteorology, 2016, 216: 203-214.

DOI

[30]
段群滔, 罗立辉. 1990—2015年青藏高原人类足迹数据集[J]. 中国科学数据(中英文网络版), 2020, 5(3): 303-312.

Duan Quntao, Luo Lihui. A dataset of human footprint over the Qinghai-Tibet Plateau during 1990—2015. China Scientific Data, 2020, 5(3): 303-312.

[31]
Liu J, Kuang W, Zhang Z et al. Spatiotemporal characteristics, patterns, and causes of land-use changes in China since the late 1980s[J]. Journal of Geographical Sciences, 2014, 24(2): 195-210.

DOI

[32]
Zhang X, Liu L, Chen X et al. GLC_FCS30: global land-cover product with fine classification system at 30 m using time-series Landsat imagery[J]. Earth System Science Data, 2021, 13(6): 2753-2776.

DOI

[33]
Chen Z, Yu B, Yang C et al. An extended time series (2000–2018) of global NPP-VIIRS-like nighttime light data from a cross-sensor calibration[J]. Earth System Science Data, 2021, 13(3): 889-906.

DOI

[34]
Li L, Zhang Y, Liu L et al. Mapping changing population distribution on the Qinghai-Tibet Plateau since 2000 with multi-temporal remote sensing and point-of-interest data[J]. Remote Sensing, 2020, 12(24): 4059.

DOI

[35]
He J, Yang K, Tang W et al. The first high-resolution meteorological forcing dataset for land process studies over China[J]. Scientific Data, 2020, 7(1): 1-11.

DOI

[36]
李青丰, 乌云其木格. 草畜平衡管理系列研究(6)——家畜数量动态变化[J]. 草业科学, 2012, 29(6): 997-1001.

Li Qingfeng, Wu Yunqi, Mu Ge. Animal demanding-feed availability balance control series(6)—Animal number monthly change. Pratacultural Science, 2012, 29(6): 997-1001.

[37]
Nicolas G, Robinson T P, Wint G R W et al. Using random forest to improve the downscaling of global livestock census data[J]. PLOS ONE, 2016, 11(3): e150424

[38]
Leyk S, Gaughan A E, Adamo S B et al. The spatial allocation of population: A review of large-scale gridded population data products and their fitness for use[J]. Earth System Science Data, 2019, 11(3): 1385-1409.

DOI

[39]
Breiman L. Random forests[J]. Machine Learning, 2001(4545): 1-5.

[40]
Gromping U. Variable importance assessment in regression: linear regression versus random forest[J]. The American Statistician, 2009, 63(4): 308-319.

DOI

[41]
Fotheringham A S, Brunsdon C, Charlton M. Geographically weighted regression: The analysis of spatially varying relationships[M].Hoboken: Wiley, 2002.

[42]
Georganos S, Kalogirou S. Package 'Spatialml'[DB/OL]. 2019. https://mran.microsoft.com/snapshot/2020-04-25/web/pack-ages/SpatialML/index.html

[43]
Ye T, Zhao N, Yang X et al. Improved population mapping for China using remotely sensed and points-of-interest data within a random forests model[J]. Science of The Total Environment, 2019, 658: 936-946.

DOI

[44]
Reed F, Gaughan A, Stevens F et al. Gridded population maps informed by different built settlement products[J]. Data, 2018, 3(3): 33.

DOI

[45]
林珲, 张鸿生, 林殷怡, 等. 基于城市不透水面–人口关联的粤港澳大湾区人口密度时空分异规律与特征[J]. 地理科学进展, 2018, 37(12): 1644-1652.

DOI

Lin Hun, Zhang Hongsheng, Lin Yinyi et al. Spatiotemporal changes of gridded urban population in the Guangdong-Hong Kong-Macao Greater Bay Area based on impervious surface-population correlation. Progress in Geography, 2018, 37(12): 1644-1652.

DOI

[46]
Tan M, Li X, Li S et al. Modeling population density based on nighttime light images and land use data in China[J]. Applied Geography, 2018, 90: 239-247.

DOI

[47]
Zhuang H, Liu X, Yan Y et al. Mapping multi-temporal population distribution in China from 1985 to 2010 using Landsat images via deep learning[J]. Remote Sensing, 2021, 13(17): 3533.

DOI

[48]
Wu S, Wang Z, Du Z et al. Geographically and temporally neural network weighted regression for modeling spatiotemporal non-stationary relationships[J]. International Journal of Geographical Information Science, 2021, 35(3): 582-608.

DOI

Options
文章导航

/

〈 〉