焦思佳，吳田軍，董世英，王長鵬

(長安大學(xué) 理學(xué)院，陜西 西安 710064)

0 引 言

20世紀(jì)60年代，以衛(wèi)星定位系統(tǒng)、地理信息系統(tǒng)和遙感為支撐的空間信息技術(shù)逐漸發(fā)展起來，大量帶有空間位置的數(shù)據(jù)被采集、處理并應(yīng)用[1-2].相比于其他數(shù)據(jù)，空間數(shù)據(jù)因其時(shí)空相關(guān)性的特性，難以使用變量獨(dú)立性假設(shè)的經(jīng)典統(tǒng)計(jì)學(xué)方法，這使得幾何空間中的牛頓推測(cè)等方法并不適用.1970年，Toblert[3]提出“地理學(xué)第一定律”，為空間數(shù)據(jù)的分析和應(yīng)用提供了理論基礎(chǔ).在該領(lǐng)域，將離散點(diǎn)的測(cè)量數(shù)據(jù)轉(zhuǎn)化為連續(xù)數(shù)據(jù)表面的空間推測(cè)是最為重要研究問題之一.經(jīng)過幾十年的持續(xù)發(fā)展，空間推測(cè)算法不斷完善，并逐漸被應(yīng)用到土壤水質(zhì)、海洋環(huán)境、地質(zhì)勘探、空氣質(zhì)量等諸多領(lǐng)域.但隨著生產(chǎn)力的進(jìn)步，社會(huì)以及經(jīng)濟(jì)生活對(duì)地理空間信息的精細(xì)性、時(shí)效性需求進(jìn)一步凸顯，這也倒逼各領(lǐng)域?qū)臻g推測(cè)算法提出更高的要求.因此，進(jìn)一步發(fā)展空間推測(cè)模型，提升專題制圖水平，具有重要現(xiàn)實(shí)意義.

目前，典型的空間推測(cè)方法大致可分為以下四類：(1)以反距離加權(quán)(Inverse Distance Weighted, IDW)為代表的確定性推測(cè)方法.IDW[4]是一種以距離作為權(quán)重的滑動(dòng)平均加權(quán)推測(cè)方法，伴隨著實(shí)際問題數(shù)據(jù)集的復(fù)雜性，基本的IDW滿足不了空間推測(cè)需求，因此，在之后的研究中其經(jīng)過不斷改進(jìn)發(fā)展，例如，王可偉等[5]在IDW中引入圓形窗口與夾角權(quán)因子，有效地提高建模的效率與精度.(2)以克里金(Kriging)為代表的地統(tǒng)計(jì)推測(cè)方法.克里金方法是1951年南非地質(zhì)學(xué)家克里金(Krige)首次提出，后經(jīng)法國著名數(shù)學(xué)家Matheron發(fā)展深化[6].由于克里金將空間相關(guān)性考慮在內(nèi)以及使用克里金標(biāo)準(zhǔn)偏差量化推測(cè)誤差這一優(yōu)點(diǎn)，成為主流方法，隨后也有一定的擴(kuò)展，例如泛克里金(Universal Kriging, UK)[7]、具有外部漂移的克里金(Kriging with External Drift, KED)[8].劉婕[9]運(yùn)用UK推測(cè)北京市六城區(qū)預(yù)估點(diǎn)的PM2.5，并驗(yàn)證統(tǒng)計(jì)值通過F檢驗(yàn)及t檢驗(yàn).鄔春明等[10]提出基于線性動(dòng)態(tài)變化因子結(jié)合柯西變異粒子群算法對(duì)變異函數(shù)的擬合模型參數(shù)進(jìn)行最優(yōu)化估計(jì)，同時(shí)在適應(yīng)度函數(shù)中引入克里金地理權(quán)重來增強(qiáng)變量的空間相關(guān)性，有效地提高推測(cè)精度并改善變異函數(shù)擬合曲線誤差過大的問題.(3)以回歸克里金(Regression Kriging, RK)代表的組合方法.Mohanasundaram等[11]運(yùn)用RK推測(cè)預(yù)估點(diǎn)的地下水位，證明推測(cè)結(jié)果優(yōu)于其他克里金方法.當(dāng)然，機(jī)器學(xué)習(xí)(Machine Learning, ML)的不斷發(fā)展同時(shí)促進(jìn)RF與克里金組合，例如，Li等[12]將RF、廣義線性模型與地統(tǒng)計(jì)方法組合，證明這些組合方法比傳統(tǒng)模型精度更高.(4)以ML為代表的推測(cè)方法.2001年，Breiman[13]提出RF，并且說明RF適用于回歸問題，同年，Rigol等[14]首次提出在運(yùn)用神經(jīng)網(wǎng)絡(luò)(Neural Network, NN)推測(cè)時(shí)，將回歸趨勢(shì)與空間關(guān)聯(lián)性一同考慮.2011年，Li等[15]提出將RF應(yīng)用到環(huán)境變量的空間推測(cè)中，并與普通克里金(Ordinary Kriging, OK)、IDW組合，表明提出方法的有效性以及對(duì)輸入變量的敏感性.WU等[16]提出基于地理圖斑的RF空間推測(cè)方法，相較于傳統(tǒng)的基于規(guī)則網(wǎng)格的方法，該方法在推測(cè)精度方面有一定的提高.盡管ML在空間推測(cè)方面非常成功，但在直接使用該類技術(shù)時(shí)大多沒有考慮到觀測(cè)值是具有地理空間自相關(guān)的.因此，在之后的研究中，經(jīng)度、緯度等地理背景相關(guān)的推測(cè)因子被引入到模型構(gòu)建中.Behrens等[17]提出將地理空間自相關(guān)的歐式距離與ML組合，并證明比RK、地理加權(quán)回歸(Geographically Weighted Regression, GWR)等方法更具優(yōu)勢(shì).Hengl等[18]在2018年提出“Random Forest for spatial prediction (RFsp)”模型，其以預(yù)估點(diǎn)到樣本點(diǎn)的緩沖距離作為推測(cè)因子，證明其相較于線性地統(tǒng)計(jì)建模與克里金等傳統(tǒng)方法，提高了推測(cè)精度.2020年，Sekulic等[19]提出“Random Forest Spatial Interpolation (RFSI)”模型，其將鄰近點(diǎn)的觀測(cè)值以及到預(yù)估點(diǎn)的距離作為推測(cè)因子引入模型中，并驗(yàn)證RFSI的推測(cè)結(jié)果優(yōu)于克里金以及RFsp.

RFsp、RFSI分別以預(yù)估點(diǎn)與所有樣本點(diǎn)的緩沖距離、鄰近點(diǎn)的觀測(cè)值與其到預(yù)估點(diǎn)的距離作為推測(cè)因子來彌補(bǔ)RF在空間推測(cè)方面的不足，但RFSI對(duì)于距離的應(yīng)用仍存在潛在問題，且模型中運(yùn)用的鄰近點(diǎn)被考慮在同一等級(jí)水平中，這并未充分體現(xiàn)地理學(xué)第一定律的空間相關(guān)性原則.有鑒于此，針對(duì)RFSI的上述不足，本文提出基于位置距離的反距離加權(quán)隨機(jī)森林(Random Forest with Inverse Distance Weighted based on location distance, RFIdw)模型，主要針對(duì)樣本點(diǎn)的觀測(cè)值與到預(yù)估點(diǎn)的距離實(shí)施反距離加權(quán)策略，將距離因素的遠(yuǎn)近考慮在內(nèi)，離預(yù)估點(diǎn)越近的樣本點(diǎn)將賦予更高的權(quán)重.另外，由于反距離加權(quán)組合之后，建模過程中的推測(cè)因子減少，隨機(jī)森林的mtry等參數(shù)設(shè)置將會(huì)在更小的范圍，從而減少模型擬合時(shí)間.本文通過Spatial Interpolation Comparison 97 (SIC97)數(shù)據(jù)對(duì)RFIdw與RK、RFsp、RFSI加以比較，驗(yàn)證RFIdw在空間推測(cè)方面的有效性.

1 研究區(qū)及數(shù)據(jù)集介紹

降水量由于受地區(qū)、海拔等各種因素影響，往往呈現(xiàn)出復(fù)雜的空間分布趨勢(shì)，因此，在空間推測(cè)研究方面被廣泛應(yīng)用.本文通過SIC97數(shù)據(jù)集所對(duì)應(yīng)的研究區(qū)以及包含的具體數(shù)值加以介紹，更好地闡明RFIdw模型的應(yīng)對(duì)問題.

1.1 研究區(qū)概況

本文選取瑞士作為研究區(qū)域如圖1所示，該區(qū)域地處歐洲中南部，位于北緯45°49′～47°48′，東經(jīng)5°57′～10°29′之間，國土面積約為4.1萬km2，地域雖小，但各地氣候差異很大.阿爾卑斯山由東向西伸展，形成了瑞士氣候的分界線，以北地區(qū)受溫和潮濕的西歐海洋性氣候和冬季寒冷夏季溫?zé)岬臇|歐大陸性氣候的交替影響，變化較大；以南地區(qū)則屬地中海氣候，全年氣候宜人.全國年降水量在 1 000～2 000 mm 之間，3/4地區(qū)平均年降水量超過 1 000 mm.該區(qū)域的降水深受地形的影響，高山峻嶺處降水量遠(yuǎn)遠(yuǎn)超過中部高原一些地區(qū)及河谷地帶.

圖1 瑞士DEM與站點(diǎn)圖Fig.1 Station locations in Swiss on top of DEM of study area

1.2 實(shí)驗(yàn)數(shù)據(jù)集

本文數(shù)據(jù)集包括了站點(diǎn)觀測(cè)的降水量數(shù)據(jù)集、DEM、CHELSA(Climatologies at high resolution for the earth’s land surface areas)降水量數(shù)據(jù)等，具體說明如下.

1) 降水量數(shù)據(jù).研究采用的降水量數(shù)據(jù)是1997年4月在環(huán)境研究所(Joint Research Centre, EC, Ispra)放射性環(huán)境監(jiān)測(cè)機(jī)構(gòu)下組織的一項(xiàng)活動(dòng)中所收集，該數(shù)據(jù)集包括1986年5月8日測(cè)量的100次降水量以及估計(jì)的367個(gè)站點(diǎn)降水量，單位為 0.1 mm[20].具體信息如表1所示.

2) DEM數(shù)據(jù).研究采用的DEM數(shù)據(jù)是從https://www.usgs.gov/獲取，空間分辨率為 1 km，具體信息如圖1所示.

3) CHELSA降水量數(shù)據(jù).由于向上的氣流加劇山頂斜坡位置的云和降水形成，而局部環(huán)流系統(tǒng)沿山谷軸線的下沉分支導(dǎo)致云溶解相應(yīng)地降低谷底的降水量這種特殊的地形降水效應(yīng)，阿爾卑斯山山頂可能會(huì)有較高的降水量.CHELSA降水量數(shù)據(jù)則是對(duì)其降水效應(yīng)進(jìn)行近似，并將其運(yùn)用到ERA-Interim氣候再分析降尺度模型中輸出的結(jié)果[21].

2 空間推測(cè)方法與模型

2.1 空間推測(cè)方法

目前，空間推測(cè)技術(shù)的發(fā)展主要分為兩個(gè)階段，前一階段主要是克里金等傳統(tǒng)方法的發(fā)展，但其理論性高，有諸多假設(shè)條件，并且由于數(shù)據(jù)集的復(fù)雜性，很難滿足.因此，后一階段二十一世紀(jì)初興起的ML彌補(bǔ)了傳統(tǒng)方法的部分缺點(diǎn)，引發(fā)了空間推測(cè)方法的進(jìn)一步提升，并被廣泛應(yīng)用到各領(lǐng)域.而本文就是在RFSI模型的基礎(chǔ)上加以改進(jìn)，提出RFIdw模型.

2.1.1 RFSI方法與模型

由于RF中忽略了樣本點(diǎn)之間的空間自相關(guān)性，可能會(huì)導(dǎo)致推測(cè)結(jié)果不準(zhǔn)確，為彌補(bǔ)這項(xiàng)不足，構(gòu)建了RFSI模型，其是在RF的基礎(chǔ)上引入鄰近點(diǎn)的觀測(cè)值以及到預(yù)估點(diǎn)的水平位置距離，公式表達(dá)式如下：

(1)

式中：covj(s0)(j=1,…,m)為預(yù)估點(diǎn)s0類似海拔、溫度、NDVI等的推測(cè)因子，z(si)為第i個(gè)鄰近點(diǎn)si的觀測(cè)值，dloci(i=1,…,n)是第i個(gè)鄰近點(diǎn)si與預(yù)估點(diǎn)s0之間的水平位置距離.RFSI將鄰近點(diǎn)的信息考慮在內(nèi)，相比RF、RFsp更加接近空間推測(cè)原理.

2.1.2 RFIdw方法與模型

RFSI采用鄰近點(diǎn)的觀測(cè)值以及到預(yù)估點(diǎn)的水平位置距離反映推測(cè)位置的信息，但是由于模型構(gòu)建的最終目的是實(shí)現(xiàn)降水量的精準(zhǔn)推測(cè)，而在RFSI模型訓(xùn)練的過程中，當(dāng)鄰近點(diǎn)到預(yù)估點(diǎn)的距離小于或者大于一定范圍時(shí)，推測(cè)的過程以相同的方式進(jìn)行，這將導(dǎo)致推測(cè)結(jié)果出現(xiàn)偏差.因此，為減小距離對(duì)模型訓(xùn)練過程中的影響以及充分體現(xiàn)地理學(xué)第一定律，本文在RFSI的基礎(chǔ)上加以改進(jìn)提出RFIdw模型，其大致可以分為反距離加權(quán)和模型構(gòu)建兩部分，核心思想是針對(duì)RFSI中選取的鄰近點(diǎn)的觀測(cè)值以及到預(yù)估點(diǎn)的距離反距離加權(quán)，其組合值與原有的環(huán)境推測(cè)因子構(gòu)建形成RFIdw模型.模型的表達(dá)式為：

(2)

(3)

對(duì)于ωi(s0)，其表達(dá)式為：

(4)

(5)

式中：(xi,yi)為第i個(gè)鄰近點(diǎn)si的位置，(x0,y0)為預(yù)估點(diǎn)s0的位置，l為反距離的指數(shù).結(jié)合SIC97降水量數(shù)據(jù)，本文基于RFIdw模型執(zhí)行空間推測(cè)過程的算法偽代碼如表2所示.

2.2 模型評(píng)價(jià)標(biāo)準(zhǔn)

為驗(yàn)證RFIdw模型推測(cè)結(jié)果的效性以及準(zhǔn)確性，選取RK、RFsp、RFSI三種方法與其進(jìn)行比較，并利用平均絕對(duì)值誤差(Mean Absolute Error, MAE)、均方根誤差(Root Mean Square Error, RMSE)、判定系數(shù)(Coefficient of Determination,R2)、一致相關(guān)系數(shù)(Concordance Correlation Coefficient, CCC)這四個(gè)評(píng)價(jià)標(biāo)準(zhǔn)加以比較，公式如下：

(6)

(7)

(8)

(9)

表2 基于RFIdw模型的空間推測(cè)算法偽代碼

3 SIC97數(shù)據(jù)空間推測(cè)結(jié)果分析與討論

3.1 空間推測(cè)制圖結(jié)果分析

基于RK、RFsp、RFSI以及RFIdw模型的空間推測(cè)結(jié)果以及不確定性如圖2、圖3所示，其中(d)為RFIdw的推測(cè)結(jié)果圖與不確定性圖.由圖2(d)推測(cè)結(jié)果可知，降水量呈由西南角到東北角帶狀分布趨勢(shì)，其中，瑞士西部區(qū)域降水量較多，中部區(qū)域降水量較少.由圖3(d)的不確定性結(jié)果可知，絕大部分區(qū)域推測(cè)標(biāo)準(zhǔn)差保持在較小的水平，較大的區(qū)域主要在瑞士的東部區(qū)域，結(jié)合圖1的DEM數(shù)據(jù)可知，瑞士的東南部區(qū)域海拔比其他區(qū)域高，導(dǎo)致站點(diǎn)數(shù)據(jù)稀缺，以致空間推測(cè)方法難以在這一區(qū)域捕捉到有效信息，如需更準(zhǔn)確地掌握該區(qū)域的降水量情況，還需要進(jìn)一步獲取樣本點(diǎn)的信息.

圖2 RK (a)、RFsp (b)、RFSI (c)、RFIdw (d)瑞士降水量推測(cè)圖Fig.2 RK (a)，RFsp (b)，RFSI (c)，RFIdw (d) spatial prediction results of Swiss rainfall

圖3 RK(a)、RFsp(b)、RFSI(c)、RFIdw(d)瑞士降水量推測(cè)標(biāo)準(zhǔn)差圖Fig.3 RK(a)，RFsp(b)，RFSI(c)，RFIdw(d) Swiss rainfall prediction standard error

結(jié)合RK、RFsp、RFSI空間推測(cè)以及不確定性圖與RFIdw相比較，降水量的推測(cè)結(jié)果總體趨勢(shì)大致相同，但在局部個(gè)別區(qū)域中差異明顯，主要集中在瑞士阿爾卑斯山以南海拔較高的區(qū)域，相比于RFsp、RFSI，RFIdw推測(cè)結(jié)果更加精確.在不確定性方面，本文以標(biāo)準(zhǔn)差為評(píng)價(jià)指標(biāo)，其中，RK標(biāo)準(zhǔn)差的表達(dá)式為：

(10)

式中：C0、C1是變異函數(shù)的參數(shù)，c0是預(yù)估點(diǎn)與樣本點(diǎn)之間的協(xié)方差向量，q是推測(cè)因子的n×(p+1)維矩陣，C是樣本點(diǎn)之間n×n維的協(xié)方差矩陣，q0是預(yù)估點(diǎn)s0處的p+1維推測(cè)因子向量，對(duì)于RFsp、RFSI、RFIdw的標(biāo)準(zhǔn)差則為：

(11)

相對(duì)RFsp、RFSI、RFIdw三種方法的標(biāo)準(zhǔn)差圖，RK標(biāo)準(zhǔn)差相對(duì)較小，但是其標(biāo)準(zhǔn)差呈現(xiàn)均勻分布的趨勢(shì)，對(duì)不確定信息的衡量包含的信息相對(duì)較少，即特殊點(diǎn)(預(yù)估點(diǎn)周圍樣本點(diǎn)相對(duì)較少)的標(biāo)準(zhǔn)差無法更加標(biāo)準(zhǔn)地度量.基于RF的三種空間推測(cè)方法RFIdw、 RFsp、RFSI在推測(cè)標(biāo)準(zhǔn)差方面不同之處主要集中在瑞士偏東南的區(qū)域中，RFsp的標(biāo)準(zhǔn)差基本保持在150相對(duì)較高的水平，RFSI雖然減小部分區(qū)域的標(biāo)準(zhǔn)差，但仍有小部分區(qū)域保持在較高水平，相對(duì)于RFsp、RFSI，RFIdw的誤差對(duì)于東南的區(qū)域明顯降低.結(jié)合四種空間推測(cè)方法的推測(cè)結(jié)果與不確定信息來看，RFIdw推測(cè)結(jié)果保持著較高的精度，并且在不確定性方面，RFIdw與RFsp、RFSI相比，標(biāo)準(zhǔn)差更小，與RK相比，不確定性更加具有信息性，因此，RFIdw對(duì)于空間推測(cè)的結(jié)果更加合理.

3.2 交叉驗(yàn)證結(jié)果分析

本文提出的RFIdw模型與RK、RFsp、RFSI交叉驗(yàn)證結(jié)果如表3所示.對(duì)比發(fā)現(xiàn)，在推測(cè)精度方面，RK最大，RFsp最小；在推測(cè)標(biāo)準(zhǔn)差方面，則相反.此外，從圖4觀測(cè)值與推測(cè)值的相關(guān)圖可得知，相較于RFIdw，RFsp、RFSI是相對(duì)分散的，同時(shí)證實(shí)了表3中RFsp、RFSI方法較高的RMSE，較低的R2、CCC.

表3 基于五折交叉驗(yàn)證四種推測(cè)方法的精確度

(a) RK相關(guān)圖 (b) RFsp相關(guān)圖 (c) RFSI相關(guān)圖 (d) RFIdw相關(guān)圖圖4 RK (a)、RFsp (b)、RFSI (c)、RFIdw (d)基于觀測(cè)值與推測(cè)值的相關(guān)圖Fig.4 RK (a),RFsp (b),RFSI (c),RFIdw (d) correlation plots based on observations and predictions

進(jìn)一步分析，由于RK是克里金與多元線性回歸的組合，其具有克里金平穩(wěn)性、殘差服從正態(tài)分布等假設(shè)條件，雖然在SIC79數(shù)據(jù)推測(cè)結(jié)果方面，其表現(xiàn)出更高的推測(cè)精度，但是，在復(fù)雜的實(shí)際問題中，這些理想化條件往往很難滿足，以致于結(jié)果可能會(huì)出現(xiàn)偏差；RFsp雖然將空間位置關(guān)系考慮在內(nèi)，但緩沖距離計(jì)算的過程往往緩慢，并且在本次實(shí)驗(yàn)中，并沒有很高的推測(cè)精度；RFSI將樣本點(diǎn)之間的空間自相關(guān)性考慮在內(nèi)，但是模型訓(xùn)練過程中距離的應(yīng)用可能導(dǎo)致推測(cè)結(jié)果的偏差.因此，相較于RK、RFsp、RFSI，RFIdw在空間推測(cè)方面不失為一種好的選擇.

4 總結(jié)與展望

為獲得精確的空間推測(cè)結(jié)果，本文發(fā)展了一種基于水平位置距離的反距離加權(quán)隨機(jī)森林RFIdw模型，不僅考慮了鄰近點(diǎn)的觀測(cè)值以及到預(yù)估點(diǎn)的距離，并對(duì)每個(gè)鄰近點(diǎn)賦以權(quán)重，從而更好地體現(xiàn)了地理學(xué)第一定律的思想.為了驗(yàn)證RFIdw在空間推測(cè)準(zhǔn)確性與不確定性等方面，本文通過SIC97數(shù)據(jù)進(jìn)行了對(duì)比實(shí)驗(yàn)，將RFIdw與RK、RFsp、RFSI這三種空間推測(cè)模型加以比較，從推測(cè)制圖效果、不確定性以及交叉驗(yàn)證精度分析，RFIdw相較于RK，減少了例如克里金模擬變異函數(shù)等的過程；相較于RFsp，減少計(jì)算緩沖距離的過程，提高模型訓(xùn)練的速度；相較于RFSI，有效地解決模型訓(xùn)練過程中應(yīng)用距離的問題，并且得出RFIdw在推測(cè)結(jié)果方面優(yōu)于RFsp、RFSI這兩種方法，在不確定性表達(dá)方面更加具有信息性.

雖然本文中的RFIdw模型在空間推測(cè)方面有一定的有效性與準(zhǔn)確性，但仍存在問題亟待解決：首先，針對(duì)本文的推測(cè)結(jié)果，RFIdw的結(jié)果稍遜于RK，可能由于在RFIdw模型中考慮的只是簡單的反距離加權(quán)，因此，未來應(yīng)該對(duì)權(quán)重設(shè)計(jì)開展更加深入的研究，例如引入鄰近點(diǎn)之間的距離以及對(duì)權(quán)重指數(shù)的復(fù)雜化[22]，使權(quán)重更加合理化；其次，本文對(duì)距離的刻畫只是兩點(diǎn)之間的水平位置距離，沒有考慮到海拔等距離的因素，導(dǎo)致可能兩點(diǎn)之間雖然幾何空間中距離相近，但在地理空間中兩點(diǎn)并不相似，以致選取的鄰近點(diǎn)可能有失偏頗，今后可以考慮地理空間中的測(cè)地距離[23]；最后，對(duì)RFIdw模型的評(píng)價(jià)方面，目前只考慮了推測(cè)結(jié)果及其不確定性，今后還可考慮其他方面的因子，以便更加全面地衡量模型性能.