【精选】【MODIS数据处理#11】例六：绘制NDVI多年变化趋势空间分布图_ndvi变化趋势图怎么做_Salierib的博客

二、操作步骤

注：需要利用ArcGIS和python，需要用到的代码见下文，python需要安装依赖库(scipy、numpy、pandas)。

2.1 处理NDVI数据

在这里插入图片描述
处理过程可以参考添加链接描述，我这边准备了近5年的NDVI栅格数据，分辨率为250m。
注意：要保证各年的NDVI栅格的分辨率和坐标系是相同的，这里我用的坐标系是WGS_1984_UTM_Zone_49N

2.2 栅格转点

2.3 多值提取至点

这一步的耗时可能较长

运行结束后打开 输入点要素 的属性表可以看到匹配好的点对

2.4 添加经纬度坐标

工具箱\系统工具箱\data management tools.tbx\要素\添加 xy 坐标

运行结束后，要素的属性表会多出来两列，分别对应x、y坐标，由于我用的是UTM分带投影坐标系，单位是m。

下面将改属性表导出为txt文件

2.5 利用python计算斜率和显著性系数

import pandas as pd
import scipy.stats as stats
读取存储有多年NDVI信息的表格，对于位置[x,y]来说有NDVI构成的时间序列
如：[NDVI_2015,NDVI_2016,NDVI_2017,NDVI_2018,NDVI_2019]
用一元线性回归分析这个时间序列的趋势和显著性
txt_file:str
    待分析的表格的完整路径
vi_names:List[str]
    NDVI时间序列对应的列名组成列表,如['2015', '2016', '2017', '2018', '2019']
out_file:str
    斜率和显著性系数结果的输出路径
txt_file = r"E:\GIS教程\趋势分析\shanxi15-19.txt"
vi_names = ['2015', '2016', '2017', '2018', '2019']
out_file = r"E:\GIS教程\趋势分析\slope and alpha.csv"
vi_years = pd.read_table(txt_file, sep=",",thousands=",")
virows = vi_years[vi_names].values
years = range(len(vi_names))
slopes = []
alphas = []
for virow in virows:
    slope, itc, r, p, *useless = stats.linregress(years,virow)
    slopes.append(slope)
    alphas.append(p)
vi_years.insert(len(vi_years.columns), "slope", slopes)
vi_years.insert(len(vi_years.columns), "alpha", alphas)
vi_years[["POINT_X","POINT_Y","slope","alpha"]].to_csv(out_file)
2.6 转为栅格
 
导入表格文件到csv
 
显示XY数据，注意要和栅格的坐标系设置为一致
 
 
保存得到的点要素
 
 
用点转栅格工具分别将要素的slope和alpha字段转为栅格，参数设置中需要将像元大小和原始的NDVI栅格大小保持一致
 
 至此我们就得到了斜率和显著性系数的栅格影像。
  
三、变化趋势分级
 
可以按照斜率slope和显著性系数alpha的大小将多年NDVI的变化趋势分为以下三类： 
alpha >=0.05;无显著变化,用数字1表示
slope>0,alpha<0.05;显著增加,用数字2表示
slope<0,alpha<0.05;显著降低,用数字3表示 
用栅格计算器完成分级操作，地图代数表达式为： 
Con("alpha.tif">=0.05,1,Con("slope.tif">0,2,3))
                    文章目录一、功能介绍二、脚本代码三、工具参数四、工具界面系列文章目录：一、功能介绍Con((“A2000M06.NDVI.tif”&lt;=0.05),0,Con((“A2000M06.NDVI.tif”&gt;=0.7),1,(“A2000M06.NDVI.tif”-0.05)/0.65))二、脚本代码#!/usr/bin/python# -*- coding: UTF-8 -*-三、工具参数四、工具界面交流群...
				我们这里展示一个简单的NDVI和EVI区域统计情况，主要是展示长时序的变化情况，但是代码会出现一些错误，我们这里看看代码的错误提示：
Line 51: No features contain non-null values of "system:time_start".
Line 51: No images in the collection intersect the specified regions.
				该植被指数数据来源于美国NASA定期发布的MODIS数据系列，原始数据集为MOD13Q1。MOD13Q1的空间分辨率为250m，时间分辨率为16d。经过提取子数据集、拼接、投影栅格、换算单位、裁剪得到16d 250m的NDVI数据后，再利用最大合成法（MVC，Maximum Value Compositing）生成对应年份的NDVI数据。
地区：中国
时间分辨率：年
空间分辨率：250m
投影坐标系：Albers conical equal area
椭球：WGS84
变形比例：1.0
中央经线：105
标准纬线：25
标准纬线：47
引用：Didan, K. (2015). MOD13Q1 MODIS/Terra Vegetation Indices 16-Day L3 Global 250m SIN Grid V006 [Data set]. NASA EOSDIS Land Processes DAAC. Accessed 2022-04-13 from https://doi.org/10.5067/MODIS/MOD13Q1.006
（二）MRT利用python批量拼接、投影、转换、提取NDVI
（三）MRT利用python批量拼接、投影、转换、提取NDVI质量控制波段
（四）利用ArcGIS、安徽省的shp文件进行批量裁剪
（五）Matlab：批量进行SG滤波、NDVI质量控制
（六）排除异常值，使NDVI值在0-1之间
（七）MVC最大化合成月、季、年数据
后续教程看浏览量，有需求再更新
.........
趋势分析法是一种通过对随时间变化的变量进行线性回归分析，从而来预测其变化趋势的方法。其计算方法如下：
式中，Slope为像元回归方程的斜率，NDVI i为第i年的NDVI的平均值，n为研究的时间长度,视自身情况而定。当Slope>0时，表示该像元NDVI为增加趋势；当Slope=0，表示该像元NDVI基本不变；当Slope<0时，表示该像元NDVI为减少趋势。
二、基于ArcGIS进行趋势分析
1、趋势分析前提
趋势分析的前提是要有NDVI年均值图像，可以直接通过arcma
				该植被指数数据来源于美国NASA定期发布的MODIS数据系列，原始数据集为MOD13A3。MOD13A3的空间分辨率为1km，时间分辨率为monthly。经过提取子数据集、拼接、投影栅格、换算单位、裁剪得到逐月1km的NDVI数据后，再利用最大合成法得到对应年份的NDVI数据。
地区：中国
时间分辨率：年
空间分辨率：1km
投影坐标系：Albers conical equal area
椭球：WGS84
变形比例：1.0
中央经线：105
标准纬线：25
标准纬线：47
引用：Didan, K. (2015). MODIS/Terra Vegetation Indices Monthly L3 Global 1km SIN Grid V006, NASA EOSDIS Land Processes DAAC.https://doi.org/10.5067/MODIS/MOD13A3.006该植被指数数据来源于美国NASA定期发布的MODIS数据系列，原始数据集为MOD13A3。MOD13A3的空间分辨率为1km，时间分辨率为monthly。经过提取子数据集、拼接、投影栅格、换算单位、
				该植被指数数据来源于美国NASA定期发布的MODIS数据系列，原始数据集为MOD13A3。MOD13A3的空间分辨率为1km，时间分辨率为monthly。经过提取子数据集、拼接、投影栅格、换算单位、裁剪得到逐月1km的NDVI数据。
地区：中国
时间分辨率：逐月
空间分辨率：1km
投影坐标系：Albers conical equal area
椭球：WGS84
变形比例：1.0
中央经线：105
标准纬线：25
标准纬线：47
引用：Didan, K. (2015). MODIS/Terra Vegetation Indices Monthly L3 Global 1km SIN Grid V006, NASA EOSDIS Land Processes DAAC.https://doi.org/10.5067/MODIS/MOD13A3.006
				该植被指数数据来源于美国NASA定期发布的MODIS数据系列，原始数据集为MOD13A3。MOD13A3的空间分辨率为1km，时间分辨率为monthly。经过提取子数据集、拼接、投影栅格、换算单位、裁剪得到逐月1km的NDVI数据。
地区：中国
时间分辨率：逐月
空间分辨率：1km
投影坐标系：Albers conical equal area
椭球：WGS84
变形比例：1.0
中央经线：105
标准纬线：25
标准纬线：47
引用：Didan, K. (2015). MODIS/Terra Vegetation Indices Monthly L3 Global 1km SIN Grid V006, NASA EOSDIS Land Processes DAAC.https://doi.org/10.5067/MODIS/MOD13A3.006
定义：NDVI（Normalized Difference Vegetation Index,归一化差分植被指数,标准差异植被指数)，也称为生物量指标变化，可使植被从水和土中分离出来。
NDVI基本原理：
植物叶片组织对蓝光(470nm)和红光(650nm)有强烈吸收，而对绿光和红外光强烈反射。叶片中心的海绵组织和叶片背面组织对近红外辐射(NIR, 700-1000nm)反射较强。从红光(Red)到红外光
from scipy.stats.stats import pearsonr
RasterFile1 = 'D:/2ndcmpttions/遥感小分队/相关系数图/NDVI波段/NDVIAPR.tif'
RasterFile2 =