康銀紅 王嘉馳 宋 鑫,2 武劍飛 賀 帥
(1.四川農業大學水利水電學院, 雅安 625014; 2.四川省綿陽市涪城區農業農村局, 綿陽 621000;3.中國農業科學院農田灌溉研究所, 新鄉 453002)
干旱是一個地區由于長期降水與蒸散收支不平衡導致水分異常短缺而產生的自然災害,具有發生頻率高、影響范圍廣、持續時間長等特點,不僅會給當地農業等部門帶來直接經濟損失,也會對水資源、土地資源及社會資源造成間接影響[1]。因此,加強干旱研究有助于減輕干旱造成的危害。目前,干旱研究大多基于干旱指數進行[2-4]。干旱指數能夠有效反映地表的干濕情況,是監測、評估干旱的重要參數,也是旱災研究的基礎。因此,開展干旱指數適用性評估,對干旱預警及資源保護具有重大意義[5]。
目前干旱評價指標主要有帕默爾干旱指數(PDSI)[6]、標準化降水指數(SPI)[7]、標準化降水蒸散指數(SPEI)[8]等。其中PDSI指數雖然考慮了溫度對干旱的影響,適用全球變暖背景下的干旱研究,但其時間尺度固定,難以表征短期干旱[9]。相較于PDSI指數,SPI指數能夠在不同空間和時間尺度上評估干旱,但其僅考慮了降水作用,并未考慮蒸散的影響。SPEI指數綜合了前兩種指數的優點,不僅可以反映多時間尺度的干濕情況,還通過引入參考作物騰發量(ET0)評估了其他氣象因素對干旱的影響,能更加客觀地描述地表干濕變化,適用于氣候變暖背景下干旱特征的分析,在國內外應用廣泛[10-12]。
由于ET0的計算精度會影響SPEI指數的干旱監測能力,因此在利用SPEI監測干旱前,需要考慮ET0的計算精度。Penman-Monteith(PM)法[13]是聯合國糧農組織(FAO)推薦使用的ET0計算方法,具有較高的權威性和準確性。由于PM法計算ET0時,需要大量氣象參數,對站點的要求較高,所以在某些資料匱乏的地區難以得到有效應用[14-16]。為此,學者們提出了許多ET0簡化算法,但這些方法都有特定的適用條件,計算ET0時,需要根據不同地區實際情況選用[17]。目前,針對不同ET0計算方法適用性的研究較多,然而考慮不同ET0算法對SPEI影響的研究還相對較少。因此,有必要厘清不同ET0計算方法的適用性及其對SPEI的影響。
四川省水資源分布差異顯著,區域降水不均、地勢復雜,容易發生區域性、季節性干旱。近年來,由于受到全球氣候持續變暖的影響,四川省干旱強度及頻率有增長趨勢[18]。特別是2022年7月28日—8月22日,持續高溫事件綜合強度為1961年有完整氣象觀測記錄以來最強,部分地區旱情嚴重,農作物減產,群眾飲水困難。因此需要選擇適宜的SPEI指數對區域降水不均、地勢復雜、容易發生干旱的四川省進行準確的干旱監測與評估,并較為準確地研究當地的干旱時空分布情況。本文以PM法計算的ET0為標準,通過誤差分析、空間插值等方法比較7種不同方法的ET0計算結果,并選取其中表現較好的3種方法計算得到其對應的SPEI指數,通過理論驗證、泰勒圖等方法找出適合四川省不同時間尺度和不同區域的SPEI指數,為全球氣候變暖背景下四川省的干旱監測和防控提供理論依據。
研究區為四川省(26°03′~34°19′N,97°21′~108°33′E),地形地貌復雜多樣,東部為盆地,西部為青藏高原,地勢西高東低。根據四川省地形特征將研究區分為川西高原、川西南山地和川中盆地(圖1)。該區地處亞熱帶,受地理位置及地形條件影響,區域氣候差異顯著。川中盆地屬于亞熱帶季風氣候,氣溫高、濕度大;川西南山地為亞熱帶半濕潤氣候區,氣溫高、雨量少;川西高原為典型的高原氣候,氣候立體變化明顯。省內水資源豐富,年平均降水量在1 000 mm以上,但區域水資源分布差異顯著,降水多集中在夏季,冬春季雨量稀少。

圖1 四川省地形圖Fig.1 Topographic map of Sichuan Province
所需的氣象數據主要包括四川省34個氣象站點1967—2016年的平均相對濕度、日照時數、日最高氣溫、日最低氣溫、日平均氣溫、日降水量、高度10 m處風速等。數據來自中國氣象數據網中國地面氣候資料日值數據集(V3.0)(http:∥data.cma.cn),選取的氣象站點資料質量控制良好,對資料中部分氣象站點的缺失數據通過同時期的臨近站點插補處理。四川省歷年來的實際旱情數據和相關描述來源于中國氣象局國家氣候中心(https:∥cmdp.ncc-cma.net/)、《中國氣象災害大典》(四川卷)[19]、《四川水旱災害》[20]等相關統計資料。
2.2.1ET0計算方法
本文采用8種常用的ET0計算方法,分別為標準綜合算法PM法[13];輻射法中的FAO-24Radiation(FAO-Ra)法[21]、Priestley-Taylor(PT)法[22]以及Makkink(MK)法[23];溫度法中的Hargreaves-Samani(HS)法[24]和Blaney-Criddle(BC)法[25];物質傳輸法中的WMO(World Meteorological Organization)法[26]和Rohwer(Ro)法[26]。各計算方法見表1。通過各方法計算得到1967—2016年的日ET0,然后將1—12月的日ET0求和作為當年ET0,各區域的ET0由相應區域內各站點的平均值來表示。

表1 ET0計算方法Tab.1 Calculation methods of reference crop evapotranspiration
2.2.2SPEI計算方法
利用四川省34個氣象站點的逐月降水量和ET0計算1、3、12個月時間尺度的SPEI(SPEI1、SPEI3、SPEI12),12個月時間尺度的SPEI用來進行年際分析,其中第12個月的SPEI被用來監測當年的干旱狀況。以PM公式為例,所對應的SPEI表示為SPEI_PM。
SPEI計算過程如下[8,27]:
首先計算逐月降水量與參考作物騰發量的差值Di,計算式為
Di=Pi-ET0
(1)
式中i——月數Pi——逐月降水量,mm
ET0——月參考作物騰發量,mm
再采用三參數的log-logistic概率密度函數對Di數據序列進行擬合,計算式為
(2)
式中α——尺度參數
β——形狀參數
γ——位置參數
3個參數可以通過線性矩方法擬合獲得。通過概率密度函數得到Di數據序列的累計概率函數,計算式為
(3)
對Di數據序列的累計概率函數進行正態標準化處理,并求得SPEI指數。計算式為
P=1-F(x)
(4)
(5)

(6)

式中系數c0、c1、c2、d1、d2、d3分別為2.515 517、0.802 853、0.010 328、1.432 788、0.189 269、0.001 308。
2.2.3精度評價指標
為了探究ET0方法及對應的SPEI在各區域的差異,以PM/SPEI_PM為標準,采用均方根誤差(RMSE)、相對誤差(RE)評價其他方法的精度,兩個指標的值越小,方法的精度越高。為了驗證水分虧缺量是否滿足log-logistic概率分布,選用Kolmogorov-Smirno(K-S)檢驗比較經驗分布與假設分布的擬合優度,選用納什效率系數(Nash-Sutcliffe efficiency coefficient,NSE)比較經驗分布與假設理論分布的接近程度。NSE越接近于1,表示二者越接近。
四川省各區域ET0計算結果見圖2。由圖2可知,在川西高原,PM法多年平均ET0為931.89 mm,最大值為977.82 mm,最小值為882.21 mm。FAO-Ra法、PT法、BC法、HS法均高估了ET0,其中FAO-Ra法偏差較大,多年ET0變化范圍為1 067.21~1 195.46 mm,相比PM法偏高166.33~234.92 mm,PT法與PM法表現出較好的吻合度,多年均值偏差為72.92 mm,偏差較小;HS法、BC法均值則分別偏大104.26、71.30 mm;MK法、WMO法均低估了ET0,其中WMO法偏低較多,偏差范圍為418.57~483.04 mm,MK法則偏低88.78~151.09 mm;Ro法年際變化不穩定,多年均值比PM法偏低103.66 mm。

圖2 四川省各區域ET0計算結果Fig.2 Calculation results of ET0 in each region of Sichuan Province
在川西南山地,PM法多年平均ET0為1 036.46 mm,最大值為1 132.88 mm,最小值為963.02 mm。BC法在川西南山地計算值最高,多年ET0變化范圍為1 245.12~1 322.58 mm,相比PM法偏高158.46~322.44 mm,偏差較大;FAO-Ra法、HS法、PT法與PM法年際變化趨勢一致,但FAO-Ra法出現了較大的偏差,多年均值比PM法偏高155.01 mm,HS法、PT法多年均值則分別偏高121.71、86.44 mm;MK法、WMO法同樣低估了ET0,其中WMO法偏低最多,為392.82~546.74 mm,MK法則偏低 98.25~198.18 mm;Ro法年際變化不穩定,多年均值比PM法偏低75.23 mm。
在川中盆地,各方法計算得到的ET0差異顯著。PM法計算得到的多年平均ET0為866.38 mm,最大值為953.31 mm,最小值為795.47 mm。PT法、HS法、BC法均比PM法計算值高,尤其BC法多年ET0變化范圍在1 369.22~1 451.76 mm之間,相比PM法偏高469.30~594.48 mm,精度最差。PT法、HS法與PM法年際變化基本一致,多年平均ET0較PM法分別偏高123.73、203.68 mm;輻射法中的FAO-Ra法與PM法計算值最為接近,變化趨勢一致,其多年平均ET0為876.75 mm,偏差范圍在-19.86~60.15 mm之間,計算結果最優;MK法、WMO法的ET0均比PM法偏低,其中WMO法低估程度明顯較大,偏差范圍為363.35~504.46 mm,在該地區適用性最差,MK法則偏低155.10 mm;Ro法年際變化不穩定,多年均值比PM法偏低93.39 mm。
利用RMSE來評價四川省各簡便方法年ET0差異。各方法在不同區域誤差變化見圖3。從圖3來看,在川西高原,PT法和BC法表現出較好的計算精度,RMSE均在76.93 mm以下;HS法、Ro法、MK法的RMSE在102.09~116.62 mm之間,FAO-Ra法的RMSE高達200.59 mm;WMO法計算精度最差,RMSE為448.25 mm。川西南山地PT法和Ro法的RMSE低于99.11 mm,計算結果較好,其次為HS法、FAO-Ra法、MK法,RMSE變化范圍在125.13~157.09 mm之間,BC法、WMO法均表現出較大的誤差,其中WMO法的RMSE為479.90 mm,精度最差。在川中盆地,FAO-Ra法RMSE最小,僅為18.66 mm,其次為Ro法、PT法、MK法、HS法,RMSE變化范圍為123.62~204.74 mm,WMO法和BC法誤差較大,RMSE均在454.08 mm以上。

圖3 四川省各區域ET0計算方法的RMSEFig.3 RMSE of ET0 calculation methods in each region of Sichuan Province
利用ArcGIS中的克里金插值法對四川省34個氣象站點年ET0的RE進行空間插值(圖4),分析不同方法在四川省的空間適用性。

圖4 四川省各區域ET0的RE空間分布Fig.4 RE spatial distributions of ET0 in each region of Sichuan Province
由圖4可以看出,在川西高原,PT法的適用性最好,大部分氣象站點的RE介于-3.8%~14.2%之間,只有高原北部的色達、若爾蓋及紅原一帶的RE相對較大,在16.8%~19.3%間變化。而WMO法的RE最差,在整個川西高原地區都低估較多,各站點低估區間為-65.1%~-25.7%,其中北部的若爾蓋、紅原、色達及石渠地區更是低估-58.4%以上,適用性最差。在川西南山地,PT法的適用性仍然最好,除了昭覺和雷波的RE在15.2%左右外,其余站點RE僅在-1.3%~13.3%變化;而WMO法在整個川西南地區都低估較多,各站點相對低估-64.5%~-28.7%,計算精度最差。在川中盆地,FAO-Ra法的適用性最佳,除了巴中的RE高于6.4%外,其余各站點的RE都相對較小,介于-2.2%~3.3%之間,而WMO法在整個川中盆地都低估較多,除了萬源和雅安外,其余各站點均低估-62.6%~-48.1%,評價效果最差。
由于四川省不同地區存在氣候差異,導致各方法在不同地區應用時差異較大。HS法在緯度較低、平均溫度高、輻射強的川西南山地ET0最大,這與HS法主要考慮平均溫度和晝夜溫差影響,同時利用太陽天頂輻射進行計算有關。BC法隨著區域海拔降低,其數值逐漸偏大,這與區域海拔降低,氣溫逐漸升高有關。PT法忽略了空氣動力學項,僅考慮了輻射項,在四川省3個區域均高估了ET0,這與符娜等[28]對云南地區進行研究得出PT法高估ET0的結論一致。由于忽略了土壤熱量,MK法在四川省3個區域均低估了ET0。RAHIMIKHOOB等[29]對伊朗北部8個氣象站點研究也發現,MK法在所有站點均低估了ET0,左德鵬等[30]、李志[31]也得出類似結論。在濕潤的川中盆地,FAO-Ra法的ET0計算值最接近PM法,這是因為FAO-Ra法在濕潤區的計算效果較優。NANDAGIRI等[32]評價了印度4個不同氣候區的氣象站點逐日ET0,發現FAO-Ra法在濕潤區的計算效果優于干旱區。樊軍等[33]的研究也表明FAO-Ra法在半干旱地區估算ET0較差。Ro法與WMO法同屬于物質傳輸法,其計算精度對使用的參數變化非常敏感[26]。盡管Ro、WMO法在3個區域計算結果一致偏小,但Ro法偏差較小,WMO明顯偏差較大,這說明Ro法的初始參數在四川省具有較好的適用性,而WMO法的參數需要根據天氣情況進行調整。
從各項誤差指標來看,MK、Ro法在四川省3個區域的表現較為穩定,PT法在川西南山地和川西高原均表現出較低的誤差,WMO法在各區域適用性最差,綜合考慮,選取PT、MK、Ro 3種簡便方法計算的ET0進行SPEI分析。
3.3.1SPEI指數理論驗證
計算SPEI指數的前提是降水量和蒸散量之差(Di)服從三參數對數邏輯分布。因此,在計算指數之前,本文通過K-S檢驗和NSE檢驗探究四川省各氣象站點的Di是否服從以上概率分布,結果如表2所示。

表2 Di分布檢驗結果Tab.2 Distribution test results of Di
由表2可知,采用三參數的對數邏輯分布函數對Di進行擬合,通過K-S檢驗和NSE評估擬合效果,結果發現34個站點的Di均通過K-S檢驗,且NSE基本都在0.95以上,表明三參數對數邏輯分布函數在四川省各月份不同時間尺度下的Di擬合效果都較好,可以用于SPEI指數的多時間尺度計算。
3.3.2時間變化
對四川省不同區域1967—2016年SPEI_PM、SPEI_PT、SPEI_MK、SPEI_Ro的年值進行比較(圖5),其中紅色圖例表示實際旱情。

圖5 四川省1967—2016年12個月時間尺度SPEI_PM、SPEI_PT、SPEI_MK、SPEI_Ro年際變化曲線Fig.5 Temporal variations of SPEI_PM, SPEI_PT, SPEI_MK and SPEI_Ro at 12-month timescale in Sichuan Province during 1967—2016
由圖5可知,在四川省的3個區域中,SPEI_PT、SPEI_MK、SPEI_Ro具有與SPEI_PM相似的變化情況,其變化值均處于-2~2之間。從四川省各區域的SPEI計算結果來看,基于不同ET0方法計算的年SPEI差異較小,變化趨勢基本一致,4種方法均顯示1967—2016年,川西高原、川西南山地的SPEI有增加的趨勢,川中盆地則存在SPEI減小的趨勢。雖然SPEI_PM、SPEI_PT、SPEI_MK、SPEI_Ro最小值不同,但時間分布一致。在有實際旱情的年份,最小值均低于0,因此4種SPEI指數都能夠識別歷史干旱事件。總的來說,在同一區域,不同方法之間會有差異,但SPEI幾乎都維持在同一水平,干濕情況交替出現。
3.3.3誤差分析
圖6為四川省1967—2016年間SPEI_PM、SPEI_PT、SPEI_MK、SPEI_Ro的泰勒圖,其中點A為觀測點,代表SPEI_PM。由圖6可知,在四川省各區域,12個月時間尺度下的SPEI_PT、SPEI_MK、SPEI_Ro均與SPEI_PM有較好的相關系數,這是因為SPEI是通過歷史同期數據的異常來監測干旱,主要與ET0變化情況有關,而受ET0本身的影響較小。從標準偏差來看,SPEI_PT、SPEI_MK、SPEI_Ro在川西高原、川西南山地差別不大,都較接近1,而在川中盆地,SPEI_Ro的標準偏差明顯大于SPEI_PT、SPEI_MK,達到1.14。從RMSE來看,SPEI_MK與SPEI_PM的RMSE最小,其次依次為SPEI_PT、SPEI_Ro。綜合來看,在四川省全域,SPEI_MK與SPEI_PM有更好的相關性。

圖6 四川省12個月時間尺度SPEI_PM、SPEI_PT、SPEI_MK、SPEI_Ro的泰勒圖Fig.6 Taylor diagrams of SPEI_PM, SPEI_PT, SPEI_MK and SPEI_Ro at 12-month timescale in Sichuan Province
相較于川西高原和川西南山地,川中盆地各方法計算的SPEI更具敏感性。對川中盆地不同時間尺度下SPEI_PT、SPEI_MK、SPEI_Ro的誤差進行分析,如圖7所示。由圖7可知,1、3、12個月時間尺度下,SPEI_PT、SPEI_MK與SPEI_PM都具有較好的相關性,而SPEI_Ro與SPEI_PM的相關性最差。從RMSE來看,1、3、12個月時間尺度SPEI_PT、SPEI_MK與SPEI_PM的RMSE大都集中在0.05~0.20區間內,而12個月時間尺度SPEI_Ro與SPEI_PM的RMSE達到0.44。綜合來看,大部分情況下,SPEI1相關性比SPEI3、SPEI12好。1個月時間尺度下,SPEI_MK與SPEI_PM有更好的相關性,相關系數達到0.99,RMSE僅為0.15,SPEI_Ro與SPEI_PM相關性最差。

圖7 川中盆地1、3、12個月時間尺度SPEI_PM、SPEI_PT、SPEI_MK、SPEI_Ro的泰勒圖Fig.7 Taylor diagrams of SPEI_PM, SPEI_PT, SPEI_MK and SPEI_Ro at 1-, 3- and 12-month timescales in central Sichuan basin
由于川中盆地各方法計算的SPEI更具敏感性,因此對川中盆地不同時間尺度下SPEI_PT、SPEI_MK、SPEI_Ro的小波實部及方差進行分析,其對比結果如圖8、9所示。圖8中等值線的正負代表偏濕潤和偏干旱時期,圖9中小波方差反映SPEI時間序列的能量波動隨時間尺度的變化特征,峰值對應的時間尺度為SPEI變化的主要周期。

圖8 川中盆地1、3、12個月時間尺度SPEI小波分析Fig.8 Wavelet analysis of SPEI at 1-, 3- and 12-month timescales in central Sichuan basin

圖9 川中盆地1、3、12個月時間尺度SPEI小波方差Fig.9 Wavelet variances of SPEI at 1-, 3- and 12-month timescales in central Sichuan basin
由圖8可以看出,在時間尺度10~15 a上,各方法的SPEI1周期振蕩變化明顯,貫穿整個研究時段,經歷了5次“偏干旱-偏濕潤”的交替變化,2016年等值線小于0且向負值中心變化,說明未來一段時間將處于偏干旱時期。在時間尺度5~10 a、15~20 a上,各方法的SPEI3周期振蕩變化明顯。結合圖9發現,SPEI3的小波方差峰值出現在9 a的3.5和 19 a的5.5,分別對應第2主周期、第1主周期。在時間尺度10~15 a上,各方法的SPEI12周期振蕩變化明顯,經歷了5次“偏干旱-偏濕潤”的交替變化。在1個月時間尺度上,SPEI_PT表現出與SPEI_PM較好的一致性,兩者第1主周期都位于13 a,小波方差差距最小,SPEI_Ro則表現出與SPEI_PM較為不同的周期性,其第1主周期位于11 a,與SPEI_PM差距最大。在3、12個月時間尺度上,各方法的SPEI都具有相似的第1主周期??傮w而言,SPEI_PT和SPEI_PM具有最相似的周期振蕩變化,SPEI_Ro和SPEI_PM的周期差距最大。
(1)PT法等7種ET0算法的空間變異性十分明顯。在川西高原和川西南山地,PT法與PM法有較好的相關性,適用性較強;在川中盆地,輻射法中的FAO-Ra法適用性最好。
(2)結合歷史干旱事件分析,SPEI_PM、SPEI_PT、SPEI_MK、SPEI_Ro 4種指數都能夠識別出四川省大多數干旱事件,具有較好的干旱監測能力,其中SPEI_MK監測的干濕情況最接近SPEI_PM。
(3)大多數情況下,SPEI1相關性比SPEI3、SPEI12好。1個月時間尺度下,SPEI_MK與SPEI_PM有更好的相關性,而SPEI_Ro與SPEI_PM相關性最差。
(4)SPEI_PT和SPEI_PM具有最相似的周期振蕩變化,SPEI_Ro和SPEI_PM的周期差距最大。
(5)四川省SPEI_MK監測的干濕情況最接近SPEI_PM,數據缺失的條件下,SPEI_MK可以作為SPEI_PM的替代指數,這為四川省的干旱監測和防控提供理論依據,能夠減少和預防干旱帶來的災害。由于本文僅研究了ET0計算方法及相應SPEI指數在四川省的適用性,下一步可深入探究不同SPEI指數監測的氣象干旱是否會在干旱傳播過程中對水文干旱、農業干旱造成影響。