随着“一带一路”倡议的提出及丝绸之路经济带的建设, 作为其重要组成部分-天山北坡的生态问题日益得到关注。天山北坡深居内陆, 干旱的大陆性气候显著, 具有特色鲜明、较为脆弱的山地-绿洲-荒漠生态系统。作为陆地生态系统的重要组成部成分, 植被是连接土壤、大气、水分的自然纽带^[1-2], 在全球能量传输、生物地球化学和水文循环等方面起着全球变化“指示器”的重要作用^[3-4]。对干旱区而言, 植被防风固沙, 涵养水源, 改善区域小气候, 防治大气污染等, 是区域重要的生态屏障, 对区域的生态安全和可持续发展起着不可替代的作用。

植被覆盖度是指单位面积内植被冠层垂直投影所占百分比^[5-6], 是量化区域植被长势的综合指标, 其动态变化反映了地球上的植被对全球变化的响应程度。对植被覆盖度动态变化进行监测不仅是评价区域生态环境质量变化的有效途径, 也是研究气候变化和地理物质循环的基础^[7-8]。20世纪以来, 随着城市化进程的不断推进, 人类活动对植被的影响日益显著^[9], 工农业生产从自然界获取大量的资源, 导致全球森林与草原植被的不断减少, 由此所引发的生态危机日益加重。因此, 对于植被覆盖度的研究逐渐变得广泛而深入^[10]。

研究植被覆盖度主要有两种方法：一是地面调查法, 二是遥感监测法^[11]。相对于地面调查法, 遥感监测法具有实时、经济、便捷等优点, 正在得到越来越广泛的使用^[12]。其中像元二分法是遥感监测法中应用范围比较广、相对比较成熟的植被覆盖度提取方法^[13]。它是利用归一化植被指数(NDVI)对植被覆盖度进行提取的一种方法。前人基于NDVI监测植被覆盖变化已取得大量成果：Tucker等^[14]在对非洲大陆土地覆盖分类的研究中, 首次对NDVI数据进行了主成分分析；袁丽华等使用MODIS NDVI产品对2000—2010年黄河流域植被覆盖度进行监测, 并分析其时空变化特征^[15]；许青云等利用MODIS NDVI产品, 对陕西省农作物进行识别研究^[16]；何宝忠等利用MODIS NDVI数据, 对新疆2005—2015年植被覆盖度进行分析^[11]。

基于此, 本文利用RS、GIS技术, 结合NDVI、数字高程模型(DEM)数据, 采用像元二分法、线性回归分析、稳定性分析等对天山北坡的植被覆盖度进行研究, 以期为区域在一带一路建设中的生态安全格局构建及可持续发展提供理论支撑。

1 研究区概况

本文选取的天山北坡研究区(图 1)位于准噶尔盆地南缘, 天山北麓中段, 北接古尔班通古特沙漠, 南部为天山山区。区域总体地势南高北低。区域属温带荒漠干旱区, 年均温约7—9℃, 年降水量约160—300 mm, 气温和降水的日较差和年较差均较大^[17]。包括乌鲁木齐市、克拉玛依市、石河子市、昌吉市、阜康市、奎屯市、乌苏市、五家渠市、呼图壁县、玛纳斯县和沙湾县11个县(市)。天山北坡区位优势明显, 能源丰富, 发展潜力巨大, 但也存在水资源短缺, 资源环境承载力低等限制性因素^[18], 随着国家“一带一路”建设的推进, 天山北坡在兴疆固边、辐射带动沿线地区中发挥着不可替代的重要作用^[19]。

2 数据源与方法 2.1 数据源

NDVI数据来源于中国科学院计算机网络信息中心地理空间数据云平台(http://www.gscloud.cn)的MODND1M数据产品, 空间分辨率为500 m。为了最大可能减少冰雪因素的干扰, 结合研究区物候特征, 选取植被覆盖比较高的月份(6—9月)作为研究时段, 并对这4个月份的NDVI数据进行合成平均。

地理数据主要包括DEM数据。DEM数据来自于资源环境数据云平台, 空间分辨率90 m, 为与NDVI数据精度一致, 将其重采样为500 m提取研究区的高程信息, 以200 m为间隔划分, 共划分为20级, 并将海拔4000 m以上的点合并为一类进行研究。

2.2 植被覆盖度计算

NDVI是反映地表植被生长状态的重要指标因子, 定义为反射率之差与反射率之和的比值：

(1)

式中, NIR是传感器的近红外波段, R是红光波段。

利用植被覆盖度和NDVI之间极显著的线性相关关系, 通过像元二分法模型对植被覆盖度进行提取。

(2)

NDVI_veg表示完全由绿色植被贡献的信息, 而NDVI_soil则表示完全由裸地提供的信息。

并在像元二分模型的基础上进行植被覆盖度提取：

(3)

式中, VFC代表植被覆盖度。截取置信度99.5%和0.5%分别作为NDVI最大值和最小值, 并参考中华人民共和国水利部2008年颁布的《土壤侵蚀分类分级标准》, 结合研究区实际情况及相关研究^[20], 将植被覆盖度划分为5个等级：< 30%(低覆盖度)、30%—45%(中低覆盖度)、45%—60%(中覆盖度)、60%—75%(中高覆盖度)、>75%(高覆盖度)。

2.3 植被覆盖度年际变化趋势

采用一元线性回归方法分析植被覆盖度的年际变化。计算公式如下^[7]：

(4)

式中:θ表示回归方程的斜率, i表示年份；C_i表示第i年的植被覆盖度。趋势的显著性用F检验, 其计算公式为：

(5)

(6)

(7)

式中, U和Q分别指误差平方和及回归平方和, y_i是第i年的植被覆盖度像元值, 为其回归值, y为2001—2015年植被覆盖度平均值, n为研究年数。根据李卓等^[7]的研究, 将其变化趋势共分为6个等级：极显著退化(θ < 0, P < 0.01)；显著退化(θ < 0, 0.01 < P < 0.05)；不显著退化(θ < 0, P>0.05)；不显著改善(θ>0, P>0.05)；显著改善(θ>0, 0.01 < P < 0.05)；极显著改善(θ>0, P < 0.01)。

2.4 植被覆盖度稳定性

引入变异系数对植被覆盖度稳定性进行研究。变异系数的大小反映了数据系列的波动程度^[21]。其计算公式为：

(8)

σ为标准差, 计算公式为

(9)

C为研究时段内多年平均植被覆盖度, 按照变异系数的划分等级, 将变异系数(C)划分为弱变异(C < 0.1), 中等变异(0.1≤0.1 < 1)以及强变异(C≥1)。

2.5 植被覆盖度景观等级变化图谱信息

为了更深入研究2001—2015年研究区植被覆盖度的时间变化趋势, 以时间序列为轴, 对植被覆盖度等级空间分布图进行空间叠加计算, 生成天山北坡植被覆盖度等级空间-属性-过程图谱^[22]。将变化类型分为6类：(1)前期变化型：只在2001—2005年变化；(2)后期变化型：只在2010—2015年变化；(3)中间过渡型：2005—2010年仅发生一次变化；(4)反复变化型：在2001—2015年至少发生两次变化且2001和2015年等级相同；(5)持续变化型：在2001—2015年持续变化；(6)稳定型：在2001—2015年没有变化。图谱的计算公式为：

(10)

式中：T为植被覆盖度等级信息变化图谱编码的时空复合体数据；n为参与计算的植被覆盖度等级空间分布栅格数据集期数, A₁, A₂, …为不同年份植被覆盖度等级栅格数据。

2.6 植被覆盖度随高程的变化

在本研究中, 引入DEM研究每个格点15a的平均植被覆盖度随海拔变化, 不同植被覆盖度在不同海拔面积占比是指各海拔中的同一植被覆盖度格点数与相同海拔内所有格点数的比值。

2.7 植被覆盖度面积加权重心转移分析

重心是地理学中描述地理要素或空间对象转变的重要空间指标。重心的动态迁移反映了地理要素空间分布的迁移轨迹, 加权重心则是通过对重心坐标赋予权重指示地理现象分布的不均匀性, 本文选用面积加权重心模型计算不同等级植被覆盖度的重心, 以反映其空间差异。可用以下公式计算^[23]：

(11)

(12)

式中：X和Y分别代表该等级植被覆盖度面积加权重心的经纬度坐标。n为图斑数量, C_t代表第t个图斑的面积, X_i和Y_i分别表示第i个图斑的重心坐标。

3 结果与分析 3.1 植被覆盖度时空变化特征

对研究区2001—2015年6—9月的平均植被覆盖度进行统计分析(图 2)。2001—2015年, 天山北坡平均植被覆盖度在0.4和0.5之间变化, 但存在年际差异：2001年最低, 约为0.4；2011年达到峰值, 约为0.49；此后基本稳定在0.45—0.49。大致可分为三个阶段：2001—2011年, 波动上升阶段；2011—2012年, 平稳波动阶段；2012—2015年, 第二次波动阶段。总体而言, 15年间平均植被覆盖度整体呈上升趋势。

对研究区不同等级的植被覆盖度进行统计分析, 结果发现：2000—2015年间, 研究区的低植被覆盖度面积大致呈减少趋势, 中低、中和中高植被覆盖度面积呈波浪式变化, 高植被覆盖度面积大致呈增长态势(图 3)。其中：(1)低植被覆盖度所占比例在32%—47%；2001年的低植被覆盖度占比最大, 为47%；此后低植被覆盖度所占比例稳中下降, 在2010年跌至最低, 约为32%；之后大致在35%—40%之间波动。(2)中低、中和中高植被覆盖度所占比例基本稳定在33%—46%, 2008年占比最低, 约为33%；2010年最高, 约为46%。(3)低植被覆盖度和中低植被覆盖度所占比例基本在49%—63%, 2001年最高, 此后波动下降, 2010年达到最低, 约为49%。区域以低和中低植被覆盖度为主。(4)研究区高植被覆盖度所占比重在13%—26%之间变化, 最低值、最高值分别出现在2001和2013年。总体看来, 高植被覆盖度所占比例有逐渐增加的趋势。

对2001—2015年间不同等级的植被覆盖度空间分布进行分析(图 4)。发现不同年份的植被覆盖度分布表现出一定的规律性：高值区基本呈西北-东南走向两个条带状贯穿区域东西, 南缘的高值区依天山北麓分布, 中北部的高值条带区贯穿阜康市南部、乌鲁木齐市中部、昌吉市中部、呼图壁县中部、玛纳斯县中北部、沙湾县中北部、乌苏市东北部以及奎屯市、石河子市、五家渠市大部分区域。这可能与山前地区、盆地边缘的降水相对丰沛, 内陆河流在山区汇集流向盆地为植被生长提供了丰富的水源有关。低值区与高值区相间分布, 区域的北部、东北部以低植被覆盖度为主, 可能和与古尔班通古特沙漠临近, 水源稀少, 蒸发强烈, 且内陆河流大多流程短、流量小, 加上沿线区域工农业用水, 植物生长所需水源不足不利于植被生长有关。总体而言植被覆盖逐渐转好, 植被覆盖度逐渐增高, 高值区逐渐扩大, 而低值区逐渐缩小。

3.2 植被覆盖度变化趋势及变化类型分析 3.2.1 植被覆盖度变化趋势分析

随着“退耕还林”、“退耕还草”等政策和生态工程的实施, 天山北坡植被覆盖度呈现改善趋势(图 5, 表 1)。具体表现为：(1)植被覆盖度改善区域所占比重(54.44%)大于退化区域(45.58%)；(2)植被覆盖度不显著变化区域约占区域总面积的70%, 多分布在区域的东部、南部以及东北部地区, 而显著变化区域占比虽小, 但以改善为主；(3)植被显著增加和极显著增加区域面积远大于显著退化和极显著退化区域面积, 区域整体植被覆盖度呈良性发展态势；(4)植被覆盖度增加区域主要在区域中部及中北部偏西地区, 大致呈两个东西向贯通的横轴分布在区域农耕区北部及以北地区；而植被覆盖度退化区域则呈环状分布于植被覆盖度改善区四周, 极显著退化区域呈零星状分布。这在一定程度上说明了人类活动对区域植被覆盖度的改善产生了积极影响。

3.2.2 植被覆盖度稳定性分析

经异常值剔除及统计分析(图 6), 发现：(1)区域99%以上的植被覆盖度稳定性介于弱变异和中等变异, 区域植被覆盖度相对稳定；弱变异区域主要分布在研究区的中部和南部, 中等变异区域主要分布在研究区中部及偏北地区；(2)强变异区域占比极小, 呈零星状主要散布在区域南部高海拔地区, 且远离城市建成区, 说明植被覆盖度强变异的主要原因是自然因素。

3.2.3 植被覆盖度等级变化图谱

植被覆盖度等级变化图谱显示, 2001—2015年间天山北坡不同等级的植被覆盖度空间变化表现出明显的差异(图 7)。整体看来, 植被覆盖度稳定型居多, 占比达50%, 且分布范围最广, 呈不规则环状分布在研究区四周；其次是反复变化型, 约占16.5%, 呈西北-东南走向、断断续续的轴状分布, 比稳定型更靠近研究区中心；接着依次为持续变化型、后期变化型、中间过渡型和前期变化型, 这四种类型虽然随着比例降低渐趋分散, 但基本呈不连续的西北-东南向带状分布于中部区域。

对比图 4、图 7, 发现低植被覆盖度和植被变化稳定型分布较为一致, 可能是低植被覆盖度主要分布于荒漠戈壁及高山高原等人类活动较少干扰地区, 因此植被覆盖度变化不大, 以稳定型为主；而越接近城市建成区, 植被覆盖度类型受到的人类活动干扰越大, 变化类型也渐趋复杂。

3.3 植被覆盖度随地形变化趋势

研究区海拔在224—5166 m之间(图 8), 80%以上的区域海拔在224—2500 m之间。北部、中部地势起伏较为和缓, 而区域南缘的天山北麓地区地势起伏较大。

不同海拔的植被覆盖度存在较大差异(图 9)：(1)海拔224—1600 m之间, 低植被覆盖度和中低植被覆盖度占绝对优势, 两者占比达50%以上；海拔800—1200 m之间低、中低植被覆盖度达90%以上；此后随着海拔升高两者占比逐渐下降, 到海拔2200—2400 m之间降至最低值不足50%；到海拔3880 m以上, 低植被覆盖度和中低植被覆盖度占绝对优势, 较高的植被覆盖度基本绝迹, 可能由于高海拔区域光热条件和水分条件较差, 土壤有机质含量低不利于植被生长；总的来看, 低植被覆盖度和中低植被覆盖度随海拔升高呈“~”型变化趋势；(2)中、中高和高植被覆盖度所占比重随海拔升高先下降后升高最后再下降, 高植被覆盖度集中于海拔1800—2600 m之间, 占比基本在50%以上；在海拔800—3800 m之间最为突出, 高达57%以上；(3)海拔2800—5166m之间, 随着海拔升高, 低植被覆盖度占比迅速上升, 几乎从0增加到100%；而高植被覆盖度所占比例从50%骤减到不足1%, 充分说明了由于海拔不同造成的光热条件悬殊对于植物生长具有极强的限制性。

3.4 不同等级植被覆盖度面积加权重心分布

总体而言, 不同时间段各等级植被覆盖度的面积加权重心位置分布既有相同之处, 也呈现出一定的差异(图 10)。除了低植被覆盖度面积加权重心位于区域北部偏东地区, 与低植被覆盖度分布基本一致, 其他各等级植被覆盖度面积加权重心大致位于区域中部地区, 且基本在沙湾县、石河子市、玛纳斯县及呼图壁县分布。低植被覆盖度面积加权重心和其他等级植被覆盖度面积加权重心相隔较远, 其他四个等级植被覆盖度面积加权重心随时间呈现由集中到分散的趋势。2001—2015年间, 低植被覆盖度和高植被覆盖度面积加权重心向西北有一定程度的移动, 而中低、中、中高植被覆盖度面积加权重心则大致向东南移动。

4 讨论

作为衡量地区生态的重要因子, NDVI从正面反映了植被覆盖度的变化, 但不同研究采取的置信度水平存在差异^{[7, 24]}, 由此可能造成结果的差异, 但基本都是利用NDVI和植被覆盖度的正相关关系进行后者的提取。在过去几十年间, 已有大量研究运用不同数据源(MODIS、Landsat、SPOT等)从不同时间尺度和空间尺度对植被覆盖度动态变化进行了探讨^[25-26]。已有研究表明, 我国总体植被覆盖有逐渐绿化的趋势^[27]。本研究得出的天山北坡植被覆盖度有逐渐增长的趋势的结论, 与何宝忠^[11]、位宏等^[20]得出的新疆、玛纳斯河流域植被覆盖度变化趋势基本一致。对生态环境极为脆弱的西北干旱区而言, 监测其生态环境变化尤为重要^[28]。随着“三北防护林”工程的深入推进, 各项退耕还林还草、退牧还草等政策的积极实施, 天山北坡的植被覆盖度可能进一步改善。

天山北坡地区深居内陆, 降水稀少, 蒸发旺盛, 植被稀疏, 植被覆盖度以低、中低为主。受到地形因子、河流流程流向、绿洲分布及人类活动的影响, 不同等级植被覆盖度大致呈西北-东南向的条带状分布。区域人口集中于绿洲, 区域的高山、荒漠地带鲜有人类活动, 植被覆盖度等级相对比较稳定, 变化图谱的稳定型多与低植被覆盖度区域重合。

对山区而言, 温度随着海拔升高逐渐下降, 而降水量随海拔升高先增多后减少。因此, 不同的海拔高度产生了不同的水热组合。一般而言, 最有利于植被生长的海拔组合往往出现在海拔中段^[29], 研究区域植被覆盖度比较高的区域集中在1800—2600 m左右。此外, 研究采用的NDVI及植被覆盖度数据都是500 m空间分辨率, 中低空间分辨率对于植被覆盖度的反映比较宏观, 缺少微观尺度的考量。因此, 对于时间和空间尺度的进一步把握是今后工作中所必须要考虑的重要问题。未来天山北坡的发展建设要更加注意中低海拔地区的生态修复, 为区域社会发展奠定坚实生态基础。

5 结论

结合多种数据源, 利用RS/GIS方法, 探究了天山北坡15年间的植被生长变化趋势, 得到如下结论：

(1) 天山北坡多年平均植被覆盖度介于0.4和0.5, 以低、中低植被覆盖度为主, 低植被覆盖度面积大致呈减少趋势, 中低、中和中高植被覆盖度面积呈波浪式变化, 高植被覆盖度面积大致呈增长态势。不同等级植被覆盖度大致呈西北-东南向相间分布, 植被覆盖度有逐渐上升的趋势, 西部较东部改善更为明显。

(2) 区域植被覆盖度稳定性多介于弱变异和中等变异之间, 植被覆盖度变化类型也以稳定型为主, 区域植被覆盖相对稳定, 波动较小。

(3) 天山北坡植被覆盖度变化受海拔高度影响明显：随海拔上升, 较低植被覆盖度比例呈现先上升后下降再次上升趋势, 中、中高和高植被覆盖度与之相反；中高、高植被覆盖度在海拔800—3800 m之间可达57%以上；海拔3880 m以上, 低植被覆盖度占绝对优势。

(4) 各等级植被覆盖度面积加权重心集中在沙湾县、石河子市、玛纳斯县及呼图壁县3县1市内分布；低植被覆盖度和高植被覆盖度面积加权重心大致向西北移动, 而中低、中、中高植被覆盖度面积加权重心则大致向东南移动, 呈现出由集中到相对分散的趋势。