基于CyGNSS數(shù)據(jù)的土壤水分與植被光學(xué)厚度反演研究

嚴清赟 金雙根,3 黃為民 賈燕 魏思遠

0 引言

深入探索水文和碳循環(huán)是對生物地球物理學(xué)的特征要素(例如土壤水分、地表生物多樣性等)進行研究的重點[1-2]．衛(wèi)星遙感因其能有效地提供大空間尺度的數(shù)據(jù),成為解決此類問題的一種有效方案．

由于土壤水分(Soil Moisture,SM)影響介電常數(shù),并且微波波段信號對其敏感[3],因此可以利用微波遙感技術(shù)獲取大范圍的SM數(shù)據(jù)．土壤水分主動-被動探測衛(wèi)星(SMAP)[4]和土壤水分和海洋鹽度探測衛(wèi)星(SMOS)[5]分別搭載了L波段輻射計,能夠提供基于亮度溫度數(shù)據(jù)的SM估計,空間分辨率約為40 km,每2～3天實現(xiàn)全球覆蓋．L波段信號對地表SM敏感,且受地表粗糙度和植被的影響較小,因此適用于SM的遙感[6]．除了被動式微波遙感接收器,主動遙感系統(tǒng),例如Sentinel-1[7]和ERRASAR-X[8]也能夠提供SM的測量．植被光學(xué)深度(VOD,或τ)作為生物多樣性的近似量,可以通過L波段SMAP和SMOS任務(wù)以及C波段平臺,例如高級微波輻射掃描計(AMSR-E)和高級甚高分辨率輻射計(AVHRR)獲得[9]．

近幾十年來,全球?qū)Ш叫l(wèi)星系統(tǒng)反射測量(GNSS-R)技術(shù)作為一種有效的遙感工具得到了廣泛的應(yīng)用．由于GNSS-R搭載的L波段信號源能夠輕易穿透云層和大氣,所以能夠?qū)崿F(xiàn)全時及全天候的遙感監(jiān)測[10-11]．GNSS-R已經(jīng)成功應(yīng)用于海面風監(jiān)測[12-14]、高程測量[15-16]、海冰遙感[17-19]、濕地分類[20]、 海嘯預(yù)警[21]等方面．隨著海量CyGNSS數(shù)據(jù)的公開,利用CyGNSS數(shù)據(jù)估算大尺度的SM的研究已成為學(xué)者們關(guān)注的熱點[22-31]．Clarizia等[24]提出了一種反射率-植被-粗糙度的算法，利用線性回歸函數(shù)實現(xiàn)SM的反演.Yan等[28]采用了類似方法,但計算了CyGNSS的觀測變量來解析表面粗糙度效應(yīng)．文獻[32-34]研究和評估了從CyGNSS數(shù)據(jù)提取植被相關(guān)參數(shù)的可行性．

值得注意的是,目前基于CyGNSS的研究,無論是估計SM還是推導(dǎo)τ,通常都依賴于其他輔助數(shù)據(jù),例如表面粗糙度、海拔高度等．本文提出了一種通過采用神經(jīng)網(wǎng)絡(luò)(Neural Network,NN)算法和暴力窮舉法,利用CyGNSS產(chǎn)品進行獨立的SM/τ反演的有效方案．基于此方案,后續(xù)的SM和τ數(shù)據(jù)能夠同時并直接地從CyGNSS數(shù)據(jù)估計出,且不依賴任何其他的輔助數(shù)據(jù)．本項工作推進了CyGNSS作為SM和τ數(shù)據(jù)的獨立來源的研究,并向公眾打開了一扇CyGNSS數(shù)據(jù)應(yīng)用新的窗口．

1 數(shù)據(jù)描述

本節(jié)首先描述了CyGNSS數(shù)據(jù)獲取與預(yù)處理過程,隨后介紹了使用的土壤水分SM與植被光學(xué)厚度VOD參考數(shù)據(jù)及來源．

1.1 CyGNSS數(shù)據(jù)

本文使用了CyGNSS 2.1版本的一級數(shù)據(jù)(參見https:∥podaac-tools.jpl.nasa.gov/drive/files/allData/cygnss/L1/v2.1)．CyGNSS星座包含了8顆微小衛(wèi)星,能提供覆蓋泛熱帶(38°S ～ 38°N)的GNSS-R數(shù)據(jù),該數(shù)據(jù)有高時空分辨率的特點．本文使用了2018年整年以及2020年內(nèi)有可用參考數(shù)據(jù)的時間段．CyGNSS一級數(shù)據(jù)包含了有地理位置信息的時延多普勒圖 (DDM)、雷達散射截面 (BRCS或σ)以及其他測量信息與地理坐標信息,如,入射角θ、SNR(信噪比)、鏡面反射點(SP)的經(jīng)緯度(lat/lon)以及SP到發(fā)射機與接收機的距離(Rt與Rr)等．本文使用了收集于陸地上且SNR大于0的CyGNSS數(shù)據(jù)．為了抑制數(shù)據(jù)中SP估計值的誤差,本文保留了BRCS峰值在時延軸第4與15位的數(shù)據(jù),類似操作也在文獻[26,28]中采用．另外,θ大于60°的數(shù)據(jù)由于質(zhì)量得不到保證而被剔除[25]．

1.2 參考數(shù)據(jù)

后續(xù)基于CyGNSS數(shù)據(jù)獲得的SM和τ結(jié)果會與SMAP第3級數(shù)據(jù)進行對比．SMAP的SM數(shù)據(jù)可以從https:∥nsidc.org/data/SPL3SMP/versions/6獲?。摂?shù)據(jù)的空間分辨率為36 km,大致每3天能實現(xiàn)全球覆蓋．SMAP數(shù)據(jù)提供了SM、VOD數(shù)據(jù)以及質(zhì)量標記等數(shù)據(jù)．本文使用了2018年全年以及2020年第240至291天質(zhì)量標記為“反演成功”的數(shù)據(jù)．為了方便后續(xù)的結(jié)果驗證與比對,對CyGNSS數(shù)據(jù)與SMAP產(chǎn)品都進行了3 d平均的操作,這樣能有效地提高空間覆蓋率．另外,本文將CyGNSS數(shù)據(jù)重新處理并投射到SMAP數(shù)據(jù)采用的EASE-Grid網(wǎng)格中．

2 估計方法

本節(jié)詳細敘述了從CyGNSS數(shù)據(jù)獨立反演SM和τ的過程,包括計算CyGNSS觀測變量、神經(jīng)網(wǎng)絡(luò)的構(gòu)建,以及通過暴力窮舉法進行結(jié)果計算3個步驟．

2.1 CyGNSS觀測量計算

通過假設(shè)植被覆蓋的光滑面的反射信號是相干的,CyGNSS的地表反射率ΓRL可由以下公式表示[26,35],即:

ΓRL(s,θ)=RRL(s,θ)γ2exp(-4k2s2cos2θ),

(1)

式中:θ是入射角;s表示土壤的介電常數(shù),通常取決于SM[3];透射率γ包含了信號由植被的衰減,通常是VOD的函數(shù);指數(shù)項代表了表面粗糙度效應(yīng),k為波數(shù),s是表面的均方根高度;表面的菲涅爾反射系數(shù)RRL在CyGNSS情況下可以表示為水平極化與垂直極化反射系數(shù)分量(Rhh與Rvv)的線性疊加[36]:

(2)

(3)

(4)

因此,在植被覆蓋的平整表面上的CyGNSSΓRL可認為是SM、θ、VOD以及表面粗糙度效應(yīng)的函數(shù)．實際中,基于相干反射的假設(shè),CyGNSSΓRL可直接從CyGNSS的BRCS(即σ)獲取[22,24,28,37]:

(5)

θ、σ、Rt與Rr均可從CyGNSS數(shù)據(jù)獲?。疚氖褂孟喔尚詷俗R(CF)來表征表面粗糙度效應(yīng),CF描述了CyGNSS數(shù)據(jù)的功率發(fā)散程度[38]．除此之外,經(jīng)緯度信息(lat/lon)也被用來作為輸入[32]．將CyGNSSΓRL、θ、CF以及l(fā)at/lon都基于EASE-Grid網(wǎng)格取逐3 d平均．作為例子,圖1展示了這5項以及SMAP的SM與VOD在2018年的平均值．

圖1 2018年CYGNSS與SMAP數(shù)據(jù)的平均值

2.2 神經(jīng)網(wǎng)絡(luò)的構(gòu)建

正如前文所提,本文假設(shè)了ΓRL、θ、CF(εCF)、lat/lon(εlat/εlon)、SM(ζSM)以及VOD(τ)之間的相互作用．自然而然地,它們之間的關(guān)系可以表達如下:

ζSM=f(ΓRL,θ,εCF,εlat,εlon,τ),

(6)

τ=g(ΓRL,θ,εCF,εlat,εlon,ζSM).

(7)

獲取這兩個公式嚴格的數(shù)學(xué)表達式有較大難度．這里,借助機器學(xué)習(xí)算法來構(gòu)建函數(shù)f與g．神經(jīng)網(wǎng)絡(luò)算法已成功地在GNSS-R遙感領(lǐng)域證明了其有效性與準確性,例如海冰監(jiān)測[39]、風速反演[14]以及土壤/植被遙感[26-27,32]．本文也采用神經(jīng)網(wǎng)絡(luò)算法．

考慮有n個樣本的訓(xùn)練集T={(x1,y1),…,(xn,yn)},其中xj是輸入向量(包含ΓRL、θ、CF、lat/lon與τ,且j=1,…,n),而yj代表了式(6)中的SM值．在這里,神經(jīng)網(wǎng)絡(luò)算法可以看作用來構(gòu)建函數(shù)f,從而實現(xiàn)f(x)高精度地近似SM(或y)．同理,為了估算VOD,可以參照式(7)．圖2給出了神經(jīng)網(wǎng)絡(luò)的結(jié)構(gòu)示意,由1個輸入層、3個隱藏層以及1個輸出層構(gòu)成．每一層的神經(jīng)元個數(shù)可由使用者任意設(shè)定,相鄰層的神經(jīng)元都由權(quán)重相連接,除輸入層外神經(jīng)元的激勵過程總結(jié)如下:

a(k+1)=φ(W(k+1)a(k)+b(k+1)),k∈[0,1,2,3],

(8)

圖2 神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)示意

式中a(k+1)、b(k+1)與W(k+1)分別是第(k+1)層的激勵向量、偏差向量以及權(quán)重矩陣．φ是激勵函數(shù)．另一方面,輸入層的表達式如下:

a(0)=x.

(9)

式(8)與(9)總結(jié)了神經(jīng)網(wǎng)絡(luò)的前向傳播路徑．一旦網(wǎng)絡(luò)完成學(xué)習(xí)過程,即權(quán)重與偏差都通過訓(xùn)練而確定,便可以通過式(6)與(7)的形式來估計SM與VOD．在網(wǎng)絡(luò)訓(xùn)練階段,采用后向傳播學(xué)習(xí)[40]以及LM算法[41],它們的本質(zhì)在于,在給定訓(xùn)練數(shù)據(jù)集T時,通過調(diào)整權(quán)重與偏差值將以下的誤差函數(shù)最小化:

(10)

總地來說,本文使用了神經(jīng)網(wǎng)絡(luò)算法來構(gòu)建f與g,進而通過式(6)與(7)的形式實現(xiàn)對SM與VOD的估計．詳細的神經(jīng)網(wǎng)絡(luò)算法使用過程可參照文獻[14,39]．本文對神經(jīng)網(wǎng)絡(luò)的構(gòu)建、訓(xùn)練以及后續(xù)分析都基于MATLAB R2019b軟件．

本研究通過CiteSpace對桂醫(yī)2007至2016年相關(guān)文獻的作者、關(guān)鍵詞與突現(xiàn)詞進行知識圖譜分析以及對桂醫(yī)年度發(fā)文量、學(xué)科發(fā)文量進行文獻統(tǒng)計分析，得到以下結(jié)論：

2.3 獨立反演方案

當f與g構(gòu)建完成后,從后續(xù)的CyGNSS數(shù)據(jù)反演SM/VOD即可通過暴力窮舉法實現(xiàn),其詳細步驟如下．首先,將式(6)中的τ替換成式(7),可以得到:

ζSM=f[ΓRL,θ,εCF,εlat,εlon,

g(ΓRL,θ,εCF,εlat,εlon,ζSM)],

(11)

圖3 獨立反演SM與VOD的方案流程

或以一種更簡潔的形式,即:

ζSM=F(ΓRL,θ,εCF,εlat,εlon,ζSM),

(12)

或者寫成:

ζSM-F(ΓRL,θ,εCF,εlat,εlon,ζSM)=0.

(13)

3 結(jié)果與評價

本節(jié)對提出的SM和τ反演方案進行了測試,并且對2018與2020年的數(shù)據(jù)進行了評估．根據(jù)SMAP數(shù)據(jù)的可用性,本文選取了2018全年以及2020年第240天至291天的數(shù)據(jù)．基于第2節(jié)描述的數(shù)據(jù)匹配方案,分別獲得了2018年1 620多萬和2020年230多萬個的數(shù)據(jù)樣本．

3.1 神經(jīng)網(wǎng)絡(luò)的訓(xùn)練

隨機選擇了2018年數(shù)據(jù)中的1%進行神經(jīng)網(wǎng)絡(luò)的訓(xùn)練,即以f的形式(式(6))從ΓRL、θ、CF、lat、lon以及VOD獲取SM,以g的形式(式(7))獲取VOD．另一方面,剩下的數(shù)據(jù),包括2018年剩余的99%和2020年的數(shù)據(jù),都將用作最終的反演驗證．

神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)的設(shè)計(如層數(shù)或神經(jīng)元的數(shù)目,以及激勵函數(shù)的選擇)可以是任意的,并可能導(dǎo)致網(wǎng)絡(luò)的精確度不同．本文測試了大量不同結(jié)構(gòu)的神經(jīng)網(wǎng)絡(luò),來得到“最優(yōu)”的網(wǎng)絡(luò),將隱藏層數(shù)設(shè)置為1、2、3,神經(jīng)元數(shù)為1～50,并將所有被測網(wǎng)絡(luò)的激勵函數(shù)設(shè)置為“tansig”[29],其形式如下:

(14)

通過對測試的神經(jīng)網(wǎng)絡(luò)進行準確度評價,最終選擇了有3個隱藏層,每個隱藏層有26個神經(jīng)元組成的網(wǎng)絡(luò)來進行對SM的估計,即f;另外,有3個隱藏層,每個隱藏層有38個神經(jīng)元構(gòu)成的網(wǎng)絡(luò)被用作為τ的預(yù)測器．圖4中顯示了這2種網(wǎng)絡(luò)性能的密度．目前為止,只有1%的2018年數(shù)據(jù)運用到了神經(jīng)網(wǎng)絡(luò)的訓(xùn)練中．

圖4 數(shù)據(jù)密度圖與性能

3.2 SM與VOD的反演

確定好神經(jīng)網(wǎng)絡(luò)f和g,按照2.3節(jié)的步驟對SM和τ進行反演．根據(jù)訓(xùn)練數(shù)據(jù)的范圍和精度,設(shè)置SM的搜索范圍為0.02 ～ 0.88 cm3/cm3,步長為10-4cm3/cm3．值得注意的是,用于反演SM/τ的輸入樣本只需要CyGNSS的數(shù)據(jù),包括ΓRL、θ、CF、lat與lon．除去之前用于神經(jīng)網(wǎng)絡(luò)訓(xùn)練的數(shù)據(jù)之外,2018年和2020年剩余的所有測試數(shù)據(jù)都用于反演結(jié)果的評估．在實際操作中,SM是根據(jù)計算出的差值Δ=(0.02∶10-4∶0.88)-F[ΓRL,θ,εCF,εlat,εlon,(0.02∶10-4∶0.88)]進行評估的,結(jié)果中對應(yīng)min(|Δ|)的值被視作SM．τ則通過式(7)計算得出．當Δ中多次出現(xiàn)交叉零點時(說明有多解),則認為反演無效．這可能是由網(wǎng)絡(luò)固有誤差引起的．盡管如此,僅約6%的數(shù)據(jù)因為這個問題而被舍棄．

上述反演過程基于式(11)—(13)的推導(dǎo).具體來說,是將變量τ替換了式(6)中的g．在現(xiàn)實中,這個過程可以采用另一種方式,即將式(7)中的f由SM代替,那么可得到與式(13)對等的形式:

τ-G(ΓRL,θ,εCF,εlat,εlon,τ)=0．

(15)

顯然,也可以使用式(15)和上述的類似方法,對SM和τ進行新一輪的反演．為提高結(jié)果的可靠性,本文利用式(15)也進行了反演,為τ設(shè)置的搜索范圍為[0,1.297 0],步長為10-4．只有當基于式(13)和(15)得到的反演結(jié)果一致時,才認為最終的反演結(jié)果是有效的．為了檢驗并量化結(jié)果一致性的影響,我們測試了用以篩選反演結(jié)果的不同閾值,范圍為[10-5,10-1]．如圖5所示,當兩種方案的反演結(jié)果吻合度越高,所得到的最終結(jié)果與參考數(shù)據(jù)也更一致．考慮到所獲得的準確結(jié)果,本文將采用10-5作為閾值進行更深層次的評價和討論．

圖5 不同反演方案所得估計值間的一致性對最終結(jié)果的影響

3.3 驗證與評估

如圖6所示,基于訓(xùn)練數(shù)據(jù)與測試數(shù)據(jù)的反演結(jié)果均展現(xiàn)了高準確度．相關(guān)系數(shù)r和均方根誤差RMSE匯總在表1中．結(jié)果表明,本文所提出的反演方案具有較好的通用性和可靠性．特別地,盡管測試數(shù)據(jù)集(2020年)與訓(xùn)練數(shù)據(jù)集(2018年)存在一定的時間差距,但仍能利用2020年數(shù)據(jù)成功地進行反演,驗證了該方案的有效性．并且,訓(xùn)練數(shù)據(jù)和測試數(shù)據(jù)的密度之間的高度一致性也是顯而易見的．

由于2018年只有1%的數(shù)據(jù)參與了神經(jīng)網(wǎng)絡(luò)的訓(xùn)練,所以測試集的結(jié)果相對于訓(xùn)練集的結(jié)果有一個不顯著下降是合理且可接受的．2018年和2020年數(shù)據(jù)的反演性能的差異可能是由于總體的數(shù)據(jù)質(zhì)量存在著潛在差異．值得注意的是,2020年CyGNSS數(shù)據(jù)的非相干求和時間(0.5 s)比2018年的(1 s)短,因此前者的空間采樣率是后者的兩倍．這對整體的數(shù)據(jù)質(zhì)量帶來了兩種相反的影響．一方面,SNR與非相干求和時間的平方根成正比[42],這使得2018年CyGNSS的個體數(shù)據(jù)都具有更高的SNR；另一方面,2020年的數(shù)據(jù)具有更高的空間采樣率,這將導(dǎo)致在每個Ease-Grid網(wǎng)格中有更多的樣本,從而提高了精度[28]．因此,2018年和2020年網(wǎng)格化后的數(shù)據(jù)質(zhì)量可能存在一定的差異．

圖6 最終SM與VOD反演結(jié)果的密度

表1 最終反演的性能統(tǒng)計

圖7 SM結(jié)果一致性比對

下面進一步驗證SMAP和CyGNSS數(shù)據(jù)產(chǎn)品時空變化響應(yīng)的一致性．因為2018年的數(shù)據(jù)量更大,這里展示了基于2018年的數(shù)據(jù)結(jié)果．如圖7所示,以年平均的形式呈現(xiàn)了SMAP和CyGNSS結(jié)果之間的空間相關(guān)性(SM見圖7a和7b,τ見圖8a及8b)，并且通過比對SM和τ時間序列的標準差(Stds),證實了時間尺度上的一致性 (圖7、圖8c、8d)．SMAP與CyGNSS數(shù)據(jù)結(jié)果間的差異很?。?/p>

圖8 VOD結(jié)果一致性比對

4 結(jié)束語

本文提出了一種基于CyGNSS數(shù)據(jù)同時反演土壤水分和植被光學(xué)厚度的方法．這種方法結(jié)合了神經(jīng)網(wǎng)絡(luò)技術(shù)與暴力窮舉算法．設(shè)計的輸入包括由CyGNSS導(dǎo)出的反射率、相干性標識、入射角以及經(jīng)緯度．兩種不同的變量替換方式形成了兩種不同的反演模型(具體形式見式(13)與(15)．為了增強該方案的可靠性,只有當兩種模型的反演結(jié)果一致時,我們才認為反演結(jié)果是有效的,并評估了其對結(jié)果準確性的影響．通過最終反演結(jié)果與SMAP數(shù)據(jù)的對比驗證,發(fā)現(xiàn)反演結(jié)果與參考數(shù)據(jù)展現(xiàn)了良好的一致性,驗證了該反演模型的有效性和準確性．針對2018年與2020年的土壤水分比對結(jié)果,得到的相關(guān)系數(shù)分別達0.86與0.84,均方根誤差分別為0.064和0.071 cm3/cm3．對于植被光學(xué)厚度,2018年與2020年的結(jié)果均得到0.98的相關(guān)系數(shù),均方根誤差分別為0.079和0.084．本文進一步從時空變化的角度驗證了反演結(jié)果的可靠性．今后,在有同時同地參考數(shù)據(jù)的情況下,將使用統(tǒng)一非相干求和時間的CyGNSS數(shù)據(jù)以便消除整體數(shù)據(jù)質(zhì)量的不確定性．如何確保神經(jīng)網(wǎng)絡(luò)系統(tǒng)的效率、適用性和精確性是至關(guān)重要的,值得進一步研究和改進．此外,本文神經(jīng)網(wǎng)絡(luò)是基于SMAP數(shù)據(jù)進行調(diào)參的,SMAP數(shù)據(jù)也可以被其他具有更高分辨率的數(shù)據(jù)或現(xiàn)場測量數(shù)據(jù)所替代．并且,本文采用暴力窮舉法獲得的最終反演結(jié)果,其初始值的設(shè)置與步長精度與最終性能高度相關(guān)．因此,擴大搜索范圍和細化步長可以進一步改善結(jié)果．暴力窮舉法的替代方法也將在未來研究中進一步探索．