利用無人機遙感反演高潛水位礦區(qū)沉陷地玉米葉綠素含量

肖 武，陳佳樂，趙艷玲，胡振琪，呂雪嬌，張 碩

(1.中國礦業(yè)大學(xué)(北京) 土地復(fù)墾與生態(tài)重建研究所，北京 100083; 2.浙江大學(xué) 公共管理學(xué)院，浙江 杭州 310058; 3.中國礦業(yè)大學(xué) 環(huán)境與測繪學(xué)院，江蘇 徐州 221116)

煤炭是我國的主體能源，煤炭資源開采 90%以上的產(chǎn)量來自于地下開采[1]，且多采用走向長壁全部垮落法開采，地表不可避免地產(chǎn)生下沉，造成大量土地的沉陷損毀[2]。我國東部礦區(qū)多為高潛水位礦區(qū)且煤炭和糧食的復(fù)合面積較大，由于煤炭開采對地表造成擾動，使大量優(yōu)質(zhì)耕地沉陷為積水區(qū)，造成原有生態(tài)系統(tǒng)的嚴(yán)重破壞,引發(fā)植被退化、水土流失、物種減少等一系列生態(tài)效應(yīng)[3]，如何對受損生態(tài)環(huán)境進(jìn)行快速監(jiān)測、確定土地受損范圍與程度一直以來都是研究的重點。

國內(nèi)外學(xué)者早期研究中對礦區(qū)環(huán)境的監(jiān)測多利用傳統(tǒng)的衛(wèi)星遙感技術(shù)來實現(xiàn)，且多集中于礦區(qū)沉陷地積水變化[4-5]、植被擾動監(jiān)測[6]、土地復(fù)墾、規(guī)劃[7]等方面，其中，歐洲的穆拉茲等[8]采用 LANDSAT-TM中等分辨衛(wèi)星數(shù)據(jù)和SPOT5較高分辨率的衛(wèi)星數(shù)據(jù)，針對華沙西南的露天煤礦，進(jìn)行了礦區(qū)植被覆蓋動態(tài)變化、礦區(qū)環(huán)境現(xiàn)狀、土地利用現(xiàn)狀及變化狀況的遙感監(jiān)測;CHRISTIAN等[9]將遙感技術(shù)與GIS相結(jié)合應(yīng)用于采礦模擬，研究了采礦活動對地下水位變化的影響，以及與地表地類變化、植被覆蓋變化的3者相互關(guān)系;國內(nèi)的研究中，商華艷等[10]基于遙感和GIS技術(shù)，對云南省安寧磷礦區(qū)的礦山開發(fā)引起的環(huán)境問題進(jìn)行調(diào)查和監(jiān)測，采用層次分析法進(jìn)行礦山地質(zhì)環(huán)境質(zhì)量評價;吳雪茜等[11]應(yīng)用開采沉陷預(yù)測技術(shù)、地理信息及衛(wèi)星遙感技術(shù),對淮南礦區(qū)土地、水域演變趨勢進(jìn)行研究并提出治理對策;徐嘉興等[12]運用遙感、GIS技術(shù)，從景觀生態(tài)學(xué)角度構(gòu)建指標(biāo)體系和評價模型，綜合評價了礦區(qū)土地景觀生態(tài)質(zhì)量。在礦區(qū)生態(tài)擾動監(jiān)測中，遙感衛(wèi)星尺度大、分辨率低、時間周期長且遙感指標(biāo)宏觀化、單一化，無法進(jìn)行地面驗證等問題一直存在。近年來，隨著無人機遙感技術(shù)的不斷成熟，加之無人機響應(yīng)快、周期短、精度高、易操作、成本低的特點，無人機遙感技術(shù)在作物長勢監(jiān)測[13]、農(nóng)田生態(tài)環(huán)境信息監(jiān)測[14]、作物產(chǎn)量估算[15]等方面得到廣泛應(yīng)用，但在礦區(qū)的應(yīng)用還尚處于起步階段，大多數(shù)研究都集中于對礦區(qū)各類采礦設(shè)施與土地利用情況分類與監(jiān)測，礦區(qū)地?zé)豳Y源分布調(diào)查，非法與越界開采識別，露天礦工程量(采剝量、堆放量等)計算等幾個方面[16]，而對于礦區(qū)生態(tài)環(huán)境監(jiān)測、土地?fù)p毀程度識別上的應(yīng)用還有待于進(jìn)一步挖掘。

采煤沉陷對耕地擾動在一定程度上可通過作物生化參數(shù)的改變來反應(yīng)，其中，葉片葉綠素含量是一個連續(xù)的變量，是農(nóng)情遙感監(jiān)測的重要參數(shù)之一，成熟期葉片葉綠素含量直接影響著作物整體的生物量，進(jìn)而決定了作物的產(chǎn)量，而在高潛水位采煤沉陷影響下，區(qū)內(nèi)地形受到沉陷下沉影響導(dǎo)致潛水位相對升高，產(chǎn)生漬害進(jìn)而影響農(nóng)作物的正常生長，其直觀反映為農(nóng)作物葉綠素含量降低，作物枯死和減產(chǎn)、絕產(chǎn)。因此把葉片葉綠素含量作為識別耕地受損程度的表征是可行的。同時作物生化參數(shù)與植被指數(shù)之間也存在著明顯的相關(guān)性，利用植被指數(shù)等遙感參量反演作物生化參數(shù)是監(jiān)測小麥、玉米等作物長勢的重要方法，如裴浩杰等[17]基于多種生化參數(shù)指標(biāo)與光譜指數(shù)結(jié)合偏最小二乘回歸構(gòu)建模型判斷出研究區(qū)小麥整體的長勢差異;楊粉團(tuán)等[18]針對玉米粘蟲災(zāi)情構(gòu)建了基于重歸一化植被指數(shù)多時相的葉片生物量定量模型,實現(xiàn)了對玉米粘蟲災(zāi)情程度的有效監(jiān)測。需要指出的是，現(xiàn)階段研究主要針對某正常農(nóng)田或受某單一脅迫影響植被的生化參數(shù)反演而展開的，并且無人機所搭載傳感器多為普通數(shù)碼相機或近紅外多光譜相機，如牛慶林等[19]基于無人機遙感平臺搭載高清數(shù)碼相機構(gòu)建低成本的遙感數(shù)據(jù)獲取系統(tǒng)，通過加入玉米生長高度作為變量提高了整個生育期的玉米LAI的監(jiān)測精度;賀英等[20]利用無人機搭載數(shù)碼相機獲取RGB影像提取常見的可見光植被指數(shù)結(jié)合多種模型對玉米冠層SPAD進(jìn)行估算，模型最佳效果R2為0.824 7，RMSE為4.3。而采用紅邊波段的多光譜數(shù)據(jù)、以礦區(qū)采煤沉陷耕地作為研究區(qū)的研究甚少。筆者以無人機多光譜相機為數(shù)據(jù)源，結(jié)合地面實測樣本數(shù)據(jù)，研究成熟期植被指數(shù)與耕地?fù)p毀程度的關(guān)系，在前人的研究基礎(chǔ)上篩選相關(guān)性好的植被指數(shù)并借助紅邊波段對葉綠素響應(yīng)敏感的優(yōu)勢加以改進(jìn)，嘗試構(gòu)建最優(yōu)的采煤沉陷影響下的玉米葉片葉綠素?zé)o人機遙感監(jiān)測模型，為土地?fù)p毀監(jiān)測與評價，土地復(fù)墾與生態(tài)修復(fù)等提供基礎(chǔ)數(shù)據(jù)與理論支撐。

1 材料與方法

1.1 試驗區(qū)選擇與樣點建立

東灘煤礦位于山東省濟(jì)寧市境內(nèi)，跨兗州、鄒城、曲阜3市(縣)，地理位置116.85°～116.95°E，35.40°～35.52°N。地處魯中低山丘陵到平原洼地的過渡地帶，整體地勢由東北向西南逐漸降低，潛水埋深為2 m左右。屬于高潛水位礦區(qū)。礦區(qū)內(nèi)土壤類型多為褐土，土質(zhì)較好，土壤肥沃，耕性良好，主要實行冬小麥與夏玉米的輪作模式,其中,夏玉米一般在當(dāng)年6月種植,10月收獲,為重要的糧食產(chǎn)區(qū)。礦區(qū)內(nèi)主要河流有白馬河與泥河，向南流入南陽湖，均為季節(jié)性河流。研究區(qū)為東灘煤礦A1工作面開采影響的范圍(圖1)。其中A1工作面位于三采區(qū)南部3號煤層，地面平均標(biāo)高+50.72 m，開采平均標(biāo)高-497.8 m。所采煤層煤厚平均8.30 m，煤層傾角平均為4°，自2014年8月開始開采后，地面已經(jīng)形成了大面積的沉陷區(qū)，地面沉陷深度已達(dá)7 m，積水深度達(dá)到2 m左右。大量的優(yōu)質(zhì)耕地沉入水里，形成典型的礦區(qū)積水盆地。沉陷區(qū)內(nèi)由內(nèi)及外依次分布為蘆葦、沼澤、灘涂、水淹的農(nóng)作物。結(jié)合該工作面的開采特征，依據(jù)采煤沉陷預(yù)計軟件獲得沉陷區(qū)的下沉等值線，按照傳統(tǒng)定義上的下沉10 mm等值線為影響邊界。保險起見，采用走向、傾向、角平分線方向布設(shè)3條實地采樣線K1,K2,K3，長度分別為550,540,620 m，其中近水區(qū)域植被長勢受采煤擾動影響明顯，需要采集更細(xì)致的樣本信息來表征植被在這一區(qū)域的生長變化趨勢，因此3條樣線上樣點布設(shè)間距均隨下沉值的減小而變大，依次為5,10,20,30,60 m，記錄所有樣點GPS位置和相應(yīng)環(huán)境信息，共計54個采樣點，其中，3條樣線K1,K2,K3最遠(yuǎn)端點布設(shè)在擾動邊界之外作為參照點(圖1)。

1.2 數(shù)據(jù)獲取與處理方法

1.2.1 多光譜影像數(shù)據(jù)獲取與處理

無人機搭載的傳感器為瑞士parrot sequoia[21]多光譜相機(下文簡稱sequoia)，是專為無人機在農(nóng)業(yè)科研、調(diào)查而研發(fā)的，能適用多種飛行器。它不僅可以獲取1 600萬像素RGB三原色照片，還能獲取120萬像素的綠光(波長550 nm、帶寬40 nm)、紅光(波長660 nm、帶寬40 nm)、紅邊(波長735 nm、帶寬10 nm)、近紅外(波長790 nm、帶寬40 nm)等波段影像。依靠自帶的光照傳感器可記錄光照條件并自動校準(zhǔn)4個多光譜傳感器的獨立亮度，同時內(nèi)置GPS和IMU。搭載sequoia的遙感平臺為四旋翼無人機大疆M100。多光譜影像的采集于2017-08-12T 11:00—12:00在沉陷區(qū)上空進(jìn)行，影像獲取時，氣象數(shù)據(jù)顯示，天空晴朗少云，光輻射強度穩(wěn)定。無人機飛行高度110 m，設(shè)定航速9 m/s。傳感器鏡頭視場角15，鏡頭垂直向下，地面分辨率13 cm，航拍面積1.1 km2覆蓋整個研究區(qū)。

圖1 研究區(qū)概況和采樣點布設(shè)Fig.1 Distribution of test sample points in study area

在地面鋪設(shè)好傳感器自帶的光譜反射率校正板，每個架次起飛前先手持飛機在校正板正上方1.5 m處拍照，獲得當(dāng)時條件下的標(biāo)準(zhǔn)反射率(圖2)，選用軟件pix4 d mapper對影像進(jìn)行處理，在處理過程中利用標(biāo)準(zhǔn)反射率校正所有的航拍影像以得到理想的處理成果，并使用ENVI5.1軟件，以研究區(qū)數(shù)碼正射影像為參考影像在圖像不同位置均勻選取 30 個參考點對多光譜影像進(jìn)行幾何精校正，經(jīng)檢驗圖像幾何糾正誤差小于 0.5 個像元，根據(jù)葉綠素地面測量對應(yīng)的樣點位置構(gòu)建感興趣區(qū)域(region of interest,ROI)，以 ROI 范圍內(nèi)地物的平均反射率光譜值作為該樣點玉米葉片反射率光譜，得到各樣點的反射率光譜數(shù)據(jù)。

1.2.2 樣點地表高程獲取

在無人機獲取遙感數(shù)據(jù)的同一天(8月12日)，采用經(jīng)礦區(qū)已知點校正過的南方銀河一號RTK沿K1,K2,K3三條樣線測得樣點地表高程。

1.2.3 葉綠素測量

在進(jìn)行航測的同時，在地面樣線上同步測量玉米葉片葉綠素值。圖1分別沿走向、傾向、角平分線方向樣線由內(nèi)及外進(jìn)行測量。因地面沉陷，部分樣點處于常年積水區(qū)域和濕地中，植被樣本多為蘆葦、雜草等，而本文研究對象為耕地，故最后采集耕地上玉米樣點共39個。測量使用專業(yè)的SPAD-502葉綠素儀進(jìn)行，每個樣點選取1 m×1 m范圍內(nèi)按五點法選取具有代表性的5株玉米，選取每株的頂葉至倒二葉，對每片葉子的葉尖、葉中、葉基3個部位各測定2次，同一片葉共測得6個SPAD值并取平均作為該株整個葉片的相對葉綠素含量[22]。在此基礎(chǔ)上對樣方內(nèi)的五株樣本取其平均值作為該樣點上的葉綠素值。圖3為采樣點實測葉綠素含量與下沉對應(yīng)關(guān)系圖。

1.3 植被指數(shù)選取及改進(jìn)

基于植被光合色素與光譜反射率特征波段具有很強的相關(guān)性，特征波段的選取需要參考色素的光譜特征信息[23]。如圖4所示:在可見光的范圍內(nèi)，葉綠素出現(xiàn)兩個強烈的吸收峰、分別為以450 nm和640～680 nm為中心的藍(lán)波段和紅波段，吸收的峰值出現(xiàn)在670 nm;一個強烈的反射峰出現(xiàn)在550 nm左右的的綠光波段，此時葉綠素吸收系數(shù)最小。在700 nm附近，也為對應(yīng)的葉綠素吸收谷值波段，因此550 nm與700 nm附近常常被選作抗干擾的特征波段，用來削弱非光合作用物質(zhì)引起的光合有效輻射;紅邊區(qū)間700～780 nm為葉綠素在紅邊波段的強吸收到近紅波段多次散射形成的高反射平臺的爬升嵴，這個區(qū)間能非常敏感的響應(yīng)植被營養(yǎng)、長勢、水分、生化參數(shù)等指標(biāo)的變化，已被國內(nèi)外的大量研究結(jié)果證實。其中，當(dāng)植被長勢良好時，紅邊會向長波方向移動(紅移)，當(dāng)植被遭受蟲害、污染、水分等因素脅迫時，紅邊則會向短波方向移動(藍(lán)移)[24];在780 nm以后，近紅外波段葉綠素對電磁波的吸收特征微弱，也常常將780～800 nm附近波段選作特征波段。

結(jié)合以上對光譜曲線的特征波段的分析與本次實驗多光譜傳感器的多通道優(yōu)勢，初步選取了現(xiàn)有研究中常見的且符合本次實驗特征波段組合的植被指數(shù)16個。表1列出了指數(shù)的計算公式及出處。

與現(xiàn)有研究中針對某一單一地類或受某單一脅迫影響的植被葉綠素反演不同的是，該研究區(qū)為礦區(qū)采煤沉陷地，植被種類重疊分布，植被生長情況受脅迫因子較復(fù)雜，選用的16種傳統(tǒng)定義的植被指數(shù)與玉米葉片葉綠素含量的相關(guān)性表現(xiàn)可能與前人研究結(jié)果不一致。因此為保證獲得最優(yōu)的反演精度，提前計算16種植被指數(shù)與實測葉綠素含量的相關(guān)性并根據(jù)其高低進(jìn)行篩選，選取其中相關(guān)性好的部分植被指數(shù)作為參照，利用紅邊波段在此基礎(chǔ)上進(jìn)行改進(jìn)獲得新的植被參數(shù)。在多種新的植被指數(shù)都與葉綠素含量的相關(guān)性良好時，若忽視植被指數(shù)間的近似線性關(guān)系則可能會使得回歸方程不穩(wěn)定，有些植被指數(shù)對葉綠素影響的顯著性被隱蔽起來，某些回歸系數(shù)的符號與實際意義不相符合等[25]。

因此根據(jù)式(1)計算改進(jìn)前后植被指數(shù)兩兩之間的方差膨脹因子(variance inflation factor，VIF)，以分析和控制變量間的多重共線性;再按照相關(guān)性高且顯著的原則，并設(shè)定變量最大 VIF<10，確定改進(jìn)植被指數(shù)及其組合;最后分別將新的植被指數(shù)及其組合作為入選變量輸入葉綠素反演模型。

VIF=1/(1-R2)(1)

式中，VIF為植被指數(shù)間的方差膨脹因子;R為植被指數(shù)間的相關(guān)系數(shù)。

1.4 數(shù)據(jù)分析方法

在野外實驗中，由于實驗環(huán)境和儀器本身等客觀因素或偶然因素的影響，可能會導(dǎo)致獲取的數(shù)據(jù)中存在異常數(shù)據(jù)。為了剔除可能存在的異常數(shù)據(jù)，提高數(shù)據(jù)的可靠性及分析結(jié)果的精度，在建立回歸方程前，對樣本數(shù)據(jù)進(jìn)行回歸分析與殘差分析[26]。根據(jù)求出的殘差和標(biāo)準(zhǔn)殘差值選擇樣本，最終選取了殘差最小的36個樣本。表2為剔除粗差后的實測樣本的統(tǒng)計特征。

表2優(yōu)化后玉米葉片葉綠素含量統(tǒng)計特征
Table2Statisticalcharacteristicsofchlorophyllcontentinoptimizedmaizeleaves

樣本數(shù)平均值(Mean)最大值(Max)最小值(Min)標(biāo)準(zhǔn)差(SD)標(biāo)準(zhǔn)誤差(SE)3654.2371.322.089.91.50

在確定好樣本集后，以預(yù)處理后的無人機多光譜圖像為數(shù)據(jù)源，基于0.05 m空間分辨率數(shù)碼正射影像精確劃分的36個實地樣方范圍，提取各個樣方的植被指數(shù)，在建立經(jīng)驗?zāi)Ｐ椭?，先分別討論初步選取的植被指數(shù)與葉綠素含量的相關(guān)性，根據(jù)結(jié)果進(jìn)行初步排序，篩選出相關(guān)性好的植被指數(shù)并加以改進(jìn)。利用SPSS statistics 19.0軟件對成熟期的36個玉米葉綠素實測值進(jìn)行隨機抽樣，選出26個做為建模樣本，10個做為檢驗樣本[27]。針對改進(jìn)后的植被指數(shù)分別使用線性、對數(shù)、冪指數(shù)等模型進(jìn)行擬合。從中篩選出與玉米葉綠素含量高度相關(guān)且精度相對較高的回歸模型。

1.5 模型評價方法

對上述構(gòu)建的各種模型選取均方根誤差、決定系數(shù)和估測精度進(jìn)行模型分析檢驗，其計算公式見表3。其中決定系數(shù)(R2)表示模擬值與實測值的擬合程度，其值越趨近于1，擬合程度越高;均方根誤差(RMSE)主要用于模型驗證，反映了模擬值與實測值的偏離度，其值越小，模型精度越高;估測精度(EA)為對反演結(jié)果進(jìn)行的精度檢驗及綜合評估，其值越趨近于100%，表明反演模型的估測精度越高。

表3類模型檢驗指標(biāo)及計算公式
Table3Threeassessmentindexesanditsformulas

指標(biāo)公式?jīng)Q定系數(shù)(R2)R2=∑ni=1(xi-x)2(yi-y)2∑ni=1(xi-x)2∑hni=1(yi-y)2均方根誤差(RMSE)RMSE=1n∑ni=1(yi-xi)2估測精度(EA)EA=1-RMSEm ×100%

2 結(jié)果與分析

2.1 采煤擾動與遙感參數(shù)關(guān)系

圖5 各距離玉米生長及對應(yīng)植被指數(shù)特征Fig.5 Corn growth maps and corresponding vegetation index maps at various distances

采煤沉陷造成采空區(qū)上覆巖層的拉伸、斷裂與彎曲等移動變形，致使高潛水位的地區(qū)地下水位相對上升，沉陷中心區(qū)域耕地變?yōu)檎訚?，玉米在?nèi)的農(nóng)作物無法繼續(xù)種植，周邊影響區(qū)域地表變形減弱了土壤持水能力和通氣狀況，影響有機物和礦物質(zhì)的分解、淋溶和沉積，以及由此而引起的土壤侵蝕，使土壤保水能力變差、養(yǎng)分流失嚴(yán)重、土質(zhì)惡化，造成耕地上玉米不同程度的減產(chǎn)[28]。筆者選取樣線K1上3處耕地?fù)p毀程度差異明顯的典型樣方(圖5(a),(b),(c))與該樣線上玉米葉片葉綠素含量的變化趨勢(圖3)相結(jié)合對比分析，可以得知樣線上玉米葉綠素呈現(xiàn)隨地表高程值升高、與積水中心距離增加而增加至最大值后,略微下降后便穩(wěn)定在某一值的現(xiàn)象。同時提取3個典型樣方內(nèi)3種相關(guān)性高的植被指數(shù)值進(jìn)行對比分析(圖5(d))，野外實測樣本點土地不同損毀程度之間植被指數(shù)之間差異明顯，表征到樣線上植被指數(shù)值上，植被指數(shù)變化趨勢與玉米葉綠素一致。這說明在采煤沉陷耕地的復(fù)雜背景中，地下的采煤擾動導(dǎo)致地表作物長勢差異進(jìn)而反應(yīng)到遙感參數(shù)上是一個連續(xù)的過程，通過選取合適植被指數(shù)構(gòu)建最優(yōu)的模型反演不同沉陷深度的玉米葉綠素，進(jìn)而表征來采煤沉陷區(qū)耕地?fù)p毀程度是可行的。

通過分析傳統(tǒng)的16種植被指數(shù)與葉綠素的相關(guān)性結(jié)果(表4)，發(fā)現(xiàn)16種植被指數(shù)與葉片葉綠素均呈正相關(guān)關(guān)系，即植被指數(shù)的降低均對應(yīng)著葉片葉綠素含量的降低。其中，結(jié)合表3中各指數(shù)的計算方式，可以得到初步的規(guī)律:本次試驗中近紅波段和紅邊的組合與葉綠素含量的相關(guān)系數(shù)最高，可達(dá)到0.7以上，其次是近紅波段與綠波段的組合也有較好效果，一般穩(wěn)定在0.6左右。近紅與紅波段的組合最差，但也達(dá)到了0.5。這也佐證了現(xiàn)有學(xué)者認(rèn)為在礦區(qū)復(fù)雜的生態(tài)環(huán)境中，紅邊波段和近紅外波段的組合對葉綠素響應(yīng)更敏感的結(jié)論[29]。因此，改進(jìn)植被指數(shù)的方向為盡可能用紅邊波段替換紅或綠波段，以期獲得更高的擬合精度。

2.2 改進(jìn)植被指數(shù)通過分析

表4中植被指數(shù)與葉綠素含量的相關(guān)性結(jié)果，發(fā)現(xiàn)RVI,NDVI,MVI,GNDVI,NLI五個指數(shù)的相關(guān)性均在0.50以上，其中已經(jīng)存在的NDVI與GNDVI改造的的紅邊指數(shù)NDVIred-edge相關(guān)性達(dá)到了0.65，為了增強其他植被指數(shù)與葉綠素的相關(guān)性，引入對葉綠素響應(yīng)更敏感的紅邊波段，通過利用紅邊波段替換紅或綠波段來提高反演精度。改進(jìn)后的植被指數(shù)與葉綠素含量的相關(guān)系數(shù)普遍高于傳統(tǒng)植被指數(shù)，表5列出了選取改進(jìn)的植被指數(shù)與葉綠素含量的相關(guān)系數(shù)。擴(kuò)展后植被指數(shù)與葉綠素含量的相關(guān)系數(shù)值提高 0.1～0.2，可見，引入紅邊波段對傳統(tǒng)植被指數(shù)進(jìn)行改進(jìn)可增強其與葉片葉綠素的相關(guān)性。

表416種植被指數(shù)和葉片葉綠素含量的相關(guān)系數(shù)(n=26)
Table4Correlationcoefficientbetweenvegetationindexandleafchlorophyllcontent(n=26)

植被指數(shù)相關(guān)系數(shù)RVI0.50**RVI10.58**NDVI0.50**SAVI0.20**DVI0.18**TVI0.18**MVI0.49**MSAVI0.30**CIred-edge0.63**CIgreen0.50**GNDVI0.60**OSAVI0.31**RDVI0.20**NLI0.41**NDGI0.21**NDVIred-edge0.65**

注:**表示在 0.01 水平上顯著相關(guān)。

表5基于植被指數(shù)與葉片葉綠素含量的相關(guān)性(n=26)
Table5Correlationbetweenvegetationindexandleafchlorophyllcontent(n=26)

植被指數(shù)相關(guān)系數(shù)改進(jìn)植被指數(shù)計算公式相關(guān)系數(shù)RVI/ RVI10.50/0.58**RVIred-edgeρNIR/ρRed-edge0.65**MVI0.49**MVIred-edge[(ρNIR-ρRed-edge)/(ρNIR+ρRed-edge)+0.5]1/20.70**NLI0.41**NLIred-edge(ρ2NIR-ρRed-edge)/(ρ2NIR+ρRed-edge)0.50**

注:**表示在 0.01 水平上顯著相關(guān)。

綜合以上的分析結(jié)果，對相關(guān)系數(shù)在0.6以上的所有的植被指數(shù)進(jìn)行雙變量相關(guān)性分析，根據(jù)相關(guān)系數(shù)計算方差膨脹因子，對選取的植被指數(shù)的多重共線性分析(表6)，分析可知:單變量中改進(jìn)后的紅邊植被指數(shù)與葉綠素的相關(guān)性最高，但是由于紅邊指數(shù)之間都存在著嚴(yán)重的共線性，以此構(gòu)建多元回歸模型會出現(xiàn)過度擬和的現(xiàn)象[30]。因此，為進(jìn)一步提高模型的相關(guān)性，考慮以最優(yōu)紅邊指數(shù)為基礎(chǔ)，構(gòu)建帶有綠、紅波段植被指數(shù)的多元回歸模型。

由表3,6分析可知:NDVI,GNDVI兩植被指數(shù)自身與葉綠素含量相關(guān)系數(shù)在0.50以上，且與各指數(shù)之間的的方差膨脹因子均低于10，基本不存在共線性。根據(jù)改進(jìn)植被指數(shù)確定原則，考慮模型計算的簡潔性，選取在MVI指數(shù)上改進(jìn)的植被指數(shù)MVIred-edge和共線性最低的多變量組合(MVIred-edge,GNDVI,NDVI);進(jìn)而分別將兩個變量作為葉綠素反演模型入選變量。

2.3 葉綠素反演模型

基于模型入選變量，分別采用線性、對數(shù)、冪、指數(shù)等方法建模并選取其中最佳的反演模型，并基于驗證樣本進(jìn)行模型的檢驗和比較(表7)。

表7顯示了引入綠、紅波段后的多元組合與最優(yōu)的單變量紅邊指數(shù)在不同回歸模型下對玉米葉綠素的最佳反演能力。在整個單變量回歸模型中，由改進(jìn)后的紅邊MVI經(jīng)驗構(gòu)建的冪模型建模精度最高，R2=0.70,RMSE=1.056 SPAD。以MVIred-edge為基礎(chǔ)構(gòu)建的多元線性回歸模型MVIred-edge,GNDVI,NDVI模型精度提高到R2=0.73，RMSE=0.938 SPAD。

表6相關(guān)性好的植被指數(shù)間的方差膨脹因子
Table6Varianceinflationfactorbetweenbettervegetationindices

參數(shù)NDVIred-edgeCIred-edgeGNDVINDVIRVIred-edgeMVIred-edgeNDVIred-edge893.212.5776.0086.00CIred-edge2.822.4292.0088.00GNDVI7.452.823.30NDVI2.422.58RVIred-edge88.00MVIred-edge

表7優(yōu)選后植被指數(shù)構(gòu)建的最佳反演模型及模型精度檢驗(n=26)
Table7Optimuminversionmodelforconstructingpost-vegetationindexandmodelaccuracytest(n=26)

植被指數(shù)模型表達(dá)式?jīng)Q定系數(shù)均方根誤差MVIred-edge冪Y=76.533X4.1760.70**1.056 SPADMVIred-edge,GNDVI,NDVI多元Y=-146.010-22.237X1+12.515X2+228.542X30.73**0.938 SPAD

經(jīng)驗回歸模型預(yù)測精度的高低與建模樣本數(shù)量關(guān)系密切，在成圖前使用未參與建模的10個樣本數(shù)據(jù)對MVIred-edge構(gòu)建的回歸模型和MVIred-edge,GNDVI,NDVI構(gòu)建的多元回歸模型進(jìn)行精度驗證，綜合評價2種模型對玉米葉綠素含量的預(yù)測能力，分析發(fā)現(xiàn)，MVIred-edge,GNDVI,NDVI3種指數(shù)間的VIF最大不超過13，說明選取的3個植被指數(shù)之間基本不存在共線性。3者構(gòu)造的多變量模型建模精度和預(yù)測精度都優(yōu)于單一的MVIred-edge模型(圖6)，即多變量模型是多種指被指數(shù)中反演玉米葉綠素的最佳植被指數(shù)組合，其構(gòu)建的線性模型預(yù)測玉米葉綠素含量的效果最理想，模型決定系數(shù)達(dá)到0.59，其他指數(shù)模型最高為0.56，整體提高了0.03～0.12，RMSE為3.43，較傳統(tǒng)指數(shù)降低1.0～1.5。其他模型中估測精度EA最高81.3%，整體提高了2.1%～4.7%，最終達(dá)到83.4%。說明利用上述模型顯著提高了采煤沉陷耕地玉米葉綠素含量的反演精度。

圖6 兩種反演模型下葉綠素含量預(yù)測值與實測值的比Fig.6 Comparison between predicted and measured values of chlorophyll content in two inversion models

2.4 反演模型的應(yīng)用

利用上述模型結(jié)合研究區(qū)玉米種植范圍，進(jìn)行采煤沉陷耕地玉米成熟期葉綠素含量反演和制圖(圖7)。獲得研究區(qū)玉米葉綠素含量介于9～79 SPAD，平均值為56 SPAD,標(biāo)準(zhǔn)差為8.7 SPAD，與研究樣本的描述性統(tǒng)計結(jié)果(表2)較為一致。利用自然間斷分類法基于數(shù)據(jù)中固有的自然分組，對分類間隔加以識別，在數(shù)據(jù)值差異相對較大的位置處設(shè)置邊界進(jìn)行分組，將研究區(qū)玉米葉綠素含量分為 5 個等級，統(tǒng)計每個等級的像元數(shù)量及所占比例,見表8。

圖7 基于多元回歸模型采煤沉陷耕地葉綠素反演Fig.7 Chlorophyll inversion map of coal mining subsidence based on multiple regression modle

遙感識別結(jié)果表明，研究區(qū)內(nèi)玉米葉綠素含量主要集中于9～61 SPAD，其面積占整體80.8%，玉米葉綠素含量高于61 SPAD的玉米面積占19.2%。研究區(qū)內(nèi)玉米葉綠素含量整體偏低，反映出耕地整體受采煤沉陷擾動較嚴(yán)重;玉米葉綠素含量空間分布基本表現(xiàn)為沿3條樣線自擾動邊界向積水區(qū)域逐漸降低，至近水區(qū)域均達(dá)到最低值。其中樣線K1最遠(yuǎn)端附近、K3覆蓋區(qū)域玉米葉綠素含量分布在61～79 SPAD，表明此處耕地大多未受影響。中部區(qū)域耕地多屬于中度和輕度損毀，反演結(jié)果與實地調(diào)查樣本基本一致。擾動邊界外部分耕地也呈現(xiàn)輕度擾動的現(xiàn)象，究其原因為研究區(qū)內(nèi)玉米均為農(nóng)民自然播種，玉米品種、施肥的質(zhì)量等其他因素的差異可能導(dǎo)致了這一現(xiàn)象。

表8基于多元模型反演的研究區(qū)玉米葉綠素含量分等統(tǒng)計
Table8Statisticaltableofchlorophyllcontentofcorninresearchareabasedonmultivariatemodelinversion

等級像元數(shù)所占比例/%Ⅰ(9～44 SPAD)1 063 5842.3Ⅱ(44～52 SPAD)1 236 8953.4Ⅲ(52～61 SPAD)24 378 15675.1Ⅳ(61～68 SPAD)4 362 53414.2Ⅴ(68～79 SPAD)1 487 6525.0

3 結(jié) 論

(1)采煤沉陷區(qū)內(nèi)玉米葉綠素含量對紅邊波段的響應(yīng)同樣敏感，致使單變量模型中紅邊和近紅波段的組合指數(shù)表現(xiàn)最佳。

(2)通過對傳統(tǒng)指數(shù)計算中引入紅邊波段進(jìn)行改造，植被指數(shù)與葉片葉綠素含量得相關(guān)系數(shù)提高0.1～0.2，RMSE降低了0.11～1.98，篩選得到相關(guān)性最好的單變量模型紅邊指數(shù)MVIred-edge，并在此基礎(chǔ)上，排除多重共線性，構(gòu)建多元線性回歸模型MVIred-edge,GNDVI,NDVI。其中模型精度為0.73、均方根誤差0.938 SPAD、估測精度為83.4%，3個指標(biāo)均為實驗?zāi)Ｐ椭凶顑?yōu)。表明該方法可有效提高葉片葉綠素含量遙感反演模型的精度。

(3)基于最佳模型對研究區(qū)玉米葉綠素含量空間分布進(jìn)行反演和分析，結(jié)果表明:遙感指標(biāo)識別獲得的研究區(qū)玉米葉綠素含量空間分布與地面樣點描述性統(tǒng)計結(jié)果一致，研究區(qū)內(nèi)80%以上的玉米長勢受到了不同程度的擾動，75.1%的玉米葉綠素含量分布在52～62 SPAD，5.7%的玉米葉綠素含量分布在9～52 SPAD。玉米葉綠素含量整體低于擾動邊界外的平均值，空間分布上葉綠素含量與地表沉陷值變化整體上呈負(fù)相關(guān)，即下沉值越大，玉米葉綠素含量越低。該研究為礦區(qū)土地?fù)p毀監(jiān)測與評價，土地復(fù)墾與生態(tài)修復(fù)等提供基礎(chǔ)數(shù)據(jù)與理論支撐。

4 討 論

研究區(qū)受采煤擾動地表已基本達(dá)到穩(wěn)沉且受損耕地上的玉米屬于成熟期，葉片發(fā)育成熟，葉片葉綠素含量可以很好的反映作物的產(chǎn)量情況。應(yīng)用成熟期的遙感影像和地面實測數(shù)據(jù)可以對作物最終的整體長勢進(jìn)行反演。在以后的研究中，應(yīng)在玉米生長的各個時期同步或準(zhǔn)同步進(jìn)行無人機多光譜影像和地面數(shù)據(jù)的采集進(jìn)一步提高災(zāi)害過程實時監(jiān)測的準(zhǔn)確性和時效性。同時還可以設(shè)置地面模擬試驗進(jìn)行多光譜的采集和分析，加強基于葉片生物量、作物產(chǎn)量模型的反演機理研究和精度提高。

高潛水位采煤沉陷耕地的特點是地下潛水位埋深高，采煤沉陷后地下水位容易上升到地表標(biāo)高以上，導(dǎo)致地表積水，積水程度不同，作物受影響程度不同，葉綠素是植被光合作用中最為重要的色素，是植被光合作用能力強弱、生理脅迫狀況、固碳能力及氮利用效率的良好指示器。因此監(jiān)測作物葉綠素含量作為表征耕地?fù)p毀程度的一個指標(biāo)是可行的。

筆者選用改造后的紅邊指數(shù)MVIred-edge在單變量指標(biāo)中相關(guān)系數(shù)最高，究其原因主要是礦區(qū)采煤沉陷地本身情況復(fù)雜，礦區(qū)植被整體長勢不佳、葉綠素含量偏低，紅邊波段向左移動，發(fā)生“藍(lán)移”現(xiàn)象，致使紅邊和近紅外波段組合指數(shù)表現(xiàn)最佳。在充分利用紅邊優(yōu)勢的情況下，排除多重共線性的影響。最終選取多元線性回歸模型MVIred-edge,GNDVI,NDVI有效提高了反演精度。在后期的研究中可以考慮利用PLSR,BPNN,SVM等模型進(jìn)一步提高模型精度。

在野外實地采集樣本的過程中，發(fā)現(xiàn)研究區(qū)內(nèi)植被的分布趨于多樣化，由沉陷中心的沉水植物到浮葉植物再到近水岸邊的挺水植物，最后才是長勢各異的玉米，此次試驗只研究了現(xiàn)有耕地中各樣觀測線上玉米葉片葉綠素含量的變化趨勢，導(dǎo)致可用的樣本較少，限制了反演精度的提高。因此如何將觀測線上所有植被的生長參數(shù)如生物量、作物產(chǎn)量等指標(biāo)標(biāo)準(zhǔn)化、統(tǒng)一化，構(gòu)建多指標(biāo)、多植被、時序性的反演模型是下一步研究工作的重點。