CN1924611A - 土地退化(荒漠化)评价参数遥感反演与监测技术方法 - Google Patents

Abstract

本发明公开一种土地退化(荒漠化)遥感评价参数指标反演与监测的技术方法。所述方法包括：对所获取遥感数据进行土地退化(荒漠化)遥感参数的选择与评价，建立土地退化(荒漠化)监测评价指标体系；全国荒漠化遥感参数的反演；荒漠化分类算法的评价与实现；利用植被萌芽事件分析植被与非植被区界线变化；利用NDVI变异系数(CoV)分析荒漠化动态变化。本发明可应用于各种土地退化类型以及生态环境的监测评价与变化监测。

Description

土地退化(荒漠化)评价参数遥感反演与监测技术方法

技术领域

本发明涉及信息领域的应用技术，尤其涉及利用遥感技术所取得的信息为土地退化(荒漠化)评价和监测及资源环境可持续发展规划提供决策信息的应用技术。

背景技术

在土地退化研究中需要获取多种地表自然条件的参数包括：土地覆盖、地表温度、土壤湿度和植被参数等。这些参数可由常规气象、水文和农田等观测站观测获取和实地调查取得。但通过实际观测所获取的只是点上的数据，虽然利用一些聚合方法可以将点的观测结果扩展到区域但其精度无法保证，不能代表这些参数大区域的空间分布特征。

遥感技术的发展为解决这一问题提供了重要手段，使点上观测的结果推广到面上成为可能。根据不同地表参数对遥感波段的反应，可以用不同的数学模型和算法获取特定的地表参数。

尽管国际上利用遥感技术做了大量土地退化评价参数的研究，但是采用遥感反演的参数只有1-2个。随着多光谱遥感技术的发展，尽可能利用遥感技术反演土地退化评价参数，并针对性中国自然特征，建立适用的模型，实现全国大范围的土地退化(荒漠化)评价参数反演，对快速进行全国土地退化的监测具有重要的科学意义和实用价值。

发明内容

因此，本发明的目的在于提供一种遥感技术的应用方法，可以通过遥感技术所获取的数据选择合适的光谱参数，建立数据模型，对遥感数据进行处理，以提取土地退化和荒漠化的信息。

为实现上述目的，本发明提出了一种土地退化和荒漠化遥感监测方法，包括：荒漠化遥感监测综合指标体系的建立；荒漠化遥感监测参数的反演；荒漠化程度遥感评价指标的确定；荒漠化遥感分类；荒漠化动态变化监测。

附图说明

图1是不同荒漠化监测指标组合的分类图；

图2是中国2001年生长季MSAVI累计值分布图；

图3是中国2001年生长季平均陆面温度分布图；

图4是NDVI-Ts特征空间简化图；

图5是中国2001年生长季平均TVDI分布图；

图6是中国2001年生长季平均反照率分布图；

图7是中国2001年最大植被覆盖度分布图；

图8是1995年中国土地荒漠化现状图；

图9是2001年中国土地荒漠化现状图；

图10是中国荒漠化土地动态变化图；

图11是中亚及中国干旱半干旱区不同年份的沙漠边界图；

图12是中亚植被非植被分布图；

图13是1982-2000年中亚及中国干旱半干旱区植被萌芽事件分布频率图；

图14是中亚及中国干旱半干旱区不同年份的cov分布图。

具体实施方式

在本发明的优选实施案例中，首先建立了荒漠化遥感监测综合指标体系，包括土地退化(荒漠化)遥感参数的选择与评价以及荒漠化程度指标体系的确定。

1)土地退化(荒漠化)遥感参数的选择与评价

由于土地退化基本上发生在干旱、半干旱和亚湿润干旱区，这些地区面积广大，自然条件恶劣、交通不便，因此，在国际上早在80年代就利用遥感技术开展了对土地退化的研究工作。但自从80年代J.Tucker建立了NOAA AVHRR植被指数(NDVI)对非洲的沙漠化进行研究以来，国际上利用遥感技术对荒漠化的研究基本上还还停留在采用NDVI等一、两个指标进行评价。

随着多光谱遥感技术的发展，有些荒漠特征参数可以利用遥感数据反演出来，本发明选取的指标有：反映荒漠化土地自然属性和生态状况的植被覆盖度(FVC)和改进型土壤调整植被指数(MSAVI)以及反映荒漠化土地物理属性的反照率(Albedo)、陆地表面温度(LST)和土壤湿度(TVDI)共5个指标。

为确保所用遥感监测指标的监测准确性和不可替代性，对荒漠化监测指标的监测精度进行了评价。以科尔沁沙地为试验对象，在选择的5个指标中，使用相同的训练样本，分别用1995年的分辨率为1公里的NOAA-AVHRR数据反演得到的单个指标和多个指标的组合，来进行科尔沁沙地荒漠化监测精度的检验，选取最合适的指标或指标组合。

荒漠化监测各个指标所用的具体数据分别是科尔沁地区的年最大植被覆盖度、MSAVI的生长季(4-10月)累积值、生长季的反照率平均值、生长季的陆面温度平均值、生长季的土壤湿度(TVDI)平均值。

分类中使用相同的训练样本数据，得到各荒漠化遥感监测指标在不同荒漠化程度等级的数值范围，使用决策树分类器分类。但分别使用不同数目的荒漠化监测指标和指标组合进行分类。分类结果的评价以总体分类精度为指标进行。总体分类精度是指分类正确的像元数据除以总像元数所得的百分比，总像元是指训练样本中的所有像元，包括分类正确的像元。所用荒漠化遥感监测指标的具体组合和分类精度见表1。

表1不同监测指标组合的荒漠化程度监测的总体分类精度表

一个指标	MSAVI	Albedo	LST		FVC	TVDI
							分类精度	68.91％	56.24％	43.56％		66.84％	30.57％
两个指标	MSAVI+Albedo	MSAVI+LST	MSAVI+TVDI		FVC+Albedo	FVC+LST
							分类精度	86.59％	84.63％	77.65％		85.31％	82.49％
三个指标	MSAVI+Albedo+LST	MSAVI+TVDI+Albedo	MSAVI+LST+TVDI		Albedo+TVDI+FVC	Albedo+LST+FVC
							分类精度	90.23％	89.18％	86.84％		88.52％	87.67％
四个指标	MSAVI+Albedo+LST+TVDI			FVC+Albedo+LST+TVDI
							分类精度	93.83％			92.05％
五个指标	MSAVI+Albedo+LST+TVDI+FVC
							分类精度	95.21％

从表1中可以看出，单独使用一个荒漠化监测指标进行分类时，MSAVI的分类精度最高，达到68.91％(图1)，其次是植被覆盖度和反照率，分类精度分别为66.84％和56.24％，分类精度最低的是TVDI，精度仅为30.57％。当使用两个监测指标分类时，MSAVI+Albedo的组合分类精度最高，达到86.59％(图1)，比仅使用一个指标时的精度有明显的提高。当使用三个指标分类时，MSAVI+Albedo+LST的组合分类精度最高，为90.23％(图1)。而使用四个分类指标分类时，虽然MSAVI+Albedo+LST+TVDI组合比FVC+Albedo+LST+TVDI组合的分类精度高(图1)，但是两者相差不大。当五个监测指标全部参与分类时，分类精度达到95.21％，在所用组合中精度最高(图1)。可见，在荒漠化监测中，每增加一个指标，反映的荒漠化信息随之增加，使荒漠化监测的精度也大大提高。从使用一个指标分类的结果也可以看出，MSAVI的分类精度最高，说明MSAVI所反映的荒漠化信息最多，对荒漠化遥感监测的贡献最大。

2)荒漠化程度指标体系的确定

本发明实施案例中，利用遥感技术确定了荒漠化程度的评价指标，其步骤如下：

所用资料是1995年中国防治荒漠化协调小组组织林业部调查规划设计院、中国林科院、林业部西北调查规划设计院、中科院综考委、北京林业大学等单位以全国沙漠、戈壁和沙化土地普查为主要依据，编制的1∶100万中国荒漠化土地分布现状图，以及由1995年NOAA-AVHRR数据反演出的5个荒漠化遥感监测指标，即1995年生长季MSAVI累积值、生长季平均反照率、生长季平均陆面温度、年最大植被覆盖度和生长季平均TVDI值。

指标体系的确定方法是，首先将选取的荒漠化程度样本图，即中国林科院荒漠化监测中心绘制的1995年的中国土地荒漠化图数字化，然后将数字化的样本图与NOAA-AVHRR反演出的各个荒漠化遥感监测指标的栅格图进行GIS叠加分析，统计不同荒漠化程度类型的各个指标值，根据统计的不同荒漠化类型大部分像元的数值，确定NOAA-AVHRR荒漠化监测指标体系。根据确定后的指标体系，使用NOAA-AVHRR数据，就可以对1995年的中国及中亚地区的荒漠化状况进行监测。

NOAA-AVHRR荒漠化监测指标体系见表2、3、4和5。

                     表2亚湿润干旱区的NOAA荒漠化监测指标体系

荒漠化程度	MSAVI	植被覆盖度	反照率	陆面温度	TVDI
						非荒漠化	＞3.0	＞0.8	＜180	＜29	＜0.39
轻度	2.7-3.0	0.7-08	180-220	29-32	0.39-0.51
						中度	2.2-2.7	0.53-0.7	220-250	32-35	0.51-0.60
重度	1.9-2.2	0.4-0.53	250-280	35-38	0.60-0.70
						极重度	＜1.9	＜0.40	＞280	＞38	＞0.70

                     表3半干旱区的NOAA荒漠化监测指标体系

荒漠化程度	MSAVI	植被覆盖度	反照率	陆面温度	TVDI
						非荒漠化	＞2.2	＞0.6	＜200	＜30	＜0.35
轻度	1.9-2.2	0.3-0.6	200-220	30-36	0.35-0.45
						中度	1.2-1.9	0.23-0.3	220-260	36-40	0.45-0.56
重度	0.8-1.2	0.15-0.23	260-300	40-46	0.56-0.65
						极重度	＜0.8	＜0.15	＞300	＞46	＞0.65

                     表4干旱区的NOAA荒漠化监测指标体系

荒漠化程度	MSAVI	植被覆盖度	反照率	陆面温度	TVDI
						非荒漠化	＞1.2	＞0.54	＜210	＜37	＜0.5
轻度	0.8-1.2	0.15-0.54	210-230	37-41	05-0.6
						中度	0.4-0.8	0.08-0.15	230-250	41-49	0.6-0.68
重度	0.1-0.4	0.03-0.08	250-300	49-53	0.68-0.78
						极重度	＜0.1	＜0.03	＞300	＞53	＞0.78

                    表5高寒区的NOAA荒漠化监测指标体系

荒漠化程度	MSAVI	植被覆盖度	反照率	陆面温度	TVDI
						非荒漠化	＞1.3	＞0.24	＜230	＜30	＜0.49
轻度	1.1-1.3	0.12-0.24	230-255	30-33	0.49-0.56
						中度	0.7-1.1	0.09-0.12	255-290	33-37	0.56-0.65
重度	0.4-0.7	0.04-0.09	290-320	37-40	0.65-0.71
						极重度	＜0.4	＜0.04	＞320	＞45	＞0.71

2001年本发明使用了MODIS数据进行荒漠化遥感监测。但2001年没有适合的样本数据用以MODIS监测指标体系的确定，并且由于MODIS和NOAA-AVHRR数据的影像特征存在较大差异，而且各个荒漠化监测指标在计算过程中，所用方法不同，因此使用与NOAA-AVHRR相同的荒漠化监测指标体系，也是不合理的。在发明中，对MODIS和NOAA-AVHRR在同一空间位置的同一指标进行回归分析，可以看出，1995年的NOAA-AVHRR和2001年的MODIS数据具有较强的相关关系(表6)。因此，利用这种相关关系可以依据NOAA-AVHRR建立起的荒漠化监测指标体系计算出MODIS所对应的荒漠化指标体系(表7、8、9和10)。根据所建立的荒漠化监测指标体系，利用2001年的MODIS数据，即可以进行中国和中亚地区的荒漠化土地监测。

                        表6NOAA与MODIS荒漠化监测指标的相关关系

监测指标	亚湿润干旱区		半干旱区		干旱区		高寒区
									关系式	相关系数	关系式	相关系数	关系式	相关系数	关系式	相关系数
	MSAVI	Y＝0.46*X+1.06	0.94	Y＝0.66*X-0.12	0.943	Y＝0.83*X-0.003	0.96	Y＝0.59*X+0.02									0.93
FVC									Y＝0.77*X+0.06	0.927	Y＝0.84*X-0.004	0.92	Y＝0.24*X+1.02	0.96	Y＝0.92*X+0.05	0.85
	Albedo	Y＝1.03*X-19.01	0.93	Y＝1.03*X-23.89	0.91	Y＝0.51*X+83.78	0.87	Y＝0.81*X+2.63									0.87
LST									Y＝0.35*X+16.19	0.77	Y＝0.54*X+9.52	0.83	Y＝0.32*X+17.87	0.78	Y＝0.91*X-12.91	0.85
	TVDI	Y＝0.65*X+0.48	0.78	Y＝0.77*X+0.25	0.82	Y＝0.56*X+0.23	0.8	Y＝1.19*X-0.15									0.93

                 表7亚湿润干旱区的MODIS荒漠化监测指标体系

荒漠化程度	MSAVI	植被覆盖度	反照率	陆面温度	TVDI
						非荒漠化	＞1.9	＞0.70	＜170	＜29	＜0.43
轻度	1.65-1.9	0.55-0.70	170-200	29-32	0.43-0.56
						中度	1.45-1.65	0.45-0.55	200-220	32-34	0.56-0.67
重度	1.25-1.45	0.35-0.45	220-240	34-36	0.67-0.76
						极重度	＜1.25	＜0.35	＞240	＞36	＞0.76

                   表8半干旱区的MODIS荒漠化监测指标体系

荒漠化程度	MSAVI	植被覆盖度	反照率	陆面温度	TVDI
						非荒漠化	＞1.50	＞0.6	＜180	＜29	＜0.54
轻度	1.20-1.50	0.43-0.6	180-205	29-32	0.54-0.60
						中度	0.90-1.20	0.31-0.43	205-230	32-34	0.60-0.69
重度	0.60-0.90	0.1-0.31	230-250	34-37	0.69-0.72
						极重度	＜0.60	＜0.1	＞250	＞37	＞0.72

                   表3-9干旱区的MODIS荒漠化监测指标体系

荒漠化程度	MSAVI	植被覆盖度	反照率	陆面温度	TVDI
						非荒漠化	＞1.10	＞0.40	＜180	＜30	＜0.53
轻度	0.80-1.10	0.20-0.40	180-200	30-35	0.53-0.58
						中度	0.50-0.80	0.20-0.32	200-250	33-35	0.58-0.61
重度	0.40-0.50	0.08-0.20	250-265	35-37	0.61-0.68
						极重度	＜0.40	＜0.08	＞265	＞37	＞0.68

                 表3-10高寒区的MODIS荒漠化监测指标体系

荒漠化程度	MSAVI	植被覆盖度	反照率	陆面温度	TVDI
						非荒漠化	＞0.80	＞0.45	＜180	＜21	＜0.4
轻度	0.70-0.80	0.34-0.45	180-200	21-25	0.4-0.55
						中度	0.50-0.70	0.21-0.34	200-220	25-29	0.55-0.67
重度	0.30-0.50	0.07-0.21	220-260	29-33	0.67
						极重度	＜0.30	＜0.07	＞260	＞33	＞0.71

优选地，本发明建立了五种土地退化(荒漠化)评价参数遥感反演算法模型。

1)改进型土壤调整植被指数(MSAVI)

本发明选用MSAVI作为荒漠化监测指标，是因为MSAVI既能尽量消除土壤背景的影响，也不需要人为确定参数，比较方便客观。

MSAVI(Modified Soil-Adjusted Vegetation Index)为改进型土壤调节植被指数，该指数的计算公式为：

$\mathrm{MSAVI}=(2\mathrm{NIR}+1-\sqrt{{(2\mathrm{NIR}+1)}^2-8(\mathrm{NIR}-R)})/2---\left(1\right)$

        SAVI＝[(NIR-R)/(NIR+R+L)]×(L+1)            (2)

R是指可见光波段，NIR是指近红外波段。

经计算取得了中国2001年生长季MSAVI累计值分布图(图2)。

2)陆地表面温度

应用Becker and Li(1990)的算法计算中国每旬的陆面温度。具体步骤为：

●计算归一化植被指数NDVI

        NDVI＝(Ch2-Ch1)/(Ch2+Ch1)                   (3)

其中，Ch1为可见光通道，Ch2为近红外通道。NDVI的值在-1到+1范围之间。

●计算地表比辐射率ε

Josef et al.(1997)在已有工作(Griend and Owe，1993；Salisbury，1994)的基础上，计算出AVHRR的ε₄和Δε，方程为：

            ε₄＝0.9897+0.029ln(NDVI)                  (4)

            Δε＝ε₄-ε₅＝0.01019+0.01344ln(NDVI)    (5)

本发明使用公式(4)和(5)计算分波段的比辐射率ε₄和ε₅。

●计算陆面温度Ts

本发表采用Becker and Li(1990)的***窗算法计算中国及中亚地区的每旬陆面温度：

T_s＝1.274+(T₄+T₅)/2{1+[0.15616(1-ε)/ε]-0.482(Δε/ε²)}    (6)

    +(T₄-T₅)/2{6.26+[3.98(1-ε)/ε]+38.33(Δε/ε²)}

其中，ε＝(ε₄+ε₅)/2，Δε＝(ε₄-ε₅)。

计算出每旬的陆面温度后，每三旬取平均求得生长季(4-10月)每月的温度，然后对生长季共7个月的陆面温度求平均，得到中国生长季的平均陆面温度值，见图3。

3)土壤湿度与旱情

根据Sandholt等(2002)利用简化的NDVI-TS特征空间提出水分胁迫指标，即温度植被旱情指数，在简化的特征空间，将湿边(Ts-min)处理为与NDVI轴平行的直线，旱边(Ts-max)与NDVI呈线性关系。该简化的NDVI-Ts特征空间见(图4)。

将NDVI-Ts特征空间简化处理为三角形，把中国分为3个农业气候区分别进行了旱情监测。根据计算出的每8天的TVDI值，求生长季(4-10月)的平均值，分别得到中国的生长季的平均TVDI值作为荒漠化监测的土壤水分指标。

$\mathrm{TVDI}=\frac{[T_s-(a_1+b_1*\mathrm{NDVI})]}{[(a_2+b_2*\mathrm{NDVI})-(a_1+b_1*\mathrm{NDVI})]}---\left(7\right)$

式中，a₁，b₁，a₂，和b₂分别是旱边和湿边拟合方程的系数。

从中国的TVDI分布图(5)可以看出，荒漠地区的TVDI值较大，青藏高原的TVDI值相对较低。TVDI值越大，说明该区的土壤水分含量越少，反之，说明土壤水分含量较大。

4)反照率的计算

采用Saunders(1990)方法计算的中国2001年的反照率。首先计算出每10天合成的反照率后，取每三旬的反照率平均值得到月反照率，然后对生长季(4-10月)共7个月的的反照率求平均，分别得到中国生长季平均反照率分布图(图6)。

5)植被覆盖度

MODIS(NOAA-AVHRR)归一化植被指数NDVI模型运行中，一个栅格的信息是由裸土与植被按面积的加权平均所组成的。同样MODIS(NOAA-AVHRR)归一化植被指数NDVI通过卫星传感器所测到的每个像元的信息φ，就可以表达为由植被区所贡献的信息和由裸土区域所贡献的信息的加权和。即图像中每个像元的NDVI值可以看成是有植被覆盖部分的NDVI与无植被覆盖部分的NDVI的加权平均，其中有植被覆盖部分的NDVI的权重即为此像元的植被覆盖度σ_v，而无植被覆盖部分的NDVI的权重即为(1-σ_v)。

             Φ＝Φ_vσ_v+(1-σ_v)Φ_s                (8)

式(8)中，下标v和s分别表示完全有植被覆盖区和裸土区的值。

将式(8)直接应用到NDVI中，即可得到植被覆盖度的最简单的表达式(9)为：

${\mathrm{\σ}}_v=\frac{\mathrm{NDVI}-{\mathrm{NDVI}}_s}{{\mathrm{NDVI}}_v-{\mathrm{NDVI}}_s}---\left(9\right)$

由于可见光R和近红外波段NIR对于能量流的贡献是呈比例的，因此将式(9)和NDVI算法相结合，得到式(10)：

${\mathrm{\σ}}_v=\frac{\mathrm{NDVI}-{\mathrm{NDVI}}_s}{{\mathrm{NDVI}}_v-{\mathrm{NDVI}}_s+(1-a)(\mathrm{NDVI}-{\mathrm{NDVI}}_v)}---\left(10\right)$

式(10)中，a＝(R+NIR)v/(R+NIR)s。当a＝1时，则式(10)变为线性的式(9)。

另一个非线性的公式(Choudhury et al.1994；Baret et al.1995；Wittich 1997)是

${\mathrm{\σ}}_v=1-{\left(\frac{{\mathrm{NDVI}}_v-\mathrm{NDVI}}{{\mathrm{NDVI}}_v-{\mathrm{NDVI}}_s}\right)}^b---\left(11\right)$

在估算大尺度的植被覆盖度时，由于地面数据的缺乏，因此用基于线性关系的方法比其它较为复杂的方法更适合。本发明估算中国的植被覆盖度的公式如下：

${\mathrm{\σ}}_v=\frac{N_{p,\max }-N_s}{N_{c,v}-N_s}---\left(12\right)$

式中，N_p，max是指该像元的年最大NDVI值；N_c，v是指每类土地覆盖类型的植被覆盖度为100％时相对应的像元NDVI值；N_s是指每类土地覆盖类型的NDVI最小值。

由于图像中不可避免的存在着噪声的影响，它可能产生过低或过高的NDVI值，如果用这样的值来计算植被覆盖度当然只会得到错误的结果。为了避免这样的错误发生，在确定N_c，v和N_s时，原则上应确定一个置信度，做出图像中每种土地覆盖类型的每个像元集合的NDVI概率分布，计算置信区间内的最大值N_c，v与最小值N_s。N_s应该是不随时间改变的，对于大多数类型的裸地表面，理论上应该接近零，一般取值0.05。对于IGBP其它土地覆盖类型的N_s，在冬天可能有更小的NDVI值，而在夏天，由于云污染和大气的影响，N_s的取值具有更大的不确定性。因此，对每一种IGBP土地覆盖类型分别确定一个N_s值是不实际的，本发明在估算植被覆盖度时，把所有土地覆盖类型的N_s都取值为0.05。

Zeng等人估算过1992年的全球植被覆盖度。

本发明根据Zeng确定的置信度，利用MODIS数据分别对不同土地覆盖类型做NDVI的累积频率图，然后确定不同土地覆盖类型的N_c，v(表11)。将2001年的NDVI最大值图像和2001年的土地覆盖分类图合成到一个文件，然后对文件的每个像元按照其所属土地覆盖类别，找到与表(11)中相对应的N_c，v值，把像元的NDVI值，N_c，v和N_s代入式(12)计算中亚及中国2001年的植被覆盖度。图7是中国2001年的最大植被覆盖度分布图。从图7可以看出，由西北向东南，植被覆盖度呈现逐渐增大的趋势，以农牧交错带为界线，以西的植被覆盖度明显小于以东区域。

本发明主要利用决策树分类器对中国的荒漠化和动态变化进行了研究，详细步骤论述如下：

为了比较，本发明分别采用了非监督分类，监督分类中的最大似然法(MLC)和决策树三种分类器，以科尔沁沙地为例，使用2001年的MODIS数据反演的五个荒漠化遥感监测指标，在选取相同的荒漠化程度训练样本情况下，进行了荒漠化程度分类，以便选取分类精度最好的分类器实现中国及中亚地区的荒漠化遥感监测。

分类完成后，随机选取500个样点，结合基于TM图像的科尔沁沙地90年代末的荒漠化监测图，对三种分类结果进行精度评价，分类精度结果见表(12)。

从表12中可以看出，在三种分类器的总体分类精度中，决策树分类结果最高，分类精度为94.63％；其次是最大似然分类法，分类精度是88.59％；分类精度最小的是非监督法，精度为80.88％。决策树分类法比最大似然和非监督分类的分类精度分别高出6.04％和13.75％。

为更科学地评价三种分类器的性能，本发明同时采用Kappa系数进行分类精度的比较。Kappa统计在评价不同分类器的分类精度上更具有统计上的辨析力。从表12的最后一行可以看出，三种分类器的Kappa系数比较结果与分类精度比较结果是相同的，即分类精度最高的是决策树分类器。因此，根据上述结果，本发明将采用决策树分类器进行荒漠化程度的分类。

表11用于MODIS植被覆盖度监测的各个土地覆盖类型的Nc，v

IGBP land cover		％	Nc，v
				1	常绿针叶林	75％	0.89
2	常绿阔叶林	75％	0.90
				3	落叶针叶林	75％	0.86
4	落叶阔叶林	75％	0.88
				5	混交林	75％	0.89
6	密灌	90％	0.88
				7	疏灌	90％	0.62
8	木质稀树大草原	75％	0.85
				9	稀树大草原	75％	0.56
10	草地	75％	0.48
				11	永久湿地	75％	0.78
12	农田	75％	0.84
				13	城市用地	90％	0.88
14	农田自然植被混合	75％	0.81
				15	冰雪
16	裸地	90％	0.60
				17	水体

             表12不同分类器的荒漠化程度分类精度表

	非监督分类	最大似然分类	决策树分类
				非荒漠化	81.03％	87.64％	96.15％
轻度	80.67％	89.75％	91.47％
				中度	75.01％	71.63％	89.30％
重度	69.41％	85.72％	93.69％
				极重度	85.56％	90.63％	98.03％
总体分类精度	80.88％	88.59％	94.63％
				Kappa系数	0.735	0.802	0.885

              表13中国荒漠化土地监测结果比较

	国家***	本发明
			研究时段	1994-1999(5年)	1995-2001(6年)
荒漠化土地增长面积(万km2)	5.2	5.68
			荒漠化土地年均增长(万km2)	1.04	0.95
轻度荒漠化土地所占比例	20％	18％*
			中度荒漠化土地所占比例	33％	31％*
重度荒漠化土地所占比例	21％	21％*
			极重度荒漠化土地所占比例	26％	30％*

^*数据以2001年为准。

在荒漠化遥感监测指标体系的基础上，对中国的1995年和2001年的荒漠化分布范围和荒漠化程度分别进行了评价。

根据本发明中国荒漠化土地遥感监测的结果，1995年，我国的荒漠化土地面积为250.87×10⁴km²，占国土面积的25.86％。具体分布见图(8)。

根据本发明荒漠化土地遥感监测的结果，2001年，本发明取得的中国荒漠化土地面积为256.55×10⁴km²，占国土面积的26.45％。具体分布见图(9)。

在中国1995年和2001年荒漠化程度分布图的基础上，得到中国的荒漠化的动态变化图。通过对荒漠化动态变化图的统计分析，可以得出不同荒漠化类型土地面积的变化以及荒漠化程度的演变情况(图10)。

国家林业部曾于1994年和1999年先后进行了两次全国性荒漠化土地普查和监测工作，按可比口径的荒漠化土地变化分析，我们将本发明的研究成果与其进行了比较(表13)。

从表13中，本发明荒漠化土地监测结果中荒漠化土地的年均增长速度小于国家林业部的监测结果，其中主要原因是两者研究范围的差异，本发明的研究范围没有包括亚湿润干旱区的岛状区域，面积相对减小。除此之外，进行荒漠化监测时所用遥感数据源的不同，遥感数据获取时间的不同等也是产生差异的原因。

在本发明的优选实施例中，利用植被萌芽事件分析植被与非植被区界线变化；以及利用NDVI变异系数(CoV)分析荒漠化动态变化。

本发明所用的方法包括两个方面，一方面是基于年际NDVI时间序列数据的沙漠界线(或植被区与非植被区的界线)提取，另一方面是基于CoV指数的多年荒漠化遥感监测。前者用于提取沙漠边际的多年变化趋势；后者在前者工作的基础上，用于监测沙漠界线外的荒漠化动态变化情况。

1)沙漠界线提取

前人的许多研究表明，NDVI的季节变化与植被物候期密切相关。植被的萌芽期是指植被受水分或温度的驱使，从静止状态开始进行活跃的光合作用的时期，在这一时期，NDVI值骤然增大。在干旱半干旱地区，大部分植被都是季节性的，从植被开始生长，NDVI值逐渐增加，直到植被生长到最茂盛时NDVI值达到最大，然后NDVI逐渐减小。本发明根据区域内是否能监测到植被的萌芽期来提取沙漠界线。

Lloyd和Reed等认为NDVI的突然增长是植被的萌芽时期。在沙漠中植被分布稀少，因此用AVHRR的NDVI数据难以测到植被的萌芽期(Loiyd，1990；Reed et al.，1994)。F.Yu等(2004)认为可以根据一个地区是否有植被萌芽事件的发生来监测沙漠的界线。

本发明使用是否有植被萌芽事件发生这一指标来判定沙漠界线。为确定NDVI的突然增长，即确定植被的萌芽期，对每年每个像元的NDVI季节变化进行了监测。一个地区若有植被的萌芽发生，则应该符合下面几个条件：①NDVI必须在一段时期内连续增长，这段时期至少是一个半月的时间，因此在本发明中NDVI至少是在5旬是连续增长的；②在这段时间内，NDVI能达到最大值；③NDVI值大于0.05；④以上几个条件必须发生在一个合理的时间段内，即发生在植被的正常生长期：4月到8月之间。在PCI图像处理软件的EASI编程环境下，实现沙漠界线的提取。本发明所指的沙漠界线即是植被区与无植被区的界线。

2)基于CoV方法的多年荒漠化遥感监测

地表植被是荒漠化评价的一个重要指标，而NDVI是用于监测植被变化的最经典植被指数(QiJ，2000)。NDVI在使用遥感图像进行植被研究以及植物物候研究中得到广泛应用，它是植物生长状态以及植被空间分布密度的最佳指示因子，与植被分布密度呈线性相关。NDVI的变异系数(CoV)是监测植被生长周期的一个比较稳定方法。CoV一般用来比较不同样本之间的变化，这些样本是时间序列数据，因此CoV表示的是一段时间内每个像元的NDVI变化情况。CoV可以用来监测一段时间内植被的动态变化，并已经被用来评价植被的变化和生态***界线的提取(Tucker，1991；L.Milich，and E.Weiss，1997)。CoV是多时相NDVI数据的综合表现，因此比单时相的NDVI变化更稳定。

在不同的时间段(如周，月，年)，每个像元的NDVI序列值都可以生成一幅CoV图像。在发明中，每年的CoV图像是根据每个像元的月NDVI最大值算出的。首先计算图像中每个像元的12个月的NDVI的标准差(σ)，然后用标准差除以该像元12个月的NDVI平均值(μ)，即得到该像元的年CoV值。具体方法如式(13)所示：

           cov_ij＝σ_ij/μ_ij                      (13)

$\mathrm{\μ}=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{p}_i/n---\left(14\right)$

$\mathrm{\σ}=\sqrt{\frac{1}{n-1}[\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{p_i}^2-\frac{{\left(\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{p}_i\right)}^2}{n}}---\left(15\right)$

式(13)中，i和j分别是每个像元的行列数，式(14)和式(15)中，pi代表像元值，n为时间序列个数。CoV值处于同一尺度内，因此更有利于不同像元间的比较。CoV可以用来量测植被的生长周期，因此CoV的变化能够监测一个区域内植被生长周期的变化。在干旱半干旱区有植被覆盖的区域，CoV值的下降与降雨的减少或生物量的降低有关。

E.Weiss等人曾用NDVI的CoV以及CoV坡度(即CoV多年变化趋势拟合曲线的斜率)来评价沙特***的草原变化情况，结果证明此方法能够成功用于干旱半干旱区(E.Weiss，S.E.Marsh，and E.S.Pfirman，2001)。

本发明在提取沙漠界线的基础上，使用CoV方法来监测沙漠界线外干旱半干旱区的荒漠化情况。首先，对每个像元应该确定多年CoV值线性拟合结果最好的方程式，然后得到CoV的回归线斜率即CoV坡度；CoV坡度反映了每个像元多年的整体变化趋势。CoV坡度的具体算法采用最小乘方OLS估计法，如式(16)所示：

$b=\frac{\mathrm{n\Σ}{x}_i{y}_i-\mathrm{\Σ}{x}_i\mathrm{\Σ}{y}_i}{\mathrm{n\Σ}{x}_i^2-{\left(\mathrm{\Σ}{x}_i\right)}^2}---\left(16\right)$

式(16)中，x为年份，本发明中x是1982，1983，…，2000；y为每年的每个像元的CoV值，i为1，2，…n，n为年总数，本发明中n为18。

如果在研究时段内CoV的变化趋势是下降的，则可推断出该区的NDVI的年内变化幅度减小，即该区植被状况变差，正处于荒漠化过程中。反之，CoV的变化趋势是上升的，则可认为该区的植被增加或生长状态变好。

利用每年36旬的NDVI数据，在PCI的EASI编程环境下，实现每年的植被区与非植被区界线的提取。图11是列举的中亚及中国干旱半干旱区六年的沙漠边界分布图，它们分别是1982，1986，1990，1993，1997和2000年。

首先提取中亚及中国干旱半干旱区从1982-2000(1994年除外)共18年的植被分布区域，然后通过GIS叠加，获得中亚及中国干旱半干旱区植被分布图(图12)。图中的非植被区域是指在1982年至2000年(1994年除外)18年间从来没有植被分布的区域。

图13是18年的植被分布区域通过GIS叠加的另一个结果，即中亚及中国干旱半干旱区1982年至2000年期间，植被萌芽事件的发生频率分布图。从图中可以看出，干旱半干旱区的沙漠界线具有很强的变化性，并且在沙漠的中心地带从来没有监测到有植被的分布。从沙漠经沙漠草原，然后过渡到典型草原带，在18年中，植被发生过的次数从2次过渡到18次。

对这个地区沙漠界线的多变性的一个可能解释是，干旱半干旱区的年降雨量的变异性比较大的缘故。在年降雨量充沛时，植被生长茂盛，沙漠界线向沙漠中心方向移动；而在年降雨量较低时，沙漠界线向草原方向移动。因此，沙漠草原在湿润的年份表现为典型草原的物候特征，而在干旱的年份则表现为沙漠的特征。

图14即是中亚及中国西部干旱半干旱区的CoV坡度分布图。图中，绿色区域表示CoV的坡度几乎为零，即植被状况稳定或多年变化不明显的区域；黄色至红色代表CoV坡度为正的区域，表示植被生长多年来是趋于变好的，并且，由黄色向红色的变化是表示CoV坡度的变化幅度是逐渐增大，即植被生长变好的程度更明显；蓝色是CoV坡度为负的区域，表示植被生长状况变坏，该区处于荒漠化中，而由浅蓝色向深蓝色的过渡说明荒漠化程度趋于严重。CoV坡度之所以分为7个等级，是为更清楚地显示荒漠化程度的强弱，没有更具体的物理意义。图14中CoV坡度为正的区域，中国主要分布在新疆天山地区，西藏南部，和内蒙古东部地区，研究结果大体与我国的防护林分布相吻合。而发生荒漠化的地区则主要分布在沙漠和沙地的边缘以及农牧交错带。中亚地区CoV坡度为负，发生荒漠化的区域主要分布在沙漠边缘，湖滨周围，以及农牧交错带等区域。

Claims (6)

1.一种土地退化(荒漠化)评价参数遥感反演与监测技术方法，包括：

获取遥感数据并进行预处理；

对预处理过的遥感数据，通过建立数学模型，反演土地退化(荒漠化)评价参数；

利用遥感反演的参数建立荒漠化监测指标体系；

利用荒漠监测指标体系取得1995和2001年中国荒漠化变化信息；

利用植被萌芽事件分析植被与非植被区界线变化；

利用NDVI变异系数(Cov)分析中亚及中国干旱、半干旱区荒漠化动态变化。

2.根据权利要求1所述的方法，其特征在于，所述建立数学模型反演的土地退化评价参数包括：植被指数(MSAVI)、植被覆盖度(FVC)、反照率(ALBEDO)、陆地温度(Ts)、土壤湿度(TVDI)。

3.根据权利要求1所述的方法，其特征在于，所述荒漠化遥感监测指标体系建立的步骤包括：

荒漠化遥感监测指标选取原则；

荒漠化遥感监测综合指标的确定；

提出适合于应用遥感进行大尺度荒漠化评价与监测的参数指标，并通过监测参数组合进行荒漠化监测精度评价，结果表明采用改进型土壤调节植被指数、植被覆盖度、反照率、陆面温度和土壤温度指数的组合在荒漠化遥感监测中分类精度最高，可以作为监测综合指标；

基于遥感数据的荒漠化程度指标体系的确定；

选取荒漠化程度图与遥感各荒漠化指标叠加分析；并根据荒漠气候类型不同，把研究区分为亚湿润干旱区、半干旱区、干旱区和高寒区四个区域，针对每一个区域分别建立荒漠化程度遥感监测指标体系。

4.根据权利要求1所建立的方法，其特征在于，所述的中国及中亚地区荒漠化变化信息获取，通过利用基于中国土地荒漠化样本建立的荒漠化遥感监测指标体系，对中国相应的四个分区1995年及2001年荒漠化状况进行了监测评价。

5.根据权利要求1所述的方法，其特征在于，所述的植被与非植被区界线变化是通过利用1982-2000年8km分辨率的NOAA AVHRR NDVI 10天合成时间序列数据集，根据区域内是否发生过植被萌芽事件，确定每年中亚及中国干旱、半干旱区植被与非植被界线，并对18年植被界线进行GIS叠加，分析多年的界线变化情况。

6.根据权利要求1所述方法，其特征在于，所述的中亚及中国干旱、半干旱区荒漠化动态变化规律是利用1982-2000年NOAA-AVHRR 10天合成时间序列数据，分别计算每年的NDVI的变异系数，然后用最小二乘法求得18年的NDVI变异系数坡度。通过分析变异系数坡度变化，实现中亚及中国干旱、半干旱区多年荒漠化状况的评价。

Images