摘 要
针对传统植被覆盖度估算方法只能估算高密度覆盖的绿色光合植被的问题，该文提出了一种可进行绿色光合植被、非绿色光合植被和裸地覆盖度估算的方法，选取博尔塔拉蒙古自治州2010年和2016年的Landsat TM/OLI影像，通过构建NDVI-DFI特征空间提取端元特征值，运用像元三分模型估算，并与像元二分模型估算结果和实地采样估算结果对比。结果表明：像元三分模型估算与实际情况相符，且精度较好。2016年博州地区植被覆盖度较2010年增长显著，生态环境得到明显改善。像元三分模型能够较好地估算光合/非光合植被覆盖度，提高遥感获取植被信息的能力，为科学评估生态环境质量提供参考。
引用格式
唐梦迎，丁建丽，夏楠，等．博尔塔拉蒙古自治州植被覆盖度估算［Ｊ］．测绘科学，２０１９，４４（７）：７４-８１．
正文
引言
植被覆盖度是干旱区荒漠化过程发展或逆转的重要指示因子。对干旱区生态系统来说，植被在物质能量交换、大气及水文循环等方面扮演的作用尤为显著[1-2]。如何掌握植被覆盖信息，对区域生态环境保护具有重要意义。遥感技术具有易操作，成本低，可周期性重复的特点能够宏观地提取植被覆盖信息。常用的像元二分法[3-4]在提取混合像元植被盖度信息时仅考虑光合植被（PV）和裸地（BS）两种组分，而自然界植被不仅包含绿色健康的光合植被（PV）也包含干枯的非光和植被（NPV）[5-6]。在干旱区存在大量NPV[7]，忽略NPV而估算植被覆盖度是不严谨的。
文献[8]首次提出运用像元三分模型估算植被覆盖度，并构建基于高光谱和多光谱图像的像元三分模型，较好地估算了澳大利亚热带稀疏草原区的fPV和fNPV的时空分布。文献[9]利用NDVI-DFI特征空间构建像元三分模型，估算出fPV和fNPV时空动态变化与牧草物候发育特征相吻合，说明该方法能较好地提取区域的fPV和fNPV。文献[10-11]运用线性模型可以较好估算区域局部的NPV。文献[12]基于多源数据计算博州地区的植被覆盖度，分析了不同月份植被覆盖度的空间差异，但对于植被变化的驱动因素考虑较少。上述研究表明了像元三分模型对于估算植被覆盖度具有一定的实用性，但相关研究较少，目前仅应用于草原植被覆盖度的估算及草原火灾预测，而干旱区同样存在较多NPV和BS，对此需进行进一步探讨。
本文研究区设定为“一带一路”沿线新疆地区发展的重要窗口之一—博尔塔拉蒙古自治州，借助不同时间尺度的landsat TM/OLI数据，利用像元三分模型估算植被覆盖度，通过高清影像和实地采样分析评价估算结果，探讨影响该区植被覆盖度变化的驱动因素，为推广以NDVI-DFI特征空间为技术条件的像元三分模型在新疆干旱半干旱地区的适用性提供研究基础，以及为该区生态环境监测和保护提供数据支撑。
1 研究区概况
博尔塔拉蒙古自治州，以下简称博州，下辖博乐市、精河县、温泉县和阿拉山口市，地理坐标范围是79°56′～83°52′E， 44°01′～45°24′N，北部和西部与哈萨克斯坦国接壤，东部和南部分别与塔城地区和伊犁地区相毗邻。地貌特征三面环山，东南部玛依力山，西北部阿拉套山，南部婆罗科努山。地形西高东低，最大高程差近1 300 m。该地属于温带大陆性气候，年均气温5.6℃，极端最高气温44℃，极端最低气温—36℃，年均降水量181 mm，年均潜在蒸发量1 500~2 000 mm。主要土壤类型为盐碱土、沙质土、钙积土和石膏土等。主要植被类型为耐旱耐盐碱的灌木和小半灌木，如梭梭(Haloxy lon ammo dendron)、柽柳（Tamarix laxa willd）、盐节木（halocnemum strobilaceum）和芦苇（Phragimites australis）等。农作物主要以棉花为支柱产业，此外还有玉米、小麦和枸杞等高产量作物。
图1 研究区及采样点分布图
Fig.1 Distribution Map of Study Area and Sampling Points
2 数据与方法
2.1数据来源和处理
本研究中植被覆盖度估算所使用的遥感影像数据为Landsat TM/OLI数据集（空间分辨率为30 m，重访周期16 d），下载于http://glovis.usgs.gov。选取时间为2010年8月9日和2016年8月16日拍摄于当地植被生长较好的时期，且云量较低。首先对每景遥感影像进行辐射定标，将DN（digital number）值转化为反射率值。然后利用艾比湖流域的矢量边界，分别对遥感数据进行镶嵌和裁剪后得到研究区域影像数据。最后用阈值法对影像中的云量和水体作掩膜处理。
2.2 采样点布设
一些研究学者提出[13-14]，垂直照相法结合数字图像处理计算植被覆盖度的方法便捷客观，计算精度高，而考虑到目估法测量植被覆盖度的不确定性，因此需采点验证。采样点的布设不仅要遵循均匀性原则，还需要考虑道路的可达性，区域代表性和研究区范围。考虑上述因素，以艾比湖、主要河流周边、绿洲农场以及沙漠边缘为重点区域，布设49个采样点，获取植被覆盖度的垂直拍摄照片，主要围绕0，0.25，0.50，0.75和1.00共5个等级范围进行验证，并用GPS记录采样点坐标。
2.3
1）归一化植被指数。归一化植被指数（normalized difference vegetation index，NDVI）是遥感观测PV的一个重要的指示因子，其特点是操作简单，植被指示性和地域适用性都较强，公式见式（1）。
式中：BandNIR是近红外波段，BandR是红光波段，其中Landsat TM影像的BandNIR都是第4波段，BandR是第3波段，而Landsat OLI的BandNIR和BandR分别是第5波段和第4波段。
2）干枯燃料指数DFI（dead fuel index）。文献[15]研究发现NPV和BS在650 nm和850 nm处的反射率高于PV；在1 650 nm NPV和BS的反射率低于PV，在2 100 nm，NPV的反射率在PV和BS之间。这些特征波段对应于Landsat TM波段由低到高依次为BandR（红光）、BandNIR（近红外）、BandSWIR1（短波红外1）和BandSWIR2（短波红外2），即Band3、Band4、Band5和Band7；对应Landsat 8_OLI的波段依次为Band4、Band5、Band6和Band7。其计算公式如式（2）所示，式中将结果乘以100，目的是为突显DFI值的差异性。
3） NDVI-DFI特征空间。假定待测像元有NPV、PV和BS 3个组分，其NDVI和DFI指数符合线性关系，其各种比例组合像元的NDVI-DFI特征空间会表现为三角形（图2）。NPV的NDVI低、DFI高，位于三角形的左上角；BS的NDVI、DFI均很低，位于三角形的左下角；PV的NDVI高、DFI几乎接近0，位于三角形的右侧。
3） NDVI-DFI特征空间。假定待测像元有NPV、PV和BS 3个组分，其NDVI和DFI指数符合线性关系，其各种比例组合像元的NDVI-DFI特征空间会表现为三角形（图2）。NPV的NDVI低、DFI高，位于三角形的左上角；BS的NDVI、DFI均很低，位于三角形的左下角；PV的NDVI高、DFI几乎接近0，位于三角形的右侧。
图2 NDVI-DFI特征空间示意图
Fig.2 Sketch map of NDVI-DFI Characteristic Space
4） 提取端元特征值。建立以NDVI和DFI为特征因子的特征空间，随后提取位于3个端元的特征值。利用遥感方法提取端元值通常包括①几何顶点提取法；②纯净像元指数法（pixel purity index method，PPI）；③连续最大角凸锥选取法。由于多光谱遥感影像各波段具有严重的噪声[16-17]，采用最小噪声分离法（minimum noise fraction，MNF）消除其影响。在提取端元值方法中PPI算法较为成熟，且能依据像元自身特点选取端元值，可极大地减少人为选择造成的误差。因此，为消除Landsat波段间的相关性，运用MNF法，选取PPI提取端元特征值。
5） 植被覆盖度等级表。结合干旱半干旱区绿洲植被特有的生态特征，为更精细研究的需要，采用文献[18]划分的植被覆盖等级，共分为5级，分别为极低覆盖度、低覆盖度、中覆盖度、中高覆盖度和高覆盖度。分类标准见表1。
表1 植被覆盖度等级划分
Tab.1 The Grades of Vegetation Cover
3结果与分析
3.1 端元值确定及NDVI-DFI特征空间分析
端元值确定要通过PPI计算获得[19-20]。首先运用ENVI 5.1软件MNF Rotation工具实现数据的变换降维，生成大量穿过数据集合内部的随机测试向量，之后仅对其中6个主要波段进行投影，设定迭代次数为2 000，阈值系数为3。将PPI指数大于3的像元视为端元，取各个顶点纯净端元的平均值作为相应端元的特征值，得到表2。表中PV端的NDVI最大，说明NDVI是表征PV的主要因子。NPV端对应的DFI最大，说明DFI是表征NPV的主要因子。
表2 端元特征值信息统计
Tab.2 The Statistics of End Member Eigenvalue Information
由2010年NDVI-DFI特征空间散点图（图3）得出，该特征空间散点分布较2016年呈现出不太理想的三角形结构，有大量的散点分布于NDVI值在0~0.25范围内，且无明显分布规律，表明植被生长盖状况不佳，大量的裸地和干枯植被是造成不规则分布的原因之一，而二者主要分布于NDVI<0.25的区域。当NDVI值在0.25~0.8时，随着NDVI值增大，DFI值整体相应地随之减小，即植被覆盖度越高，干枯植被覆盖越少。2016年该特征空间的散点分布呈现出较为规则的三角形，说明2016年植被生长状况整体较好。NDVI值在0~0.2阈值范围内，NDVI值增大，DFI随之增大；当NDVI值在0.2~1之间时，随NDVI值增大，DFI呈递减趋势。
图3 NDVI-DFI特征空间散点图
Fig.3 Scatter Plots of NDVI-DFI Characteristic Space
3.2 植被覆盖度估算
1)植被覆盖度估算结果分析。将上述端元特征值带入像元三分计算模型，得到的fPV、fNPV和fBS的RGB合成图（图4）。图中亮色表示有PV覆盖更多的区域，灰色表示有NPV覆盖更多的区域，深灰色表示裸地覆盖更多的区域，其间混合的颜色表示混合3种类型的区域，黑色是将水体和积雪进行掩膜处理的部分。
由图4可知，2010年fPV主要集中在人类活动区，包括农田、林地等；fNPV主要分布在山地和山前洪积冲击扇；fBS分布在博州北部和东部的山区及山前平原，以及东南部的沙漠。相比2010年，2016年fPV、fNPV和fBS变化显著，其中fPV沿博尔塔拉河向南北方向扩展，东南部和西南部山区也出现大量绿色植被，覆盖面积增长了15.37%。fNPV随着fPV面积的增加而面积发生了一定变化，相比2010年略微减少了4.81%。fBS主要分布在艾比湖周边以及绿洲边缘，面积显著减少了15.37%。（该数据由通过ENVI软件的统计工具得出）
图4 fPV、fNPV和fBS估算结果图
Fig.4 The Map of fPV, fNPVand fBS Evaluating Result
绿洲荒漠交错带生态关系错综复杂，生态系统相对敏感。选择绿洲荒漠交错带典型区域进行局部分析（图5）。从中看出，2016年影像中部绿洲范围植被较2010年有所增多，具体表现为PV平均增加了18.68%，NPV平均增加10.61%，BS平均减少29.28%。绿洲面积扩大，生态系统得到一定改善。
图5 绿洲荒漠交错带局部植被覆盖状况图
Fig.5 Local Vegetation Cover Map of Oasis Desert Ecotone
2）与像元二分模型估算结果的比较。为探究像元三分模型对于光合植被、非光合植被和裸地的区别能力，将上述植被估算结果与像元二分模型反演结果进行对比，参考相关反演方法[3]得到2016年像元二分模型估算植被覆盖度空间分布图。从图6中看出该反演结果与图4中2016年植被覆盖状况近似，说明二者均能一定程度地表述研究区植被覆盖度。为进行详细对比，另选取一绿洲荒漠交错区域，并添加谷歌地球高清影像进行对比验证。该区域为裸露山地，植被覆盖度极低，周围有农田包裹，经实地考察发现其西侧为撂荒的旱地，早期种有打瓜和枸杞，因此存在一些植物残体。像元二分估算结果在该区域未显示出植被分布差异，均表现为极低的植被覆盖，无法判断是否为裸地或非光合植被。而像元三分模型结果显示出更多的裸地和非光合植被信息，比较符合实际情况，表明像元三分模型能较好地分辨光合植被、非光合植被和裸地。
图6 基于两种模型的绿洲荒漠交错带局部对比分析
Fig.6 Comparison of local oasis desert ecotone based on the two models
3） 估算结果验证。将室外采样结果与对应点的fPV像元信息进行比较（通过ArcGIS 10.2 Extraction工具实现），得出遥感估算的fPV最小值、最大值、均值和标准差分别是0.000、1.000、0.221和0.302；采样的fPV最小值、最大值、均值和标准差分别是0.000、1.000、0.177和0.199。二者相关系数为0.767（P<0.01），表明二者存在较好的相关性。再将该两种结果进行线性拟合（见图7）。图中可以看出遥感反演结果与实地采样结果趋于y=0.6869x+0.112直线，其决定系数R2为0.7515。另外可发现遥感方法对于低植被覆盖区（fPV<0.05）估算结果与实际存在一定偏差，可能与复杂的地表覆盖特征有关。但对于博州地区而言，像元三分模型的估算效果是符合精度要求的。
图7 fPV与采样估测值的线性拟合图
Fig.7 Linear Fitting of fPVValue and Sampling Estimation Value
3.3数据分析
1） 植被覆盖度等级变化。通过对fPV反演结果的统计得出，博州2010年和2016 年两个时相的平均植被覆盖度分别为13.16%和28.52%。再按照植被覆盖等级将各年fPV分为5级，得到表3。由表可知2010年I级植被所占百分比最大，为75.42%（面积16884.81 km2），2016年同样是I级占比较大，为38.47%（9011.22 km2），但相比2010年显著减少了36.95%（7873.59 km2）。除I级植被面积有所减少，2016年其余各级植被覆盖度均有所增加，表现为II级植被增加最明显，增加了17.65%，III级植被增加了8.33%，四级和五级植被分别增加了5.57%和5.39%。以上说明自2010年到2016年博州地区植被覆盖度增长显著。高等级植被的显著增加，说明该地区生态环境得到了一定的改善。2010年~2016年6年间各级植被覆盖度年均增长率分别为-9.94%、23.74%、17.14%、16.93%和10.19%。
表3 植被覆盖度分级统计
Tab.3 The Grades of PVStatistics
2） 植被覆盖度的年际转移状况。利用 ENVI 软件转移矩阵功能模块，计算出2010年与2016年研究区不同等级植被覆盖度的面积转移状况进行统计和分析（表4）。与2010年相比，2016年共有61 938.15 km2的植被未发生变化。在恢复过程中，IV向V的转换面积最大，为4 011.63 km2，占总面积的26.67%，III向V的转换面积次之，为2 671.60 km2占总面积的17.76%，II向V的转换面积，为1 467.39 km2占总面积的9.76%；在退化过程中，II向I的转变面积最大，为575.21 km2占总面积的43.43%， V向IV的转变面积次之，为544.46 km2占总面积的19.88%，IV向II的转换面积，为457.17 km2占总面积的0.8%。总的来说，2010—2016年，由低级转向高级植被覆盖面积>高级转向低级植被覆盖面积，即恢复面积>退化面积。
3.4 驱动力分析
1)自然因素。降水和气温同样影响植被覆盖度的大小和植被类型的差异[21]。从博州各县市的气象数据中获取2010~2016年各年1~7月的降水量和气温数据，并绘制成图8。由图得，各年7个月的降水总量差异明显，但整体呈递增趋势，尤其是2016年降水量突破了这六年以来的最大值，在近十年也属于降水较多的年份。这是由于在2016年6月和7月博州出现爆发式降雨，强度大，范围广，持续时间长[22]。而当年博州各县市中，精河县和阿拉山口市降水量较常年多约60%，博乐市和温泉县较常年也偏多近10%。降水作为干旱区的气象指标对该地区植被生长具有至关重要的作用，故一定程度地解释了2010和2016年植被分布状况。2015年7月，新疆出现历史罕见的大范围、长时间、高强度的高温天气去哪将平原地区出现35℃以上的高温。而2015年7月博州地区出现最高气温突破历史极值[23]。
2）人为因素。植被覆盖度受自然因素影响，同时也受人类活动的影响。一方面政策驱动，国家提出的“18亿亩耕地红线”的严格保持，生态文明建设的有效推行以及当地部门对艾比湖国家级湿地保护区的有效管理，有利地保护了研究区的耕地、林草地和艾比湖流域的生态环境，为自然生态恢复提供了有力保障。另一方面，近年来全民环保意识不断增强，人类活动产生的废弃物得到有效处理，加之对水资源的合理利用也极大促进植被健康生长。但在干旱区生态环境极为脆弱，为维持区域生态可持续仍需加强环境保护力度。
4结束语
通过构建NDVI-DFI特征空间的方法，估算博州地区光合/非光合植被覆盖度取得较好结果，尤其在绿洲和荒漠交错区域，光合植被、非光合植被和裸地三者间的覆盖变化表现得较为突出。相比像元二分模型，该方法同样适合研究区植被覆盖度估算。
2016年较2010年光合植被覆盖度增加，非光合植被和裸地减少，生态环境得到改善。将估算结果与实地采样结果进行对比，二者具有较好的线性关系。植被动态变化监测中，2016年较2010年II级、III级、VI级、V级植被有所增加，I级植被减少。降水量的增多导致植被覆盖度的增加。
本文系统地分析了博州地区的植被覆盖状况，对评估该区域的生态环境综合治理效果具有重要意义。然而目前NDVI-DFI特征空间中的DFI多用于草原干枯植被的遥感估算，本研究针对干旱区绿洲和荒漠地区的DFI估算具有一定的适用性，但其作用机制是否与草原植被相同，仍需进一步探讨。
参考文献（略）
关于《测绘科学》
主管：中华人民共和国自然资源部
主办：中国测绘科学研究院
邮箱：chkx@casm.ac.cn
网站：http://chkd.cbpt.cnki.net