植被指數(shù)最佳地理軸線研究

姜欣彤，周祖昊，黎 曙，龐 亮，嚴(yán)子奇

（1.中水珠江規(guī)劃勘測設(shè)計(jì)有限公司，廣東 廣州 510610；2.中國海洋大學(xué) 工程學(xué)院，山東 青島 266100；3.中國水利水電科學(xué)研究院 流域水循環(huán)模擬與調(diào)控國家重點(diǎn)實(shí)驗(yàn)室，北京 100038；4.黃河水利委員會 西寧水文水資源勘測局，青海 西寧 810008）

1 引 言

在地形復(fù)雜、高程差異顯著的高原山區(qū)，植被的空間分布格局具有明顯的特點(diǎn)。有學(xué)者發(fā)現(xiàn)高原山區(qū)的植被分布受高程影響，植被指數(shù)（NDVI）空間分布具有明顯的垂向條帶性特征。如陳菁［1］對植被垂直分布空間格局進(jìn)行研究，并以位于四川省西南部的邛崍市為例探究植被與高程帶之間的關(guān)系，結(jié)果表明不同高程帶內(nèi)植被分布具有差異性；Li等［2］針對喜馬拉雅山東部雅魯藏布江流域內(nèi)NDVI沿海拔梯度的分布模式進(jìn)行分析，發(fā)現(xiàn)NDVI在海拔3 000 m附近區(qū)域內(nèi)達(dá)到最大；劉曉婉等［3］對雅魯藏布江流域NDVI的空間分布特征進(jìn)行分析，發(fā)現(xiàn)該流域內(nèi)NDVI與高程具有相依性，總體呈現(xiàn)隨高程提高而減小的趨勢。另有學(xué)者研究發(fā)現(xiàn)高原山區(qū)的植被分布除受高程影響外，地理要素也是其變化的影響因子，NDVI空間分布具有明顯的水平地帶性特征。如趙建軍［4］將阿帕拉契亞山NDVI分別與高程和地理因子的相關(guān)關(guān)系進(jìn)行分析，研究表明高程和經(jīng)緯度均影響NDVI空間分布，該區(qū)域內(nèi)高程每升高100 m，其NDVI下降0.006 2個單位，緯度每向北移動1°NDVI下降0.01個單位；張文江等［5］在對青藏高原中東部NDVI空間分布格局進(jìn)行分析時發(fā)現(xiàn)，青藏高原呈現(xiàn)獨(dú)特的三維地帶性特征，在3 500 m高程以下NDVI與緯度相關(guān)關(guān)系較強(qiáng)，而3 500 m以上NDVI與經(jīng)度相關(guān)關(guān)系更強(qiáng)；趙健赟等［6］對青海高原植被空間分布特征進(jìn)行研究發(fā)現(xiàn)，研究區(qū)內(nèi)NDVI與經(jīng)度有明顯的正相關(guān)關(guān)系，與緯度有明顯的負(fù)相關(guān)關(guān)系，NDVI分布有明顯的水平地帶性特征。

縱觀前人針對高原山區(qū)植被空間分布特征的研究，大都聚焦于NDVI與高程、經(jīng)度和緯度的相關(guān)關(guān)系上，以此說明NDVI具有垂向條帶性和經(jīng)緯向地帶性。但是，在地形條件復(fù)雜的高原山區(qū)，植被指數(shù)并非嚴(yán)格按照經(jīng)度或緯度的方向呈現(xiàn)條帶性變化趨勢，或許是沿著某一方向表現(xiàn)出顯著的條帶性分布特征，因此單用高程或經(jīng)緯度作為指標(biāo)難以定量描述。筆者在前人研究的基礎(chǔ)上，將經(jīng)緯度作為組合因子，提出植被指數(shù)空間變化最佳地理軸線的概念，即在此方向上，植被指數(shù)的空間分布具有最顯著的條帶性特征，以此深入探究高原山區(qū)的植被指數(shù)空間分布條帶性變化規(guī)律。

位于黃河干流唐乃亥水文斷面以上的黃河源區(qū)，集水面積僅占黃河流域面積的16.2%，但徑流量卻占黃河流域總徑流量的三分之一，故有黃河的“水塔”之稱［7］。黃河源區(qū)位于青藏高原東北部，因特殊的地理位置和氣候特征，其生態(tài)環(huán)境質(zhì)量對下游地區(qū)的生態(tài)環(huán)境有重要的影響，因此對黃河源區(qū)植被空間分布規(guī)律進(jìn)行研究具有重要意義。幅員遼闊、地形復(fù)雜的黃河源區(qū)植被空間分布格局遵循怎樣的規(guī)律？地理因子如何影響黃河源區(qū)植被空間分布格局？造成黃河源區(qū)植被空間分布格局的機(jī)制是什么？本文將從最佳地理軸線的角度出發(fā)，揭示黃河源區(qū)植被指數(shù)的空間分布條帶性變化規(guī)律，并對其形成如此空間格局的原因作出合理解釋。

2 研究區(qū)概況及數(shù)據(jù)

黃河源區(qū)地理坐標(biāo)范圍為北緯32°09′—36°34′、東經(jīng)95°54′—103°24′，總面積約為13萬km2?？傮w地形特點(diǎn)主要為起伏和緩的丘狀高原，區(qū)域內(nèi)地勢西北高、東南低。其中分布有峽谷、高山、盆地、冰川等地貌，西部和中部地勢高，最高處位于阿尼瑪卿山，最低處位于龍羊峽，海拔在2 409～6 254 m之間（見圖1）［8］。源區(qū)內(nèi)自然環(huán)境類型多樣，高寒植被分布廣泛，約占源區(qū)面積的3／4。主要的植被生態(tài)類型包括高寒草甸、高寒草原、高寒沼澤、高寒灌叢、常綠針闊葉林等［9］。

圖1 黃河源區(qū)高程及氣象站

本研究采用的數(shù)據(jù)如下：①美國NASA航天中心發(fā)布的MOD13Q1歸一化植被指數(shù)，時間分辨率為16 d，空間分辨率為1 km，數(shù)據(jù)選用年份為2003—2017年；②中科院地理空間數(shù)據(jù)云網(wǎng)站下載的數(shù)字高程（DEM）數(shù)據(jù)，空間分辨率為30 m；③基于研究區(qū)內(nèi)12個氣象站通過反距離權(quán)重法空間插值后再進(jìn)行高程修正的氣溫、降水量逐月數(shù)據(jù)，將2003—2017年每年逐月降水量求和作為年降水量，逐月氣溫求平均作為年氣溫。

結(jié)合MRT技術(shù)、ArcGIS及MVC法對MOD13Q1遙感數(shù)據(jù)產(chǎn)品進(jìn)行處理，得到NDVI值。通過MRT技術(shù)可以對遙感影像進(jìn)行拼接和重投影處理，并提取NDVI波段。利用ArcGIS掩膜工具提取研究區(qū)范圍，利用柵格計(jì)算工具將原始NDVI資料從-10 000～10 000的浮點(diǎn)數(shù)轉(zhuǎn)化為-1～1之間的標(biāo)準(zhǔn)NDVI值。本文在計(jì)算時只考慮有植被的區(qū)域，因此將柵格中河流、湖泊、巖石和裸土等NDVI值小于0的非植被區(qū)域像元設(shè)置為空值以剔除。MVC法是目前國際上通用的最大化合成法［10］，即將每月內(nèi)的NDVI數(shù)據(jù)取最大值進(jìn)行像元合成，以此消除云、大氣和太陽高度角的部分干擾［11］。將通過MVC法進(jìn)行像元合成后的NDVI作為月NDVI值。植被生長季為4—10月，本文將每年4—10月NDVI值計(jì)算平均值作為該年植被生長季NDVI，得到研究區(qū)2003—2017年15幅生長季NDVI空間分布圖，再計(jì)算15 a植被生長季NDVI多年平均值，以此為數(shù)據(jù)基礎(chǔ)研究黃河源區(qū)植被指數(shù)空間分布規(guī)律（見圖2）。采用ArcGIS工具在研究區(qū)內(nèi)等間距選取1 990個樣本點(diǎn)，提取樣本點(diǎn)位置所在經(jīng)緯度和15 a植被生長季平均NDVI值。

圖2 黃河源區(qū)植被指數(shù)空間分布規(guī)律

3 最佳地理軸線計(jì)算方法及原理

本文提出最佳地理軸線的概念，即將經(jīng)度和緯度作為組合因子，定義一種新的坐標(biāo)系，定義方法為：從地理平面上選定一個原點(diǎn)，從原點(diǎn)向不同角度作極軸，計(jì)算樣本點(diǎn)xi在軸線上的投影距離Di。對于每一角度的軸線，可建立樣本點(diǎn)NDVI值與樣本點(diǎn)在此軸線上投影距離D的回歸關(guān)系并計(jì)算決定系數(shù)R2，R2最大值對應(yīng)的軸線即為“最佳地理軸線”，沿著此軸線方向，植被指數(shù)空間分布呈現(xiàn)顯著的條帶性變化趨勢。

將經(jīng)線、緯線分別作為橫縱坐標(biāo)，單位為度。以橫坐標(biāo)大于所有樣本點(diǎn)經(jīng)度最大值、縱坐標(biāo)小于所有樣本點(diǎn)緯度值為原則選取一個原點(diǎn)。從原點(diǎn)向研究區(qū)內(nèi)做不同角度射線，作為地理軸線，將經(jīng)線記作K1，緯線記作K2，穿過研究區(qū)內(nèi)部的軸線記作Ki。將研究區(qū)內(nèi)所有樣本點(diǎn)NDVI值與樣本點(diǎn)在軸線Ki上的投影距離D建立回歸關(guān)系，并計(jì)算決定系數(shù)R2。

決定系數(shù)R2也稱為擬合優(yōu)度，它表征回歸方程在多大程度上解釋了因變量的變化，即表征因變量NDVI的變異中有多少可由控制的自變量投影距離D來解釋。其計(jì)算公式為

式中：SSE為殘差平方和；SST為總離差平方和；y為樣本點(diǎn)觀測值；ˉy為樣本點(diǎn)均值；＾y為根據(jù)擬合回歸關(guān)系得到的樣本點(diǎn)擬合值。

R2的取值范圍為0～1，R2越大，自變量對因變量的解釋程度越高，自變量引起的變動占總變動的百分比越高，觀察點(diǎn)在回歸直線附近越密集。對于樣本點(diǎn)NDVI值與投影距離D的回歸關(guān)系，R2越高，表明對于同一投影距離D，NDVI分布越密集。過原點(diǎn)做半徑為D的圓，過不同角度的軸線（Ki）與圓的交點(diǎn)做切線（ti）。由于該切線上所有樣本點(diǎn)在軸線上的投影距離都為D，因此若此條切線上NDVI分布較散亂、大小不一，擬合NDVI值與樣本點(diǎn)NDVI值殘差平方和就會大，導(dǎo)致回歸模型決定系數(shù)R2變??；若此條切線上樣本點(diǎn)NDVI值分布較集中，樣本點(diǎn)NDVI值在回歸線附近密集，則決定系數(shù)R2就大。選擇R2最大值對應(yīng)角度的軸線作為最佳地理軸線，表明在此軸線的垂線方向上NDVI分布范圍較為集中，沿此軸線NDVI空間分布具有明顯條帶性變化規(guī)律。

4 分析與討論

以黃河源區(qū)NDVI為例，計(jì)算其最佳地理軸線，如圖3所示。把研究區(qū)NDVI值劃分為＜0.2（藍(lán)）、0.2～0.4（紅）、0.4～0.6（綠）和＞0.6（紫）4個區(qū)間，將樣本點(diǎn)所在經(jīng)緯度和NDVI值做氣泡圖。由圖3可知，研究區(qū)NDVI值從東南向西北方向遞減。以（103.5°，31.5°）為原點(diǎn)，從該原點(diǎn)出發(fā)分別把斜率k在tan 90°～tan 180°范圍內(nèi)的每一度軸線作為地理軸線，建立樣本點(diǎn)NDVI值與每條軸線上投影距離D的回歸關(guān)系并計(jì)算決定系數(shù)R2和均方誤差MSE，并對決定系數(shù)的顯著性進(jìn)行F檢驗(yàn)。

圖3 黃河源區(qū)植被指數(shù)最佳地理軸線計(jì)算

不同角度軸線樣本點(diǎn)NDVI值與軸線上投影距離D回歸關(guān)系如圖4所示。R2、F統(tǒng)計(jì)值和回歸模型的均方誤差MSE隨角度變化規(guī)律如圖5所示。

圖4 不同角度軸線樣本點(diǎn)NDVI值與投影距離D回歸關(guān)系

從與x軸正方向夾角90°（經(jīng)線）到夾角180°（緯線）內(nèi)的每一度軸線，NDVI與投影距離的回歸模型均通過了0.05顯著性水平F檢驗(yàn)，這說明地理因子對黃河源區(qū)NDVI確實(shí)有顯著影響。由圖5可知，當(dāng)軸線與x軸正方向夾角為128°，即西偏北52°時，NDVI與此軸線上的投影距離回歸關(guān)系最佳，R2達(dá)到0.649 6，此時回歸模型的F統(tǒng)計(jì)值最大，MSE最小，在此方向上，NDVI在擬合回歸線附近分布較密集，同一投影距離上NDVI分布較為集中，即在此方向上NDVI具有明顯的條帶性，因此西偏北52°為研究區(qū)NDVI空間分布的最佳地理軸線。

圖5 NDVI與投影距離回歸模型統(tǒng)計(jì)指標(biāo)隨軸線角度變化

為深入探究NDVI空間分布格局形成的原因，分析影響NDVI的氣象因子（降水量和氣溫）與NDVI的相關(guān)關(guān)系。對樣本點(diǎn)多年平均NDVI值分別與多年平均降水量和多年平均氣溫進(jìn)行皮爾遜相關(guān)分析，得到降水量與NDVI的相關(guān)系數(shù)為0.72，氣溫與NDVI的相關(guān)系數(shù)為0.40。結(jié)果表明，降水量與NDVI有強(qiáng)相關(guān)性。鑒于此，用最佳地理軸線計(jì)算方法對NDVI主要影響因子降水量進(jìn)行分析，以探析研究區(qū)內(nèi)降水的空間分布格局是否具有同樣的顯著條帶性分布規(guī)律。提取研究區(qū)內(nèi)1 990個樣本點(diǎn)的15 a平均降水量，將其與在不同軸線上的投影距離分別建立回歸關(guān)系，R2、F統(tǒng)計(jì)值和回歸模型的均方誤差MSE隨軸線角度變化規(guī)律如圖6所示。對于降水量，每一條軸線對應(yīng)的回歸模型也均通過了0.05顯著性水平F檢驗(yàn)，說明降水也受地理因子影響。由圖6可知，降水的最佳地理軸線角度為西偏北53°，此時回歸模型F統(tǒng)計(jì)值最大，MSE最小，R2最大為0.791 4。

圖6 降水量與投影距離回歸模型統(tǒng)計(jì)指標(biāo)隨軸線角度變化

對比圖5和圖6可以看出，隨著軸線角度在經(jīng)線和緯線之間轉(zhuǎn)移，NDVI與投影距離回歸模型和降水與投影距離回歸模型基本展現(xiàn)出了較為一致的變化規(guī)律，均能找到R2最大、F統(tǒng)計(jì)值最大、均方誤差MSE最小的一個角度，且NDVI的最佳地理軸線為西偏北52°，降水的最佳地理軸線為西偏北53°，二者僅偏差1°。同時也可看出，在最佳地理軸線附近一定區(qū)域內(nèi)的角度，其R2與最佳地理軸線對應(yīng)的R2相差不大，為量化最佳地理軸線與其他角度上的R2差異，根據(jù)斜率計(jì)算公式做出R2變化斜率圖（見圖7），由圖7可知，無論是降水量還是NDVI，最佳地理軸線附近的軸線R2變化平緩，距離最佳地理軸線遠(yuǎn)的軸線R2變化相對明顯，角度偏向緯線方向的軸線比角度偏向經(jīng)線方向的軸線R2變化平緩。若以斜率的絕對值0.001作為閾值，則可以確定一個最佳地理軸線范圍，即對于NDVI，123°～137°是其具有顯著條帶變化的方向，對于降水量，123°～133°是其具有顯著條帶變化的方向。

圖7 R2變化斜率隨軸線角度變化

通過對黃河源區(qū)植被指數(shù)空間分布最佳地理軸線角度的計(jì)算，以及影響植被指數(shù)主要?dú)庀笠蜃咏邓康目臻g分布最佳地理軸線角度計(jì)算發(fā)現(xiàn)，研究區(qū)內(nèi)NDVI的最佳地理軸線角度為西偏北52°，降水的最佳地理軸線角度為西偏北53°，且NDVI與降水量均在最佳地理軸線附近一定區(qū)域內(nèi)，呈現(xiàn)出較為顯著的條帶性分布特征，二者的最佳地理軸線角度僅偏差1°，這也很好地解釋了為何研究區(qū)內(nèi)NDVI與降水量強(qiáng)相關(guān)。由此可見，在空間分布上，黃河源區(qū)降水量與NDVI的條帶變化規(guī)律具有較高的一致性，但存在微小差異，這是因?yàn)橹脖环植纪瑫r還受氣溫、土壤等因素影響。綜上所述，黃河源區(qū)的植被空間分布受地理因子影響，但地理因子只是影響植被分布的間接因子，造成黃河源區(qū)植被空間分布格局的機(jī)制主要是降水的分布規(guī)律。

5 結(jié) 論

以經(jīng)緯度為組合因子，提出植被指數(shù)空間變化最佳地理軸線的概念，深入探究黃河源區(qū)植被指數(shù)空間分布條帶性變化規(guī)律。選取研究區(qū)1 990個樣本點(diǎn)，將樣本點(diǎn)NDVI值與不同角度地理軸線投影距離建立回歸關(guān)系，找到研究區(qū)植被指數(shù)條帶性變化的方向，并對呈現(xiàn)此種空間變化格局的原因做出合理解釋，得出的主要結(jié)論如下。

（1）在地形復(fù)雜的高原山區(qū)植被指數(shù)往往不是嚴(yán)格按照經(jīng)緯向呈現(xiàn)水平地帶性，而可能是沿著某一特定軸線角度呈現(xiàn)出顯著的條帶性，因此將經(jīng)緯度作為組合因子，定義一種新的計(jì)算方法，可以得到植被指數(shù)空間分布最佳地理軸線。

（2）植被指數(shù)與降水量的相關(guān)性較強(qiáng)，相關(guān)系數(shù)為0.72，與氣溫的相關(guān)性較弱，相關(guān)系數(shù)為0.40。利用最佳地理軸線計(jì)算方法分別對黃河源區(qū)的植被指數(shù)和降水的空間分布規(guī)律進(jìn)行探析，結(jié)果表明植被指數(shù)的最佳地理軸線為西偏北52°，降水量的最佳地理軸線為西偏北53°。

（3）黃河源區(qū)植被分布主要受降水影響，二者在空間上的分布規(guī)律基本一致，但植被分布同時還受氣溫、土壤等因素影響，因此NDVI與降水的“最佳地理軸線”存在微小差異。