岷江上游小流域NDVI与环境要素的多尺度耦合关系研究

熊亚兰，张科利，朱雅晗，夏 丽

（1.四川农业大学资源学院，成都 611130；2.北京师范大学地理科学学部地理学院，北京 100875）

植被是陆地生态系统的重要组成部分，对全球的气候和环境变化具有重要的指示作用，此外植被能有效地减少水土流失、山体滑坡和泥石流等。归一化植被指数（normalized difference vegetation index，NDVI）与生物量、叶面积指数、植被覆盖度和土地利用等密切相关，是目前最为常用的表征植被状况的指标[1]。在山区，由于地形控制了太阳辐射和降水的空间再分配过程，因此地形因子是植被格局的重要塑造力量[2]。

对植被覆盖度与地形因子关系的研究随着小波分析技术的应用已从静态、单一尺度逐渐发展到动态、多尺度[3-5]。基于一维小波变换的研究中样带的选择是关键，而二维小波变换能克服样带代表性等限制，较一维小波变换能更全面深入地反映植被和环境要素空间分布格局的多尺度关系。二维小波变换可将栅格数据分解为低频和高频部分。低频系数反映实测值经小波变换滤去细节后所剩余的本征信息，能反映全局变化态势，分解尺度越大，其变化规律也趋于清晰和简单化[6-8]。

岷江上游茂县段是长江流域重要的生态屏障和水源涵养地。自20世纪以来，受人类活动和气候变化的影响，该区脆弱的生态环境受到破坏，生态环境破坏已严重地制约了当地人民生产生活水平的提高，削弱了岷江上游地区的生态屏障主体功能。本研究拟采用二维小波变换对该区典型小流域NDVI与环境要素的多尺度耦合关系进行研究，研究结果可为生态植被恢复、水土流失治理和农业产业规划提供基础资料。

1 材料和方法

1.1 研究区概况

渭门小流域（103°48′42″～103°54′52″E，31°44′52″～31°51′18″N）位于岷江上游干旱河谷茂县北侧的渭门乡、永和乡境内，流域面积74.04 km2（图1）。该区属高原型季风气候，垂直分带气候显著，多年平均降水量484.1 mm，蒸发量大；多年平均气温11.1℃，气温年较差小，日较差大。全区地势西高东低，海拔约为1 550～3 394 m，山顶形态类型多尖顶，沟谷深切，地形崎岖。土地利用以林地、草地为主，耕地面积小。全区自然植被类型为典型的干旱河谷灌丛植被，以矮灌木与半灌木占优势，几乎无大面积森林存在。

图1 研究区域及位置示意图Figure 1 Map showing the scope and location of the study area

1.2 数据来源与处理

NDVI数据基于2015年1月Landsat-8 OLI影像（http://www.gscloud.cn/）经几何校正和辐射校正后计算得到，空间分辨率为30 m。NDVI的计算公式为[9]：

式中：NIR为近红外波段反射率；R为红光波段反射率；在TM影像中分别为波段5和波段4。

数字高程模型（digital elevation model，DEM）来源于美国太空总署（NASA）和日本经济产业省（METI）联合推出的ASTER GDEM数据集，空间分辨率为30 m，坐标系为WGS84。利用ArcGIS软件从DEM中提取高程、坡度和坡向。将坡向按西北-东南方向重分类为 5 级，其范围分别为 292.5°～337.5°、337.5°～360°和 0°～22.5°以及 247.5°～292.5°、22.5°～67.5°和202.5°～247.5°、67.5°～112.5°和 157.5°～202.5°、112.5°～157.5°，并分别赋值为 1、2、3、4、5。NDVI与坡向相关系数的计算，以及回归方程的拟合均采用坡向的赋值。

居民点和河流图层根据Google卫星地图目视解译得到，在解译之前需对Google卫星地图进行几何校正和投影变换。利用ArcGIS的欧氏距离分析工具生成距最近居民点和距最近河流的栅格图，空间分辨率为30 m。

1.3 小波多尺度分析的基本原理

设任意一个采样点坐标为（x，y），点坐标（x，y）变化时各要素的变化就形成了一个离散变化的二维信号（函数）f（x，y）。小波变换就是通过基本小波（母小波，mother wavelet）与待分析函数内积，以达到分解原函数的目的。二维小波变换的公式如下[10-12]：

式中，WTf为小波系数，a为尺度因子，b1和b2分别为x和y方向的平移因子，ψ（x）为母小波，f（x，y）为待分析函数（栅格值）。原函数经过母小波与待分析函数作内积之后，可分解得到不同尺度下的小波系数。本研究中小波母函数选择在小波分析中常用的db6，算法采用Mallat塔式算法。根据塔式算法，实测值与分解后小波系数的关系为：

式中，j为分解的级数；A为低频系数；H、V和D分别为水平、垂直和对角方向的高频系数；G为多级分解尺度下3个方向高频系数之和。本研究主要对低频系数进行分析，低频系数的提取在MATLAB中完成。根据小流域范围和栅格大小可对NDVI和环境要素进行6级分解，分解后每一级所对应的栅格大小分别为60 m×60 m、120 m×120 m、240 m×240 m、480 m×480 m、960 m×960 m 和 1 920 m×1 920 m。通过二维小波重构可以确保不同分解尺度下小波系数的个数一致。

2 结果和分析

2.1 NDVI与环境要素的基本统计特征分析

由表 1 可得，NDVI的范围为-0.11～0.42，平均值为0.11；海拔变化较大，其范围为 1 550～3 394 m，均值为 2 442 m；坡度范围为 0°～70°，均值为 30°。距最近居民点和距最近河流的距离分别为0～7 078 m和0～7 445 m，均值分别为2 082 m和1 526 m。因受本地区地形地貌等条件的控制，居民点及其工程-经济活动主要分布于地理位置相对较低的河谷及支沟地带。

表1 NDVI及环境要素的基本统计特征Table 1 Descriptive statistical characteristics of NDVI and environmental factors

2.2 NDVI与环境要素的多尺度空间分布特征

由图2可知，渭门小流域植被覆盖偏低，结合该区的实际植被覆盖情况，将植被覆盖度划分为6个等级：≤0（非植被）、0～0.05、0.05～0.10、0.10～0.15、0.15～0.25 和＞0.25。小波变换对要素具有一定的平滑效应，每进行一级分解，细节信息逐渐被忽略，要素的空间分布规律趋于清晰和简单化。二维小波变换对各要素的多尺度平滑效应不同。NDVI和坡度、坡向在30～120 m尺度的平滑效应不明显，在240～1 920 m尺度平滑效应逐渐明显。高程、距最近河流距离和距最近居民点距离的多尺度平滑效应不明显，表明这3个要素在研究区本身具有一定的空间连续性。

图2 NDVI和环境要素的多尺度空间分布图Figure 2 Multi-scale spatial patterns of NDVI and environmental factors

综合7个尺度看，在5个环境要素中，坡向和NDVI的空间分布格局最为相似，其次是高程。从水系的分布特征来看，在靠近河流处，NDVI通常较低。坡度、距最近居民点距离的多尺度空间分布特征与NDVI差异较大。将NDVI的多尺度空间分布图与高清遥感影像以及实地调研进行对比分析后发现NDVI在 0～0.10 的区域主要为农用地，0.10～0.15 主要为灌丛林地和农用地交错的区域，＞0.15主要为灌木林地。总的来看，NDVI的空间分布特征呈现出东南低西北高、河谷低山地高的特征。

2.3 NDVI与环境要素的多尺度小波相关分析

多尺度相关分析与单一尺度的相关分析相比，能更好地厘清不同要素之间的耦合特性。由表2可知，NDVI与各环境要素均呈极显著相关，相关性大小受尺度变化的影响较明显。多尺度下坡向与NDVI的相关性最强，其次是高程。坡向与NDVI的相关系数为 0.405～0.543，并且随尺度的增加表现为先减小后增加。高程与 NDVI的相关系数为 0.226～0.362，随分解尺度的增加而增大。坡度与NDVI呈极显著负相关，相关性在30～480 m尺度呈逐渐增加的趋势，随后在960 m尺度略有降低，最后在1920 m尺度出现陡增，相关系数达-0.358。距最近河流距离、距最近居民点距离与NDVI的相关性较小，且随尺度的变化也较小，说明该区这两个因素对NDVI的影响较小。

表2 不同尺度下NDVI与环境要素的相关性Table 2 Correlation coefficients between NDVI and environmental factors at different scales

为了进一步厘清NDVI与环境要素间的定量关系，构建了多尺度下NDVI的回归方程（表3）。为了消除建模因子共线性的问题，在SPSS软件中采用逐步回归法对NDVI进行了拟合，结果在7个尺度的拟合方程中均没有环境要素被剔除，表明在研究尺度范围内NDVI的空间分布特征是由高程、坡度和坡向等5个环境要素共同决定。回归方程的决定系数随分解尺度的增大表现为先降低后增加，表明在480 m尺度以下，分解尺度越小各要素对NDVI的控制作用越大；而在480 m尺度以上，分解尺度越大，各要素对NDVI的控制作用越大，480 m是否是该小流域的特征尺度有待进一步研究。

表3 不同尺度下NDVI的拟合方程Table 3 Regression prediction model and R2for NDVI at different scales

3 讨论与结论

在对渭门小流域NDVI与环境要素的多尺度相关研究中发现，将坡向按西北-东南方向赋值所得的相关性优于按南-北或东-西方向赋值。NDVI在东南坡表现出明显的优势，这与张学玲等[13-15]在武功山、连江流域和赤水河流域的研究结果一致，分析其原因主要是渭门小流域位于横断山区岷江上游，其山脉走向大体垂直于东南季风，加之明显的焚风效应使得NDVI沿迎风坡至背风坡逐渐减少[16]。坡向与NDVI的相关性在小尺度上表现更为明显，除了上述风向的原因外，还与该地山顶形态类型多尖顶，谷深坡陡，地形崎岖有关。

高程与NDVI在各尺度下均呈极显著正相关，其原因是海拔较低处靠近河谷且地势平缓，农业垦殖率高，而随着海拔的升高，人类活动较少，有利于天然植被生长。坡度与NDVI在各尺度下均呈极显著负相关，其原因坡度越大越容易发生滑坡和泥石流。高程、坡度与NDVI不论是正相关还是负相关，在较小尺度上与NDVI的相关系数都比较小，强相关性均出现在较大尺度上，这与李双成等[17-20]在青藏高原、武夷山的研究结果一致，说明高程和坡度趋向于在宏观尺度上制约NDVI的空间分布。距最近河流距离、距最近居民点距离对NDVI的影响较小，说明地形是该区导致NDVI产生空间分异的主要原因。

岷江上游渭门小流域NDVI的空间分布特征为东南低西北高、河谷低山地高。多尺度下NDVI和坡向的空间分布相似度最高，其次是高程。NDVI与各环境要素均呈极显著相关，与地形因子的相关性大小受尺度变化的影响较明显。多尺度下坡向与NDVI的相关性最强，其次是高程。坡向与NDVI的最强相关性出现在较小尺度上，高程、坡度与NDVI的最强相关性出现在较大尺度上。在30～1 920 m尺度NDVI的空间分布特征是由高程、坡度和坡向等5个环境要素共同决定，但是各要素对NDVI的控制作用大小随尺度发生改变。

小波变换可以揭示出NDVI和环境要素的多尺度特征，但就小波变换方法来说还需考虑以下一些问题。首先，数学领域的研究已经表明，选择不同的小波基函数，会导致小波分解后结果发生变化，因此需进一步探索小波基的筛选方法；其次，目前所采用的小波分析方法仅能给出断续尺度（即尺度必须为2n，n为整数）上的小波系数，无法在连续尺度上开展计算和分析，这也是在今后的研究中需解决的问题之一。