1982—2015年永定河山區(qū)植被變化及對天然徑流的影響

杜 勇，李建柱，牛凱杰，馮 平，胡慶芳，鄭彥辰

（1. 天津大學(xué) 水利工程仿真與安全國家重點實驗室，天津 300072；2. 永定河流域投資有限公司，北京 100193；3. 南京水利科學(xué)研究院 水文水資源與水利工程科學(xué)國家重點實驗室，江蘇 南京 210029）

1 研究背景

植被作為陸地生態(tài)系統(tǒng)的重要組成部分，具有連接大氣、土壤、水文等生態(tài)要素的功能，在自然生態(tài)系統(tǒng)和人類生產(chǎn)生活中發(fā)揮著不可替代的作用［1-3］。歸一化植被指數(shù)（Normalized Difference Vegeta?tion Index，NDVI）是公認(rèn)的反映大尺度地表植被覆蓋和生長狀況的有效指標(biāo)，與光合有效輻射、綠葉密度、植被生產(chǎn)力以及累積生物量等均具有密切的關(guān)系，通常用來表征植被活動的強弱［4-6］。近年來，由于氣候變化和人類活動的影響，區(qū)域乃至全球尺度的植被均呈現(xiàn)顯著增加趨勢［7-10］，植被變化不僅能影響氣候，而且可以通過改變植被蒸騰、冠層截流和土壤蒸發(fā)等水文過程而改變河川徑流，故能顯著影響流域水文循環(huán)過程［11-14］。因此，分析植被的時空演替特征，揭示在不同時空尺度下導(dǎo)致植被變化的驅(qū)動因子的空間異質(zhì)規(guī)律，并解析以植被變化為重要特征的下墊面變化對徑流變化的貢獻及影響，可為應(yīng)對氣候變化和人類活動下生態(tài)系統(tǒng)的適應(yīng)能力提供理論依據(jù)，同時有助于深刻理解植被活動對水循環(huán)過程的影響，并對完善流域水資源管理措施具有重要理論和實踐意義。

作為植物生長發(fā)育所必需的外界環(huán)境因子，氣候要素（水分和溫度等）對植被活動起重要作用，同時人類活動對植被活動也具有不同程度的促進或抑制作用［15-18］。以黃土高原為例，已有研究表明過去30年氣候變化對黃土高原的植被恢復(fù)有一定促進，同時植被恢復(fù)工程等人類活動也有利于NDVI的增加，因此氣候變化和人類活動共同驅(qū)動黃土高原的植被變化，但兩者的驅(qū)動作用呈現(xiàn)出顯著的空間異質(zhì)性［2，19-20］。同時大量研究證實，植被覆蓋作為下墊面的重要性質(zhì)，其變化將導(dǎo)致水循環(huán)過程變化，影響河川徑流過程［21-23］。目前有關(guān)植被變化等研究已在黃土高原等生態(tài)脆弱區(qū)和人類活動密集區(qū)全面開展，但這些研究多采用相關(guān)性分析等方法對植被變化的動因進行定性研究，關(guān)于不同時間尺度下各種驅(qū)動因素影響的空間異質(zhì)性的定量解析還不夠深入［24-25］，且當(dāng)前關(guān)于NDVI 變化對天然徑流量影響的文獻較少，不能清晰地從源頭理解植被變化對河川徑流過程的影響。

本文以黃土高原上的典型流域永定河山區(qū)為代表，基于氣溫、降水、徑流以及長時序NDVI 數(shù)據(jù)集，通過構(gòu)建地理加權(quán)回歸模型，研究流域內(nèi)NDVI 和降水、平均氣溫的回歸參數(shù)、殘差均值及其趨勢的空間異質(zhì)性，并進一步探討以植被變化為突出特征的下墊面變化對流域典型斷面天然徑流的影響。

2 研究方法

2.1 NDVI 與影響因素的空間趨勢分析方法NDVI 及其影響因素的空間趨勢分析可采用基于柵格尺度的普通最小二乘法（OLS），通過模擬每個柵格NDVI 及其影響因素在相應(yīng)時段的變化趨勢，從而反映研究區(qū)在不同研究時段內(nèi)的NDVI 及其影響因素空間特征?？臻g趨勢分析的計算公式為：

式中：θslope為線性回歸的斜率，表示不同時間尺度下NDVI、回歸參數(shù)、殘差等研究對象在流域或像元內(nèi)的年際變化趨勢；n 為研究時段年份數(shù)；Yi為第i年研究對象數(shù)值。

2.2 NDVI 與影響因素的空間異質(zhì)性表征方法地理加權(quán)回歸（Geographically Weighted Regression，GWR）是由Brunsdon 等［26］提出的一種空間分析方法，有助于揭示一定的研究區(qū)域內(nèi)自變量與因變量之間的非平穩(wěn)空間相關(guān)性。GWR 模型是對普通線性回歸的重要拓展，該方法的參數(shù)是空間位置的函數(shù)，該模型的表達形式為：

式中： yi、 xik、 εi分別為空間上i 點的因變量、自變量和模型殘差；(ui，vi)為i 點的空間位置；m為樣本量i 點處自變量的個數(shù)，本文取2，分別表示降水量和氣溫； βk為i 點上的回歸系數(shù)； β0為截距。在本文中，年均NDVI 為因變量，降水、年平均氣溫為相應(yīng)自變量，用以分析NDVI 對氣候要素的空間響應(yīng)關(guān)系。

GWR 模型采用局部加權(quán)最小二乘法求解，其中的權(quán)重是待估點到各觀測點的空間距離的函數(shù)，具有距離衰減效應(yīng)。參數(shù)公式為：

式中：β( ui，vi)為回歸系數(shù)的無偏估計；W (ui，vi)為權(quán)重矩陣，一般采用高斯模型作為權(quán)重函數(shù)，具體參見文獻［27］；X 和Y 分別為自變量和因變量的矩陣。

本文通過建立GWR 模型得到年降水量和年平均氣溫對應(yīng)的回歸參數(shù)，并計算得到NDVI 的模擬值（NDVIcc），用來表征氣候因素對NDVI 的影響；然后計算NDVI 實際觀察值（NDVIobs）與模擬值（NDVIcc）之間的差值，表示人類活動對NDVI 的影響：

2.3 植被變化對徑流影響的研究方法流域水量平衡方程為：

式中：R 為徑流深，mm；P 為降水量，mm；ET 為實際蒸散發(fā)量，mm； ΔS 為流域儲水量的變化，mm。當(dāng)研究時間尺度較大時， ΔS ≈0，本文為了減少永定河山區(qū)流域ΔS的影響，取各變量時間尺度為10年的滑動平均值。

在流域尺度上，降雨量可直接通過雨量站點觀測數(shù)據(jù)進行空間插值獲得，實際蒸散發(fā)可基于Bu?dyko 假設(shè)的水熱耦合平衡方程計算：

式中：ET0為年均潛在蒸散發(fā)量，mm；n 為反映流域下墊面特征的參數(shù)（簡稱下墊面參數(shù)），與地形、土壤和植被等一系列因素有關(guān)。對于流域下墊面參數(shù)n 的估計，諸多學(xué)者提出了不同的估算方法，這些方法均認(rèn)為植被是影響n 的重要因素。Li 等［27］針對全球典型大流域研究分析發(fā)現(xiàn)，參數(shù)n 與NDVI存在較好的擬合關(guān)系：

式中：a、b 為回歸參數(shù)，對特定流域為常數(shù)，a 的大小表征了n 對NDVI 變化的敏感性程度，b 則反映了NDVI 以外的其它因素對n 的影響。

3 研究區(qū)域和數(shù)據(jù)

3.1 研究區(qū)域永定河流域位于海河流域，地理位置介于東經(jīng)111.95°—116°22′，北緯38°90′—41°16′，流經(jīng)內(nèi)蒙古、山西、河北、北京、天津五個省級行政區(qū)，集水面積約4.7 萬km2，占海河水系的14.7%。2016年流域總?cè)丝诩s813 萬，城鎮(zhèn)化率為58.3%。永定河流域上游一級支流為桑干河和洋河，兩支合并后始稱為永定河。永定河流域?qū)侔霛駶櫚敫珊导撅L(fēng)氣候，多年平均降水量介于360 ~650 mm，6—8月降水量約占年均降水量的80%。干流官廳水庫以上流域?qū)儆谏絽^(qū)，面積為4.51 萬km2，占全流域總面積95.8%，是主要產(chǎn)流區(qū)。冊田水庫、響水堡水庫分別為永定河一級支流桑干河、洋河的主要水庫（圖1），集水面積分別為1.71 萬km2和1.45 萬km2。

圖1 永定河山區(qū)氣象水文站點和河道斷面分布

3.2 數(shù)據(jù)來源與處理

3.2.1 NDVI 數(shù)據(jù) 本文收集了1982—2015年的永定河山區(qū)逐月NDVI 數(shù)據(jù)，數(shù)據(jù)來源于GIMMS ND?VI3g 數(shù)據(jù)集（https：//ecocast.arc.nasa.gov/data/pub/gimms/），空間分辨率為8 km，根據(jù)研究區(qū)域空間范圍，提取了相應(yīng)網(wǎng)格單元逐月NDVI 數(shù)據(jù)，經(jīng)進一步處理后得到季尺度和年尺度的NDVI 數(shù)據(jù)。

3.2.2 氣象數(shù)據(jù) 本文收集了1982—2015年的MSWEP 逐日降水?dāng)?shù)據(jù)和中國區(qū)域地面氣象要素驅(qū)動數(shù)據(jù)集中的逐日平均氣溫數(shù)據(jù)。兩者通過網(wǎng)絡(luò)共享獲?。ňW(wǎng)址分別為：http：//www.gloh2o.org/和http：//www.tpdc.ac.cn/zh-hans/），空間分辨率均為0.1°×0.1°。同時，在流域內(nèi)可獲取國家氣象科學(xué)數(shù)據(jù)中心（http：//cdc.cma.gov.cn/home.do）提供的7 個國家基本氣象站點1982—2015年逐日常規(guī)地表氣象數(shù)據(jù)（含雨量和氣溫）。

對永定河山區(qū)0.1°×0.1°的逐日MSWEP 網(wǎng)格降水?dāng)?shù)據(jù)和中國區(qū)域地面氣象要素驅(qū)動數(shù)據(jù)集中的平均氣溫數(shù)據(jù)，采用IDW 方法，分別采用國家基本氣象站的降水與氣溫數(shù)據(jù)進行校正，之后對校正后的逐日降水、平均氣溫網(wǎng)格數(shù)據(jù)進行重采樣處理，得到與NDVI數(shù)據(jù)空間分辨率一致的逐日降水、平均氣溫格網(wǎng)數(shù)據(jù)，最后通過時間尺度上的聚合處理，得到月、季、年等尺度的降水、平均氣溫數(shù)據(jù)。

3.2.3 地表徑流數(shù)據(jù) 根據(jù)海河流域水文年鑒和各省市水資源公報，收集了冊田水庫、響水堡水庫和官廳水庫所在斷面1956—2016年的實測徑流量，采用分項還原法，對天然年徑流量計算有影響的水量進行還原，將被人類開發(fā)利用的水量（如農(nóng)業(yè)用水、工業(yè)用水、生活用水、環(huán)境用水耗損量、跨流域引調(diào)水量等）逐項還原到天然年徑流量中，得到永定河典型斷面的天然徑流量。

4 結(jié)果分析與討論

4.1 NDVI 的時空變化特點圖2 給出了1980年和2015年研究區(qū)域土地利用類型分布，從中可知永定河山區(qū)主要土地利用類型是耕地、林地、草地、水域、建設(shè)用地及未利用土地，其中又以耕地、林地和草地為主。2015年較1980年林地面積增加了2364.4 km2，相對增幅34.6%；耕地、草地面積分別減少1610.1 和1056.8 km2，相對降幅分別為7.7%和7.5%；建設(shè)用地增加了926.9 km2，相對增幅116.2%。

圖2 永定河山區(qū)1980年和2015年土地利用類型

表1 為永定河山區(qū)土地利用面積轉(zhuǎn)移矩陣，反映出1980—2015年土地利用變化的主要特征是耕地、林地和草地之間的相互轉(zhuǎn)移以及部分耕地、草地轉(zhuǎn)化為建設(shè)用地；其結(jié)果是，流域內(nèi)耕地占總面積的比重由1980年的47.3%下降至2015年的43.7%，林地和草地面積之和占總面積的比重由47.3%上升至50.2%，而建設(shè)用地占總面積的比重由1.8%上升至3.9%，不同類型土地利用變化量占流域總面積的比重較小且變化范圍介于0.4% ~ 5.3%之間，因此，總體而言，1980—2015年期間永定河山區(qū)土地利用格局變化不大。

表1 1980—2015年永定河山區(qū)土地利用變化轉(zhuǎn)移矩陣 （單位：km2）

圖3 給出了永定河山區(qū)1982—2015年年均NDVI 變化，從3 可知，年均NDVI 呈波動增加趨勢。在研究時段內(nèi)，永定河山區(qū)年均NDVI 的波動范圍是0.2722 ~ 0.3394，最小值和最大值分別出現(xiàn)在1985年和2013年?？傮w上，1982—2015年年均NDVI 線性趨勢率為0.0013/a（p<0.01），這說明表明永定河山區(qū)植被恢復(fù)較明顯。圖4 表明，1982—2015年永定河山區(qū)NDVI 的年際變化趨勢具有一定的空間異質(zhì)性；年均NDVI 呈增加和減少趨勢的區(qū)域面積分別占總面積的95.36%和4.64%。其中，年均NDVI 增長較快（slope≥2×10-3/a）的面積約占研究區(qū)域總面積的12.32%，主要分布在桑干河上游南部林地和永定河官廳水庫附近的濕地；年均NDVI 呈降低趨勢的區(qū)域分布比較零散，主要在城市擴張的地區(qū)，且降低速度率較小，多在-1×10-3/a 以內(nèi)。

圖3 1982—2015年永定河山區(qū)年均NDVI年際變化

圖4 1982—2015年永定河山區(qū)年均NDVI 趨勢空間分布

圖5 給出了永定河山區(qū)1982—2015年春、夏、秋季的NDVI 變化情況，從中可知春、夏、秋季平均NDVI 均呈波動增加趨勢。在研究時段內(nèi)，永定河山區(qū)春季平均NDVI 的變化范圍是0.2196 ~0.2746，最小值和最大值分別出現(xiàn)在1985年和2014年，其趨勢率為0.0013/a（p<0.01）；夏季平均NDVI的變化范圍是0.4317 ~ 0.5512，最小值和最大值分別出現(xiàn)在1985年和2013年，趨勢率為0.0018/a（p<0.01）；秋季平均NDVI 的變化范圍是0.2626 ~ 0.3669，最小值和最大值分別出現(xiàn)在1984年和2013年趨勢率為0.0020/a（p<0.01）。這表明在春、夏、秋季永定河山區(qū)植被均有明顯恢復(fù)。由圖6 可知，1982—2015年永定河山區(qū)NDVI 變化趨勢在春、夏、秋季均具有一定的空間異質(zhì)性。其中，春季平均NDVI呈增加和減少趨勢的區(qū)域面積分別占總面積的92.03%和7.97%， NDVI增長較快（slope≥2×10-3/a）的面積約占18.26%，主要分布在桑干河上游南部林地和永定河官廳水庫附近的濕地地區(qū)；夏季平均NDVI呈增加和減少趨勢的區(qū)域面積分別站總面積的90.58%和9.42%，NDVI增長較快（slope≥2×10-3/a）的面積所占比例高達46.40%，主要分布在桑干河中上游地區(qū)林地、洋河流域北部的林草地和官廳水庫附近的濕地地區(qū)；秋季平均NDVI 呈增加和減少趨勢的區(qū)域面積分別站總面積的98.41%和1.59%，NDVI 增長較快（slope≥2×10-3/a）的面積約占51.30%，主要分布情況與夏季類似。

圖5 1982—2015年永定河山區(qū)春、夏、秋季平均NDVI年際變化率年際變化

圖6 1982—2015年永定河山區(qū)春、夏、秋季平均NDVI 變化率的空間分布

因此，在土地利用格局變化不大的情況下，永定河山區(qū)NDVI 在年度和春、夏、秋季均有較大幅度的增加，可能由于受降水、氣溫的氣象因素或流域內(nèi)植被恢復(fù)工程、農(nóng)林草地管理水平的提高等人類活動的影響，使年平均和各季度的NDVI 均有較大提高。

4.2 NDVI 與相關(guān)驅(qū)動因子的空間非平穩(wěn)性

4.2.1 驅(qū)動因子回歸參數(shù)與殘差的均值 為驗證植被NDVI 與降水、氣溫的空間相關(guān)關(guān)系是否具有非平穩(wěn)性，首先對回歸參數(shù)均值進行空間自相關(guān)分析，得到各要素GWR 回歸參數(shù)的Moran 指數(shù)與Z 值如表2 所示。表2 結(jié)果表明，NDVI 與降水、平均氣溫的GWR 模型回歸參數(shù)均值的Moran 指數(shù)均大于0，說明參數(shù)具有正的空間自相關(guān)，且Z 值均大于2.58，說明計算參數(shù)在0.01 顯著水平下具有統(tǒng)計意義。因此采用GWR 方法分析降水、平均氣溫等要素在不同地理位置對NDVI 的定量影響是合理的。通過GWR 模型殘差項可用來近似表征除降水、平均氣溫以外的其他因素對NDVI 的定量影響，包括地理位置、太陽輻射、人類活動等因素，結(jié)合流域的特點和前人研究［28］，人類活動影響最大，本文用殘差體現(xiàn)人類活動的影響作用；另外，殘差作為NDVI 實測值與模擬值之差具有一定的不確定性，這種不確定性主要來源于輸入信息（如遙感數(shù)據(jù)、實測數(shù)據(jù)）的不確定性、模型結(jié)構(gòu)的不確定性等方面，對結(jié)果影響不大。

表2 各要素GWR 回歸參數(shù)的Moran 指數(shù)與Z 值

永定河山區(qū)1982—2015年年均NDVI 與降水、平均氣溫的GWR 空間回歸參數(shù)及殘差的均值分布如圖7 所示。從圖7 可知，在1982—2015年期間，總體上NDVI 與降水、平均氣溫表現(xiàn)為正相關(guān)關(guān)系，殘差也為正值，這說明降水、平均溫度、人類活動均對植被活動有一定的促進作用。同時，回歸參數(shù)及殘差的均值表現(xiàn)出顯著的空間異質(zhì)性。NDVI 與年降水量的GWR 空間回歸參數(shù)在桑干河、洋河中上游地區(qū)介于5.9×10-4~ 7×10-4，表明這些地區(qū)年均降水量每升高100 mm 會導(dǎo)致NDVI 增加5.9×10-2~ 7×10-2；同時NDVI 對降水變化的敏感性從上游至下游逐漸加大，下游地區(qū)的GWR 空間回歸參數(shù)提升至1.2×10-3~ 1.5×10-3。NDVI 與年均氣溫的GWR 空間回歸參數(shù)介于0 ~ 7.2×10-3，其中約45%的區(qū)域低于1×10-3，主要分布在桑干河和洋河下游以及官廳水庫以下地區(qū)。這說明流域內(nèi)約一半的范圍NDVI 對平均溫度的變化不甚敏感，溫度每升高1℃導(dǎo)致植被NDVI 的增加值低于1×10-3。GWR模型殘差值在各像元內(nèi)均為正，表明永定河山區(qū)人類活動對NDVI 為促進作用，且在流域周邊地區(qū)殘差值較高，這可能由于永定河山區(qū)周邊地區(qū)主要為林地、草地，而中心區(qū)域主要為耕地，近年來人工植被恢復(fù)措施導(dǎo)致林草地NDVI 增長較耕地更快。

圖7 1982—2015年永定河山區(qū)年均NDVI 與降水、平均溫度的GWR 回歸參數(shù)和殘差均值空間分布

1982—2015年永定河山區(qū)不同季度的NDVI 與各影響因素的GWR 回歸參數(shù)均值分布如圖8 所示。從圖8 可以看出，各季度GWR 回歸參數(shù)分布與年平均情況下具有相似空間分布特征，且回歸參數(shù)均為正值。

圖8 1982—2015年永定河山區(qū)各季度NDVI 與降水、平均溫度的GWR 回歸參數(shù)和殘差均值空間分布

NDVI 與降水GWR 的空間回歸參數(shù)值在春、秋季大致介于2.5×10-3~ 7.5×10-3之間，NDVI 變化對降水敏感性最強的區(qū)域分布在永定河下游地區(qū)；而在夏季，GWR 回歸參數(shù)值介于7.6×10-4~ 3.68×10-3，表明這一季節(jié)NDVI 對降水的敏感性程度相對較弱，洋河流域北部地區(qū)敏感性程度相對較高。春、夏、秋季NDVI 與年平均氣溫的GWR 回歸參數(shù)大于2×10-3的面積占總面積的比重分別為27.4%、51.2%和61.2%，可以看出春季NDVI受溫度影響的作用較夏、秋季明顯偏弱；在空間上，各季度NDVI對溫度變化較敏感的區(qū)域主要分布在桑干河、洋河流域上游。不同季節(jié)GWR 模型的殘差項均為正，其中春、夏、秋季殘差項大于0.05 的范圍分別占流域面積的10.4%、22.5%、13.0%，主要分布在流域周邊林草地和濕地，這說明人類活動在春、夏、秋季對NDVI 均為促進作用，且對夏季NDVI 受人類活動影響作用最大，秋季次之、春季最少。

4.2.2 驅(qū)動因子回歸參數(shù)與殘差的趨勢 永定河山區(qū)1982—2015年年均NDVI 與降水、平均氣溫的GWR 回歸參數(shù)及殘差的時間變化趨勢如圖9 所示。從圖9 可知，永定河山區(qū)NDVI 與年降水的GWR空間回歸參數(shù)呈下降趨勢的面積占研究區(qū)域總面積的比例接近17%，主要分布在洋河流域，而其余83%的面積對應(yīng)的像元回歸參數(shù)均呈上升趨勢，這表明1982—2015年流域內(nèi)多數(shù)地區(qū)NDVI 受降水的影響程度呈增加趨勢；NDVI 與年平均氣溫的GWR 回歸參數(shù)呈增加和下降趨勢的像元比例為62.9%和37.1%，其中洋河流域和官廳水庫以下區(qū)域主要呈下降趨勢，而桑干河流域主要呈上升趨勢，這說明洋河流域和官廳水庫以下流域NDVI 受氣溫的影響程度在減弱，而桑干河流域NDVI 受氣溫的影響程度增強；殘差呈增加和下降趨勢的像元比例分別為39.1%和60.9%，表明人類活動對大多數(shù)區(qū)域NDVI的促進有所減弱。

圖9 1982—2015年永定河山區(qū)年均NDVI 與降水、平均溫度的GWR 回歸參數(shù)和殘差趨勢空間分布

1982—2015年永定河山區(qū)不同季度下NDVI 與各影響因素的GWR 空間回歸參數(shù)趨勢分布如圖10所示。降水方面，在春、夏季流域大部分地區(qū)（除東部外）NDVI 與降水的GWR 的回歸參數(shù)斜率為正，而秋季NDVI 與降水的GWR 回歸參數(shù)斜率在全流域內(nèi)為負(fù)，表明在永定河山區(qū)大部分區(qū)域春、夏季NDVI 受降水的促進作用不斷增加，秋季NDVI 受降水的促進作用在全流域內(nèi)呈減小趨勢；平均氣溫方面，春季流域內(nèi)大多數(shù)區(qū)域平均氣溫與NDVI 的GWR 空間回歸參數(shù)斜率為負(fù)，其中洋河流域回歸參數(shù)斜率小于-20×10-5/a，表明春季流域總體上平均溫度對NDVI 的促進作用呈減緩趨勢，且以洋河流域降低程度最為顯著；夏季平均氣溫與NDVI 的GWR 空間回歸參數(shù)斜率為負(fù)的區(qū)域主要分布在洋河流域中下游地區(qū)，占總面積比例為31.3%，其余占流域68.7%的地區(qū)NDVI 受平均氣溫的促進作用呈增加趨勢；秋季NDVI 受溫度的促進作用呈減緩趨勢的面積占44.2%，分布區(qū)域較夏季有所擴大。殘差方面，春、夏、秋季殘差趨勢呈增加的面積占流域面積的比重分別為45.9%、40.7%和45.5%，且主要分布在永定河山區(qū)周邊的林草地地區(qū)，而中部大部分區(qū)域NDVI 受人類活動的促進作用呈減弱趨勢。

圖10 1982—2015年永定河山區(qū)春、夏、秋季NDVI 與降水、平均溫度的GWR 回歸參數(shù)及殘差的趨勢空間分布

4.3 NDVI 變化對流域天然徑流的影響人類活動對河川徑流的影響主要包括兩部分內(nèi)容，其一主要為流域內(nèi)人類取用水和水利工程調(diào)節(jié)的影響，其二主要為下墊面變化（土地利用和覆被變化）的影響［29］。對于天然徑流，鑒于已對社會經(jīng)濟取用水和水利工程調(diào)節(jié)進行還原，天然徑流變化可不考慮前一種影響，僅考慮下墊面變化的的影響。圖11 為永定河山區(qū)1956—2016年典型斷面（冊田水庫斷面、響水堡水庫斷面、官廳水庫斷面）天然徑流量，可以看出流域天然徑流量的年際衰減趨勢十分明顯，這是氣候變化、土地利用和覆被變化共同作用的結(jié)果。

為了進一步探究NDVI 變化對永定河山區(qū)天然徑流的影響，首先對天然徑流變異性進行了診斷，采用滑動T 檢驗法對各斷面1956—2016年天然徑流量進行了突變診斷（采用的顯著性水平α=0.05，統(tǒng)計量臨界值為1.67），診斷結(jié)果如圖12 所示。由該圖可知，冊田水庫、響水堡水庫、官廳水庫斷面的滑動T 檢驗統(tǒng)計量分別在1982年、1983年、1983年取得最大值，且均顯著高于臨界值，這說明桑干河、洋河和永定河上游年際天然徑流量在1982年、1983年、1983年發(fā)生了突變。為便于分析，本文統(tǒng)一將1983年作為3 個斷面天然年徑流量的變異點，同時下文將1956—1982年記作基準(zhǔn)期，將1983—2016年記作變異期。結(jié)合圖11分析，官廳水庫斷面基準(zhǔn)期年均天然徑流量43.4 mm，變異期為21.2 mm，相對減少51.1%；冊田水庫斷面基準(zhǔn)期、變異期年均天然徑流量分別為37.8 和22.5 mm，相對減少40.5%；響水堡水庫斷面基準(zhǔn)期、變異期年均天然徑流量分別為51.0 和21.9 mm，相對減少57.1%。

圖11 1956—2016年永定河山區(qū)主要斷面地表天然徑流量

圖12 1956—2016年永定河山區(qū)控制斷面天然地表徑流量滑動T 檢驗統(tǒng)計量

圖13 給出了永定河山區(qū)干支流斷面在不同階段的年降水量-天然徑流量散點圖及線性回歸關(guān)系。從圖13 可知，基準(zhǔn)期內(nèi)降水量與天然徑流量有較好相關(guān)性，冊田水庫、響水堡水庫和官廳水庫斷面降水-天然徑流量線性回歸方程的決定系數(shù)分別為0.57、0.63 和0.59，表明基準(zhǔn)期降雨量對天然徑流有較強的解釋能力；而變異期內(nèi)降水-徑流量線性回歸方程的決定系數(shù)變?yōu)?.17、0.07 和0.10，降水量對天然徑流的解釋能力明顯變差。變異期年降水-天然徑流量散點基本位于基準(zhǔn)期散點下方，同等年降水量條件下永定河山區(qū)產(chǎn)流能力大幅度下降。根據(jù)兩個時期的年降水-天然徑流回歸關(guān)系，計算得到變異期內(nèi)由降水減少導(dǎo)致各斷面天然徑流累計縮減量分別為3.64、4.26 和3.34 mm，對天然徑流減少的貢獻率分別為16.4%、14.6%和21.8%。相關(guān)研究［30］表明，海河流域潛在蒸散發(fā)變化對徑流的變化影響微弱，降水減少和由于人類活動導(dǎo)致的下墊面變化是流域徑流量減少的主導(dǎo)因素。由于降水減少對永定河徑流減少的貢獻率僅為20%左右，因此導(dǎo)致永定河流域天然徑流減少的主要原因是土地利用和覆被等下墊面條件的變化，由前述分析可知，永定河山區(qū)土地利用格局總體變化不大，而NDVI 呈顯著變大的趨勢，這進一步表明永定河山區(qū)天然徑流衰減是由NDVI 的增加所引起的。

圖13 永定河山區(qū)干支流斷面年降雨量-天然徑流量關(guān)系

通過式（5）與式（6）得到冊田水庫、響水堡水庫、官廳水庫控制流域內(nèi)在不同年份的下墊面參數(shù)n，發(fā)現(xiàn)永定河山區(qū)三個水庫控制流域n 的多年平均值分別為3.18、2.72、2.74，均大于傅抱璞模型推薦的n 的缺省值［31］，因此永定河山區(qū)的植被水分利用效率要高于平均水平，這是由于永定河山區(qū)處于中國半干旱地區(qū)，實際蒸散發(fā)量占降水量的比重高于90%，屬于典型的水分限制區(qū)，因此參數(shù)n 的均值大于平均值符合水熱平衡規(guī)律，且在各控制流域中，冊田水庫控制流域的持水能力最大，而響水堡水庫控制流域的持水能力最小。同時發(fā)現(xiàn)，與天然徑流的變化趨勢相反，永定河山區(qū)下墊面參數(shù)n和年均NDVI 隨時間呈增加趨勢，其中冊田水庫、響水堡水庫、官廳水庫下墊面參數(shù)n 的年均增長率分別為0.0359、0.0426、0.0433，對應(yīng)各流域年均NDVI 的增長率分別為0.0011、0.0013、0.0013，其中，響水堡水庫與官廳水庫的下墊面參數(shù)n 和NDVI 增長情況表現(xiàn)出高度一致性，均大于冊田水庫控制流域。繪制了不同年份滑動平均的下墊面參數(shù)n 與NDVI 點據(jù)分布，對下墊面參數(shù)n 與NDVI 進行擬合，如圖14 和表3 所示，可以看出下墊面參數(shù)n 與NDVI 有較好的正相關(guān)關(guān)系，線性回歸系數(shù)a 分別為17.90、22.14、20.97，截距b 分別為-2.19、-2.84、-3.45，表明NDVI 每增加1%時，下墊面參數(shù)n將增加20%左右，這意味著永定河流域產(chǎn)流性能對植被變化有極強的敏感性。回歸方程的決定系數(shù)分別為0.44、0.45、0.46，對應(yīng)相關(guān)系數(shù)擬合效果均達到顯著性水平（P<0.01）。

表3 永定河山區(qū)Budyko 模型中下墊面參數(shù)n 與植被NDVI 擬合結(jié)果

圖14 永定河山區(qū)典型流域內(nèi)年均NDVI 與下墊面參數(shù)的關(guān)系

永定河山區(qū)NDVI 的增加極大促進了流域河川天然徑流量的衰減作用，結(jié)合前述NDVI 與相關(guān)驅(qū)動因素的GWR 模型計算結(jié)果，NDVI 增加與降水量、平均氣溫作用下的自然生長有關(guān)，也與植被恢復(fù)相關(guān)的退耕還林、風(fēng)沙源治理等大規(guī)模人類活動有關(guān)。根據(jù)調(diào)查，自21 世紀(jì)以來永定河流域內(nèi)開展實施了多項植被恢復(fù)工程，主要有2000年開始實施的京津風(fēng)沙源治理一期工程、2002年實施的退耕還林工程以及2013年實施的京津唐風(fēng)沙源治理二期工程等，而2000年以前流域內(nèi)無相關(guān)大型植被恢復(fù)工程，在一定程度上說明，2000年前植被的生長主要受氣候條件的影響，而2000年后植被生長同時受氣候條件和人類活動的共同影響作用。結(jié)合圖3 永定河山區(qū)年均NDVI 變化圖可知，1982—2000年期間年均NDVI 增加值為0.012，而2000 至2015年年均NDVI 增加0.043，后者遠大于前者，表明2000年后大型植被恢復(fù)工程極大促進了永定河山區(qū)NDVI 的增加，同時也引起流域生態(tài)系統(tǒng)蒸散作用增強，導(dǎo)致永定河山區(qū)在同等降水條件下產(chǎn)流能力顯著降低。有研究已指出，大型植被措施對濕潤地區(qū)的水資源影響較小，而在干旱半干旱地區(qū)往往會顯著減少水資源［32］，顯然永定河山區(qū)植被與河川天然徑流的變化情況與后者是相符的。1982—2015年永定河山區(qū)冊田水庫、響水堡水庫、官廳水庫控制流域的平均濕潤指數(shù)（即降水量與潛在蒸散發(fā)的比值）分別為0.39、0.40、0.41，相應(yīng)下墊面參數(shù)n 的平均值分別為3.18、2.71、2.74，表明永定河山區(qū)屬于持水能力高的半干旱地區(qū)，結(jié)合文獻［33］給出的中國人工植被恢復(fù)的原則，建議未來永定河流域植被恢復(fù)應(yīng)以自然修復(fù)為主，對人工植被恢復(fù)工程應(yīng)進行適當(dāng)?shù)目刂?，同時遵守植被分布規(guī)律來界定植被分區(qū)［34］。

5 結(jié)論

通過構(gòu)建地理加權(quán)回歸模型探究了1982—2015年永定河山區(qū)NDVI 和降水、平均氣溫的空間非平穩(wěn)關(guān)系，并進一步剖析了NDVI 變化對流域典型斷面天然徑流的影響，結(jié)果如下：（1）1982—2015年永定河山區(qū)土地利用格局變化不大，但流域尺度上植被恢復(fù)明顯，年均和春、夏、秋季NDVI平均趨勢率為0.0013/a、0.0013/a、0.0018/a、0.0020/a，且流域內(nèi)絕大多數(shù)空間像元NDVI 的年際增加趨勢均顯著。（2）1982—2015年，永定河山區(qū)年、季尺度下的NDVI 和降水、平均氣溫的回歸參數(shù)、殘差均值以及兩者的年際變化率均具有明顯空間異質(zhì)性。其中回歸參數(shù)及殘差均值在空間像元內(nèi)均為正，表明流域內(nèi)降水、平均氣溫等氣候因素和人類活動因素對NDVI均表現(xiàn)為增加作用，但在不同空間位置的差異亦較為顯著。（3）1956—2016年永定河山區(qū)天然徑流量具有顯著的衰減趨勢，結(jié)合水熱耦合平衡方程發(fā)現(xiàn)流域天然徑流量對植被變化具有高度敏感性，以植被變化為突出特征的下墊面變化對永定河天然徑流衰減具有主導(dǎo)作用。（4）永定河山區(qū)屬于持水能力高的半濕潤半干旱地區(qū)，流域產(chǎn)流能力對植被覆蓋狀況的變化較為敏感。大規(guī)模的恢復(fù)工程的建設(shè)對河川天然徑流衰減的可能影響應(yīng)給予高度重視，建議未來永定河流域植被恢復(fù)應(yīng)以自然修復(fù)為主，對人工植被恢復(fù)工程應(yīng)進行適當(dāng)控制。