肖星星，呂偉才，仲臣，曹進(jìn)進(jìn)，李延杰，張燦

（ 1. 安徽理工大學(xué) 空間信息與測(cè)繪工程學(xué)院, 安徽 淮南 232001；2. 安徽理工大學(xué) 礦區(qū)環(huán)境與災(zāi)害協(xié)同監(jiān)測(cè)煤炭行業(yè)工程研究中心, 安徽 淮南 232001；3. 安徽理工大學(xué) 礦山采動(dòng)災(zāi)害空天地協(xié)同監(jiān)測(cè)與預(yù)警安徽普通高校重點(diǎn)實(shí)驗(yàn)室, 安徽 淮南 232001；4. 安徽廣航建設(shè)工程有限公司, 安徽 蕪湖 241002 ）

0 引 言

歸一化植被指數(shù)(NDVI)被用作對(duì)植物覆蓋面積生長(zhǎng)狀況的長(zhǎng)期監(jiān)測(cè)，是衡量植被覆蓋的標(biāo)準(zhǔn)化方法，對(duì)因植被覆蓋變化導(dǎo)致的區(qū)域氣候變化，農(nóng)作物種植面積等研究提供定量的判斷依據(jù)[1-2]. NDVI作為分析人員在遙感中最常用的指標(biāo)，形成了改進(jìn)型甚高分辨率輻射計(jì)(AVHRR)、中分辨率成像光譜儀(MODIS)、可見光紅外成像輻射儀(VIIRS)等數(shù)據(jù)采集形式，但此類方法對(duì)觀測(cè)條件要求較高，且背景噪聲消除復(fù)雜，難以完成連續(xù)觀測(cè)[3-4].

相較于上述方法， GPS對(duì)地觀測(cè)不僅能克服NDVI數(shù)據(jù)采集的缺陷，而且具有大范圍、全天時(shí)全天候、連續(xù)運(yùn)行等先天優(yōu)勢(shì). GPS是通過地面接收設(shè)備接收來自空間衛(wèi)星發(fā)射的無線電信號(hào)，進(jìn)行空間定位的空間導(dǎo)航系統(tǒng)，其中L波段微波信號(hào)有良好的大氣穿透性，會(huì)被植被含水量的細(xì)微變化產(chǎn)生干涉，導(dǎo)致載波信號(hào)振幅的改變. 這一特性為使用GPS干涉反射測(cè)量(GPS-IR)進(jìn)行NDVI反演提供了條件.LARSON等[5]利用GPS反射信號(hào)特點(diǎn)，對(duì)地表積雪厚度進(jìn)行了測(cè)量，證明了GPS進(jìn)行地表覆蓋形態(tài)觀測(cè)的可行性；CHEW等[6]同時(shí)對(duì)土壤濕度進(jìn)行了反演，進(jìn)一步說明了利用GPS反演土壤理化性質(zhì)具有較高精度；SMALL等[7]也通過計(jì)算GPS偽距多路徑均方根誤差(RMSE)進(jìn)行植被含水量的反演，證實(shí)了其與植被各階段生長(zhǎng)狀況不同含水量存在較好的相關(guān)性；隨后LARSON[8]等提出了一種歸一化微波反射指數(shù)(NMRI)，發(fā)現(xiàn)在植被指數(shù)的反演過程中具有較好的實(shí)用性. WAN等[9]建立了GPS反射信號(hào)振幅與植被含水量之間的回歸模型，發(fā)現(xiàn)反演出的植物含水量模型具有較高精度；鄭南山等[4]將BP神經(jīng)網(wǎng)絡(luò)應(yīng)用于NDVI模型反演中進(jìn)一步提升模型精度. 上述研究成果為利用GPS反射信號(hào)來進(jìn)行NDVI反演提供了理論依據(jù).

綜上，GPS-IR反演NDVI具有觀測(cè)成本低廉、實(shí)時(shí)連續(xù)、不受傳感器和大氣環(huán)境限制等優(yōu)點(diǎn). 本文基于上述特點(diǎn)對(duì)兩個(gè)不同區(qū)域GPS參考站數(shù)據(jù)進(jìn)行實(shí)驗(yàn)研究，進(jìn)行了NMRI與NDVI之間的相關(guān)性分析，實(shí)現(xiàn)了NDVI植被指數(shù)的反演.

1 原理與方法

1.1 GPS-IR進(jìn)行NDVI反演原理

圖1為GPS接收機(jī)在對(duì)地觀測(cè)中記錄的信號(hào)強(qiáng)度變化. 其中L波段有較好地大氣穿透性并且對(duì)于植被水分含量變化比較敏感，因此GPS接收機(jī)接收到L波段的直射信號(hào)，與地面植被反射信號(hào)產(chǎn)生干涉.

圖1 GPS-IR進(jìn)行NDVI反演幾何模型

由于多路徑效應(yīng)導(dǎo)致接收到的反射信號(hào)振幅發(fā)生變化，對(duì)振幅進(jìn)行處理計(jì)算歸一化微波指數(shù). 其中GPS-IR可以采用常規(guī)單天線測(cè)量型接收機(jī)進(jìn)行地基觀測(cè). 地面上有效觀測(cè)反射區(qū)是指接收機(jī)能夠接受到反射信號(hào)區(qū)域，在測(cè)區(qū)較為平坦的條件下，此區(qū)域可以用第一菲涅耳反射區(qū)進(jìn)行較為全面的描述.

如圖2所示，Q點(diǎn)為空間衛(wèi)星在地面的鏡像點(diǎn)，由惠更斯-菲涅爾原理[10]可知，以Q點(diǎn)和接收機(jī)為焦點(diǎn)的旋轉(zhuǎn)橢球面中空間區(qū)域?yàn)榉颇鶢柗瓷鋮^(qū)，此橢圓與地面相交區(qū)域(圖2中黃色區(qū)域)即地面反射的有效區(qū)域?yàn)榈谝环颇瓷鋮^(qū)，由長(zhǎng)半軸a、短半軸b決定的橢圓，其長(zhǎng)、短半軸可根據(jù)下式計(jì)算：

圖2 菲涅爾區(qū)示意圖

式中： λ 為GPS衛(wèi)星載波波長(zhǎng)，單位為m；H為接收機(jī)天線高，單位為m；e為衛(wèi)星高度角.

由式(1)可知，當(dāng)GPS觀測(cè)條件作為固有變量保持不變時(shí)，第一菲涅爾反射區(qū)形狀影響參數(shù)僅為衛(wèi)星高度角，當(dāng)使用載波波長(zhǎng)為0.19 m的L1載波信號(hào)作為計(jì)算參數(shù)時(shí)，接收機(jī)天線架設(shè)保持1.8 m高度不變時(shí)，以自變量高度角5°，有效反射區(qū)域?yàn)?41.443 8 m2；當(dāng)高度角選擇7°時(shí)，有效觀測(cè)面積縮小至72.341 4 m2；衛(wèi)星高度角為10°，有效觀測(cè)面積為35.631 6 m2.為了方便觀察上述性質(zhì)，選用衛(wèi)星方位角為0°進(jìn)行繪圖，圖2中接收機(jī)位置為(0，0)點(diǎn)天線高為1.8 m.

由圖3可知，接收機(jī)天線架設(shè)穩(wěn)定后衛(wèi)星高度角越大，反射信號(hào)覆蓋區(qū)越小，且逐漸向接收機(jī)方向縮小[11-13].

圖3 不同高度角下的有效散射區(qū)大小

反射點(diǎn)位置是在接收反射信號(hào)的過程中，其反射信號(hào)與接收信號(hào)的最短路徑與反射面的交點(diǎn)，其中該最短路徑稱為徑向距離. 如圖1所示，反射信號(hào)在地面的反射點(diǎn)O為反射點(diǎn)位置，由接收機(jī)天線高和衛(wèi)星高度角定義

式中：r為徑向距離，單位為m；H為接收機(jī)天線高；e為衛(wèi)星高度角. 隨著GPS衛(wèi)星與接收機(jī)相對(duì)位置的變化，反射點(diǎn)軌跡也隨之變化，在不同的弧度表現(xiàn)出不同的走向和分布情況. 圖4給出了徑向距離與衛(wèi)星高度角、接收機(jī)天線高度的關(guān)系圖.

圖4 徑向距離與衛(wèi)星高度角、天線高度的關(guān)系圖

由圖4可知，隨著GPS衛(wèi)星與接收機(jī)相對(duì)位置的變化，反射點(diǎn)軌跡也隨之變化，在不同的弧度表現(xiàn)出不同的走向和分布情況. 徑向距離與衛(wèi)星高度角的變化成正相關(guān)，與天線高度成負(fù)相關(guān). 綜合多顆衛(wèi)星反射點(diǎn)軌跡和菲涅爾反射區(qū)計(jì)算，可判定反射信號(hào)充分?jǐn)y帶參考站周圍的植被信息，所以準(zhǔn)確確定菲涅爾反射區(qū)，直接關(guān)系到NDVI遙感的準(zhǔn)確性和精度.

1.2 NMRI

NMRI是衡量GPS反射信號(hào)振幅變化的一項(xiàng)重要指標(biāo)，其中對(duì)于L1載波上偽距多路徑RMSE用MP1的計(jì)算可表示為[14]

式中：P1是以m為單位在L1載波上偽距觀測(cè)值；f1=1 575.42 MHz，f2=1 227.60 MHz；λ1=0.19 m，λ2=0.24 m；

φ1和φ2是以周為單位的L1、L2載波相位觀測(cè)值. 首先對(duì)接收機(jī)觀測(cè)數(shù)據(jù)進(jìn)行去均值化處理并整合，獲得全部MP1的 RMSE 數(shù)值；然后基于衛(wèi)星尋優(yōu)進(jìn)行分類加權(quán)平均處理，以此得到更為準(zhǔn)確的MP1的 RMSE數(shù)值. NMRI的計(jì)算方法為

式中：MP1RMSE是單日MP1的RMSE值；max(MP1RMSE)是將MP1RMSE序列數(shù)值進(jìn)行由大到小排列后前5%的MP1RMSE的平均值，其結(jié)果均在0～1，不排除少數(shù)負(fù)值出現(xiàn)的可能. 與NDV1保持一致.

1.3 NDVI

NDVI是通過近紅外波段和紅光波段二者反射率的差值與和值進(jìn)行作比產(chǎn)生的用于反映地面綠色植被量的測(cè)度. 其計(jì)算公式為

式中， ρch2和 ρch1分別表示近紅外波段和紅光波段經(jīng)過誤差校正后的地面反射率. 紅光波段是植物將入射輻射吸收的波段，近紅外波段是植物將入射輻射強(qiáng)烈反射波段. NDVI的取值范圍是[-1，1]，負(fù)值部分表示觀測(cè)區(qū)域內(nèi)存在對(duì)可見光高反射性地物，如湖泊、冰雪層等，0的出現(xiàn)表示土地巖石直接裸露，正值部分?jǐn)?shù)值越大表示地面綠色植被覆蓋率越高.

2 實(shí)驗(yàn)與分析

2.1 數(shù)據(jù)來源

本次用于實(shí)驗(yàn)的數(shù)據(jù)來自美國(guó)大陸板塊觀測(cè)網(wǎng)(PBO)中P041站點(diǎn)2012—2016 年和P049站點(diǎn)2008—2012 年的GPS觀測(cè)數(shù)據(jù)，站點(diǎn)位于美國(guó)中西部草原地區(qū)，站點(diǎn)周圍地勢(shì)平坦，無大型障礙物遮擋和人為活動(dòng)因素干擾，利于GPS數(shù)據(jù)信息采集. 且站點(diǎn)周圍的植被有明顯的年際變化和物候特征，參考站的觀測(cè)環(huán)境如圖5所示.

圖5 GNSS 參考站

為保證實(shí)驗(yàn)條件的一致性，選擇使用相同接收機(jī)和天線的參考站進(jìn)行數(shù)據(jù)采集，其中P041和P049使用接收機(jī)和天線型號(hào)分別為TRIMBLE NETRS和TRM29659.00，觀測(cè)姿態(tài)為天線指向天頂(U)方向，采樣間隔設(shè)置為15 s. 參考站概略信息如表1 所示.

表1 GPS 參考站概況 (°)

本文NDVI數(shù)據(jù)來源于美國(guó)地質(zhì)調(diào)查局陸地過程分布式數(shù)據(jù)檔案中心，由于MODIS 產(chǎn)品數(shù)據(jù)MOD13Q1準(zhǔn)確度較高作為對(duì)比數(shù)據(jù)使用，此數(shù)據(jù)空間分辨率為250 m，時(shí)間分辨率為16 天.

2.2 數(shù)據(jù)處理

首先對(duì)測(cè)站GPS數(shù)據(jù)進(jìn)行處理，經(jīng)式(4)計(jì)算得出反射信號(hào)的歸一化微波指數(shù)NMRI，得出的2個(gè)參考站的NMRI時(shí)間序列繪制成圖，如圖6所示.

由圖6可知，2個(gè)站點(diǎn)的NMRI呈現(xiàn)明顯的周期性且數(shù)據(jù)波峰和波谷出現(xiàn)時(shí)間相似，這與植被生長(zhǎng)隨季節(jié)性變化表現(xiàn)出相同的特性. 同時(shí)2個(gè)站點(diǎn)均存在明顯的粗差點(diǎn)，這是由于GPS參考站在觀測(cè)過程中易受到降雨、降雪等天氣影響，導(dǎo)致觀測(cè)值與實(shí)際值產(chǎn)生偏差. 為削弱這種偏差，通過參照GPS氣象文件進(jìn)行手動(dòng)剔除降雨時(shí)間段觀測(cè)數(shù)據(jù)；為減少離群噪聲點(diǎn)影響，將觀測(cè)值與平均值之差的數(shù)據(jù)中大于2倍中誤差的觀測(cè)數(shù)據(jù)剔除，其結(jié)果如圖7所示.

圖6 參考站 NMRI

圖7 剔除粗差后的NMRI（藍(lán)色）與NDVI（綠色）

與圖4對(duì)比可知，NMRI由于環(huán)境因素引起的誤差離群點(diǎn)明顯減少，NDVI散點(diǎn)圖均出現(xiàn)在NMRI數(shù)據(jù)周圍且波峰和波谷出現(xiàn)時(shí)間吻合度較高，因此兩變量在時(shí)間上具有相似的周期規(guī)律. 由上述結(jié)果可初步判斷NMRI與NDVI間存在相關(guān)性，為實(shí)現(xiàn)NMRI反演NDVI增加了可能性，進(jìn)而進(jìn)行模型驗(yàn)證實(shí)驗(yàn).

2.3 模型驗(yàn)證

為驗(yàn)證該理論，選取80%的實(shí)驗(yàn)數(shù)據(jù)在 NMRI與 NDVI間分別建立一元線性回歸模型和BP神經(jīng)網(wǎng)絡(luò)反演模型，并使用剩余數(shù)據(jù)進(jìn)行模型驗(yàn)證. 回歸分析結(jié)果如表2所示，建模結(jié)果如圖8所示.

圖8 建模結(jié)果

表2 回歸分析結(jié)果

相關(guān)研究表明，來自數(shù)據(jù)量和參考站觀測(cè)系統(tǒng)處理反射信號(hào)的方式不同，因此P041和P049兩測(cè)站回歸系數(shù)不同，且接收機(jī)類型甚至固件的版本都會(huì)引起MP1RMSE數(shù)值的發(fā)生尺度變化[15-16]. 兩站點(diǎn)數(shù)據(jù)通過建立一元線性回歸模型得到的相關(guān)系數(shù)分別是0.744 28和0.734 96. 在顯著性檢驗(yàn)中，顯著性水平選擇0.05時(shí)，T檢驗(yàn)全部滿足精度要求. 在模型驗(yàn)證中，通過表2所示線性回歸模型將實(shí)驗(yàn)預(yù)留剩余NMRI時(shí)序數(shù)據(jù)進(jìn)行擬合得出相應(yīng)的NDVI值，再將遙感衛(wèi)星計(jì)算獲取的NDVI值作為地面驗(yàn)證的真實(shí)值使用，最終獲得每站NDVI的反演結(jié)果的RMSE分別為0.051 29、0.055 08.

為進(jìn)一步提高反演效果，引入3-8-1(三個(gè)輸入層八個(gè)隱藏層一個(gè)輸出層)的拓?fù)浣Y(jié)構(gòu)BP神經(jīng)網(wǎng)絡(luò)模型. 如圖9所示，增加植被背景濕度信息作為一個(gè)輸入項(xiàng)與NDVI、NMRI共同組成三個(gè)輸入層. 采用八個(gè)隱藏層，主要負(fù)責(zé)通過激活函數(shù)進(jìn)行信息處理與變換，后傳向輸出層各神經(jīng)元；最后將NDVI數(shù)值作為輸出層.

圖9 三層 BP神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)圖

其中實(shí)驗(yàn)數(shù)據(jù)與上述數(shù)據(jù)分配一致，將該實(shí)驗(yàn)結(jié)果與一元線性回歸模型對(duì)比得到P041和P049站點(diǎn)的相關(guān)系數(shù)分別提升2%和6%，RMSE分別為0.049 2和0.054 01，提升了反演的準(zhǔn)確性.

最后將NDVI 反演結(jié)果的相對(duì)誤差區(qū)間分布規(guī)律進(jìn)行統(tǒng)計(jì)分析，圖10為繪制的反演NDVI模型檢驗(yàn)的相對(duì)誤差概率統(tǒng)計(jì)圖.

由圖10可知，經(jīng)過模型預(yù)測(cè)的NDVI相對(duì)誤差值峰值部分均落在-0.1～0.1，精度較高；統(tǒng)計(jì)得出測(cè)站的反演相對(duì)誤差落在±0.2區(qū)間內(nèi)的數(shù)據(jù)量占比均達(dá)到70%以上，相對(duì)誤差集中于0附近，伴隨著相對(duì)誤差增大后，其所占比例越來越小，圖中顯示誤差分布總體上接近于正態(tài)分布，證明建立的NMRI反演歸一化植被指數(shù)模型是有實(shí)際意義的.

圖10 模型檢驗(yàn)的誤差統(tǒng)計(jì)

3 總 結(jié)

為克服AVHRR、MODIS等遙感數(shù)據(jù)產(chǎn)品存在的土壤濕度噪聲、觀測(cè)成本高及不易獲取高時(shí)間分辨率歸一化植被指數(shù)等缺陷，本文利用GPS-IR作為遙感信息源進(jìn)行NDVI反演，通過對(duì)大量連續(xù)時(shí)序?qū)崪y(cè)數(shù)據(jù)進(jìn)行分析，結(jié)果表明：

1) GPS多路徑數(shù)據(jù)產(chǎn)生的NMRI與NDVI時(shí)序數(shù)據(jù)具有相同的起伏形態(tài)且與實(shí)際氣候周期相關(guān)聯(lián).

2) 測(cè)區(qū)及周圍進(jìn)入雨水充足的季節(jié)時(shí)，植被量和水分含量增長(zhǎng)，反射信號(hào)振幅減小，NMRI增大. 同時(shí)雨雪天氣易產(chǎn)生的噪聲點(diǎn)容易分辨剔除.

3) NMRI 與 NDVI 在統(tǒng)計(jì)數(shù)據(jù)內(nèi)具有一定的相關(guān)性，在此基礎(chǔ)上建立反演NDVI的一元性回歸模型，相關(guān)系數(shù)分別是0.744 28、0.734 96，經(jīng)過模型驗(yàn)證反演出NDVI 相對(duì)誤差的RMSE為0.051 29、0.055 08；同時(shí)引入BP神經(jīng)網(wǎng)絡(luò)模型，進(jìn)一步提高反演相關(guān)性. 因此通過GPS-IR技術(shù)進(jìn)行NDVI反演結(jié)果具有可靠性和準(zhǔn)確性.

文中為減少?gòu)?fù)雜高矮植被、湖泊等外部環(huán)境因素影響，選擇開闊草地的測(cè)站點(diǎn)數(shù)據(jù)采集并進(jìn)行NDVI反演. 對(duì)于復(fù)雜環(huán)境反演效果及有效散射區(qū)進(jìn)行計(jì)算等問題還需在后續(xù)進(jìn)一步深入研究. 此研究結(jié)果為區(qū)域性NDVI數(shù)據(jù)采集提供了新的方式和思路.