唐超禮, 陶鑫華, 魏圓圓, 戴聰明, 魏合理

1. 安徽理工大學(xué)電氣與信息工程學(xué)院, 安徽 淮南 232001;

2. 中國科學(xué)院空間天氣學(xué)國家重點實驗室, 北京 100190;

3. 安徽大學(xué)互聯(lián)網(wǎng)學(xué)院, 安徽 合肥 230039;

4. 中國科學(xué)院合肥物質(zhì)科學(xué)研究院安徽光學(xué)精密機械研究所, 中國科學(xué)院大氣光學(xué)重點實驗室, 安徽 合肥 230031

地表溫度(land surface temperature, LST)是反映陸地與大氣之間能量交換的關(guān)鍵物理量, 與生態(tài)環(huán)境密切相關(guān), 通常用于地表通量估算、干旱監(jiān)測, 并且可直接或間接作為驅(qū)動模型用于全球環(huán)境監(jiān)測和氣候變化分析(Weng et al, 2009; Karnieli et al, 2010; Stephens et al, 2015)。如今, 我們可以利用地面觀測和遙感衛(wèi)星的長期觀測資料, 對所選取的地區(qū)進行LST 的研究。我們也有多種分析手段, 如基于風(fēng)云三號D 星(FY—3D)微波成像儀對中國陸地區(qū)域LST反演及驗證進行比較(Nemani et al, 1996); 又如采用蒙大拿大學(xué)地表反演算法(University of Montana land surface retrieval algorithm, UMT)中的物理模型反演中國區(qū)域LST(王博 等, 2022), 并利用該衛(wèi)星與地面站點的數(shù)據(jù)進行比較; 以及利用單窗算法對數(shù)據(jù)進行反演, 隨后用歸一化指數(shù)研究區(qū)域范圍內(nèi)的時空變換相關(guān)性(Becker et al, 1990; 王艷慧 等, 2014); 這些都是對局部數(shù)據(jù)進行分析研究。本文利用2003—2020 年衛(wèi)星反演的LST 數(shù)據(jù), 通過M-K 突變檢驗、線性回歸、EOF 等方法, 分析東亞及沿海區(qū)域的LST 時空模態(tài)特征, 并運用SARIMA 模型預(yù)測LST 的變化趨勢。

1 數(shù)據(jù)來源與研究方法

1.1 研究區(qū)域及數(shù)據(jù)來源

本文研究區(qū)域為東亞及西太平洋地區(qū), 東亞地理位置特殊, 靠近赤道, 此處地域廣闊, 東西跨度大, 其LST 主要受太陽輻射和氣候的雙重影響。又因為東亞及西太平洋地區(qū)屬于典型的季風(fēng)氣候, 該地區(qū)地貌類型多樣, 四季分明。因此, 研究該地區(qū)的LST 變化規(guī)律對東亞及西太平洋地區(qū)的氣候及災(zāi)害預(yù)警有著重要的參考意義。

本文利用Aqua 衛(wèi)星自2003 年1 月—2020 年12月 LST 數(shù)據(jù)及第五代歐洲中期天氣預(yù)報中心(European centre for medium range weather forecasting,ECMWF)再分析資料集(The fifth generation ECMWF reanalyses, ERA5)對LST 進行分析研究。大氣紅外探測儀(atmospheric infrared sounder, AIRS) 提供了溫度、濕度等數(shù)據(jù), 將大氣參數(shù)劃分成1°×1°的網(wǎng)格單元, 經(jīng)度從–180°到+180°, 緯度從–90°到+90°, 運用升軌、降軌數(shù)據(jù)即白天(A)和夜晚(D)數(shù)據(jù)進行處理。除此之外, 還使用了2018—2020年的ERA5再分析數(shù)據(jù),它是ECMWF 在陸地范圍的模擬數(shù)據(jù), 數(shù)據(jù)包含了月、日、小時數(shù)據(jù)集。在前人的文章中, 利用伊比利亞半島分析衛(wèi)星的應(yīng)用設(shè)施, 對ERA5 的LST 數(shù)據(jù)進行評估, 發(fā)現(xiàn)ERA5 的日最大LST 的標(biāo)準(zhǔn)偏差約為2～3 ℃(Johannsen et al, 2019)。

1.2 研究方法

本文采用了M-K(Mann-Kendall)突變檢驗, 線性回歸分析, 季節(jié)自回歸移動平均模型(seasonal autoregression integrated moving average, SARIMA),經(jīng)驗正交分解(empirical orthogonal function, EOF)等方法。其中, M-K 突變檢驗是一種非參數(shù)統(tǒng)計檢驗方法, 分析突變的LST 時間點可以準(zhǔn)確知其突變情況; SARIMA 模型對未來年份做季節(jié)時間序列預(yù)測,把LST 變成隨機序列。當(dāng)隨機序列有很強的季節(jié)周期性變化時, 可以對其進行季節(jié)差分, 再利用2003—2020 年月平均溫度做一階季節(jié)性差分序列的自相關(guān)和偏相關(guān)函數(shù)圖以及用差分次數(shù)確定模型參數(shù)的數(shù)值, 建立模型進行預(yù)測。此模型可以從時間序列的過去值及現(xiàn)在值來預(yù)測未來值。

2 東亞及西太平洋LST 時空變化

2.1 LST 年際變化分析

圖1a 展示了東亞及西太平洋地區(qū)LST 在2003—2020 年的年際變化曲線。如圖, 東亞及西太平洋地區(qū)的LST 在6～30℃之間, 分別在2011 年、2017 年達到最小值和最大值。

圖1 LST 年際變化圖Fig. 1 LST interannual variation chart. (a) Monthly trend chart of East Asia and various regions; (b) anomaly map

整體看, 該地區(qū)大部分是沿海地區(qū), 溫度變化有限; 局部看, 蒙古國 LST 溫差大, 最高值超過40 ℃, 最低值低于–10 ℃; 日本則相反, LST 溫差不大, 差值不足13℃。蒙古國與日本的LST 變化幅度差異, 主要原因是太陽輻射。蒙古國緯度高遠離赤道而日本緯度低靠近赤道, 且夏季風(fēng)與西太平洋漂流區(qū)的海表溫度為正相關(guān), 與赤道中東部海表溫度為負相關(guān)(沈柏竹, 2007), 整體呈自西向東走勢, 冬季相反。東亞各地區(qū)LST 的距平圖如圖1b 所示, 為進一步給出該地區(qū)是否存在LST 突變的變化過程,我們運用M-K 突變檢驗法對東亞及西太平洋LST進行突變分析(陳錦年 等, 2008; 王宏娜 等, 2009)。

由圖2 可見,UFk和UBk曲線在2007 年和2016年及2018 年附近存在突變點。2003—2011 年UFk統(tǒng)計值均大于0, 說明該地區(qū)LST 在此期間呈遞增趨勢, 并且統(tǒng)計值均未超出臨界值, 說明上升趨勢不顯著; 2011—2014 年末UFk統(tǒng)計量小于0, 說明LST 在此期間呈遞減趨勢; 2015 年以后呈逐年上升趨勢且在2019 年后超出臨界值, 說明此后上升趨勢顯著。這些年份出現(xiàn)突變點離不開當(dāng)年全球發(fā)生的自然災(zāi)害或大規(guī)模人為事件。

圖2 M-K 突變檢測結(jié)果Fig. 2 M-K mutation detection results

2.2 LST 季節(jié)變化

將2003—2020 年AIRS 的LST 數(shù)據(jù), 按季節(jié)劃分為春季(3—5 月)、夏季(6—8 月)、秋季(9—11 月),和冬季(12—2 月)。圖3a 展示了東亞及西太平洋區(qū)域自2003—2020 年春季的月平均LST 的空間分布。青藏高原等西南地區(qū)LST 較低, 中國新疆及周邊地區(qū)LST 較高, 從而出現(xiàn)了明顯的溫度差異。在平均海拔都超過3000m 的青藏高原上, 它的地勢明顯西側(cè)高于東側(cè), 中間低而南北側(cè)高, 山地走勢非常復(fù)雜, 且80%的青藏高原都是高原地貌, 全年平均氣溫均低于10℃, 這也導(dǎo)致了氣溫差異大(曹曉云 等,2021)。

圖3 2003—2020 年東亞及西太平洋LST 四季空間分布a. 春季; b. 夏季; c. 秋季; d. 冬季; e. 四季變化折線圖; f. 四季距平圖。a—d 基于國家測繪地理信息局標(biāo)準(zhǔn)地圖服務(wù)網(wǎng)站下載的審圖號為GS(2016)2938 號的標(biāo)準(zhǔn)地圖制作Fig. 3 Spatial distribution of LST seasons in East Asia and West Pacific from 2003 to 2020. (a) Spring; (b) summer; (c)autumn; (d) winter; (e) line chart of four seasons; (f) four seasons anomaly map

圖3b 為東亞及西太平洋地區(qū)自2003—2020 年夏季的月平均LST 的空間分布。其中, 夏季普遍高溫, 主要原因是太陽輻射。而夏日輻射量高, 冬日輻射量少, 是因為季節(jié)變化導(dǎo)致太陽高度角的變化,進而引起的太陽輻射的差異(黃盼 等, 2017)。然而,夏季陸地的LST 比海洋的LST 高, 是因為在相同的太陽輻射條件下, 陸地的升溫或者降溫都比海洋更加劇烈, 利于陸地LST 高于海洋。圖3c 是東亞及西太平洋區(qū)域自2003—2020 年秋季的月平均LST的空間分布, 變化趨勢與春季大體相似, 東亞北部地區(qū)LST 局部較高于周圍(高弋斌 等, 2019)。圖3d為東亞及西太平洋區(qū)域自2003—2020 年冬季的月平均LST 的空間分布。如圖, 沿海地區(qū)LST 是明顯高于內(nèi)陸的。除了太陽輻射的影響, 地溫少雨的冬季風(fēng)也是導(dǎo)致陸地溫度低的重要原因。海洋比熱容大, 地面比熱容小, 所以冬季陸地相對于海洋溫度低(王靜 等, 2009; 郭甜甜 等, 2017)。綜上可知, 本文選取的數(shù)據(jù)存在顯著的季節(jié)周期變化。圖3e 為該地區(qū)四季折線圖, 可見春、夏、秋三季幅度變化不明顯。其中, 春季和秋季溫度相近; 夏季均溫超過30 ℃; 冬季溫度較低, 最低約為8℃。圖3f 所示折現(xiàn)圖為四季圖3e 的距平圖。

2.3 LST 空間變化分析

圖4a 為AIRS 探測器2003—2020 年平均LST分布圖。沿海地區(qū)的LST 高于東亞內(nèi)陸及青藏高原區(qū)域。圖4b 中曲線表示2003—2020 年東亞及西太平洋地區(qū)LST 隨緯度的變化趨勢, 可見LST 隨著緯度的升高而變小。其中, LST 最低溫度為2℃的主要原因是緯度的增加, 沿海地區(qū)的季風(fēng)氣候以及陸地的散熱和吸熱能力小于海洋。圖4c 是ERA5 再分析數(shù)據(jù)2020 年的東亞及西太平洋的平均溫度圖, 其溫度范圍與AIRS 的結(jié)果相似。

圖4 2003—2020 年東亞及西太平洋LST 空間分布圖a. AIRS 空間分布; b. 緯度分布; c. ERA5 再分析空間分布。a、c 基于國家測繪地理信息局標(biāo)準(zhǔn)地圖服務(wù)網(wǎng)站下載的審圖號為GS(2016)2938 號的標(biāo)準(zhǔn)地圖制作Fig. 4 Spatial distribution of LST in east Asia and west Pacific from 2003 to 2020. (a) Spatial distribution of AIRS; (b)latitude distribution; (c) ERA5 reanalysis spatial distribution

EOF 是一種分析矩陣數(shù)據(jù)中的結(jié)構(gòu)特征和提取主要數(shù)據(jù)特征量方法, 也是氣象分析時空特征的一種重要方法。將LST 組合排列成m空間點數(shù)(根據(jù)文章中的經(jīng)緯度10°N—60°N 和70°E—140°E, 即50×70 個點數(shù)據(jù))和n時間點數(shù)(根據(jù)本文選用的18年數(shù)據(jù)), 得到Xm×n。

理論上, 對LST 的分析應(yīng)比濕度更直觀。但由于LST 受眾多因素的影響, 使得影響LST 變化趨勢的因素較為復(fù)雜(劉振元 等, 2018)。所以我們也可以先查看濕度在時間、空間上的變化規(guī)律, 再考查LST。

圖5a、c、e 為濕度EOF 分解第一、第二、第三模態(tài)的空間分布, 其中模態(tài)解釋方差貢獻率分別為26.18%, 18.70%和10.12%, 累計貢獻率達55%,能較好地代表濕度的變化。圖5b、d、f 為濕度EOF分解第一、第二、第三模態(tài)的時間系數(shù)。系數(shù)符號決定了模態(tài)的方向, 正號表示與模態(tài)方向相同, 負號則相反, 時間系數(shù)數(shù)值絕對值越大, 表明該年分布型式越典型(魏鳳英, 1999)。圖5a 為第一模態(tài)的空間分布圖, 貢獻率為26.18%, 該模態(tài)表征了2003—2020 年東亞及西太平洋地區(qū)濕度的最主要的分布特征。可以看出, 在空間分布上, 東亞及西太平洋大部分地區(qū)呈負值分布, 極少數(shù)地區(qū)呈正值, 說明該模態(tài)下濕度在大部分地區(qū)空間變化趨勢具有一致性,即同時上升(下降)。結(jié)合圖5b 的時間系數(shù)來看, 時間系數(shù)在2004—2009 年均為正位, 在2010—2020年負位居多, 即在2010 年前后基本存在正負相位,并呈下降趨勢, 說明2003—2020 年濕度整體是下降的。

圖5 2003—2020 年東亞及西太平洋地區(qū)濕度EOF 空間分布圖和時間系數(shù)圖a、c、e 分別為第一、第二、第三模態(tài)的空間分布圖, 基于國家測繪地理信息局標(biāo)準(zhǔn)地圖服務(wù)網(wǎng)站下載的審圖號為GS(2016)2938 號的標(biāo)準(zhǔn)地圖制作; b、d、f 為第一、第二、第三模態(tài)的時間系數(shù)圖Fig. 5 Spatial distribution map and time coefficient map of humidity EOF in east Asia and west Pacific from 2003 to 2020.(a), (c), (e) spatial distribution map; (b), (d), (f) time coefficient diagram

圖6a、c、e 表示東亞及西太平洋地區(qū)LST 年際變化的EOF 分解第一、第二、第三模態(tài)的空間分布,其中模態(tài)解釋方差貢獻率分別為29.58%, 13.35%和10.92%, 累計貢獻率達53.85%, 能較好地揭示LST的分布特征。圖6b、d、f 表示LST 數(shù)據(jù)的EOF 分解第一、第二、第三模態(tài)的時間系數(shù)。圖6a 為第一模態(tài)的空間分布圖, 貢獻率為29.58%, 高于其他模態(tài), 所以選取第一模態(tài)來分析(韓雪 等, 2014)。在空間分布上, 第一模態(tài)以中國昆侖山脈、秦嶺為分界線, 把東亞及西太平洋地區(qū)劃分成兩個明顯的上下區(qū)域, 負值主要分布在35°N 以南, 110°N 以西區(qū)域,其余地區(qū)大部分為正值, 正值顯著區(qū)主要在東北地區(qū), 表明該地區(qū)LST 變化趨勢更為明顯。結(jié)合圖6b時間系數(shù)來看, 時間系數(shù)在2009—2013 年均為負位,2014—2020 年正位居多, 時間系數(shù)呈上升趨勢, 說明2003—2020 年LST 整體是上升趨勢。圖6d 從2013—2018 年均正位, 其余年份正、負相位交替出現(xiàn), 基本沒有連續(xù)正(負)位。圖6f 和圖6d 變化相似。從時間系數(shù)來看, LST 在2012 年之前逐漸變高, 總體為正位, 之后就逐漸變低, 總體為負位(圖6f)。圖6b 中, 2003 年的時間系數(shù)為負位, 與模態(tài)空間變化方向相反, 此年, 東亞及西太平洋地區(qū)的模態(tài)空間分布表現(xiàn)為以昆侖山脈與秦嶺為分界線的東北區(qū)域溫度降低, 余下區(qū)域溫度升高的變化趨勢。反之, 時間系數(shù)為正, 呈相反的分布模式。

圖6 2003—2020 年東亞及西太平洋地區(qū)LST 的EOF 空間分布圖和時間系數(shù)圖a、c、e 分別為第一、第二、第三模態(tài)的空間分布圖, 基于國家測繪地理信息局標(biāo)準(zhǔn)地圖服務(wù)網(wǎng)站下載的審圖號為GS(2016)2938 號的標(biāo)準(zhǔn)地圖制作; b、d、f 為第一、第二、第三模態(tài)的時間系數(shù)圖Fig. 6 Spatial distribution map and time coefficient map of EOF in east Asia and west Pacific of LST from 2003 to 2020. (a),(c), (e) spatial distribution map; (b), (d), (f) time coefficient diagram

3 東亞及西太平洋LST 影響因素分析及預(yù)測

3.1 LST 影響因素分析

圖7 為LST 與濕度趨勢變化圖, 可見LST 一直呈上升趨勢, 這與本文關(guān)于東亞及西太平洋地區(qū)LST 的M-K 分析結(jié)果一致。其次, 濕度是影響LST的關(guān)鍵因素, 而溫度可以改變濕度。如圖7a, 最高溫度為21 ℃, 整體濕度隨溫度波動, 除2007 年因當(dāng)年發(fā)生的熱帶風(fēng)暴溫濕度驟變, 其余年份基本符合趨勢。圖7b 是夜晚的溫、濕度, 兩者變化情況基本一致。測量LST 和濕度, LST 反演精度約為1.63, 濕度反演精度約為 0.063%, 且溫度升高濕度減少。

圖7 LST 與濕度趨勢變化圖a. 白天; b. 夜晚Fig. 7 Trend of surface temperature and humidity. (a) Day; (b) night

選用了2003—2020 年的東亞及西太平洋地區(qū)的相關(guān)AIRS 溫濕度數(shù)據(jù), 圖8a 的相關(guān)性分析結(jié)果顯示, 各地區(qū)的白天溫、濕度之間相關(guān)性良好, 最佳擬合線的截距為183.409, 回歸系數(shù)0.403, 相關(guān)性到達0.505(陳錦年 等, 2012); 圖8b 為溫度與濕度的折線圖, 可見溫、濕度整體波動幅度相似。

圖8 LST 與濕度相關(guān)性分析(a)及時間序列比較(b)Fig. 8 Correlation analysis and time series comparison between LST and humidity. (a) Correlation analysis diagram between LST and humidity; (b) comparison diagram of LST and humidity time series

3.2 預(yù)測

運用SARIMA 進行一階季節(jié)性差分時間序列分析, 推斷東亞及西太平洋地區(qū)LST 變化趨勢和未來預(yù)測。

如前文所述, 在圖9 原始圖中, 2011 年和2012年LST 驟降, 存在異常, 且2012—2020 年LST 為明顯上升趨勢。經(jīng)過多方資料查證, 這與當(dāng)時的厄爾尼諾、火山爆發(fā)等事件相關(guān)(李剛 等, 2015;Simmons et al, 2017), 也是陸地溫度和海表溫度共同影響的結(jié)果。

圖9 原始圖Fig. 9 Original drawing

圖10 為年平均溫度的預(yù)測溫度值, 假設(shè)LST 的時間序列為w(t),v(t)=w(t)–w(t–1),z(t)=v(t)–v(t–12),z(t)=a1*z(t–1)+m(t)+b1*m(t–1)+c1*z(t–12)+d1*m(t–24),其中w是誤差項。預(yù)測時套用以上公式, 只是先預(yù)測z(t), 再預(yù)測v(t)和w(t)。根據(jù)過去值, 現(xiàn)在值以及預(yù)測的函數(shù)方程, 預(yù)測出未來值。2020 年之后, 東亞及西太平洋地區(qū)LST 變化平穩(wěn), 若當(dāng)年或者前后出現(xiàn)重大自然問題, 可能溫度有浮動。

圖10 預(yù)測圖Fig. 10 Prediction chart

4 結(jié)論

本文利用2003—2020 年衛(wèi)星觀測反演的LST數(shù)據(jù), 通過M-K 突變檢驗、線性回歸、EOF 模態(tài)分析等方法, 分析東亞及沿海區(qū)域的LST 時空模態(tài)特征, 并運用SARIMA 模型預(yù)測LST 的變化趨勢。得出以下結(jié)論。

1) 利用LST 年月平均值法, 可知東亞陸地LST變化幅度高于沿海。例如, 蒙古國身處東亞內(nèi)陸, 年平均最高值超 40 ℃, 最小值低于–10 ℃, 差值達50 ℃; 而東亞其他地區(qū)溫差不足20℃。

2) 利用M-K 突變檢驗法分析18 年LST 數(shù)據(jù)。整體看, 東亞及西太平洋地區(qū)LST 于2003—2011年, 2015 年直至以后呈遞增趨勢, 且2019 年超出臨界值, 上升趨勢顯著; 局部看, 東亞及西太平洋地區(qū)LST 在2003—2020 年春季、夏季及秋季波動平穩(wěn), 而冬季波動大。這與人類活動及火山爆發(fā)等事件相關(guān)。

3) 東亞及西太平洋區(qū)域LST 季節(jié)變化特征明顯。除青藏高原常年處于明顯低溫(主要是青藏高原80%以上的高原地貌, 年均不超過10 ℃), 而夏、冬季LST 恰好呈相反模式, 這是受太陽輻射和季風(fēng)氣候的雙重作用。夏季輻射量最高, 冬季最少, 是因為季節(jié)變化導(dǎo)致太陽高度角的變化, 進而引起太陽輻射的差異; 海洋比熱容大, 陸地比熱容小, 所以冬季陸地相對于海洋溫度低。秋季與春季變化相似。

4) 利用EOF 分解法對東亞及西太平洋區(qū)域18年LST 分析, 第一模態(tài)與LST 變化趨勢在空間上相符, 以昆侖山脈、秦嶺為分界線, 形成東亞內(nèi)陸與沿海地溫走勢不同的總態(tài)。