分布式屬性相似空間插值算法及其對流域土壤屬性插值的計算性能評估*

周腳根，雷秋良，李 勇，張?zhí)禊i

（1.江蘇省塘庫智能觀測與水環(huán)境生態(tài)管控工程研究中心，江蘇 淮安 223300；2.淮陰師范學(xué)院城市與環(huán)境學(xué)院，江蘇 淮安 223300；3.中國科學(xué)院亞熱帶農(nóng)業(yè)生態(tài)研究所，湖南 長沙 410125；4.農(nóng)業(yè)農(nóng)村部面源污染控制重點實驗室，北京 100081；5.中國農(nóng)業(yè)科學(xué)院農(nóng)業(yè)資源與農(nóng)業(yè)區(qū)劃研究所，北京 100081）

0 引言

數(shù)字土壤制圖是以空間分析和數(shù)學(xué)方法為技術(shù)手段，高效表達土壤屬性空間分布的土壤調(diào)查和制圖技術(shù)[1,2]?；跀?shù)字土壤制圖技術(shù)產(chǎn)生的各類土壤專題圖可為生態(tài)環(huán)境、土壤、水資源、農(nóng)業(yè)規(guī)劃等基礎(chǔ)性學(xué)科研究以及國土資源調(diào)查和保護提供基礎(chǔ)數(shù)據(jù)支撐[3-5]，具有非常重要的意義。數(shù)字土壤制圖技術(shù)可分為空間自相關(guān)理論和土壤景觀模型理論兩大方法?？臻g自相關(guān)理論假設(shè)土壤屬性的空間分布服從空間自相關(guān)，故未知位點土壤屬性可用與其鄰近的已知樣點估算[2]。基于該理論的空間插值技術(shù)有很多，如反距離加權(quán)(IDW)[6]、普通克里格(OK)[7]、協(xié)同克里格[8]、回歸克里格[2]、地理加權(quán)回歸(GWR)[9]等。土壤景觀模型理論的核心是經(jīng)典的土壤發(fā)生和形成理論。直觀上，土壤景觀模型理論可表述為：相似地理環(huán)境要素組合下有相似的土壤；環(huán)境組合越相似，其土壤越相似[10,11]。無論數(shù)字土壤制圖技術(shù)采用何種支撐技術(shù)，其核心過程依舊是實現(xiàn)土壤屬性的點源數(shù)據(jù)向面源數(shù)據(jù)的空間拓展。該過程被稱為空間插值或空間估值。傳統(tǒng)空間插值方法采用與未知的待插值位點空間上鄰近的一定數(shù)量觀測樣點，估算未知點的目標(biāo)屬性值。這種思路遵循了地理學(xué)第一定律：空間對象越鄰近，其屬性越相似[12]。不過，由于地形地貌、氣候、植被、成土母質(zhì)以及人類活動等成土因素的復(fù)雜性，客觀地理世界中土壤屬性的空間異質(zhì)非常強。復(fù)雜地理環(huán)境下空間對象空間鄰近，但對象屬性未必相似[13]。Zhou 等[14,15]研究結(jié)果表明：復(fù)雜地理環(huán)境下用非空間屬性比空間屬性更利于表征土壤樣點之間的相似性。隨著數(shù)字土壤制圖應(yīng)用對區(qū)域空間分辨率及研究尺度的需求不斷提升[16]，中大尺度上數(shù)字土壤制圖所面臨的海量數(shù)據(jù)，使得傳統(tǒng)集中式空間插值的計算模式逐漸遇到瓶頸[17]。以輸出空間分辨率30m*30m 的土壤養(yǎng)分分布圖為例：輸出一個市（如長沙）的柵格圖需要計算0.13 億個網(wǎng)格，輸出一個省（湖南）的柵格圖則需要計算2.35 億個網(wǎng)格，當(dāng)范圍擴大到全中國和整個全球（陸地）時，計算量則分別急劇擴增到107億和1 656億。如用傳統(tǒng)集中式IDW算法，計算一個網(wǎng)格用時0.253ms（假定CPU為Intel 2.27GHz、單核），理論上完成上述市、省及全國的柵格計算分別耗時約1h、16.5h 及30d。可見，隨著實際數(shù)據(jù)規(guī)模擴大，傳統(tǒng)的集中計算模式已經(jīng)很難適應(yīng)海量的數(shù)據(jù)處理要求，尤其是對計算有實時性要求時。

為了高效處理海量數(shù)據(jù)，目前信息技術(shù)領(lǐng)域的主流方法是引入分布式計算技術(shù)在多臺計算機上同時執(zhí)行計算，減少程序的運行時間[18][21]。但是傳統(tǒng)的分布式計算技術(shù)有3 個明顯缺點：（1）需專用并行編程接口對原串行算法進行并行化處理，應(yīng)用靈活性小；（2）對集群計算節(jié)點的故障問題處理欠佳；（3）缺乏節(jié)點計算資源動態(tài)平衡的調(diào)控能力，尤其是當(dāng)各計算節(jié)點性能差別較大時，更容易出現(xiàn)節(jié)點交互的通訊瓶頸[17]。云計算平臺Hadoop是一種面向云計算的分布式計算模型框架，其核心框架是MapReduce[20]。和傳統(tǒng)分布式計算技術(shù)一樣，Hadoop 平臺讓程序自動分布到一個由多個計算節(jié)點組成的計算集群上并行執(zhí)行。Hadoop 平臺的優(yōu)勢在于：Hadoop 會跨越機器集群的程序執(zhí)行調(diào)度和處理機器失效，以及管理機器之間的通訊請求[21]。這種模式大大降低用戶對并發(fā)處理或者分布式系統(tǒng)的處理難度。由于Hadoop 云計算平臺具備算法移植簡單、系統(tǒng)可靠性高及各節(jié)點負載動態(tài)平衡的優(yōu)勢，當(dāng)前Hadoop 平臺成為當(dāng)前云計算的主流技術(shù)，廣泛應(yīng)用于各個研究領(lǐng)域[22-24]。

鑒于當(dāng)前數(shù)字土壤制圖領(lǐng)域的空間插值方法依舊以傳統(tǒng)集中式計算為主，文章旨在將LASWR 空間插值方法[14]從集中式計算拓展為分布式計算，建立一種基于Hadoop 平臺的分布式空間插值算法，并評估測試數(shù)據(jù)規(guī)模對分布式空間插值算法計算性能的影響。該算法可為數(shù)字土壤制圖領(lǐng)域當(dāng)前應(yīng)用的空間插值方法，由傳統(tǒng)集中式計算拓展為分布式計算提供技術(shù)參考。

1 材料與方法

1.1 局部屬性相似空間插值算法（LASWR）

對流域尺度土壤屬性的空間插值案例表明：LASWR 算法的插值精度優(yōu)于經(jīng)典的IDW、OK 及GWR，且空間插值的輸出圖更能精細表征土壤目標(biāo)屬性與其影響因子之間的空間互作關(guān)系[14]。LASWR 算法的基本原理是：用非空間屬性而不是空間屬性評估未知土壤樣點與土壤觀測樣點之間的相似性；鄰近樣點對未知樣點的影響權(quán)重由樣點之間非空間屬性的相似性衡量。這種思路本質(zhì)上符合Zhu等[11]提出的地理相似理論。LASWR 算法的關(guān)鍵點在于：一是，在非空間屬性上搜索與未知樣點的鄰近觀測樣本；二是，基于非空間屬性相似加權(quán)計算鄰近樣點對未知樣點的影響權(quán)重，并由此估算未知樣點的目標(biāo)屬性值。未知點目標(biāo)屬性值的估算過程可具體參考文獻?？傮w上，LASWR算法的主要實現(xiàn)步驟如下。

（1）將數(shù)據(jù)集劃分為訓(xùn)練樣本集（觀測樣點）及待估值樣點集；給定土壤目標(biāo)屬性，提取與目標(biāo)屬性相關(guān)聯(lián)的土壤景觀環(huán)境因子；

（2）對待估值樣點集中任何一點，基于非空間屬性相似度量于訓(xùn)練樣本集中搜索k個鄰近觀測樣點；

（3）計算鄰近樣本對未知樣點的影響權(quán)重，并輸出估算值；重復(fù)過程（2）和（3）直至遍歷所有未知點。

1.2 分布式局部屬性相似加權(quán)空間插值算法

該文選擇將LASWR 算法進行分布式計算拓展，構(gòu)建出分布式局部屬性相似加權(quán)空間插值算法，簡稱分布式空間插值算法（DLASWR）。DLASWR 是基于Hadoop 的MapReduce 框架實現(xiàn)。MapReduce 是負責(zé)Hadoop平臺的海量數(shù)據(jù)計算處理框架[20-23]。MapReduce是一個將大規(guī)模分布式計算表示成一個對鍵/值對數(shù)據(jù)集作分布式操作的序列化編程范例。MapReduce 計算分為Map 和Reduce 兩個階段，輸入數(shù)據(jù)是一個鍵/值對數(shù)據(jù)集。在Map 階段，框架根據(jù)用戶指定的方式將輸入的數(shù)據(jù)集拆分成片段，并且將每個片段分配給一個Map 任務(wù)，每個Map 任務(wù)由集群中的一個節(jié)點執(zhí)行。對于每個輸入的鍵/值對(K，V)，Map 任務(wù)調(diào)用用戶定義的映射函數(shù)，將輸入轉(zhuǎn)變成一個不同的鍵/值(K'，V')。Map 階段之后，框架對中間數(shù)據(jù)(K'，V')鍵值對排序，使所有與每個特定鍵值有關(guān)的值出現(xiàn)在一起，生成一個(K',V'*)元組集合；然后將這個元組分解，片段數(shù)量和Reduce 任務(wù)的數(shù)量相等。和Map 任務(wù)一樣，在Reduce 階段，框架會將每一個(K'，V'*)元組片段分配給一個Reduce任務(wù)執(zhí)行，每個Reduce任務(wù)由集群中的某個節(jié)點執(zhí)行。每個Reduce 任務(wù)讀取分配給它的(K'，V'*)元組片段，并調(diào)用用戶定義的Reduce 函數(shù)，將其轉(zhuǎn)變成一個輸出(K，V)鍵/值對。

概括起來，MapReduce計算任務(wù)分為Map和Reduce兩個處理階段。在Map階段，MapReduce框架將計算任務(wù)拆分成大量的Map 任務(wù)并自動分配給各個節(jié)點執(zhí)行；Reduce 階段將各個Map 任務(wù)的輸出結(jié)果分配給Reduce任務(wù)，最終由Reduce任務(wù)輸出結(jié)果。

基于MapReduce 框架的邏輯計算流程，DLASWR 算法的核心計算流程是：將整個空間插值任務(wù)拆分成無數(shù)個子任務(wù)，每個子任務(wù)遵循集中式LASWR 算法計算流程。簡單講，在MapReduce 框架下DLASWR算法的實現(xiàn)流程為：為待估值樣點分配Map任務(wù)，每個Map任務(wù)執(zhí)行一次LASWR運算；而在Reduce階段將每個Map任務(wù)的結(jié)果進行匯總，作為插值結(jié)果輸出。

需要注意的是，Map 任務(wù)對每個未知點執(zhí)行LASWR 算法時都會讀取觀測樣本數(shù)據(jù)，并將結(jié)果輸出到執(zhí)行Reduce 任務(wù)的計算節(jié)點。該過程所涉及的讀寫磁盤數(shù)據(jù)及網(wǎng)絡(luò)數(shù)據(jù)傳輸操作耗時比較大。同時，Map任務(wù)的啟動和結(jié)束，也需要耗時。因此，Map任務(wù)的多少影響DLASWR算法的計算效率，需要給Map任務(wù)分配的待估值點數(shù)量設(shè)置一個參數(shù)Nmap，以調(diào)控DLASWR 算法的計算效率。Nmap決定給每個map任務(wù)分配的待估值樣點數(shù)量。給定一個空間插值任務(wù)，可以通過評估隨Nmap變化DLASWR 的計算效率，選擇合理的Nmap值。另外，受人力和物力的影響，實際所獲取的觀測樣點數(shù)量遠小于未知點數(shù)量。因此，無需將觀測樣本數(shù)據(jù)分配給每個Map任務(wù)，而是直接存儲到各計算節(jié)點上，以減少數(shù)據(jù)訪問的耗時。

總體上，MapReduce框架下DLASWR算法的主要步驟如下。

（1）將訓(xùn)練樣本集分發(fā)到各個節(jié)點上，設(shè)置初始參數(shù)Nmap。

（2）依據(jù)參數(shù)Nmap，將待估值點集分割成一定數(shù)量的子數(shù)據(jù)集。每個子集數(shù)量固定為Nmap，并被分配給各Map任務(wù)。

（4）給定任意未知點，計算節(jié)點上的Map 任務(wù)執(zhí)行一次LASWR 算法操作，此過程不斷重復(fù)直至遍歷整個待估值點子集。

（5）執(zhí)行Reduce 任務(wù)的各節(jié)點歸并源自Map 任務(wù)的輸出結(jié)果后，將其傳輸?shù)郊褐鞴?jié)點，作為最終結(jié)果輸出。

1.3 試驗平臺、評估指標(biāo)及數(shù)據(jù)源

試驗采用的云計算平臺Hadoop 集群，由22 個配置相似的計算節(jié)點(PC 機)組成。每個計算節(jié)點的主要硬件配置為：Inter 雙核處理器（2.93GHz、3GB 內(nèi)存），運行環(huán)境配置了Linux 操作系統(tǒng)及Hadoop 1.2.1 軟件包；集群節(jié)點間以千兆局域網(wǎng)連接和通訊。給定待估值點集數(shù)據(jù)總量為Tdata、集群節(jié)點數(shù)N，每個Map任務(wù)分配的待估值子集數(shù)量為Nmap，則分發(fā)給每個計算節(jié)點的Map任務(wù)Tmap為：

Nmap過大，則Tmap過小，這可能直接造成Map 任務(wù)過少，集群部分節(jié)點分配不到任務(wù)從而降低整體計算效率。Nmap過小，則Tmap過大，則導(dǎo)致Map 任務(wù)過多，同樣會增加節(jié)點間的通信開銷，降低計算效率。鑒于分布式計算廣義屬于并行處理計算，為有效反映參數(shù)Nmap取值變化對DLASWR算法性能的影響，該文采用加速比（Speedup）指標(biāo)以表征DLASWR 算法的計算效率。加速比是同一個任務(wù)在單處理器系統(tǒng)和并行處理系統(tǒng)中運行時間的比率[25,26]為：

式（2）中，Ts是單處理器系統(tǒng)中的運行時間，Tp是P個節(jié)點的運行時間。計算效率越接近1，算法的效率越高。

試驗數(shù)據(jù)源自湖南省長沙縣金井流域。金井流域?qū)儆诘湫偷膩啛釒Ъt壤丘陵地貌，土地利用類型以農(nóng)田及林地為主[27]，年降雨1 200～1 500mm 左右，面積105km2。試驗數(shù)據(jù)集分為土壤觀測樣點集及待估值樣點集。土壤觀測樣點集有1 023 個樣本，每個樣本包含地理坐標(biāo)（X、Y）、土壤有機碳含量及與土壤有機碳含量相關(guān)的4 個環(huán)境協(xié)變量[28]（高程、坡度、濕地指數(shù)、土地利用類型）共7 個屬性字段，觀測樣點及協(xié)變量空間分布見圖1 所示。待估值樣點集有20 萬個待估值點，每個點包含2 個地理坐標(biāo)字段及高程、坡度、濕地指數(shù)和土地利用4個環(huán)境協(xié)變量字段，總共大小近10M。

圖1 土壤樣點有機碳含量及環(huán)境協(xié)變量的空間分布

中、大尺度上真實的訓(xùn)練樣本數(shù)據(jù)很難獲取，同時待估值數(shù)據(jù)集的大小是影響估值計算效率的關(guān)鍵因素。因此，該文保持上述土壤觀測樣點集不變，通過多次復(fù)制待估值數(shù)據(jù)集，生成3個不同大小尺度的數(shù)據(jù)：100M、1G、10G，以評估DLASWR的計算性能。

2 結(jié)果與討論

對3 個不同規(guī)模測試數(shù)據(jù)集（100M、1G 及10G）的試驗結(jié)果表明，待估值數(shù)據(jù)量、分配給Map 任務(wù)的數(shù)據(jù)子集大小以及參與運算的節(jié)點數(shù)量影響DLASWR算法的計算性能。

在小規(guī)模測試數(shù)據(jù)集案例中(圖2)，隨分配給Map 任務(wù)的子集數(shù)量Nmap增加，DLASWR 算法的運行時間消耗逐漸減少；不過隨參與計算的節(jié)點數(shù)增加，DLASWR 算法完成插值任務(wù)的耗時也在減少（由峰值638s降到54s）。這表明Nmap過小導(dǎo)致Map任務(wù)過多，過多的Map任務(wù)增加了節(jié)點間的通訊消耗，從而拉低了分布式算法的運算性能。不過，隨計算節(jié)點的增加，整個分布式算法整體性能趨向加速。實際上，圖2a顯示：當(dāng)Nmap=5 000且計算節(jié)點數(shù)為2和4以及Nmap=10 000且計算節(jié)點數(shù)為2時，DLASWR算法的運行耗時（平均耗時511s）反而高于集中式LASWR 的。這實際上反證了Map任務(wù)過多拉低算法的計算性能。在計算效率方面（圖2b），由于100M數(shù)據(jù)量相對偏小，分布式計算的計算效率優(yōu)勢并未體現(xiàn)。DLAS‐WR 效率較低，各處理的加速比均值為0.21，變化范圍為0.1～0.51。與上述耗時分布趨勢相似，DLASWR算法的計算效率隨Map任務(wù)增加而下降。但是隨參與計算節(jié)點數(shù)增加，DLASWR效率程下降趨勢。

圖2 100M測試數(shù)據(jù)集對分布式DLAWR算法計算性能的影響

在中等規(guī)模1G 測試數(shù)據(jù)集案例中（圖3），DLASWR 算法完全凸顯對集中式LASWR 算法的計算性能優(yōu)勢。不同Nmap處理下DLASWR 平均耗時517s，變化范圍為294～1 502s，遠小于集中式LASWR 算法的2 652s。中等規(guī)模數(shù)據(jù)處理下，DLASWR 算法的計算性能與小規(guī)模數(shù)據(jù)集處理案例類似：算法的耗時總體上與參數(shù)Nmap取值成反比；節(jié)點數(shù)的增加降低DLASWR 算法的耗時(圖3a)。與小規(guī)模數(shù)據(jù)集處理不同，Nmap=1 000 000 及Nmap=2 000 000 處理下，DLASWR 算法的耗時以節(jié)點數(shù)14 為界，呈現(xiàn)先減后增加的趨勢。在計算效率方面，中等規(guī)模1G 測試數(shù)據(jù)處理下DLASWR 算法的計算效率明顯優(yōu)于小規(guī)模100M 數(shù)據(jù)處理的。6 種Nmap取值處理下DLASWR 算法的加速比平均值為0.32，高于小規(guī)模的0.21。與上文小規(guī)模數(shù)據(jù)集的實驗結(jié)果相似，中規(guī)模數(shù)據(jù)處理下DLASWR算法的計算效率整體呈現(xiàn)隨計算節(jié)點和Nmap取值增加(Nmap=1 000 000及Nmap=2 000 000處理除外)而減少的變化趨勢（圖3b）。

圖3 1G測試數(shù)據(jù)集對分布式DLAWR算法計算性能的影響

對大規(guī)模10G 測試數(shù)據(jù)集處理（圖4），DLASWR 算法的耗時及加速比的分布趨勢有異于上述小和中規(guī)模數(shù)據(jù)處理的。DLASWR 算法的計算耗時隨Nmap取值的增加而減少。并且以節(jié)點數(shù)14 為界，算法耗時先減后增。相對小和中規(guī)模數(shù)據(jù)處理，算法的計算效率進一步提升，平均值達0.47，變化范圍為0.1～0.91。并且，算法的計算效率隨Nmap取值增加而增加但隨計算節(jié)點的增加而減少。

圖4 10 G測試數(shù)據(jù)集對分布式DLAWR算法計算性能的影響

綜合3 個不同規(guī)模測試數(shù)據(jù)集的試驗處理結(jié)果可知，數(shù)據(jù)規(guī)模、節(jié)點數(shù)量及Nmap取值影響分布式DLASWR 算法的計算性能。數(shù)據(jù)規(guī)模小時（如100M 小規(guī)模案例），DLASWR 算法處理的性能不如集中式DLASWR 算法，但是隨數(shù)據(jù)規(guī)模擴大，其計算性能完全超越集中式算法，且算法的計算效率也提升明顯。

計算節(jié)點數(shù)增加整體上提升DLASWR 算法的計算性能。不過，在1G 和10G 兩個測試數(shù)據(jù)處理案例中，DLASWR 算法性能都出現(xiàn)了以節(jié)點數(shù)14 為界先減后增的現(xiàn)象。這主要是因為隨著節(jié)點增多，雖然發(fā)給每個節(jié)點的任務(wù)減少、計算時間縮短，但是節(jié)點增加會引發(fā)網(wǎng)絡(luò)通信量的顯著提升，從而拉低計算效率[20]。

Nmap取值大小決定Map 任務(wù)量，從而影響DLASWR 算法計算性能。不過，Nmap取值大小與DLAS‐WR算法計算性能并不是簡單的線性關(guān)系。例如，上述試驗結(jié)果顯示：數(shù)量為100M時，Nmap取值增加提升算法的計算性能（圖2）；但是數(shù)據(jù)量為10G時，Nmap取值增加則降低算法的計算性能（圖4）。造成這種現(xiàn)象的主要原因是：Hadoop 框架中Map 任務(wù)的切換存在額外的時間開銷。當(dāng)數(shù)據(jù)量比較小時，Map 任務(wù)的時間開銷將占總耗時的絕大部分。當(dāng)數(shù)據(jù)量足夠大時，Map任務(wù)的時間開銷在總耗時中所占比例將變得較小。通常，Nmap取值原則是以每個節(jié)點分配1～2 個Map 任務(wù)為宜；當(dāng)數(shù)據(jù)集較大時，可適當(dāng)增加Map任務(wù)數(shù)[26]。

3 結(jié)論

海量數(shù)據(jù)處理是現(xiàn)代數(shù)字土壤制圖技術(shù)向精細化、大尺度化發(fā)展所面臨的技術(shù)挑戰(zhàn)。該文構(gòu)建了一個基于云計算平臺Hadoop 的分布式空間插值算法DLASWR。對流域土壤樣點屬性的空間插值測試結(jié)果表明，DLASWR 算法的計算效率比集中式LASWR 算法顯著升。數(shù)據(jù)規(guī)模、節(jié)點數(shù)量及Map任務(wù)數(shù)量是影響DLASWR 算法計算性能的主要因素。該算法可為數(shù)字土壤制圖領(lǐng)域當(dāng)前應(yīng)用的空間插值方法由傳統(tǒng)的集中式計算拓展成分布式計算提供技術(shù)參考。受現(xiàn)有研究條件限制，未來需要在更大規(guī)模數(shù)據(jù)條件下對DLASWR算法進行計算性能驗證。