李 想,李國平
(1.成都信息工程大學(xué)大氣科學(xué)學(xué)院,四川 成都 610225;2.氣象災(zāi)害預(yù)報預(yù)警與評估省部共建協(xié)同創(chuàng)新中心,江蘇 南京 210044)
降水是天氣預(yù)警預(yù)報、氣候預(yù)測預(yù)估的最重要的氣象要素之一[1-2]。同時獲取精度更高的降水資料,可以更準(zhǔn)確認(rèn)識降水的空間分布[3-6]。
更加精確的空間降水?dāng)?shù)據(jù),除了完善觀測設(shè)備和觀測站網(wǎng)之外,也依賴于更高質(zhì)量的空間插值方法。插值技術(shù)的不斷成熟改進(jìn)了氣象臺站有限和空間分布不均的不足,但實際業(yè)務(wù)中沒有一種插值方法是普適性的[7],一般根據(jù)研究目的和區(qū)域地理特征采用適宜的插值方法。為尋求最佳插值方法,已有研究對同一區(qū)域的降水?dāng)?shù)據(jù)采用多種空間插值方法進(jìn)行比較,例如克里金插值法、樣條函數(shù)法、趨勢面法等方法[8-10]。目前國內(nèi)主要根據(jù)中國高密度氣象站點的信息采用不同插值法[11]。例如中國地面降水0.5°×0.5°格點數(shù)據(jù)集(V2.0)(1961—2013)利用薄盤樣條法,并且引入數(shù)字高程模型(digital elevation model,DEM)以盡可能地消除中國區(qū)域獨特地形條件下高程對降水空間插值精度的影響。同時由于協(xié)同克里金方法能夠引入其余影響因子,也經(jīng)常被用于業(yè)務(wù)中[12-14]。例如安徽省的降水插值中效果最佳的是引入經(jīng)緯度的協(xié)同克里金方法[15]??紤]DEM數(shù)據(jù)的協(xié)同克里金方法提高了湖南省的降水插值精度[16]。綜合來看,協(xié)同克里金方法在降水空間插值上具有較佳的表現(xiàn)。
山區(qū)的降水空間插值更為復(fù)雜。傳統(tǒng)插值依賴于樣本空間密度,缺少氣象要素的演變機理,忽略了數(shù)據(jù)的空間分布特征和地形高度等影響因素[17]。對于地形復(fù)雜且站點分布不均的山區(qū),降水空間插值的結(jié)果往往不太理想[18]。四川山區(qū)受熱帶季風(fēng)、副熱帶季風(fēng)以及青藏高原環(huán)流系統(tǒng)的影響,夏季暴雨頻發(fā),極易造成泥石流、山體滑坡等次生災(zāi)害。因此,獲取四川山區(qū)夏季精確的降水?dāng)?shù)據(jù)對防震減災(zāi)具有重要意義。本文重點針對傳統(tǒng)插值效果較差的四川地區(qū),分區(qū)域?qū)ふ腋鲄^(qū)域最優(yōu)的降水影響因子組合并且對不同方法的插值結(jié)果進(jìn)行比較,找出最優(yōu)插值方法。
四川省(97°21′E—108°12′E,26°03′N—34°19′N)位于我國西南地區(qū)內(nèi)陸,地處長江上游,總面積約4.86×105km2,地跨青藏高原、橫斷山脈、四川盆地等地貌單元,地勢西高東低。西部是青藏高原東南緣和橫斷山脈的一部分,海拔為3000~4500 m;東部由盆地底部和盆地邊緣組成,底部,海拔300~700 m,由成都平原、眉山-峨眉平原組成;邊緣地區(qū)以山地為主,海拔為1500~3000 m(圖1)。
圖1 四川省自動氣象站分布Fig.1 Distribution of automatic weather stations in Sichuan Province
所用資料:(1)四川省氣象局共157個自動氣象站2010—2019年6—8月逐小時降水資料,用于計算10 a的多年平均夏季降水量;(2)DEM數(shù)據(jù)來源于先進(jìn)星載熱發(fā)射和反射輻射儀全球數(shù)字高程模型(advanced spaceborne thermal emission and reflection radiometer global digital elevation model,ASTER GDEM)30 m空間分辨率數(shù)據(jù);(3)NDVI數(shù)據(jù)來源于LAADS(Level-1 and Atmosphere Archive & Distribution System)DAAC(Distributed Active Archive Center)提供的空間分辨率為250 m的MOD13Q1第6版產(chǎn)品,此次使用的數(shù)據(jù)是的植被指數(shù)(VI)值,有兩個主要植被層,本文使用的是第一個歸一化差異植被指數(shù)(normalized difference vegetation index,NDVI)。
首先對四川地區(qū)進(jìn)行區(qū)域劃分,并對各因子與研究區(qū)和各分區(qū)的多年平均夏季降水量進(jìn)行相關(guān)性分析,找出每個區(qū)域主要的影響因子,再將這些因子進(jìn)行組合,與降水量進(jìn)行多元線性回歸分析,找出每個區(qū)具有最佳擬合效果的因子組合;其次分別對各區(qū)進(jìn)行插值,并對結(jié)果進(jìn)行交叉檢驗,找出適合各區(qū)域的最佳插值方法,同時根據(jù)插值結(jié)果,給出各區(qū)域降水主導(dǎo)的地理影響因子。
共使用6種插值方法,分別是反距離加權(quán)(inverse distance weighted,IDW)、徑向基函數(shù)(radical basis function,RBF)、普通克里金(ordinary Kriging,OK)、協(xié)同克里金金(CoKriging,CoK)、局部多項式(local polynomial interpolation,LPI)和經(jīng)驗貝葉斯克里金(empirical Bayesian Kriging,EBK)插值。其中,IDW方法運算快、效率高,但外推能力差,適用于站點分布盡可能均勻且布滿整個插值區(qū)域的樣本數(shù)據(jù)集[19-20]。RBF工作量小且精度相對較高,適用于樣本數(shù)據(jù)集大,地形平緩的情況[21-22]。OK法又稱空間局部插值法,該方法計算速度較慢,適用于區(qū)域化變量存在空間相關(guān)性的區(qū)域[23-24]。OK是應(yīng)用最廣的克里金插值方法,而CoK是OK的擴(kuò)展形式,它將主變量的空間自相關(guān)性和主輔變量之間的交互相關(guān)性結(jié)合起來,用于無偏最優(yōu)估值中[19,25]。LPI是一種局部加權(quán)最小二乘擬合法,多用于解釋局部變異現(xiàn)象、建立平滑表面和確定變量的小范圍變異[26-27]。EBK是一種地統(tǒng)計插值方法,它與OK方法不同,是通過估計基礎(chǔ)版變異函數(shù)來說明所引入的誤差,因此大大降低了預(yù)測的標(biāo)準(zhǔn)誤差[21]。
運用ArcGIS中的交叉驗證法對各種插值法進(jìn)行誤差分析。其基本思想是將原始數(shù)據(jù)隨機地分為訓(xùn)練集和驗證集,先對訓(xùn)練集進(jìn)行訓(xùn)練,再利用驗證集來測試訓(xùn)練得到的結(jié)果,以此作為評價指標(biāo)[28-29]。選取平均誤差、平均絕對誤差、均方根誤差和綜合相對誤差作為評判標(biāo)準(zhǔn)。
聚類分析是研究多要素(多個變量)的客觀分類方法。它的聚類原則是根據(jù)某些相似性的指標(biāo)進(jìn)行聚類,把對象的個體(樣品)進(jìn)行聯(lián)合,用分裂或添加的方法進(jìn)行聚類或串組,故也稱串組分析。根據(jù)四川各縣區(qū)的經(jīng)緯度以及海拔高度,運用聚類分析的7種方法(最遠(yuǎn)鄰元素、最近鄰元素、質(zhì)心連接、組內(nèi)連接、組間連接、中位數(shù)以及快速聚類)并均采用平方歐氏距離進(jìn)行篩選比對,篩選條件包括最終分得區(qū)域數(shù)目合適、分得同一區(qū)域所在位置相對聚集、區(qū)域站點數(shù)目合適并平均、分得區(qū)域地形特征統(tǒng)一。最終選定組間連接方法將四川地區(qū)分為4個區(qū)域,分別是區(qū)域1(南部地區(qū))、區(qū)域2(東北部地區(qū))、區(qū)域3(西北部地區(qū))、區(qū)域4(中東部地區(qū))(圖2)。
圖2 四川省分區(qū)結(jié)果Fig.2 The results after the division in Sichuan Province
將四川全區(qū)及各區(qū)多年平均夏季降水量與各因子進(jìn)行相關(guān)分析,相關(guān)系數(shù)如表1所示??梢钥闯?,四川各區(qū)多年平均夏季降水量與6種因子均相關(guān),其中坡向因子與全區(qū)和各區(qū)相關(guān)性最差且未通過α=0.05的顯著性檢驗。全區(qū)多年平均夏季降水量與海拔相關(guān)系性最好,為-0.461,且通過α=0.01的顯著性檢驗。四川包含多種復(fù)雜地貌,使得夏季降水整體上受海拔影響很深,其次是坡度,如迎風(fēng)坡、背風(fēng)坡、陡坡、緩坡等對降水的影響占比較高。區(qū)域1與經(jīng)度相關(guān)性最好,且通過α=0.05的顯著性檢驗,其次是坡度、海拔、NDVI,但都沒有通過α=0.05的顯著性檢驗。區(qū)域2與海拔和坡度的相關(guān)性通過α=0.01的顯著性檢驗,這是由于區(qū)域2西側(cè)小部分區(qū)域為盆地邊緣地區(qū),地形起伏較大,地形會影響到當(dāng)?shù)氐慕邓?。區(qū)域3與NDVI的相關(guān)性最好,通過α=0.01的顯著性檢驗,說明在川西高原地區(qū),植被覆蓋率對山區(qū)降水具有一定的影響。區(qū)域4與經(jīng)度、緯度相關(guān)性最好,均通過α=0.01的顯著性檢驗。
表1 四川全區(qū)及各區(qū)多年平均夏季降水量與各因子相關(guān)系數(shù)Tab.1 The correlation coefficients between multi-year average summer precipitation in Sichuan and various factors in the whole region and each district
綜合全區(qū)和4個分區(qū)來看,坡向與全區(qū)與各分區(qū)的相關(guān)性都很低且沒有通過顯著性檢驗,所以在因子組合中去除坡向這一因子。
CoK只能加入3個輔助變量,從5個因子(經(jīng)度、緯度、海拔、NDVI、坡度)中任意選3個因子作為一個組合,總共有10種組合。以每3個因子作為自變量,多年平均夏季降水量作為因變量,進(jìn)行多元線性回歸分析,結(jié)果如表2所示。對于全區(qū)而言,R2最高的三種組合方式為AlNdSl、LaAlSl和LaAlNd組合,但擬合效果均不佳,其中AlNdSl組合與多年夏季平均降水量擬合度最高(R2=0.293),R2最高的三種組合方式中均有海拔因子,說明海拔對四川夏季降水影響很大。對于區(qū)域1來說,不同組合差距懸殊,R2最高的三種組合方式為LoNdSl、lLoLaSl和LoAlSl,其中最佳擬合組合為LoNdSl,R2=0.556。區(qū)域2,10種組合得出的擬合效果差距不大,R2最高的三種組合方式為LoLaAl、LoNdSl和LoLaSl,其中最佳擬合方式是LoLaAl的組合,R2=0.397。綜合來看,區(qū)域3擬合效果較好,R2最高的三種組合方式為LoAlNd、AlNdSl和LaAlNd,其中R2最高達(dá)0.625,為LoAlNd的組合。區(qū)域4總體而言擬合效果不如前三個區(qū)域的擬合效果好,10組組合中,R2最高的三種組合方式為LoLaNd、LoLaSl和LoLaAl,其中擬合效果較佳的是LoLaNd組合,R2=0.267,其次是Lo-LaSl和LoLaAl,3組組合都含有經(jīng)度和緯度,說明四川中東部地區(qū)夏季降水受地理位置影響比較大。
表2 全區(qū)及各區(qū)因子組合與多年平均夏季降水量相關(guān)性R2Tab.2 The correlation R2 between the combination factors and the multi-year average summer precipitation in each area and the whole area
圖3為全區(qū)多年平均夏季降水插值結(jié)果比較。從四川全區(qū)來看,插值結(jié)果整體呈現(xiàn)為西北小、東南大的空間分布,即盆地區(qū)域降水多,高原地區(qū)降水少。川中東部地區(qū)出現(xiàn)條狀的大值區(qū),為西北—東南向的橢圓形,該區(qū)域西部、西北部和南部的高山形成喇叭狀地形,使得東來的太平洋東南暖濕氣流與盆地周邊山地下沉的冷濕氣流交匯于此處,形成著名的“華西雨屏”現(xiàn)象[30]。不同的空間插值方法在全區(qū)具有一致性,但在局地存在一定差異,OK、3種組合的CoK插值結(jié)果表現(xiàn)為川東北地區(qū)的大值區(qū)與川南北部地區(qū)的部分大值區(qū)相連,而其他插值結(jié)果明顯形成了斷裂并且在盆地及其東部地區(qū)出現(xiàn)小范圍的“牛眼”現(xiàn)象,這種現(xiàn)象是IDW常出現(xiàn)的一種現(xiàn)象,這些“牛眼”能突出站點的特征,并在一定程度上提高了插值精度。OK插值以及3種組合的CoK插值差異肉眼幾乎難以分辨,必須后續(xù)對插值結(jié)果進(jìn)行更進(jìn)一步的交叉檢驗來判斷幾種方法的優(yōu)劣。
圖3 全區(qū)多年平均夏季降水插值結(jié)果比較(單位:mm)(a)OK,(b)IDW,(c)RBF,(d)AlNdSl組合CoK,(e)LaAlSl組合CoK,(f)LaAlNd組合CoK,(g)LPI,(h)EBKFig.3 Comparison of interpolated results of multi-year average summer precipitation in the whole region(Unit:mm)(a)OK,(b)IDW,(c)RBF,(d)AlNdSl CoKriging interpolation,(e)LaAlSl CoKriging interpolation,(f)LaAlNd CoKriging interpolation,(g)LPI,(h)EBK
圖4為區(qū)域1多年平均夏季降水插值結(jié)果比較。區(qū)域1大值區(qū)主要集中在其北部地區(qū),且有深入其中部地區(qū)的趨勢,東南地區(qū)也有小范圍大值區(qū),而東、西邊緣地區(qū)值較小。IDW、RBF以及EBK插值結(jié)果仍有部分“牛眼”現(xiàn)象。除OK、LPI和LoLa-Sl組合CoK外,其余插值均在東南邊緣地區(qū)有一大值區(qū)。LPI、LoLaSl組合CoK以及OK、EBK插值結(jié)果表現(xiàn)為北部大值區(qū)較為分散,其余插值則較為集中。
圖4 區(qū)域1多年平均夏季降水插值結(jié)果比較(單位:mm)(a)OK,(b)IDW,(c)RBF,(d)LoNdSl組合CoK,(e)LoLaSl組合CoK,(f)LoAlSl組合CoK,(g)LPI,(h)EBKFig.4 Comparison of interpolated results of multi-year average summer precipitation in District 1(Unit:mm)(a)OK,(b)IDW,(c)RBF,(d)LoNdSl CoKriging g interpolation,(e)LoLaSl CoKriging interpolation,(f)LoAlSl CoKriging interpolation,(g)LPI,(h)EBK
區(qū)域2插值結(jié)果表現(xiàn)為整體從西北到東南呈小、大、小分布。高值區(qū)均呈帶狀,主要位于盆地北部邊緣地區(qū),低值區(qū)不同插值方法各不相同但均在南部地區(qū)有一個大范圍的較低值區(qū)域。除3種組合的CoK以及OK插值外,其他插值在阿壩州的茂縣、綿陽的北川地區(qū)均有明顯的小值區(qū)(圖5)。
圖5 區(qū)域2多年平均夏季降水插值結(jié)果比較(單位:mm)(a)OK,(b)IDW,(c)RBF,(d)LoLaAl組合CoK,(e)LoNdSl組合CoK(f)LoLaSl組合CoK,(g)LPI,(h)EBKFig.5 Comparison of interpolated results of multi-year average summer precipitation in District 2(Unit:mm)(a)OK,(b)IDW,(c)RBF,(d)LoLaAl CoKriging interpolation,(e)LoNdSl CoKriging interpolation,(f)LoLaSl CoKriging interpolation,(g)LPI,(h)EBK
區(qū)域3海拔總體較高,降水主要集中在其東南部的九龍一帶。RBF、EBK與AlNdSl組合CoK插值結(jié)果表現(xiàn)為在阿壩州馬爾康和壤塘地區(qū)出現(xiàn)相連帶狀區(qū)域,而OK、IDW以及LPI插值則為較為分散的點狀。LoAlNd組合CoK插值在甘孜州石渠縣的值明顯高于其余幾種插值,并且LoAlNd和LaAlNd組合的CoK插值在大部分地區(qū)都為大值區(qū)域(圖6)。
圖6 區(qū)域3多年平均夏季降水插值結(jié)果比較(單位:mm)(a)OK,(b)IDW,(c)RBF,(d)LoAlNd組合CoK,(e)AlNdSl組合CoK,(f)LaAlNd組合CoK,(g)LPI,(h)EBKFig.6 Comparison of interpolated results of multi-year average summer precipitation in District 3(Unit:mm)(a)OK,(b)IDW,(c)RBF,(d)LoAlNd CoKriging interpolation,(e)AlNdSl CoKriging interpolation,(f)LaAlNd CoKriging interpolation,(g)LPI,(h)EBK
區(qū)域4各種插值方法的插值結(jié)果差異不大,高值多集中在其西北部,低值分散在東部地區(qū)。除3種組合的CoK以外,其他插值方法在其西北部插值結(jié)果呈現(xiàn)出塊狀明顯、邊緣清晰的特征。LPI在宜賓筠連附近出現(xiàn)大值區(qū),EBK、OK、IDW、RBF在筠連的大值區(qū)不太明顯,而3種組合CoK在此地沒有表現(xiàn)出大值區(qū)域(圖7)。
圖7 區(qū)域4多年平均夏季降水插值結(jié)果比較(單位:mm)(a)OK,(b)IDW,(c)RBF,(d)LoLaNd組合CoK,(e)LoLaSl組合CoK,(f)LoLaAl組合CoK,(g)LPI,(h)EBKFig.7 Comparison of interpolated results of multi-year average summer precipitation in District 3(Unit:mm)(a)OK,(b)IDW,(c)RBF,(d)LoLaNd CoKriging interpolation,(e)LoLaSl CoKriging interpolation,(f)LoLaAl CoKriging interpolation,(g)LPI,(h)EBK(Unit:mm)
總體而言,同區(qū)域不同插值方法得出的結(jié)果差異不大,插值結(jié)果大體走向與增減幾乎一致,要得到更精確的對照分析,必須對結(jié)果進(jìn)行交叉檢驗。
表3為全區(qū)與各區(qū)多年平均夏季降水交叉檢驗結(jié)果。對于全區(qū)而言最好的插值方法是EBK,綜合相對誤差(每個區(qū)域內(nèi)點的相對誤差求得的平均值)只有9.11%,其次是RBF和OK。其原因是在大量樣點數(shù)據(jù)中,RBF能減少誤差,而OK在全區(qū)更能發(fā)揮優(yōu)勢。全區(qū)3種組合的CoK均沒有達(dá)到預(yù)期效果,原因可能是全區(qū)的降水是多種因子復(fù)合影響的結(jié)果,對于僅有的3種因子的協(xié)同插值達(dá)不到預(yù)期的效果。但在3種組合CoK中,AlNdSl組合方式最優(yōu),說明對于全區(qū)來說,海拔、NDVI和坡度對于降水的影響較大。
區(qū)域1的插值結(jié)果為幾個區(qū)中最佳,其中最佳插值方法是LoNdSl組合CoK。RBF和IDW插值平均誤差雖小,但從標(biāo)準(zhǔn)均方根誤差來看略差于其余插值方法,說明這兩種插值方法的誤差極值效應(yīng)與CoK相比較大。3種組合的CoK效果均較好,說明此區(qū)域在插值方法中引入最能影響此區(qū)域的環(huán)境因子能夠很大程度上提高插值結(jié)果的精確性。區(qū)域1的降水主要受經(jīng)度和坡度影響,結(jié)合區(qū)域1地勢來看,區(qū)域1處于山地,降水受迎風(fēng)坡或背風(fēng)坡的影響較大。
區(qū)域2插值結(jié)果最佳的是EBK,其次是RBF,這是由于區(qū)域2站點多而密集,有利于RBF的插值。3種組合的CoK在區(qū)域2的插值效果不好,其原因可能是多方面的:一是該區(qū)域大部分為平原地區(qū),地勢開闊,海拔平均較低,能影響該區(qū)域的因子較為單一,用3種因子組合的CoK方法反而會降低插值精度;二是此區(qū)域站點較多,對于RBF和EBK來說更為適合。
區(qū)域3插值效果最佳的是AlNdSl組合CoK,其次是RBF。區(qū)域3處于川西高原地區(qū),海拔較高且起伏較多,另外區(qū)域3降水與NDVI相關(guān)性很高,超過了0.5且通過了α=0.01的顯著性檢驗,所以引入AlNdSl的協(xié)同克里金在區(qū)域3的插值效果最佳。
區(qū)域4插值效果最佳的是RBF和EBK,相對誤差均在10%左右,3種組合的CoK插值效果一般,但比OK插值效果佳,說明此區(qū)域的降水確實由多種因子相互影響,僅3種因子的協(xié)同插值不足以體現(xiàn)出該區(qū)域的降水分布。若將此區(qū)域更細(xì)致地劃分,會得到更好的插值結(jié)果。
綜合比較全區(qū)與4個區(qū)域所有插值方法,RBF和EBK的綜合相對誤差相對于其他插值方法更小,其原因是選用四川全區(qū)時站點密集,有利于RBF插值,并且EBK自帶減少誤差的子程序;其次為OK和Cok,在選取因子合適的區(qū)域,Cok插值可以達(dá)到更佳的插值效果,而在影響因子眾多或者地形單一的區(qū)域,RBF和EBK更勝一籌??傮w而言IDW插值效果最差,分區(qū)之后所有插值方法的插值精度要高于分區(qū)前的精度。
基于四川省157個自動氣象站點近10 a(2010—2019)夏季降水?dāng)?shù)據(jù),采用聚類分析進(jìn)行分區(qū)后通過相關(guān)性分析和多元回歸分析篩選出各區(qū)域降水量的地理影響因子。使用協(xié)同克里金插值方法的同時,采用傳統(tǒng)插值方法進(jìn)行對比并對插值結(jié)果進(jìn)行交叉檢驗,結(jié)論如下:
(1)影響全區(qū)降水的主要地理因子為海拔、坡度、經(jīng)度和緯度,主導(dǎo)的影響因子組合是海拔、NDVI和坡度;分區(qū)細(xì)化之后,對南部地區(qū)降水影響的主要因子是經(jīng)度、坡度和海拔,主導(dǎo)的影響因子組合是經(jīng)度、NDVI和坡度;對東北部地區(qū)而言,降水主要因子為坡度,其次是海拔,其主導(dǎo)的因子組合方式是經(jīng)度、緯度和海拔;西北地區(qū)的降水主要因子是NDVI和經(jīng)度,其主導(dǎo)的因子組合方式為經(jīng)度、海拔和NDVI;中東部地區(qū)的主要降水影響因子是經(jīng)度和緯度,其主導(dǎo)的因子組合方式是經(jīng)度、緯度和NDVI。
(2)全區(qū)而言,最佳的插值方式是經(jīng)驗貝葉斯克里金插值;對四川南部地區(qū)而言,最佳插值方法是經(jīng)度、坡度和NDVI因子組合的協(xié)同克里金插值;對四川東北部地區(qū)而言,最佳插值方法是經(jīng)驗貝葉斯函數(shù);對四川西北部地區(qū)而言,最佳的插值方式是海拔、NDVI和坡度為主導(dǎo)因子組合的協(xié)同克里金插值函數(shù);對四川中東部地區(qū)而言,最佳的插值方式是徑向基函數(shù)插值方法。
(3)分區(qū)后插值精度高于分區(qū)前插值精度,在所選區(qū)域降水影響因素數(shù)目適中時,選用協(xié)同克里金插值效果更佳。所選區(qū)域降水影響因素數(shù)目單一或眾多的情況下,選用經(jīng)驗貝葉斯克里金插值或徑向基函數(shù)的效果更佳。