基于混合像元分解的南方地区植被覆盖度遥感监测——以广州市为例

张志新，邓孺孺，李 灏，陈 蕾，2，陈启东，何颖清

(1.中山大学地理科学与规划学院，广州 510275;2.国家海洋局南海海洋工程勘察与环境研究院，广州 510300)

基于混合像元分解的南方地区植被覆盖度遥感监测
——以广州市为例

张志新1，邓孺孺1，李 灏1，陈 蕾1，2，陈启东1，何颖清1

(1.中山大学地理科学与规划学院，广州 510275;2.国家海洋局南海海洋工程勘察与环境研究院，广州 510300)

通过测量图像端元的地表反射率，对遥感图像进行精确大气校正;在对混合像元分解模型进行改进的基础上，建立了基于地表反射率的线性混合像元分解(Liner Spectral Unmixing，LSU)模型，有效地避免了因大气时间、空间差异所造成的多时相误差，实现了多时相对比;通过增加土壤湿度因子，消除了土壤湿度差异造成的误差，建立了适用于南方地区的植被覆盖度遥感监测模型。实地验证结果表明，该模型具有较高的精度。将该模型运用于广州市1998～2009年植被覆盖度时空变化监测，认为城市化、大型工程建设与陡坡开荒是造成广州市植被覆盖度变化的主要原因。

植被覆盖度;线性混合像元分解;地表反射率;土壤湿度;Landsat TM图像

0 引言

植被是城市生态系统中的重要组成部分。合理的植被景观空间布局，能够有效地缓解城市热岛效应，改善城市热环境。城市发展与城市化是城市植被覆盖度与植被景观变化的主要因素，尤其在快速城市化过程中，工程建设的大量实施，耕地面积的大幅减少，造成了城市内部及其周围地区植被覆盖度的大幅降低。植被覆盖度是区域内植被垂直投影面积占地表面积的百分比［1］，对其定量监测于城市生态环境的评估与统筹规划具有重要意义。

依据植被光谱信息与植被覆盖度之间的关系对大范围区域进行遥感监测，是获取区域植被覆盖度的重要手段［2］。目前较为普遍的植被覆盖度遥感监测方法有经验模型法、植被指数法、混合像元分解法和植被冠层模型法等［3－8］。经验模型法适用于小区域监测，对实测数据依赖性强;植被指数容易受到大气、土壤湿度等的影响和制约［9］，导致其与实际植被覆盖度的关系不清，因此植被指数法的定量化程度较低;植被冠层模型通常依据植被叶面积指数(LAI)计算植被覆盖度，主要分为以个体随机集合为特征的几何光学模型、以均匀散射薄层为特征的辐射传输模型及二者的混合模型等，该模型充分考虑了叶片与冠层的结构以及辐射传输等问题，但模型结构复杂、变量过多且难以测量;混合像元分解法则通过计算混合像元中各种地物占像元面积的百分比(丰度)［10］实现植被覆盖度的反演，该方法简单易行，主要针对植被指数或地物光谱特征对混合像元进行分解。

基于植被指数的混合像元分解法通过计算各种基本地物类型的植被指数(例如纯土壤覆盖像元的植被指数 NDVIsoil，纯植被覆盖像元的植被指数NDVIveg等)对混合像元的植被指数的贡献，实现对植被覆盖度的估计［11－15］。由于应用了植被指数，该方法同样存在定量化程度不高的问题。少数研究依据地物光谱特征建立混合像元分解模型［16－19］，该方法物理机制明晰，但仍存在以下不足:①此类模型一般通过直接读取或图像分析的方式［20－22］从遥感图像上获得图像端元的光谱特征，不考虑大气的时间、空间差异，导致产生时序变化监测误差;②此类模型通常不考虑土壤含水量对土壤光谱的影响，不宜用于土壤湿度较大、土壤含水量空间差异显著的我国南方地区的植被覆盖度监测。

为了更好地实现我国南方地区植被覆盖度的长期、定量遥感监测，本文对已有线性混合像元分解(Liner Spectral Unmixing，LSU)模型进行了如下改进:①测量各种基本地物(图像端元)的地表反射率，并将经过大气精校正的遥感图像的地表反射率作为模型的输入值，以消除大气时空差异的影响，提高时序变化监测结果的可对比性;②简化含水量与土壤反射率的关系，建立增加了土壤湿度因子的LSU模型。

以广州市中北部地区为实验区，借助于林地冠层分析仪测量结果对本文建立的LSU模型进行检验。然后利用该模型对广州市1998、2004和2009年3个时相的Landsat TM图像进行植被覆盖度的反演，并对结果进行了分析。

1 基于图像端元的LSU模型

1.1 增加土壤湿度因子的LSU理论

在遥感探测中，若忽略地形的影响，并将地表视为朗伯体，则行星反射率r可表示为

式中，T为电磁波入射与反射过程中的大气综合透过率;ρu为上行路径散射反照率;R为地表反射率。

大气中的可变成分(水汽、气溶胶、云雾水滴等)的分布随时间和空间的变化较大，导致T与ρu具有很强的时空不确定性，造成时空变化监测误差，因此必须对遥感图像进行大气精校正，以获得本文中LSU技术的分解对象——混合像元的地表反射率R。

线性光谱模型利用一个线性关系表达遥感系统中一个像元内各地物的类型、比例与光谱响应［23］。该模型假设在瞬时视场下，混合像元的地表反射率是由该像元内的各种图像端元的地表反射率线性组合而成，其中各种图像端元在该像元内所占的面积比为线性方程的权重［24，25］，即

式中，e为随机误差;Rv、Rs、Rw、…分别为植被、土壤、水体等图像端元的地表反射率;Av、As、Aw、…分别为植被、土壤、水体等在像元中所占的面积比。

传统的混合像元分解研究与应用，通常忽略不同湿度土壤间的反射率差异，但对于土壤湿度空间差异显著的我国南方地区，土壤湿度差异的影响应予以考虑。邓孺孺等［26，27］的研究指出，在土壤含水量明显低于土壤田间持水量时，湿润土壤反射率服从Beer定律，即

其中，Rso为干燥土壤的反射率;Vws为土壤体积含水量;aws为水分的吸收系数。

在土壤含水量较低的情况下，式(3)可改写为

式(4)表明，在土壤含水量较低的情况下，土壤体积含水量Vws与土壤反射率Rs之间线性相关，因此可以用干燥土壤与水分的反射率的线性叠加代替土壤的反射率，即将土壤作为干燥土壤和水分两种图像端元的组合。比照As，定义干燥土壤及土壤水分在混合像元中所占的面积比分别为Aso与Aws;定义水分的反射率为Rws，则式(2)可改写为

将式(5)改写为矩阵的形式，即

本文选用最小二乘法求解上述多元线性方程组，当方程组的离差平方和|E|2最小时，可以得到模型参数的估计表达式为

1.2 图像端元地表反射率的获取

利用美国ASD公司的ASD FieldSpec光谱仪测量各种地物的辐射亮度，然后计算其光谱反射率。为了与遥感图像的各个波段相对应，在Landsat TM传感器各个波段的波长范围内分别对每种地物的光谱反射率积分，并将积分结果相对于波段宽度的均值作为该地物在相应波段中的地表实测反射率。图1示出实测的(左)和TM图像上的(右)6类地物的光谱反射率曲线。

图1 图像端元实测(左)与TM(右)反射率光谱Fig.1 Spectral reflectance of basal land covers based on measurement(left)and TM(right)

原始图像端元包括图1中显示的6类地物。其中，灰白石砖、粗糙水泥地与柏油路面的实测反射率与TM影像各波段的反射率都极为相近，由于这3类地物是城市建筑与道路的主要组成部分，为了简化模型，将它们归并为同一类，即“建筑与道路”;并将上述3种地物反射率的均值作为“建筑与道路”的地表反射率。最终，本研究以植被、干燥土壤、纯净水体、建筑与道路为图像端元，将其地表反射率作为参数，用于建立LSU模型。

1.3 用于植被覆盖度反演的TM波段选择

在可见光区，红光被植被叶绿素强吸收，是光合作用的代表性波段;植被在0.55 μm附近有一较低的反射峰，近红外波段对植被差异和长势反映敏感;因此，TM 2、3和4波段是植被信息最为丰富的波段，而TM5波段对土壤含水量敏感。综上所述，本文选择上述4个波段用于植被覆盖度的反演。

对土地利用复杂的地区，采用同一种图像端元组合方式(例如植被与土壤)对全区的混合像元进行分解的方法并不符合实际，为此将地表覆盖简化为几种主要类型，定义其包含的图像端元类型，设定相应的反演波段(表1)。

表1 混合像元类型Tab.1 Types of mixed pixels

在反演植被覆盖度时，分别用各种图像端元组合类型对混合像元进行分解，然后根据式(6)与式(8)计算的各种类型反演结果计算离差平方和(|E|2)。比较并选取|E|2最小的类型作为混合像元的地表覆盖类型，再根据式(7)与式(8)计算像元的植被覆盖度。

2 精度验证

2.1 实验区与遥感数据

选择广州市北部作为实验区，区内土地利用与植被覆盖度空间差异显著:北部为丘陵地区，以林地、耕地与乡村建筑为主，植被覆盖程度较高;南部为广州市中心区，以城市建筑与道路为主，植被覆盖程度较低。首先，对实验区2009年1月2日获取的Landsat TM图像(图2(a))进行植被覆盖度反演，然后根据地表仪器测量结果，对反演精度进行定量评价。

图2 实验区 TM5(R)、TM4(G)、TM3(B)假彩色合成图像及验证点分布Fig.2 TM5(R)、TM4(G)、TM3(B)pseudo color composite image and distribution of the verifying points in the study area

2.2 数据预处理与植被覆盖度反演

2.2.1 数据预处理

根据式(7)可知，在进行植被覆盖度反演时，模型输入值为地表反射率R。地表反射率可以通过遥感图像计算获得(包括行星反射率计算与大气校正)。

(1)根据遥感图像DN值计算行星反射率r。首先对遥感图像进行几何纠正，获得具有地理坐标的图像数据，确保能够在图像上找到与地面验证点(由GPS确定其地理坐标)相对应的像元;然后对图像进行辐射校正，以消除传感器误差;最后计算行星反射率。

(2)利用大气校正方法消除大气程辐射的影响，计算地表反射率R。从行星反射率图像上选取纯净水体像元和位于阴影区的浓密植被像元作为暗像元，采用多暗像元法［28］对图像进行大气校正;然后根据式(1)计算地表反射率。

2.2.2 植被覆盖度反演

将预处理得到的地表反射率作为植被覆盖度反演模型的输入值，根据式(7)和式(8)计算实验区内各个像元的植被覆盖度。

2.3 精度检验

图3 野外测量与模型反演结果对比Fig.3 Comparison between field survey and Remote sensing retrieval results

则F的观察值为

根据F分布的临界值表，得到

当0.56＜F ＜1.79，表明遥感图像反演值与仪器测量值的方差之间无显著差异。全部验证点的平均相对误差为11.96%;位于植被覆盖度较高的白云山公园(图2(c))的验证点的平均相对误差为3.68%。结果表明，该模型能够较为准确地实现南方地区植被覆盖度的反演，对于高密度植被覆盖区的反演效果更好。

3 广州市1998～2009年植被覆盖度时空变化监测

近十年来广州市城市化进程不断加快，房地产侵占农田，丘陵地区陡坡开荒，北部原生植被减少。植被作为城市生态系统的重要组成部分，其覆盖度时空变化显著。本文按照2.2节中的数据预处理与植被覆盖度反演步骤分别对1998年12月、2004年11月和2009年1月获取的广州市Landsat TM图像进行了植被覆盖度的计算，并对结果进行了分析。

3.1 植被覆盖度时空变化

广州市各年份平均植被覆盖度分别为43.734%(1998年)、45.151%(2004 年)和 45.749%(2009年)，各辖区植被覆盖度如图4所示。

图4 广州市各辖区植被覆盖度对比Fig.4 Comparison of vegetation coverage in Guangzhou

从图4可以看出，广州市北部及东北部区域(从化市、增城市、萝岗区、花都区、白云区)平均植被覆盖度明显高于南部区域(番禺区、南沙区);市中心区(越秀区、天河区、黄埔区、海珠区、荔湾区)平均植被覆盖度明显低于非中心城区。

为了清晰地显示研究区域植被覆盖度时空变化，通过 ArcInfo软件进行栅格计算并分别生成1998～2004年与2004～2009年广州市植被覆盖度时空变化图(图5)。

图5 广州市植被覆盖度变化Fig.5 Change of vegetation coverage in Guangzhou

图5 显示，在1998～2004年和2004～2009年这2个时段中，广州市植被覆盖度空间变化特征较为一致:中部和南部平原地区植被覆盖度降低较为明显，呈斑杂状分布;北部和东部丘陵地区部分区域植被覆盖度增加;其他地区变化不明显。

3.2 植被覆盖度变化原因分析

参照原始遥感图像解译结果与实地考察发现，近十年来影响广州市植被覆盖度变化的3个主要因素为城市化、城市规划大型工程建设和丘陵地区陡坡开荒。图6～8为上述3种因素的典型区域及其表现(每幅图中的(a)与(b)分别为相应区域1998～2004年(时段1)与2004～2009年(时段2)的植被覆盖度变化空间分布图，(c)为相应区域2009年1月获取的TM5(R)、TM4(G)、TM3(B)假彩色合成图像)。

图6 城市化对植被覆盖度的影响Fig.6 Urbanization impact on the vegetation coverage

图6 选于增城市南端，为城市化对城市植被覆盖度图像的一个案例:在时段1中，西南部地区植被覆盖度整体呈大幅度降低(图6(a));在时段2中，植被覆盖度大幅度降低现象已经扩展到东北部地区(图6(b)。时段1中植被覆盖度降低的区域在TM图像中显示为紫色，已经成为城市建筑区;而时段2中植被覆盖度降低的区域呈白色版块状，为新开发区(图6(c))。截至2009年1月，该区大部分植被已经或正在被城市建筑替代。图6表明，城市化过程对植被的破坏具有一定程度的不可逆转性和时间上的延续性。

图7为城市发展进程中一类具有代表性的生态影响因素——大型工程建设。图上居中部分为番禺小谷围岛，于2003年1月开始实施广州市大学城建设工程。从图7(a)中可以看出，至2004年11月，该区域植被覆盖度较1998年整体下降，这主要是由广州市大学城的建设造成的。图7(b)中，中心区域植被覆盖度有所回升，且植被覆盖度显著增加的区域呈现较为规整的图斑，这与2006年12月大学城建设项目的完工以及后续两年的区域绿化工程密切相关;而在西北部，即位于小谷围岛南面的官洲岛的植被覆盖度变化则是因为自2008年底起广州国际生物岛工程建设的实施。与城市化相比，大型工程建设初期常常会造成较大程度的植被破坏，但工程建设中后期合理的绿地规划与建设，则能在一定程度上补救和改善建设初期对工程区生态的损害。

图7 大型工程建设对植被覆盖度的影响Fig.7 Large－scale architectural engineering impact on the vegetation coverage

图8 显示了陡坡开荒活动对植被覆盖度变化的影响。图8(a)与图8(b)中的深绿色斑块为成熟的陡坡开荒区，白色斑块为新开垦区，紫色斑块则为处于陡坡开荒过渡阶段的区域。陡坡开荒活动虽然能够带来短期、较高的经济效益，却破坏了地面植被，且又翻松了土壤，造成了产生严重土壤侵蚀的条件，导致水土流失加剧，引起局部生态环境恶化和山体滑坡、泥石流等地质灾害的发生，是威胁居民生命与财产安全的不安定因子。

图8 陡坡开荒对植被覆盖度的影响Fig.8 Reclamation impact on the vegetation coverage

4 结论

本研究在对线性混合像元分解(LSU)模型进行分析和改进的基础上，建立了具有较高精度的植被覆盖度监测模型。该模型具有以下特点:

(1)以地表反射率为LSU的作用对象，削弱了大气环境与遥感影像成像时间差异的影响，提高了时序变化监测结果的可对比性。

(2)根据土壤含水量与土壤光谱的关系，增加了土壤湿度因子，消除了土壤湿度空间差异的影响，适用于土壤湿度较大、土壤含水量空间差异显著的我国南方地区的植被覆盖度反演。

(3)地表仪器测量结果证明该模型具有较高的反演精度。利用该模型对广州市近十年来的植被覆盖度进行监测，分析了影响该地区植被覆盖度变化的主要因素，证明该模型用于我国南方地区植被覆盖度时空变化的监测可以获得良好的效果。

致谢:感谢蒋少锋在野外测量过程中给予的帮助。

［1］ 程红芳，章文波，陈 锋.植被覆盖度遥感估算方法研究进展［J］.国土资源遥感，2008(1):13 －18.

［2］ 秦 伟，朱清科，张学霞，等.植被覆盖度及其测算方法研究进展［J］.西北农林科技大学学报:自然科学版，2006，34(9):163－170.

［3］ 但尚铭，但 玻.若尔盖地区典型湿地NDVI动态特点分析［J］.四川气象，2006，26(3):115 －117.

［4］ Busetto L，Meroni M，Colombo R.Combining Medium and Coarse Spatial Resolution Satellite Data to Improve the Estimation of Sub －pixel NDVI Time Series［J］.Remote Sensing of Environment，2008，112(1):118 －131.

［5］ 阳小琼，朱文泉，潘耀忠，等.基于修正的亚像元模型的植被覆盖度估算［J］.应用生态学报，2008，19(8):1860 －1864.

［6］ 陈效逑，王 恒.1982～2003年内蒙古植被带和植被覆盖度的时空变化［J］.地理学报，2009，64(1):84 －94.

［7］ 董 印，焦 黎，杨光华，等.基于SPOT－VGT数据的新疆1998～2007年植被覆盖变化监测［J］.水土保持通报，2009，29(2):125－129.

［8］ 李 娟.兰州市南北两山植被覆盖度动态变化遥感监测［J］.矿山测量，2010(1):34－38.

［9］ 田庆久，闵祥军.植被指数研究进展［J］.地球科学进展，1998，13(4):327－333.

［10］ Iehoku C，Karnieli A.A Review of Mixture Modeling Techniques for Sub － pixel Land Cover Estimation［J］.Remote Sensing Reviews，1996，13(3 －4):161 －186.

［11］ Gutman G，Ignatov A.The Derivation of the Green Vegetation Fraction from NOAA/AVHRR Data for Use in Numerical Weather Prediction Models［J］.International Journal of Remote Sensing，1998，19(8):1533 －1543.

［12］ Qi J，Marsett R C，Moran M S，et al.Spatial and Temporal Dynamics of Vegetation in the San Pedro River Basin Area［J］.Agricultural and Forest Meteororology，2000，105(1 －3):55 －68.

［13］陈 晋，陈云浩，何春阳，等.基于土地覆盖分类的植被覆盖率估算亚像元模型与应用［J］.遥感学报，2001，5(6):416 －422.

［14］唐世浩，朱启疆，王锦地，等.三波段梯度差植被指数的理论基础及其应用［J］.中国科学(D 辑)，2003，33(11):1094－1102.

［15］阳小琼，朱文泉，潘耀忠，等.基于修正的亚像元模型的植被覆盖度估算［J］.应用生态学报，2008，19(8):1860 －1864.

［16］ Jasinski M F.Estimation of Subpixel Vegetation Density of Natural Regions Using Satellite Multispectral Imagery［J］.IEEE Transaction on Genscience and Remote Sensing，1996，34(3):804 －813.

［17］ Gilabert M A，García － Haro F J，Melià J.A Mixture Modeling Approach to Estimate Vegetation Parameters for Heterogeneous Canopies in Remote Sensing［J］.Remote Sens Environ，2000，72(3):328－345.

［18］ Paringit E C，Nadaoka K.Using Spectral Mixture Modeling Techniques to Derive Land－cover Parameters for Distributed Ediment Yield Estimation［C］//Proceedings of the 21st Asian Conference on Remote Sensing.Taipei，Taiwan，2000:52 －57.

［19］王天星，陈松林，马 娅.基于改进线性光谱分离模型的植被覆盖度反演［J］.地球信息科学，2008，10(1):114 －120.

［20］ Metternicht G I，Fermont A.Estimating Erosion Surface Features by Linear Mixture Modeling［J］.Remote Sensing of Environment，1998，64(3):254 －265.

［21］ Lobell D B，Asner G P.Cropland Distributions from Temporal Unmixing of MODIS Data［J］.Remote Sensing of Environment，2004，93(3):412 －422.

［22］ Bannari A，Pacheco A，Staenz K，et al.Estimating and Mapping Crop Residues Cover on Agricultural Lands Using Hyperspectral and IKONOS Data［J］.Remote Sensing of Environment，2006，104(4):447－459.

［23］赵英时，杨立明，陈冬梅，等.遥感应用分析原理与方法［M］.北京:科学出版社，2003.

［24］蔡 薇，郭洪海，隋学艳，等.基于混合像元分解提取小麦种植面积的技术与方法研究［J］.山东农业科学，2009(5):12－16.

［25］陶雪涛.线性模型下多通道遥感图像混合像元分解方法研究［D］.上海:复旦大学，2008:10.

［26］邓孺孺.青藏高原地表反照率反演及冷热源分析［D］.北京:中国科学院研究生院遥感应用研究所，2002:16－19.

［27］邓孺孺，田国良，柳钦火，等.粗糙地表土壤含水量遥感模型研究［J］.遥感学报，2004，8(1):75 －80.

［28］齐志新，邓孺孺.多暗像元大气校正方法［J］.国土资源遥感，2007(2):16－19.

Remote Sensing Monitoring of Vegetation Coverage in Southern China Based on Pixel Unmixing:A Case Study of Guangzhou City

ZHANG Zhi－ xin1，DENG Ru － ru1，LI Hao1，CHEN Lei1，2，CHEN Qi－ dong1，HE Ying － qing1
(1.School of Geographic Science and Planning，Sun Yat－ sen University，Guangzhou 510275，China;2.South China Sea Marine Engineering and Environment Institute，SOA，Guangzhou 510300，China)

Based on the measurement of the ground spectral reflectance of basal land covers and the accurate atmospheric correction for Landsat TM data，the authors improved the linear spectral mixture model(LSU)and developed a vegetation coverage retrieval model suitable for southern China.The effects of the atmospheric environment and the imaging time of remote sensing data were both reduced，contributing to the multi－ temporal comparison，by the utilization of the ground spectral reflectance from field survey.The soil moisture factor was considered to eliminate its remarkable spatial differentiation error in southern China.The vegetation coverage retrieval model was proved to be efficient with high precision over the in situ field verification and was applied to extract the vegetation coverage information in Guangzhou from 1998 to 2009.It is inferred that the urbanization，the large－scale architectural engineering and the reclamation activities constitute the main factors responsible for the formation of the spatio－temporal vegetation change in this area.

Vegetation coverage;Linear spectral unmixing(LSU);Ground spectral reflectance;Soil moisture;Landsat TM image

TP 79

A

1001－070X(2011)03－0088－07

2010－12－13;

2011－01－09

国家自然科学基金项目(编号:40671144)和水利部948项目(编号:200820)共同资助。

张志新(1985－)，女，硕士研究生，主要从事水质遥感和微波遥感研究。

邓孺孺(1963－)，男，博士，教授，主要从事环境遥感和 GIS应用研究。E－mail:eesdrr@mail.sysu.edu.cn。

(责任编辑:刘心季)