孫時雨,宋承運,周 露
(1.安徽理工大學 空間信息與測繪工程學院,安徽 淮南 232001;2.安徽理工大學 礦山采動災害空天地協同監測與預警安徽普通高校重點實驗室,安徽 淮南 232001)
土壤水分是能量、水和碳循環的重要組成部分,為陸地-大氣模型、水文模型和降水預測模型提供了重要信息[1]。在多尺度相互作用中,土壤水分在控制地表能量分配、調節地表徑流和土壤排水、調節冠層蒸騰和碳排放等方面發揮著重要作用[2]。長時序的土壤水分產品可用于更準確和可靠地估計深層土壤水分和蒸散發等,也可用于數據同化、天氣預報和干旱監測等。因此,迫切需要可靠、準確和完整的土壤水分信息[3]。
由于土壤水分空間異質性強,用來觀測土壤水分的設備成本高,不僅耗費人力,還耗費物力,因此很難建立區域和全球范圍的高分辨率觀測網絡。傳統的土壤水分觀測方法一般只適合范圍比較小的土壤水分的信息采集,無法應用到大范圍的地表水分的測量。隨著遙感技術的快速發展,衛星遙感是獲取地表土壤水分(通常在0~5 cm深度)最有效的工具[4]。特別是微波遙感被認為是監測土壤水分最主流的方法,因為它能夠全天工作,不受天氣的限制,可穿透云層、植被和土壤[5]。各種衛星研究和應用機構的微波傳感器可以生成全球遙感產品,如高級微波掃描輻射計EOS(AMSR-E)、土壤水分主-被動遙感衛星(SMAP)、歐洲航天局(ESA)的土壤濕度和海洋鹽度觀測衛星(SMOS)、日本的先進微波掃描輻射計2(AMSR2)和中國的微波輻射成像儀(MWRI)[6]。然而上述衛星產品的時間覆蓋范圍太短,無法滿足水文氣象的應用需求。因此, ESA CCI發布了一個完整和一致的全球土壤水分數據集。但由于衛星軌道變化、射頻干擾(RFI)和微型傳感器的物理限制,該數據集存在數據缺失[7]。數據缺失限制了ESA CCI 土壤水分產品的實際應用,特別是在需要空間完整和時間連續的土壤水分數據的領域,如氣候模擬、干旱監測和水資源管理[8]等。
國內外已經提出了一些方法來填補衛星數據集的空白。這些方法主要可分為統計、插值和機器學習等幾類。統計方法(例如多元線性回歸、非線性回歸)主要依賴輸入和目標數據的線性和非線性關系。Wang等[9]測試了一種基于三維離散余弦變換(DCT-PLS)的最小二乘回歸的統計建模方法來填補衛星產品空白。Yan等[10]運用了反距離加權(IDW)根據目標變量的空間模式填充缺失數據。其他方法已用于重建,例如Turlapaty等[11]使用改進的奇異光譜分析來填補AMSR-E土壤水分數據集中的缺失。近年來,支持向量機(SVM)和人工神經網絡(ANN)被廣泛用于填補遙感數據的缺失[12]。然而,基于機器學習算法模型的優劣很大程度上取決于輸入變量的選擇以及輸入變量和目標變量的關系,模型的結構與輸入變量、目標變量和周期高度相關[13]。收斂速度慢、易過擬合是傳統的機器學習模型及BP神經網絡的缺點,而且得到數據的精度誤差大[14]。
基于以上分析,以青藏高原那曲為研究區域,利用ESA CCI V04.2組合土壤水分產品、數字高程數據(DEM)以及基于MODIS數據重建的地表溫度數據、歸一化植被指數(NDVI)和反照率數據,采用原子搜索優化(Atom Search Optimization,ASO)算法去優化BP神經網絡模型的權值、閾值,從而建立ASO-BP土壤水分重建模型,得到空間連續土壤水分,并利用青藏高原中部土壤水分多尺度觀測網(Soil Moisture and Temperature Monitoring Network on the Central TP,TP-SMTMN)數據對結果進行驗證分析。
研究區域為青藏高原那曲地區(83°55′E~95°05′E,29°55′N~36°30′N)。那曲約有40萬km2的面積,地勢開闊,干濕季分明,屬于高原亞寒帶氣候區,月平均降水量在7月最大,在12月或1月最小,年平均降水量為476.5 mm[15]。
那曲地區在3月下旬凍土開始消融,土壤水分增加。從 5月開始表層土壤水分會繼續增加,受季風降水的影響,7月達到峰值,8月以后降水減小,土壤水分逐漸減小,10月溫度降低,土壤水分逐漸由液態變為固態,土壤水分減小,1月達到谷值,一直持續到 3月,完成一年的循環周期[16]。研究區高程及地面觀測站點示意如圖1所示。

(a)研究區高程

(b)地面觀測站點
(1)ESA CCI土壤水分數據。ESA CCI是一種長時間序列、多衛星融合的土壤水分產品,它基于主動和被動微波傳感器,生產的數據集包含主動、被動和融合數據集3種,空間分辨率為0.25°,時間分辨率為1 d[17]。它也是第一個基于遙感的長時間序列土壤水分產品,具有很好的一致性,而衛星軌道、植被、凍土、積雪和冰川等因素影響了反演后土壤水分的空間連續性與精度[18]。本文采用ECV V04.2組合土壤水分產品(每日,0.25°),作為主要數據源 (http:∥www.esa-soilmosture-cci.org/)[19]。
(2)LST、NDVI、Albedo與DEM遙感數據。LST、NDVI和Albedo數據是由MODIS MOD11C1 LST產品(0.05°)、MOD13C1 NDVI產品(0.05°)與MOD43C1反照率產品(0.05°)分別采用多時相穩健回歸方法[20]和基于時間濾波統計方法[21]進行數據重建,得到空間連續的LST、NDVI和Albedo數據,分辨率為0.05°(數據來自 http:∥reverb.echo.nasa.gov/)。美國國家航空航天局(National Aeronautics and Space Administration,NASA)航天飛機雷達地形任務(SRTM)公開的DEM(30 m)(數據來自 http:∥reverb.echo.nasa.gov/DEM)。
(3)地面觀測數據。由研究區內青藏高原中部土壤水分多尺度觀測網的數據集提供。此數據網提供3種尺度和4種深度的土壤溫度與土壤水分的觀測數據,3種尺度分別為大尺度1°、中尺度0.3°與小尺度0.1°,4種深度分別為 0~5、10、20 、40 cm,觀測間隔為0.5 h。由于被動微波土壤水分產品的空間分辨率為0.25°,選取中尺度(如圖1(b)所示) 范圍內深度為 0~5 cm的土壤水分觀測數據去驗證分析結果[22]。
BP神經網絡是一種多層前饋網絡,其核心算法為誤差反向傳播,在函數逼近、模式識別、分類等領域應用較為廣泛,具有自組織、自學習能力等特點[23]。
基本原理:基本BP算法計算誤差輸出時是按從輸入到輸出的方向執行的,而調整權值和閾值則相反。通過調整輸入節點與隱層節點、隱層節點與輸出節點的聯接強度以及閾值,從而使BP網絡的預測輸入不斷接近期望輸出。BP神經網絡還可以用來建立多元非線性關系,因此,在遙感反演估算中被廣泛應用[24]。BP神經網絡結構如圖2所示。

圖2 BP神經網絡結構Fig.2 BP neural network structure
ASO算法是一種基于分子動力學模型的新穎智能算法,它通過種群間各個原子間的相互作用力來指導群體,從而進行智能優化搜索[25]。
基本原理:根據分子動力學中原子的物理運動規律建立算法模型。由牛頓第二定律可知,如果Fi是作用在第i個原子上的相互作用力,Gi為作用在第i個原子上的約束力,原子質量為mi,第i個原子的加速度可表示為:
(1)
式中:t為當前迭代的次數,d為原子所在的維數。下面對Fi、Gi、mi三個值的求解進行闡述。
(2)
式中:Kbest為適應度函數值比較好的k個原子的集合。
(3)

(4)

(5)
(6)
式中:xi和xj分別表示原子i和原子j的位置。
(7)

(3)原子質量mi(t)。第t次迭代中第i個原子的質量mi(t)由當前種群個體的適應度值大小決定,定義如下:
(8)
(9)
式中:fi(t)為第i個原子的適應度值,fmax(t)、fmin(t)分別表示原子種群中最大適應度值和最小適應度值。
(10)
在每一次迭代過程中,根據得到的加速度更新原子i的速度和位置,更新如下:
(11)
(12)

本文應用一種基于ASO算法的反演模型,其關鍵在于通過優化原子位置來確定BP神經網絡的初始權值和閾值。首先利用ASO算法對原子位置進行多次迭代搜索,以找到最優的參數組合。再利用優化后的BP神經網絡模型進行土壤水分的反演。算法的流程如圖3所示,其算法流程分為3步:
① 根據遙感影像選取所需的輸入數據,將其進行重建和均值升尺度處理后的數據作為整體的樣本集。在確定BP神經網絡的拓撲結構之后對BP神經網絡的權值與閾值進行初始化。
② 計算ASO算法的決策變量長度,選取均方誤差作為優化的目標函數,對BP神經網絡的權值與閾值進行優化。ASO計算原子適應度,更新目前找到的最優解。根據式(8)、式(9)分別計算原子的個體加速度、速度及位置更新。設置當達到最優適應度或者最大迭代次數時終止迭代的算法停止準則,輸出原子個體的最優位置,否則繼續重復ASO。
③ 使用優化算法對神經網絡的權值和閾值參數進行優化,將優化后得到的權值和閾值參數賦給BP神經網絡。輸入樣本集與測試集進行反演,數據歸一化后最終獲得土壤水分。

圖3 ASO BP神經網絡土壤水分重建流程Fig.3 Flowchart of ASO BP neural network soil moisture reconstruction
由NDVI、LST、Albedo與SRTM DEM數據采用均值升尺度后得到空間分辨率為0.25°的數據,并與ECV V04.2作為輸入數據,優化后的土壤水分數據作為輸出,建立BP神經網絡。BP神經網絡中采用了正切S型傳遞函數tansig作為傳輸函數,最大的迭代次數為1 000,目標0.000 001。采用ASO后的BP神經網絡訓練及土壤水分反演。
為了驗證土壤水分反演的有效性,評價傳統BP神經網絡模型和ASO-BP神經網絡模型的反演精度,選取了3個指標:相關系數(R)、均方根誤差(RMSE)和無偏均方根誤差(ubRMSE)。一般來說,R值越大,RMSE、ubRMSE越小,則反演結果越好、精度越高。
由ASO-BP神經網絡模型、BP神經網絡模型重建的2015年5月1日土壤水分如圖4所示。由圖4可以看出,重建后的土壤水分與重建前ESA CCI土壤水分空間分布一致,中部和西部土壤水分值較小,東部值較大。通過重建,得到了空間連續的土壤水分。

(a)ESA CCI土壤水分圖

(b)BP重建土壤水分圖

(c)ASO-BP重建土壤水分圖
2015年1—12月ESA CCI土壤水分、重建的土壤水分和地觀測站土壤水分的時間序列如圖5所示。由圖5可以看出,重建后的土壤水分在對ESA CCI土壤水分中缺失的土壤水分進行重建生成有效值的同時,與實測土壤水分保持一致的變化趨勢。夏季土壤水分較高,冬季較低,在8月左右有一個明顯的低谷,ESA CCI土壤水分與重建土壤水分在此時間段也表現出明顯的低谷。但重建后的土壤水分值在4—6月以及10月低于實測值,主要是由于原始ESA CCI土壤水分在此期間低估了土壤水分,造成重建后的土壤水分值偏低。同時,彌補了在ESA CCI土壤水分1—3月以及11—12月缺少的數據,但與地面實測值相關較大,這主要是由于此時段為凍土,且在重建模型中,缺少相應時期的土壤水分數據作為輸入數據。因此,提高人工神經網絡輸入數據的質量是進一步提高反演精度的關鍵步驟。

圖5 2015年土壤水分的時間序列比較Fig.5 Time series comparison of soil moisture in 2015
選取地面觀測網中尺度(0.3°)站點觀測平均值與反演得到的土壤水分進行對比分析,圖6為BP重建土壤水分(圖6(a))和ASO-BP重建土壤水分(圖6(b))與地面觀測值的散點圖。 由圖6可以看出,在研究期內, ASO-BP重建土壤水分值不高于地面觀測值,土壤水分在大部分日期低于地面觀測值。而在土壤較為干燥時,尤其是在土壤水分值小于 0.2 cm3·cm-3時,土壤水分重建值則高于地面觀測值。這可能是由于土壤濕潤時,多處于季風期,降雨較多,植被生長較為茂盛,采用被動微波反演土壤水分時植被透過率估算的誤差,造成重建得到的 ESA CCI土壤水分低于地面觀測0~5 cm土壤水分觀測值,而ESA CCI土壤水分在大部分日期偏低也是造成BP重建土壤水分與ASO-BP重建土壤水分值偏低的一個重要原因。相對于BP重建土壤水分,ASO-BP重建土壤水分在土壤較為干燥時,與地面觀測值更為接近,主要原因是ASO對BP神經網絡初始值的優化,減小了干燥時反演土壤水分的誤差。
同時,ASO-BP重建土壤水分與地面觀測值的RMSE值(0.029 cm3·cm-3),ubRMSE值為0.028 cm3·cm-3,略高于BP 重建土壤水分(0.034 cm3·cm-3),ubRMSE值為0.033 cm3·cm-3,整體精度得到提高。
整體而言,重建后的土壤水分的變化趨勢與站點實測數據基本保持一致,可以較好地監測土壤水分變化情況。
為了進一步分析不同植被覆蓋度下ASO-BP土壤水分精度,選取NDVI值0.4作為不同植被覆蓋度閾值[26],分別分析NDVI≥0.4與NDVI<0.4時土壤水分精度。表1NDVI<0.4與NDVI≥0.4時,土壤水分與地面觀測值的精度比較。

表1 不同土壤水分值與地面觀測值的精度比較
由表1可以看出,在高植被覆蓋度下(NDVI≥0.4),ASO-BP土壤水分與BP土壤水分的R值較接近,RMSE值分別為0.031、0.026 cm3·cm-3, ubRMSE從0.024 cm3·cm-3降到0.022 cm3·cm-3。在植被覆蓋區域,植被在一定程度上降低了微波信號對土壤水分的敏感性,造成微波反演土壤水分時產生誤差。而BP與ASO-BP算法能夠在一定程度上對植被覆蓋度較高區域的土壤水分進行優化,提高精度。當NDVI<0.4時, ASO-BP算法相關系數R值(0.58),高于BP重建土壤水分(0.39), ubRMSE從0.033 cm3·cm-3降到0.028 cm3·cm-3。
因此,ASO-BP優化模型能夠提高土壤水分精度,并且在低植被覆蓋度區域(NDVI<0.4)時的效果比BP神經網絡模型的效果要更好。
本文使用ASO算法優化BP神經網絡,融合了NDVI、LST、DEM、Albedo、ESA CCI土壤水分數據,建立了土壤水分重建模型。實驗結果表明,本文算法能夠實現土壤水分遙感產品的重建,同時,提高了土壤水分反演精度,得出主要結論如下:
① 通過ASO算法對BP神經網絡進行優化,改善了BP神經網絡算法收斂速度慢、易過擬合的缺點。在算法中,采用ASO對BP神經網絡的初始值進行優化,使得土壤水分重建精度較傳統BP算法獲得有效提高。
② 在不同植被覆蓋度下,體現了算法的適用性。高植被覆蓋度區域(NDVI≥0.4),ASO-BP算法優于BP神經網絡算法。ASO-BP土壤水分和BP土壤水分與地面觀測值的R值從0.79提高到0.80以上, ubRMSE從0.024 cm3·cm-3提高到0.022 cm3·cm-3。在低植被覆蓋度區域(NDVI<0.4),ASO-BP優化模型土壤水分精度(R=0.58,RMSE=0.029 cm3·cm-3)比傳統BP神經網絡精度(R=0.39,RMSE=0.034 cm3·cm-3)提高更為明顯。
本文中土壤水分空間分辨率低且時間序列較短,而影響土壤水分的因素復雜多樣,在一定程度上會影響土壤水分結果精度驗證。同時,在結果驗證方面,采用像元內地面觀測站點的平均值作為真實值,但易受站點的布設位置、區域環境等因素的影響。在進一步研究中,將分析降水、土壤類型、土壤表面粗糙度等參數對模型的影響,以及土壤水分真實性驗證等。