CN112907113A - 一种考虑空间相关性的植被变化成因识别方法 - Google Patents

Abstract

本发明公开了一种考虑空间相关性的植被变化成因识别方法。包括：步骤一、收集影响植被生长的因子数据集，包括自然因子、气象因子和社会因子，以及NDVI数据；步骤二、数据处理；步骤三、利用空间分析理论对不同区域影响NDVI时间变化的主控因子进行筛选和识别；步骤四、构建地理探测器模型，分析影响因子的空间分异性对NDVI空间分布的影响，并确定适宜植被生长的各影响因子范围；步骤五、利用残差分析方法定量区分自然因素和人类活动因素对植被变化的影响。本发明方法可以在考虑不同区域因子空间相关性的条件下，准确识别影响植被时空变化的主控因素，定量描述自然因素和人类活动对植被变化的影响，为区域植被恢复和生态保护提供合理依据。

Description

一种考虑空间相关性的植被变化成因识别方法

技术领域

本发明涉及植被变化影响因素动态识别领域，具体涉及一种考虑空间相关性的植被变化成因识别方法。

背景技术

地表植被是陆地生态系统的重要组成部分，对区域生态环境具有重要作用。植被的生长过程受到多种因素的协同作用，同时植被也能够在短时间内对周围环境的变化作出响应。

遥感技术的发展，为研究植被变化的驱动因素提供了助力，归一化因子植被指数(NDVI)是最常用的指标之一，它与植被生产力和植被覆盖度有密切的正相关关系。基于连续的长系列遥感气象数据、土地利用数据等，全球气候条件、人类活动等因素与植被变化的关系也得到了深入研究。但是影响植被生长的因素众多，各因素对植被生长的影响力不同，因此在众多因素中准确识别出主控因素也存在一定难度。

不同区域气象条件、自然地理特征和社会经济等要素存在差异，但由于大气系统的连续性，气象因子在不同的空间位置的变化存在互相依赖关系，即空间自相关；同时，由于因子的协同作用，不同因子在空间上也存在相关性，即空间互相关。因子的空间相关对植被变化具有重要影响，但目前关于这方面的研究报道很少，相关研究有待进一步加强。

发明内容

本发明的目的在于提供一种识别影响植被生长主控因素的方法，旨在解决定量确定各因素及其时空分异规律对植被生长影响的问题。

在空间相关分析的基础上进行主控因子识别可以更加准确地反映植被变化的区域差异及区域相关特性。同时影响因素的空间异质性对植被生长的空间分布具有重大影响，且因素之间的协同作用对植被生长的作用也不可忽视。因此有必要在识别主控因素的基础上，综合考虑各因素的时空异质性以及多因子的协同效应对植被生长的影响。

为实现上述目的，本发明提供如下技术方案：

一种考虑空间相关性的植被变化成因识别方法，包括以下步骤：

步骤一、数据收集：收集下载反映植被覆盖状况的归一化植被指数(NDVI)数据；收集下载植被变化的影响因子数据集，包括自然因子、气象因子、社会因子；植被变化以归一化植被指数(NDVI) 反映。

步骤二、数据前处理：数据前处理包括缺失数据的插补延长、数据空间插值和重采样、数据转换；缺失数据利用线性内插或是线性外延方法进行插补延长；数据的空间插值采用反距离加权平均方法，并利用GIS软件中的Create Fishnet功能按1km分辨率在空间上生成均匀分布的采样点，以此获取采样点位置处的数据序列；数据转换是指根据研究需要对数据进行如对数化、无量纲化或比例化等处理。

利用Arcmap软件生成采样点，采用反距离加权平均方法对收集下载的网格数据集进行重采样和插值处理，获取采样点位置处的数据序列；对数据进行相应的转换如对数化、无量纲化或比例化，在后续的计算过程中根据情况选择原始数据或是转换后的数据。

根据情况选择原始数据或是转换后的数据是指相关性分析中当相关系数总体较小或不易获取时，则利用转换后的数据，否则利用原始数据。

步骤三、主控因子识别：

1)计算各影响因子和NDVI数据序列的相关系数，并利用t统计量检验相关关系是否显著，剔除不显著的的因子，初步筛选影响NDVI 变化的影响因子种类；

该步为因子初筛，利用相关系数的过滤算法，构建影响因素影响力判别模型。根据判别模型，初步确定影响植被变化的影响因子种类，其中影响力未通过显著性检验的影响因子将被剔除；通过筛选的影响因子自动组成一个新的集合，作为下一步骤的输入。

2)基于逐步回归方法对影响NDVI变化的影响因子进行分析，进一步减少影响因子的个数。

3)根据影响因子在空间上的聚集和离散情况对研究区进行分区，根据空间相关分析理论分析各子区域内影响因子的空间相关性，以确定不同区域适用的空间回归模型，分析不同区域内影响因子与NDVI 的相关性，进一步识别影响植被生长的主控因素。

步骤四、因子空间分异分析：

利用地理探测器模型对影响因子进行空间分异性分析，量化单个因子的具体影响、多个因子之间的交互作用，确定影响因子的空间分异对植被空间变化的影响，确定适合植被生长的影响因子范围。

步骤五、定量识别人类活动的影响：

分离自然因素、气象因素和人类活动因素对NDVI变化的影响。将植被变化分为未明显受到人类活动影响及明显受到人类活动影响两个时期，以未明显受到人类活动影响时期的NDVI数据及自然因子和气象因子数据逐像元建立回归模型，并用该模型计算明显受到人类活动影响时期的NDVI预测值，计算NDVI真实值和预测值之间的残差，以此量化人类活动对植被变化的贡献率。

进一步的优化方案，步骤一种所述自然因子包括高程、坡度、坡向、土壤属性、植被类型；所述气象因子包括气温、气压、相对湿度、降水、蒸发、风速风向、日照、辐射；所述社会因子包括城市建成区面积、土地利用变化、造林面积。

进一步的，步骤二中数据转换的方法包括对数化、无量纲化及比例化，其计算公式分别为：

对数化：

以10为底，对NDVI及影响因素取对数：

X_i＝logx_i (1)

以e为底，对NDVI及影响因素取对数：

X_i＝lnx_i (2)

其中x为变量变换前的数值，X为变换后的数值；i表示第i个变量，i＝1,...,n；n为因子个数；

无量纲化：

其中i、x、X的意义同前文所述；j表示变量的第j个值，j＝1,...,l， l为各自变量的序列长度；

比例化：

X_i＝αx_i (4)

其中i、x、X的意义同前文所述，α为比例系数。

进一步的，步骤三1)中t统计量计算公式为：

其中r为相关系数，l为序列长度；当P＜0.05时，则表明通过显著性检验。

进一步的，步骤三2)建立植被NDVI和影响因素之间的回归模型，通过逐步回归方法进一步剔除影响较小的因素；利用F统计量对因素影响力(回归系数)的显著性进行检验；所述的逐步回归方法具体计算步骤为：

首先，建立回归自变量X₁，...，X_n同因变量Y的一元回归模型 (其中n为自变量的个数)：

Y＝β₀+β_iX_i+ε(i＝1，…，n) (6)

其中β₀，...，β_n为回归系数，ε为残差，i表示第i个变量；计算各自变量回归系数的F统计量，记为

取其最大值记为

若

超过显著性水平临界值，则将自变量X_i1引入回归模型，记 I₁为选入变量指标集合；建立因变量Y与自变量子集 {X_i1，X₁}，{X_i1，X₂}，…，{X_i1，X_n}的二元回归模型；此时自变量回归系数的F统计量记为

最大值记为

若

超过显著性水平临界值，则将自变量X_i2引入回归模型；否则，停止引入变量；重复以上步骤，直至无变量引入为止；

进一步的，步骤三3)中所述的空间回归模型表示方式如下：

其中，Y为因变量，如NDVI等；X是自变量，如气温、降水等；ρ为空间回归系数，β为自变量回归系数，W₁和W₂分别是自变量、残差的空间邻接矩阵；ε为残差项；δ为残差项的空间回归系数；μ为随机误差。

进一步的，步骤四中利用地理探测器对影响因子的空间分异性进行分析包括因子影响力探测、因子交互作用探测、风险区探测、生态检测。

因子影响力探测：

SST＝Nσ² (12)

其中，PD为影响因子的解释力；h＝1，...，L反映了影响因子的分级情况，各级分别对应各影响因子的某个范围或类型，h表示第h级， L为总的分级数；N_h和N分别代表了不同级内及全区内的样本个数，

和σ²则分代表了不同级内及全区内对应的方差；SSW为级内方差之和，SST为全区总方差；PD的取值范围为[0，1]，PD值的简单变换满足非中心F分布，PD转换过程及非中心F分布如下：

其中，λ为非中心参数；

为级h的NDVI均值。

因子交互作用探测：

首先，分别计算两个影响因素的解释力PD(a₁)和PD(a₂)；然后，计算因素交互后的影响力PD(a₁∩a₂)；最后比较PD(a₁∩a₂)和 PD(a₁)及PD(a₂)。

风险区探测：

风险区探测通过对影响因子分级，并检测各级对应的NDVI平均值是否有显著差异，来识别出有利于植被生长的影响因子范围或类型。风险区探测基于t统计量。计算步骤如下：

首先，计算各影响因子不同级别下的NDVI的平均值；然后，根据t检验判断各影响因子不同级别下的NDVI平均值之间是否存在显著差异。最后，对应各影响因子，分别找出计算出的NDVI平均值的最大值及其对应的级别，如果该值与其他级别所对应的NDVI平均值有显著差异，则表示该级别所对应的影响因子范围/类型是最合适植被生长的；否则，找出第二大NDVI平均值及其对应级别，并比较它是否与除最大值以外的其他值有显著差异，如果是，则表示该影响因子有两个最合适植被生长的范围/类型；否则，重复此步骤，直到满足条件。t统计量的计算公式为：

其中

表示级h内的NDVI平均值；N_h为级h内的样本数量；var 表示方差。自由度的计算方法为：

在95％的置信水平下，认为两个子级间的NDVI均值存在明显差异。

生态检测：

生态检测用于确定两个影响因子(a₁和a₂)对NDVI分布的影响是否有显著差异(即a₁对NDVI分布的影响是否会比a₂更显著)；F统计量被用于检查SSW_a1和SSW_a2之间的差异，其中SSW_a1和SSW_a2分别表示由a₁和a₂形成的分级的级内方差之和。

进一步的，步骤五步骤如下：

(1)根据调查确定研究区大规模人类活动起始年份；或是利用 MK突变分析检测NDVI突变点，以突变点所在年份作为大规模人类活动的临界年份，假设突变点以前研究区植被覆盖仅受自然因子和气象因子的影响。MK突变分析计算过程为：

E(d_k)＝k(k-1)/4 (18)

var(d_k)＝k(k-1)(2k+5)/72 (19)

其中，d_k为统计量，l为序列长度，m_i表示样本y_i大于y_j(1≤j≤i) 的累积数，i、j分别表示第i、j个样本；E(d_k)和var(d_k)分别为d_k的均值和方差；UF_k为d_k标准化后的值，组成的曲线为UF，用相同的方法计算反序列UB曲线，曲线在置信区间内(＞1.96或＜-1.96)的交点即为突变点。

(2)以突变点以前的NDVI数据及自然因子和气象因子数据逐像元建立回归方程，并用该方程计算突变点后的NDVI预测值，并计算NDVI真实值和预测值之间的残差，且认为该差异是由人类活动导致：

ε＝NDVI_真实-NDVI_预测 (21)

若ε大于0，则表明人类活动对植被变化产生积极影响，否则人类活动产生消极影响。

附图说明

图1为本发明实施例的技术流程图；

图2为本发明实例提供的影响因子自相关分析及因子影响力分析结果示意图；

图3～图9分别为本发明实例提供的影响NDVI空间分布的因子识别结果示意图。图3为降水、图4为气温、图5为日照、图6为高程、图7为坡度、图8为土壤类型、图9为土地利用类型。

具体实施方式：

实施例1

以汉江上游(丹江口水库以上)流域为案例区，以2001-2018年为研究时段对本发明技术方法的具体应用作进一步说明，具体应用包括以下步骤：

步骤一、数据收集：

从中科院中国科学院资源与环境数据云平台(http://www.resdc.cn) 下载NDVI、土壤类型和土地利用数据。NDVI数据的时间序列为2001-2018年，空间分辨率为1km；土壤类型数据的空间分辨率为1km，并采用传统“土壤发生分类”系统，分为12土纲，61个土类，227个亚类；土地利用类型数据的空间分辨率为1km，土地利用分类体系为 6个一级类，25个二级类。30m分辨率的DEM数据从USGS/NASA 航天飞机雷达地形测绘任务(SRTM)数字高程数据库(http://srtm.csi)下载。降水、气温及日照时数数据从中国气象局国家气象信息中心(NMIC)(http://cdc.cma.gov.cn)下载，数据准备见表1：

表1汉江上游植被变化影响因素数据准备

步骤二、数据预处理

采用反距离加权平均法对各气象因子进行空间插值，利用GIS 软件的CreateFishnet功能按照1km的空间分辨率生成采样点，对 NDVI数据及影响因子数据进行重采样，并获取采样点处的因素数据序列。

步骤三、主控因子识别

1)因子初筛：对NDVI和影响因子进行相关分析，通过相关系数的过滤算法初步确定影响因子种类。结果表明降水、气温和日照时数均对植被的时间变化具有显著影响，其中降水与NDVI显著正相关，气温和日照时数与NDVI显著负相关。

2)因子逐步回归：利用初步筛选的因子与植被生长(NDVI)建立回归模型，如果因子的回归系数太小或不易获取，则对数据进行转换后重新建立回归模型；通过逐步回归，进一步剔除对植被生长影响较小的因子。

3)研究区分区及主控因子识别：利用GIS软件的Spatial Autocorrelation(MoransI)功能对影响因子进行自相关分析，根据 Morans’Index和Z-score判断因子的聚集和离散程度(图2)；筛选影响植被生长的主要影响因子，Morans’Index越接近1表示越聚集，越接近-1表示越离散；Z-score越大表示聚集(离散)程度越显著。自相关分析结果表明，可以根据土壤类型和土地利用类型对研究区进行分区。根据空间分析理论分别确定不同区域适用的空间回归模型，分析不同区域影响因素与植被生长的相关性，进一步识别主控因素。

步骤四、因子空间分异性分析

根据各类标准、规范及专家意见对主控因素进行分级(表2)，利用地理探测器对NDVI及分级后的因素进行地理探测，确定单个因子对植被空间变化的解释力(图2)、不同因子的空间分布对植被空间变化规律的影响(图3)、适宜植被生长的因子范围/类型(表3) 及因子交互作用对植被空间变化的影响(表4)。

表2影响因子分级

表3适宜植被生长的因子范围/类型

因子	范围/类型
		PRE	>1100mm
TEM	10-12℃
		SUN	1500-1700hours
ELEV/m	1500-2500m
		SLP/°	>21°
SOIL	淋溶土
		LU	林地

表4因子的交互作用

步骤五、定量识别人类活动影响

利用残差分析分离气候变化和人类活动对植被变化的影响，量化人类活动的影响。结果表明，汉江上游流域NDVI发生突变的年份为 2009年，平均残差为-0.014，说明人类活动使流域NDVI有所下降。人类活动导致NDVI值下降了1.74％，气候变化使NDVI值上升了7.28％；人类活动造成NDVI下降的区域约占流域面积的46％，但是人类活动的影响总体较小。

以上所述仅对本发明的实例实施而已，并不用于限制本发明，本发明中对数据的转换及空间回归模型的选取，可根据需求及具体研究区设定。凡是在本发明的权利要求限定范围内，所做的任何修改、等同替换、改进等，均应在本发明的保护范围之内。

Claims (6)

1.一种考虑空间相关性的植被变化成因识别方法，其特征在于，包括以下几个步骤：

步骤一、数据收集：收集下载NDVI数据；收集下载影响植被变化的影响因子数据集，所述影响因子包括：自然因子、气象因子、社会因子；

步骤二、数据前处理：包括缺失数据的插补延长、数据空间插值和重采样、数据转换中的一项或几项；

步骤三、主控因子识别：

1)计算各影响因子和NDVI数据序列的相关系数，并利用t统计量检验相关关系是否显著，剔除不显著的的因子，初步筛选影响NDVI变化的影响因子种类；

2)基于逐步回归方法对影响NDVI变化的影响因子进行分析，进一步减少影响因子的个数；

3)根据影响因子在空间上的聚集和离散情况对研究区进行分区，根据空间相关分析理论分析各子区域内影响因子的空间相关性，以确定不同区域适用的空间回归模型，分析不同区域内影响因子与NDVI的相关性，进一步识别影响植被生长的主控因子；

步骤四、因子空间差异分析：

利用地理探测器对影响因子的空间分异性进行分析，确定影响因子的空间分异特征对植被空间变化的影响以及适宜植被生长的影响因子范围；

步骤五、定量识别人类活动影响：

将植被变化分为未明显受到人类活动影响及明显受到人类活动影响两个时期，以未明显受到人类活动影响时期的NDVI数据及自然因子和气象因子数据逐像元建立回归模型，并用该模型计算明显受到人类活动影响时期的NDVI预测值，计算NDVI真实值和预测值之间的残差，以此量化人类活动对植被变化的贡献率。

2.根据权利要求1所述的一种考虑空间相关性的植被变化成因识别方法，其特征在于，步骤一中所述自然因子包括高程、坡度、坡向、土壤属性、植被类型；所述气象因子包括气温、气压、相对湿度、降水、蒸发、风速风向、日照、辐射；所述社会因子包括城市建成区面积、土地利用变化、造林面积。

3.根据权利要求1所述的一种考虑空间相关性的植被变化成因识别方法，其特征在于，步骤三2)中所述的逐步回归方法具体计算步骤为：

首先，建立回归自变量X₁，...，X_n同因变量Y的一元回归模型，其中n为自变量的个数：

y＝β₀+β_iX_i+ε(i＝1，…，n)

其中β₀，...，β_n为回归系数，ε为残差，i表示第i个变量；计算各自变量回归系数的F统计量，记为F_i ⁽¹⁾，取其最大值记为

若

超过显著性水平临界值，则将自变量X_i1引入回归模型，记I₁为选入变量指标集合；建立因变量Y与自变量子集{X_i1，X₁}，{X_i1，X₂}，…{X_i1，X_n}的二元回归模型；此时自变量回归系数的F统计量记为F_i ⁽²⁾，最大值记为

若

超过显著性水平临界值，则将自变量X_i2引入回归模型；否则，停止引入变量；重复以上步骤，直至无变量引入为止。

4.根据权利要求1所述的一种考虑空间相关性的植被变化成因识别方法，其特征在于，步骤三3)中所述的空间回归模型表示方式如下：

其中，Y为因变量；X是自变量；ρ为空间回归系数；β为自变量回归系数；W₁和W₂分别是自变量、残差的空间邻接矩阵；ε为残差项；δ为残差项的空间回归系数；μ为随机误差。

5.根据权利要求1所述的一种考虑空间相关性的植被变化成因识别方法，其特征在于，步骤四中利用地理探测器对影响因子的空间分异性进行分析包括因子影响力探测、因子交互作用探测、风险区探测、生态检测。

6.根据权利要求1所述的一种考虑空间相关性的植被变化成因识别方法，其特征在于，步骤五的具体步骤为：A、利用MK突变分析检测NDVI突变点，以突变点所在年份作为大规模人类活动的临界年份，假设突变点以前研究区植被覆盖仅受自然因子和气象因子的影响；B、以突变点以前的NDVI数据及自然因子和气象因子数据逐像元建立回归方程，并用该方程计算突变点后的NDVI预测值，并计算NDVI真实值和预测值之间的残差，该残差大于0，则表明人类活动对植被变化产生积极影响，否则人类活动产生消极影响。

Images