基于遥感分析20年来民勤绿洲防护林带植被变化研究

张永，杨自辉，郭树江，王强强，詹科杰，张剑挥，魏怀东

(1.甘肃农业大学林学院，甘肃 兰州730070；2.甘肃民勤荒漠草地生态系统国家野外科学观测研究站，甘肃 民勤733300)

绿洲防护带(农田防护林，防风固沙林)固定沙丘、防止流沙向绿洲蔓延、消减沙尘暴发展，是保护绿洲的生态屏障。因此，绿洲防护带植被的动态监测与评价，是探讨绿洲生态安全的重要手段，传统的监测，主要通过人工实地考察，建立定量观测样地调查等，不适合较大范围定量监测。近几年，随着国内外争相发射资源卫星，如国产高分系列，资源卫星系列，欧洲航空局哨兵卫星系列，利用卫片和航拍进行荒漠化动态监测成为首选的方法之一。

植被覆盖度(fractional vegetation cover，FVC)是指植被冠层的垂直投影面积占区域总面积的百分比,是描述防护带(包括天然林和人工林)长势与密度的重要指标。目前通过卫星数据进行植被覆盖度的研究，主要从数据源的选择、计算指数的选取、变化成因3个方面开展。传统的数据源主要是中低分辨率数据，如Modis、Landsat系列数据[1-2]。近几年，随着高分数据[3]存量的增加，以及小型无人机[4]的快速发展，高分辨率数据作为数据源成为热点。此外，利用ASD光谱仪，野外实测获取光合、非光合植被、裸土、阴影的端元光谱信息，使用线性光谱模型计算覆盖度[5]。计算指数的选择关系到覆盖度估算的精度。因此，基于NDVI[6](normalized difference vegetation index)归一化植被指数、EVI[7](enhanced vegetation index)增强型植被指数、TAVI[8](terrain adjust vegetation index)地形调节植被指数等基于像元二分模型[9-10]计算植被覆盖度被广泛采用。

土地利用变化遥感监测，能够定量反演研究区各土地利用类型在时间和空间上的动态变化。利用卫星数据进行土地类型分类，其精度成为关键。目前，主要从分类数据源、分类方法和分类器选择等方面提高分类质量。对于大尺度区域来说Modis、Landsat[11]成为首选数据。中小尺度区域利用高分数据、无人机航片[12]比较适合。基于CART决策树[13-15]、面向对象分类方法[16]、支持向量机法[17]等分类器可以提高分类精度。此外,利用纹理和地形[18]、历史土地利用数据[19]作为辅助，也成为分类过程中的重要参考指标。

本研究以民勤绿洲区及绿洲边缘防风固沙林为研究对象，通过植被覆盖度和基于CART自动决策树分类技术组合多源数据进行土地利用分类，开展防护林带监测与评价，探讨民勤绿洲防护带植被变化趋势，为绿洲防护带建设提供基础数据支撑。

1 材料与方法

1.1 研究区概况

民勤绿洲位于甘肃河西走廊东段，是石羊河延伸到沙漠腹地出现的非地带性景观，是指石羊河下游红崖山水库以下冲积成的狭长而平坦的绿洲带。甘肃民勤县东西北三面被巴丹吉林和腾格里沙漠包围，近九成面积为荒漠化土地，是中国重要的沙尘暴策源地之一。民勤绿洲像插入沙漠中的楔子，阻挡两大沙漠的合拢。20世纪50年代，民勤绿洲边缘(这里指主要受西北风危害的绿洲西北边缘一线)地下水位为1～3 m，主要植被为红柳(Tamarixramosissima)，俗称柴湾，植被生长旺盛，绿洲有石羊河多条水系的地表径流滋润，生态环境良好。1958年建红崖山水库，流域终端湖青土湖来水断流。而后，随着流域农业发展，耕地面积逐步扩大，地表径流不能满足流域农业生产需求，绿洲内来水逐渐减少，靠打井提灌满足农业生产，地下水位逐渐下降,时至今日，地下水位下降至25 m左右，植被无法利用，造成衰败死亡。在这个过程中，为了绿洲不受风沙危害，20世纪60年代，开始不间断在民勤绿洲边缘大面积治沙造林，2007年石羊河流域重点治理工程实施，促进绿洲防护带建设进展，几十年的防护带建设，绿洲植被如何变化，需要有一个量化的说明。本研究主要探讨该区20年来植被变化。研究区范围是结合卫星遥感图，通过绿洲外围分布的固沙林实地调查生长状况，包括植被的盖度、密度等来确定，见图1。

图1 民勤绿洲与固沙林区Fig.1 Minqin oasis and sand-fixing forest

1.2 数据预处理

Landsat TM数据分辨率为30 m，这种数据对于中尺度区域开展遥感定量监测，具有成本低，处理周期短，效率高等优点。因此，本研究选择1996、2006和2016年Landsat TM、OLI影像。云覆盖均小于10%，成像质量较好，时间范围为7月末至8月初，每期3景。由于选用的数据属于未经处理的原始数据，因此，先对9景数据进行辐射校正和大气校正，并采用二次多项式进行几何校正，将校正误差控制在一个像元以内。

选用12景2016年8 m多光谱、2 m全色高分1号数据，时间范围在6-9月，通过大气校正、正射校正等处理最后融合为2 m多光谱影像。对比两种分辨率下计算结果的误差大小。所有数据经过镶嵌和裁剪处理，获得研究区影像。

1.3 研究方法

1.3.1植被覆盖度的计算 在原有沙地、裸地的基础上栽植人工固沙林，以及在低覆盖度土地栽植固沙林恢复植被，在卫星遥感图上表现为植被覆盖度的增加。因此利用遥感影像反演植被覆盖度，监测和评价绿洲防护带的植被变化。本研究采用基于NDVI的像元二分模型[9]计算植被覆盖度。计算公式为：

(1)

式中：fc为植被覆盖度，NDVI是通过计算得到的3期研究区NDVI实际值，NDVIveg和NDVIsoil分别对应研究区范围内纯植被的NDVI值和土壤的NDVI值。根据研究区NDVI值的分布范围，以置信度5%分别获取3期NDVI上下阈值，作为NDVIveg和NDVIsoil。最后利用Arcgis中的重分类工具，结合民勤绿洲植被类型实际分布和荒漠化地区土地利用类型划分方法，将植被覆盖度划分为5级(包括低植被覆盖度0～5%、中低植被覆盖度5%～15%、中等植被覆盖度15%～30%、中高植被覆盖度30%～60%、高植被覆盖度60%以上)。此种分类方法，较为精确的反映干旱荒漠地区各覆盖度类型空间分布和面积变化。

1.3.2基于CART决策树分类与面积变化统计 CART决策树分类技术，是通过对测试变量和目标变量构成的训练数据集的循环二分形成二叉树形式的决策树结构[20]。本研究结合野外实测数据选取训练样本和多源数据，基于ENVI 5.3的CART计算工具，对民勤绿洲进行土地利用类型分类。并对分类后结果进行面积统计，最后建立3期影像各地物类型面积转移矩阵。多源数据[14]相应计算公式为：

图2 1996年植被覆盖分布Fig.2 Vegetation cover distribution in 1996

图3 2006年植被覆盖分布Fig.3 Vegetation cover distribution in 2006

图4 2016年植被覆盖分布Fig.4 Vegetation cover distribution in 2016

(2)

(3)

(4)

式中：NDVI为归一化植被指数；NDWI为归一化差异水体指数(normalized difference water index)；NDBI为归一化裸露指数(normalized difference bare index)。bandgreen和bandswir1、bandnir、bandred分别对应TM影像中的绿波段、中红外波段、近红外波段和红波段。NDVI可以消除地形和大气的影响，对植被信息较为敏感。NDWI可以较好地消除土壤和植被背景的影响，从而突出水体信息；NDBI可以有效地分辨裸地信息。

2 结果与分析

2.1 不同覆盖度类型的空间分布与面积变化

根据植被覆盖度计算方法，利用Arcgis进行密度分割和制图，获得3期植被覆盖度图(图2～4)。

根据3期植被覆盖度图，利用Arcgis重分类工具，对各覆盖度类型的面积计算统计。统计结果见表1。

根据图2和表1可以看出：1996年民勤绿洲防护带的低覆盖度植被主要分布在绿洲外围，占总面积的50.57%。中低覆盖度植被主要分布在中部、西南部地区，占总面积的25.90%。中等覆盖度植被主要分布在西南部、中西部、东南部地区，占总面积的6.48%。中高覆盖度植被主要分布在农田周边和西南部地区，占总面积的5.49%。高覆盖度植被主要分布在中部农田附近，占总面积的11.56%。

根据图3和表1可以看出：2006年民勤绿洲防护带的低覆盖度植被主要分布在绿洲外围且面积较大，范围较广，占总面积的67%，中低覆盖度植被主要分布在中部农田附近和东南部地区，所占比例8.73%，中等覆盖度植被分布在东南部、西南部、中部地区，所占比例5.32%，中高覆盖度植被分布在东南部、中部绿洲内部，所占比例5.35%，高覆盖度植被主要呈条带状分布在中部农田附近，所占比例3.62%。

图4和表1看出2016年各覆盖度类型：低覆盖度植被主要分布在西南部和中南部，所占比例49.39%，中低覆盖度植被主要分布在东部、中部、北部、东南部，范围较广，所占比例26.68%，中等覆盖度植被分布在东南部和中部，所占比例5.53%，中高覆盖度植被分布在东南部、西南部、中部地区，所占比例4.06%，高覆盖度植被主要分布于东南部和西南部农田附近以及中部地区，所占比例5.22%。

表1 1996、2006和2016年植被覆盖度等级面积及所占比例Table 1 1996, 2006, 2016 vegetation coverage area and proportion

表2 1996、2006和2016年不同植被覆盖度等级面积变化及所占比例Table 2 1996, 2006 and 2016 different vegetation coverage grade area changes and the proportion

结合图1～4和表2可以发现，1996-2006年10年间，低覆盖度和高覆盖度面积有所增加，分别增加32.48%和27.65%。中低覆盖度、中等覆盖度、中高覆盖度面积有所减少，分别减少66.28%、17.94%、2.61%，重点监测样区中，刘家地和老虎口由原来的中低覆盖度退化为低覆盖度，南湖乡部分区域由中低覆盖度转变为低覆盖度。2006-2016年10年间，低覆盖度和中高覆盖度面积有所减少，分别减少26.28%和24.00%。中低覆盖度、中等覆盖度、高覆盖度面积增加，分别增加了205.47%、4.12%、44.22%。重点监测样区中，青土湖区由低覆盖度改善为中低覆盖度和中等覆盖度，刘家地和老虎口由低覆盖度转变为中低覆盖度，南湖乡由低覆盖度转变为中低覆盖度。总体来看1996-2016年20年间植被覆盖度总体呈现上升趋势，整个主绿洲区生态环境朝着良性发展。

图5 1996年土地利用Fig.5 Land use in 1996

2.2 绿洲防护林带空间分布和面积变化

本研究主要监测绿洲防护带的空间分布和面积变化，因此将民勤绿洲土地利用类型分为3类：农田防护林[新疆杨(Populusalbavar.pyramidalis)]、防风固沙林[梭梭(Haloxylonammodendron)、白刺(Nitrariatangutorum)、柽柳(Tamarixchinensis)]、其他地类(农田、水域、裸岩、盐碱地、沙地)。此外，利用本研究采用的土地利用分类方法得到分类结果(图5～7)。

对各土地利用类型面积进行计算统计，结果见表3和表4。

根据图5和表3可知1996年各土地利用类型的空间分布和具体面积。其中农田防护林主要呈条带状零星分布在绿洲中部、中北部地区，占总面积的0.29%，防风固沙林主要分布在绿洲外围，分布范围较广，占总面积的63.64%，其他地类主要呈狭长形分布于绿洲区外围和绿洲中部地区，占总面积的36.07%。

图6 2006年土地利用Fig.6 Land use in 2006

图7 2016年土地利用Fig.7 Land use in 2016

结合图6和表3可以看出，2006年各土地利用类型：农田防护林主要分布于中部和中北部地区，分布范围较小，占总面积的0.72%，防风固沙林主要分布于绿洲西北部、东部边缘地带，东南部有小范围集中分布，所占比例74.10%，其他地类主要分布在绿洲中部、北部、东部边缘地带，东南部有少量分布，所占比例25.18%。

结合图7和表3可以看出，2016年各土地利用类型中，农田防护林主要分布于绿洲中部、中北部，占总面积的0.47%。防风固沙林主要分布于中部农田外围大部分地区，分布范围较广，面积较大，占总面积的76.73%。其他地类小范围分布于绿洲中部、北部、东南部、西部地区，占总面积的22.79%。

根据图1、5～7和表4看出，1996-2006年农田防护林面积增加了25.75 km2，增幅148.68%，说明在这10年间，阻止沙漠侵吞农田，开展农田防护林建设，取得一定成果。防风固沙林分布范围明显向绿洲边缘扩大，增加面积630.12 km2，增幅达16.47%。此外，老虎口和南湖乡，防风固沙林面积增加明显，青土湖和刘家地以沙地等其他地类为主，没有较大变化。2006-2016年的10年间，农田防护林面积出现减少，分布范围发生萎缩，实际减少14.64 km2，减少比例33.99%。防风固沙林面积分布范围进一步扩大，实际增加面积有158.38 km2，增幅3.55%。此外，青土湖、刘家地、老虎口，由裸地改善为防风固沙林地，南湖乡也有部分裸地转变为防风固沙林。总体上看，1996-2016年的20年间，农田防护林出现先快速增长后小幅减少趋势。防风固沙林分布范围和面积出现先快速增长后缓慢增长的趋势，说明民勤绿洲，20年来人工造林工程成果显著，有效遏制沙漠向绿洲推进。

表3 1996-2016年土地利用面积统计Table 3 Land use area statistics from 1996 to 2016

表4 1996-2016土地利用面积变化Table 4 Change of land use area from 1996 to 2016

通过表5可以发现，2006年较1996年，农田防护林面积有52.53%保持不变，有近一半的面积发生转移，其中大部分面积转化为防风固沙林和其他地类，所占比例分别为26.57%和20.90%；防风固沙林有88.97%的面积保持不变，其余主要转变为其他地类，所占比例10.36%。2016年较2006年，农田防护林只有10.87%的面积保持不变,主要转变为防风固沙林和其他地类，所占比例65.77%和23.38%；防风固沙林有82.33%面积保持不变，其余主要转变为其他地类，所占比例17.42%。

2.3 植被覆盖度和土地利用分类精度验证

高分1号融合后的2 m多光谱数据，较中低分辨率数据，具有地物纹理和轮廓更加清晰的优点，可以较为准确的判断地物类别以及更加准确的获取纯植被和纯土壤的NDVI值。因此，适合基于像元二分模型求取植被覆盖度结果和中低分辨率土地利用分类结果的精度验证。本研究将处理后的高分1号数据利用ENVI进行波段计算，获取研究区NDVI值的范围，计算植被覆盖度用于对Landsat TM/OLI数据处理结果的精度验证。同时在高分1号影像均匀选取每类地物的样点(可以弥补中尺度区域选点困难，时间和人力成本高的缺点)，结合野外实测数据建立验证样本，对基于CART决策树分类方法得到的土地利用分类结果进行精度验证。

根据表6植被覆盖度精度验证，可以看出，低覆盖度、中低覆盖度和中等覆盖度分别有83.14%、87.83%和84.47%被正确分类。中高和高覆盖度分类精度较低，有79.17%和73.42%被正确分类。总体分类精度为81.62%，Kappa系数为0.78。说明基于二分像元模型计算植被覆盖度，总体分类精度良好，但是由于TM数据分辨率较低，混合像元的存在对于提取纯植被和纯裸地NDVI值有一定的干扰。

根据表7土地利用分类精度可以发现，防风固沙林分类精度最高，有91.58%被正确分类。农田防护林分类精度最低，只有65.21%被正确分类，这是因为农田防护林为零星分布，一般宽度为1～3 m，而TM数据单位像元的边长为30 m，所以很难有纯净的农田防护林像元存在，只有范围较大的林区和经济林田可以较好地被分离出来。总体分类精度为74.95%，Kappa系数为0.71.说明总体分类精度较好。

表5 绿洲防护带面积转移矩阵Table 5 Oasis protection zone area transfer matrix (%)

表6 植被覆盖度精度验证Table 6 Verification of vegetation coverage accuracy (%)

注：总体分类精度81.62%，Kappa系数：0.78。

Note：Overall accuracy 81.62%; Kappa coefficient 0.78.

表7 土地利用分类精度验证Table 7 Land use classification accuracy verification (%)

注：总体分类精度74.95%，Kappa系数：0.71。

Note：Overall accuracy 74.95%; Kappa coefficient 0.71.

3 讨论与结论

3.1 讨论

甘肃民勤地区处于石羊河下游，三面受巴丹吉林和腾格里沙漠包围，是中国四大沙尘源区之一。特殊的地理环境，迫使民勤人民长期开展绿洲生态环境治理，2007年以来石羊河流域重点治理规划实施，治理效果、生态环境变化趋势等是人们关注的热点问题，本研究利用Landsat TM结合高分1号数据，探讨20年来防护林带面积及其植被覆盖度变化，研究结果与滑永春等[21]、韩涛等[22]、王红霞等[23]的研究结论基本一致。

利用中等分辨率遥感数据，基于光谱特征进行监督分类，在建立分类样本时容易出现“同物异谱”、“同谱异物”现象，导致各类别可分离度低，最终影响分类结果的精度。尤其是在荒漠地区，植被稀疏，覆盖度低，获取的植被NDVI值会低于一般理论值。因此，单一数据辅助分类效果不佳。CART决策树分类，通过组合多源数据和目视解译建立分类样本进行土地利用分类，较传统分类方法，精度容易把握，分类结果质量较高。

3.2 结论

1996-2016 年间，民勤绿洲及其外围防护林带各等级植被覆盖度中，低覆盖度植被的面积先增后减，整体稍有减少，实际减少面积比例为2.33%，这一类型为流动沙丘、稀疏白刺群落以及一些沙地草本植物。中低覆盖度(增加面积比例3.00%)、高覆盖度(增加面积比例84.10%)植被的面积整体处于增加趋势；中等覆盖度(减少的面积比例14.56%)、中高覆盖度(减少的面积比例25.98%)植被面积整体处于减少趋势。民勤绿洲地下水位持续降低，1996、2006、2016年地下水位分别为14、21、23 m左右，已处于民勤生态地下水位临界范围(7～12 m)[24]，靠地下水生长的植物(如柽柳、人工梭梭林)衰退死亡，一些种群(白刺)逐步适应依靠降水生长。低覆盖度植被生长主要受降雨的影响，20年来变化不大；受地下水位下降的影响，原先营造的梭梭林退化，靠地下水位生长的柽柳林死亡，是中等覆盖度和中高覆盖度植被面积减少的原因；大面积营造固沙林、封沙育林育草、退耕还林等工程措施促进了低覆盖度和高覆盖度植被面积增加。

总体来看，20年来，民勤绿洲农田防护林面积先增加后减少，整体处于增加趋势，增幅64.09%，面积增加了11.10 km2；防风固沙林面积整体在增加，先快速增长，后缓慢增长，增加了20.58%，增加面积为788.50 km2。随着石羊河流域综合治理规划的实施，节水灌溉工程减少了绿洲防护林的生态用水量，农田防护林衰败死亡，但大面积的生态经济林建设，绿洲防护林面积增加。民勤生态环境建设，大面积固沙造林，绿洲外围形成了近10 km宽的防风固沙林带，防风固沙林面积增加。

在绿洲防护带面积转移变化中：1996-2006年10年间，有47.47%农田防护林转变为防风固沙林和其他地类，所占比例分别为26.57%和20.90%。2006-2016年的10年，有89.13%农田防护林转变为防风固沙林和其他地类，所占比例为65.77%和23.38%。1996-2006年的前10年，防风固沙林主要转变为其他地类，所占比例为10.36%；2006-2016年的后10年，防风固沙林有17.42%转变为其他地类。

利用同期高分1号数据计算植被覆盖度，对Landsat TM数据结果进行精度验证。植被覆盖度精度验证: 总体分类精度为81.62%，Kappa系数为0.78。说明基于二分像元模型计算植被覆盖度，总体分类精度良好。土地利用分类精度验证：防风固沙林达91.58%，农田防护林只有65.21%。