陳明，王永前，吳錫

(1.成都信息工程大學(xué) 資源環(huán)境學(xué)院，成都 610225；2.成都信息工程大學(xué) 計(jì)算機(jī)學(xué)院，成都 610225)

0 引言

大氣中的水分僅占全球總水量的極小比例，但在水文循環(huán)過程中起著關(guān)鍵作用，與大氣環(huán)流和氣候變化關(guān)系密切[1]。水汽是影響地球能量循壞的溫室氣體之一，對(duì)全球氣候變暖有著顯著的作用[2]。青藏高原地區(qū)上空的大氣水汽分布以及變化與亞洲季風(fēng)、東亞大氣環(huán)流以及全球氣候變化息息相關(guān)[3]。因此，準(zhǔn)確地獲取青藏高原大氣水汽含量及其分布對(duì)于氣候變化的相關(guān)研究具有重要意義。

大氣水汽的獲取方法有很多，主要有探空法、地基法和衛(wèi)星遙感探測法[4]。探空法和地基探測的水汽精度普遍較高，但是只能進(jìn)行點(diǎn)觀測，且數(shù)量有限[5]。利用不同時(shí)空分辨率的衛(wèi)星遙感數(shù)據(jù)反演大氣水汽可以克服探空法和地基探測水汽的不足。目前，根據(jù)所利用遙感數(shù)據(jù)波段的不同，衛(wèi)星遙感反演大氣水汽可以分為熱紅外、可見光近紅外和微波[6]。一般在晴空條件下，利用近紅外或熱紅外觀測的遙感數(shù)據(jù)可以反演出準(zhǔn)確的大氣水汽含量[7-8]。近紅外及紅外遙感容易受云、霧和其他天氣狀況的影響，因此它的觀測能力受到限制。微波波段可以在減少大氣約束的情況下估算大氣水汽，能夠獲得有云條件下的大氣水汽[9]。因此，如何利用微波波段反演大氣水汽算法具有重要的意義。

水汽在微波波段中(1～300 GHz)有兩條轉(zhuǎn)動(dòng)吸收譜線，一條在22.235 GHz，另一條在183.31 GHz[10]。因此，利用微波進(jìn)行大氣水汽的反演，也多圍繞著這兩個(gè)波段展開。Wang等[11]提出了一種利用AMSR-E的18.7 GHz和23.7 GHz亮溫極化差反演大氣水汽的參數(shù)化模型。經(jīng)驗(yàn)反演方法是直接建立遙感信號(hào)與變量之間的關(guān)系，通過增加輸入輸出數(shù)據(jù)對(duì)，提高估算的穩(wěn)定性和精度[12]。新的發(fā)展趨勢是大量使用機(jī)器學(xué)習(xí)方法，用于訓(xùn)練的輸入輸出數(shù)據(jù)多來自于測量數(shù)據(jù)或者是基于物理模型模擬。Zabolotskikh等[13]利用AMSR2數(shù)據(jù)采用神經(jīng)網(wǎng)絡(luò)方法反演海洋上空的大氣水汽含量，具有較好的效果。Lee等[14]利用三種機(jī)器學(xué)習(xí)方法(隨機(jī)森林、神經(jīng)網(wǎng)絡(luò)、梯度提升)，基于葵花八號(hào)衛(wèi)星的熱紅外波段數(shù)據(jù)進(jìn)行大氣水汽反演，并基于歐洲中期天氣預(yù)報(bào)中心的水汽數(shù)據(jù)對(duì)反演結(jié)果進(jìn)行驗(yàn)證，三種方法效果較好。多種機(jī)器學(xué)習(xí)方法的反演結(jié)果較好，但是基于AMSR2數(shù)據(jù)采用機(jī)器學(xué)習(xí)方法反演大氣水汽的研究較少。

青藏高原地理位置特殊，是我國的水汽輸送中轉(zhuǎn)站[15]。衛(wèi)星遙感反演水汽主要在晴空條件下，對(duì)于有云條件下大氣水汽研究較少。本文基于AMSR2數(shù)據(jù)，利用隨機(jī)森林模型，在青藏高原地區(qū)反演有云和晴空條件下的大氣水汽，并利用ERA-Interim再分析數(shù)據(jù)對(duì)其進(jìn)行對(duì)比分析。本文為青藏高原及其周邊地區(qū)大氣水汽的獲取提供方法支撐。

1 研究區(qū)與數(shù)據(jù)源

1.1 研究區(qū)概況

青藏高原位于26°00′N～39°47′N，73°19′E～104°47′E之間，一般海拔在3 000～5 000 m之間，平均海拔超過4 000 m，不僅是我國長江、黃河的發(fā)源地，也是東亞、東南亞和南亞許多大河流的發(fā)源地，被稱為“世界屋脊”[16]。青藏高原總面積約為2.5×106km2，擁有獨(dú)特的地形地貌，地表類型較為多樣化。近年來青藏高原各地區(qū)上空云量較多，尤其是青藏高原東南部地區(qū)?；?016—2021年的MODIS云產(chǎn)品數(shù)據(jù)得到青藏高原地區(qū)云出現(xiàn)概率分布圖。可以發(fā)現(xiàn)，青藏高原云出現(xiàn)的概率大體呈現(xiàn)為從西北方向向東南方向遞增，青藏高原東南部地區(qū)云出現(xiàn)概率高達(dá)0.9。

1.2 數(shù)據(jù)源

1)遙感數(shù)據(jù)。

(1)AMSR2數(shù)據(jù)。2012年，日本宇宙航空研究開發(fā)機(jī)構(gòu)發(fā)射了高級(jí)微波掃描輻射計(jì)二號(hào)(advanced microwave scanning radiometer 2，AMSR2)，作為AMSR-E的替代品。AMSR2搭載在GCOM-W1衛(wèi)星上，是一種用于測量大氣表面和地球表面的微弱微波發(fā)射的工具[17]。AMSR2提供了在地球上空約700 km處對(duì)微波傳播和散射強(qiáng)度的精確測量。AMSR2的白天和夜間信息每兩天能覆蓋99%以上的地球表面。在水平和垂直極化下，AMSR2包含7個(gè)頻率通道，有6個(gè)頻率(6.9、10.6、18.7、23.8、36.5、89 GHz)用來得到亮溫，7.3 GHz頻率是用來減輕射頻干擾。傳感器分辨率隨著頻率變化而變化，范圍約為5 km(89 GHz)到60 km(6.9 GHz)，其中18.7 GHz和23.8 GHz通道的分辨率為22 km[18]。與MODIS和AIRS等其他傳感器相比，AMSR2的微波測量受多云條件和其他大氣效應(yīng)的影響較小，因此本文采用AMSR2升軌的中等分辨率L1B雙極化亮溫?cái)?shù)據(jù)。AMSR2升軌一般在世界時(shí)六點(diǎn)左右過境青藏高原地區(qū)。通過G-portal平臺(tái)下載研究區(qū)2017年1月和7月的AMSR2 L1B亮溫?cái)?shù)據(jù)。

(2)數(shù)字高程模型數(shù)據(jù)。數(shù)字高程模型(digital elevation model，DEM)是利用地形高程數(shù)據(jù)實(shí)現(xiàn)對(duì)地形曲面的數(shù)字化模擬，DEM的高程數(shù)據(jù)一般指的是絕對(duì)高程(即從標(biāo)準(zhǔn)水平面起算)。通過“地理空間數(shù)據(jù)云”平臺(tái)下載了研究區(qū)2017年SRTMDEMUTM數(shù)字高程產(chǎn)品，分辨率為90 m，投影方式為WGS 1984。

對(duì)下載的影像做如下預(yù)處理。①提取研究區(qū)域。利用IDL軟件對(duì)AMSR2 L1B亮溫?cái)?shù)據(jù)及DEM高程數(shù)據(jù)進(jìn)行重采樣、重投影(WGS1984)、拼接以及裁剪，得到青藏高原地區(qū)影像。②提取亮溫和高程信息。為了得到在晴空和有云條件下GPS站點(diǎn)對(duì)應(yīng)位置的高質(zhì)量亮溫及高程數(shù)據(jù)，首先，根據(jù)AMSR2影像的掃描時(shí)間選擇時(shí)間相近(正負(fù)半小時(shí))的葵花八號(hào)衛(wèi)星的云產(chǎn)品數(shù)據(jù)，得到GPS站點(diǎn)周圍3×3區(qū)域的像元值，若9個(gè)像元值中有超過5個(gè)像元是云像元，則該點(diǎn)為云像元，否則為晴空像元。然后，得到有云和晴空GPS站點(diǎn)的各頻率水平及垂直極化的亮溫和高程。

2) GPS水汽。通過西藏自治區(qū)獲得青藏高原GPS水汽站的觀測資料?；贗DL，首先，獲取AMSR2掃描時(shí)間的前后15分鐘時(shí)間段的平均GPS水汽值，剔除無效GPS觀測水汽值；其次，結(jié)合葵花衛(wèi)星云分類產(chǎn)品獲取青藏高原有云像元的GPS水汽值；最后，得到2017年1月和7月共 2 803個(gè)GPS水汽樣本數(shù)據(jù)。具體信息見表1。

表1 青藏高原GPS有效站點(diǎn)水汽統(tǒng)計(jì)

3)再分析數(shù)據(jù)。歐洲中期天氣預(yù)報(bào)中心(ECMWF)發(fā)布了對(duì)全球大氣數(shù)據(jù)集的ERA-Interim再分析資料，數(shù)據(jù)提供時(shí)間范圍為1979—2019年。再分析數(shù)據(jù)包含一天四次(00、06、12、18 UTC)，預(yù)報(bào)數(shù)據(jù)在世界時(shí)00和12點(diǎn)提供了四種步長(3、6、9、12)[19]。結(jié)合AMSR2升軌過境青藏高原的時(shí)間，在ECMWF 官網(wǎng)下載研究區(qū)2017年7月每天6點(diǎn)(UTC)的ERA-Interim水汽數(shù)據(jù)，分辨率為0.125°×0.125°。已有研究表明，ERA-Interim水汽數(shù)據(jù)能夠很好地展示青藏高原地區(qū)的水汽含量及其分布情況[20-21]，因此以該數(shù)據(jù)為參考，對(duì)比分析反演結(jié)果。

2 研究方法

2.1 變量選取

遺傳算法是一種受查爾斯·達(dá)爾文自然選擇思想啟發(fā)的數(shù)學(xué)模型[22]。機(jī)器學(xué)習(xí)中遺傳算法可以用來篩選變量以創(chuàng)建更優(yōu)的預(yù)測模型。其基本原理是通過適應(yīng)函數(shù)確定方向，通過選擇、交叉和突變組合，不斷地迭代搜索最優(yōu)的變量[23]?；陔S機(jī)森林回歸模型估算大氣水汽是以不同時(shí)空分辨率的遙感影像為主要數(shù)據(jù)源。以往較多研究是基于AMSR2 23.8 GHz和18.7 GHz等波段組合反演水汽。Du等[24]在考慮地形影響的基礎(chǔ)上，利用AMSR2數(shù)據(jù)，運(yùn)用參數(shù)化算法反演了大氣水汽。因此，本文選取了13個(gè)潛在變量因子參與變量篩選(表2)。

表2 參與變量篩選的13個(gè)潛在預(yù)測變量

利用python構(gòu)建遺傳算法(迭代次數(shù)10)，從13個(gè)潛在變量選取敏感變量。圖1為在迭代過程中，每個(gè)潛在變量被選取的次數(shù)。36.5 GHz和18.7 GHz亮溫(水平和垂直極化)和DEM被選中了10次；23.8 GHz的水平亮溫被選中了10次，其垂直極化下的亮溫被選中了6次；89 GHz的亮溫(水平和垂直極化)只被選中了1次。23.8 GHz位于大氣水汽吸收帶附近，23.8 GHz通道對(duì)大氣水汽反應(yīng)敏感[25]。地形強(qiáng)迫作用對(duì)溫度和水汽的分布有很大影響，水汽密度隨著對(duì)流層高度呈指數(shù)遞減趨勢[26]。地表溫度是大氣水汽反演必不可少的參數(shù)，而地表溫度的估算基本是基于36.5 GHz亮溫[27]。因此，本研究選擇23.8、36.5、18.7 GHz頻率的水平(H)和垂直(V)極化下的亮溫和高程來參與模型構(gòu)建。

圖1 潛在變量被選取的次數(shù)

2.2 樣本數(shù)據(jù)準(zhǔn)備

AMSR2 亮溫(除了7.3、10.6、89 GHz)和 DEM作為隨機(jī)森林回歸模型的輸入變量，GPS/PWV為模型的目標(biāo)變量。為了得到最優(yōu)模型，除了將像元按照晴空和有云分類，還需將樣本數(shù)據(jù)集按照季節(jié)分類(1月為干季、7月為濕季)，劃分為不同的樣本數(shù)據(jù)集。在樣本數(shù)據(jù)集中隨機(jī)選擇85%的數(shù)據(jù)作為訓(xùn)練樣本集，剩下的15%作為模型的驗(yàn)證數(shù)據(jù)集。

2.3 隨機(jī)森林回歸模型

隨機(jī)森林算法是集成學(xué)習(xí)眾多方法里的一種代表性方法[28]。隨機(jī)森林算法在以決策樹為基學(xué)習(xí)器構(gòu)建Bagging集成的基礎(chǔ)上，進(jìn)一步在決策樹的訓(xùn)練過程中引入了隨機(jī)屬性變量選擇。具體來說，隨機(jī)森林算法中，每一棵CART使用的訓(xùn)練數(shù)據(jù)集是從總數(shù)據(jù)集有放回的隨機(jī)選擇，對(duì)于基決策樹的每個(gè)結(jié)點(diǎn)，都是從所有特征集中隨機(jī)挑選特征[29]。

在python環(huán)境中，利用敏感變量構(gòu)建隨機(jī)森林回歸模型以及對(duì)變量進(jìn)行重要性排序。有兩個(gè)需要調(diào)優(yōu)的基本模型參數(shù)：樹的數(shù)量(K)和在每個(gè)節(jié)點(diǎn)拆分時(shí)隨機(jī)抽樣的變量數(shù)量(m)。除了兩個(gè)最具代表性的參數(shù)外，其他參數(shù)都是通過應(yīng)用默認(rèn)值來使用的。為了得到最優(yōu)參數(shù)，以1為間隔，對(duì)K和m的取值范圍分別設(shè)置為0～1 000、0～30。

2.4 評(píng)價(jià)指標(biāo)

機(jī)器學(xué)習(xí)回歸模型的常用評(píng)價(jià)指標(biāo)有均方根誤差(RMSE)、平均絕對(duì)誤差(MAE)、平均絕對(duì)比例誤差(MAPE)、皮爾遜相關(guān)系數(shù)(R)、方差等。本文采用R、RMSE和MAE作為評(píng)估模型預(yù)測能力的評(píng)價(jià)指標(biāo)。R可以反映出預(yù)測值與實(shí)測值之間的相關(guān)性強(qiáng)度，RMSE和MAE越小代表著模型預(yù)測結(jié)果越準(zhǔn)確。

3 結(jié)果與分析

3.1 模型實(shí)現(xiàn)

在不同大氣條件和不同季節(jié)下，以隨機(jī)森林算法為基礎(chǔ)，通過變量選取和參數(shù)優(yōu)化，得到模型預(yù)測值與實(shí)測值之間的擬合關(guān)系。圖2為GPS水汽與模型預(yù)測的水汽之間的散點(diǎn)圖。

圖2 隨機(jī)森林回歸模型精度評(píng)價(jià)

分析可得以下結(jié)論。①基于隨機(jī)森林回歸模型預(yù)測的大氣水汽和GPS水汽有顯著的相關(guān)性，RMSE和MAE整體值較小。有云像元反演的三個(gè)模型精度都很高，說明隨機(jī)森林回歸模型用來反演多云地區(qū)的大氣水汽是可行的。整體而言，相較于所有有云樣本數(shù)據(jù)參與構(gòu)建的模型，按季節(jié)構(gòu)建的模型精度更優(yōu)。干季有云像元大氣水汽預(yù)測值和GPS水汽相關(guān)性高達(dá)0.928，RMSE為1.15 mm，MAE接近0；濕季有云像元大氣水汽預(yù)測值和GPS水汽相關(guān)性為0.882，RMSE為4.08 mm，MAE為3.08 mm；所有有云像元大氣水汽預(yù)測值和GPS水汽相關(guān)性高達(dá)0.919，RMSE為4.74 mm，MAE為3.39 mm。干季有云像元水汽反演精度高于濕季，7月的青藏高原地區(qū)是雨季，呈多云，有云像元數(shù)量高于晴空，水汽含量高。②在干、濕季節(jié)中，有云像元反演精度略高于晴空像元。以GPS水汽為參考，干季有云和晴空像元的R分別為0.928和0.901，RMSE分別為0.127 cm和0.119 cm，MAE分別為0.095 cm和0.088 cm；濕季有云和晴空像元的R分別為0.882和0.845，RMSE分別為0.408和0.530，MAE分別為0.308 cm和0.439 cm。結(jié)果表明，AMSR2水汽反演精度受云層的影響較小。

3.2 與MODIS水汽產(chǎn)品比較

已有研究表明，MODIS水汽產(chǎn)品可以較好地反映青藏高原地區(qū)的水汽特征[30]，光學(xué)遙感受天氣影響較大，而微波具有一定的穿透性。因此，將GPS實(shí)測值作為參考，分別評(píng)估利用AMSR2數(shù)據(jù)采用隨機(jī)森林回歸模型反演的水汽和MODIS水汽產(chǎn)品數(shù)據(jù)在青藏高原有云條件下的精度。

基于模型反演的有云像元，選出相應(yīng)的MOD05近紅外和熱紅外有云像元，結(jié)果如圖3所示。由圖3可知，干季MOD05近紅外和熱紅外的有云像元水汽與GPS水汽R均低于0.6，普遍低于真實(shí)值。干季AMSR2有云像元水汽反演結(jié)果與GPS的觀測水汽結(jié)果符合良好，RMSE為0.11 cm，MAE為0.082 cm;濕季AMSR2有云像元水汽反演結(jié)果與GPS水汽的R為0.847，RMSE為0.274 cm。濕季MOD05近紅外和熱紅外的有云像元水汽與GPS水汽之間的相關(guān)性很低，R分別為0.38和0.308。MODIS現(xiàn)有的水汽產(chǎn)品數(shù)據(jù)在晴空條件下具有較高的精度，但是對(duì)于有云條件下的水汽反演明顯不足。因此，基于AMSR2和高程數(shù)據(jù)利用隨機(jī)森林算法反演的有云像元大氣水汽含量精度高于光學(xué)遙感有云像元水汽的反演精度。

圖3 AMSR2水汽和MOD05(近紅外/熱紅外)水汽與GPS遙感水汽之間的關(guān)系圖

3.3 青藏高原大氣水汽反演

AMSR2需要兩天遙感影像覆蓋整個(gè)青藏高原，本研究采用濕季有云條件下的隨機(jī)森林回歸模型估算2017年7月1日—2日兩天的大氣水汽，并合成完整的青藏高原地區(qū)大氣水汽。大氣水汽含量由東南向西北方向減小，高原東南地區(qū)和西南地區(qū)為水汽分布的高值區(qū)，其中高原東南地區(qū)最濕，其含量可高達(dá)5 cm。反演結(jié)果表明，青藏高原地區(qū)大氣中的水汽呈現(xiàn)南濕北干的空間分布特征，該結(jié)論與周長艷等[31]的結(jié)論一致。

為了驗(yàn)證AMSR2亮溫?cái)?shù)據(jù)反演的青藏高原大氣水汽含量及分布，本文選用了再分析ERA-Interim水汽數(shù)據(jù)與其進(jìn)行對(duì)比分析。青藏高原大氣水汽分布情況與隨機(jī)森林回歸模型反演大氣水汽分布圖相比，具有較好的一致性。為了更好地對(duì)比分析AMSR2和ERA-Interim的大氣水汽數(shù)據(jù)，將ERA-Interim大氣水汽數(shù)據(jù)進(jìn)行空間插值至AMSR2大氣水汽數(shù)據(jù)的空間分辨率，并進(jìn)行點(diǎn)對(duì)點(diǎn)的匹配。以AMSR2大氣水汽數(shù)據(jù)為參照面，利用評(píng)價(jià)指標(biāo)，分析其與ERA-Interim大氣水汽數(shù)據(jù)之間的相關(guān)關(guān)系?？梢园l(fā)現(xiàn)，AMSR2和ERA-Interim的大氣水汽數(shù)據(jù)顯著相關(guān)，兩者大氣水汽的RMSE為5.48 mm，MAE為4.35 mm。

4 結(jié)束語

青藏高原及其周邊地區(qū)的大氣水汽含量時(shí)空分布對(duì)全球和區(qū)域的氣候變化、大氣循環(huán)等有著重大影響。近年來，機(jī)器學(xué)習(xí)在遙感領(lǐng)域一直為熱點(diǎn)問題，其中，隨機(jī)森林算法憑借其高效、準(zhǔn)確的性能被廣泛應(yīng)用于遙感估算。本文基于GPS實(shí)測大氣水汽數(shù)據(jù)、AMSR2亮溫?cái)?shù)據(jù)，利用隨機(jī)森林算法構(gòu)建非線性回歸模型，對(duì)青藏高原多云地區(qū)大氣水汽含量進(jìn)行了遙感估算與分析，得到以下結(jié)論。

1)AMSR2亮溫?cái)?shù)據(jù)結(jié)合隨機(jī)森林算法適用于青藏高原有云地區(qū)的大氣水汽含量反演。無論是干季還是濕季，利用AMSR2反演有云像元的大氣水汽的精度較高，其中有云像元的水汽估算精度高于晴空像元，濕季模型擬合預(yù)測的水汽精度低于干季。與MODIS水汽產(chǎn)品相比較，AMSR2反演的有云像元水汽含量準(zhǔn)確性更高，可以用來彌補(bǔ)光學(xué)遙感反演水汽的不足。

2)23.5 GHz是大氣水汽吸收帶，水汽密度隨著對(duì)流層高度呈指數(shù)遞減趨勢。研究表明，AMSR2 23.5 GHz和36.5 GHz頻率的亮溫和DEM高程是隨機(jī)森林回歸模型中對(duì)大氣水汽含量具有重要意義的預(yù)測變量。大氣水汽估算涉及的變量較多，會(huì)受到不止一種因素的影響。

3)由于衛(wèi)星遙感數(shù)據(jù)的分辨率、大氣水汽本身特性等眾多限制因素影響，對(duì)有云像元水汽反演而言，按季節(jié)分類構(gòu)建的模型精度高于所有樣本數(shù)據(jù)參與模型的精度。濕季模型估算2017年7月1日—2日的青藏高原大氣水汽，水汽分布為南濕北干，與2017年7月1日—2日的再分析ERA-Interim水汽分布趨勢基本一致。