點(diǎn)源時(shí)間序列數(shù)據(jù)缺失值的估值方法比較——以小流域氣象和水文數(shù)據(jù)為例*

甘 蕾，周腳根，石 錦，李 希，沈健林，呂殿青，李裕元，吳金水

?

點(diǎn)源時(shí)間序列數(shù)據(jù)缺失值的估值方法比較——以小流域氣象和水文數(shù)據(jù)為例*

甘 蕾1,2，周腳根2**，石 錦2,3，李 希2，沈健林2，呂殿青1，李裕元2，吳金水2

（1.湖南師范大學(xué)資源與環(huán)境科學(xué)學(xué)院，長(zhǎng)沙 410081；2.中國(guó)科學(xué)院亞熱帶農(nóng)業(yè)生態(tài)研究所亞熱帶農(nóng)業(yè)生態(tài)過(guò)程重點(diǎn)實(shí)驗(yàn)室，長(zhǎng)沙 410125；3.湖南農(nóng)業(yè)大學(xué)工學(xué)院，長(zhǎng)沙 410128）

對(duì)點(diǎn)源時(shí)間序列數(shù)據(jù)缺失值進(jìn)行有效估值能提升其數(shù)據(jù)質(zhì)量。為探究不同估值方法對(duì)點(diǎn)源時(shí)間序列數(shù)據(jù)缺失值的估值效果及其影響因素，以亞熱帶典型小流域長(zhǎng)期定位觀測(cè)的每日氣象和水文數(shù)據(jù)（最高氣溫、最低氣溫、太陽(yáng)輻射量、降雨量及地表徑流量）為例，以均方根誤差（RMSE）、絕對(duì)平均誤差（MAE）和Pearson相關(guān)系數(shù)（r）為性能驗(yàn)證指標(biāo)，比較了線性內(nèi)插法（LIM）、K-最近鄰插值法（KNNM）、樣條插值法（SIM）、多項(xiàng)式插值法（PIM）和核密度估值法（KDEM）5種估值方法的估值性能差異及其主要影響因素。結(jié)果表明：（1）LIM、SIM和KDEM的估值性能總體上優(yōu)于其它2種方法；（2）5種估值方法對(duì)氣象數(shù)據(jù)（最高氣溫、最低氣溫和太陽(yáng)輻射量）缺失值估值的RMSE為1.81～6.35，MAE為1.30～4.20，r為0.70～0.98（P＜0.05），而對(duì)水文數(shù)據(jù)（降雨量和地表徑流量）缺失值估值的RMSE為12.54～26.28，MAE為3.60～14.21，r為0.07～0.72。可見(jiàn)，各估值方法對(duì)氣象數(shù)據(jù)的估值性能強(qiáng)于對(duì)水文數(shù)據(jù)；（3）上述數(shù)據(jù)集的變異系數(shù)（CV）與估值評(píng)估指標(biāo)（RMSE、MAE及r）線性相關(guān)（P＜0.05），是影響估值性能的重要因素。

缺失值；估值方法；變異系數(shù)；時(shí)間序列

時(shí)間序列數(shù)據(jù)是生態(tài)環(huán)境、水文及氣象等研究領(lǐng)域必不可少的基礎(chǔ)數(shù)據(jù)，這些領(lǐng)域的相關(guān)研究通常需要對(duì)環(huán)境參數(shù)進(jìn)行長(zhǎng)期定位監(jiān)測(cè)采集，但是由于儀器設(shè)備故障、環(huán)境惡劣或人為操作失誤等原因，采集到的觀測(cè)數(shù)據(jù)難免出現(xiàn)數(shù)據(jù)缺失問(wèn)題[1]，從而影響觀測(cè)數(shù)據(jù)的質(zhì)量。有效估算時(shí)間序列數(shù)據(jù)的缺失值，可以完善時(shí)間序列數(shù)據(jù)的質(zhì)量，提升數(shù)據(jù)使用效率，是空間分析與統(tǒng)計(jì)領(lǐng)域研究的熱點(diǎn)之一[2]。時(shí)間序列的估值問(wèn)題，目前的研究主要涉及兩方面：（1）面源尺度上對(duì)未觀測(cè)位點(diǎn)環(huán)境參數(shù)屬性值的估算；（2）點(diǎn)源尺度上對(duì)觀測(cè)參數(shù)缺失值的估算。

由于人力和物力的有限性，面源尺度上環(huán)境參數(shù)通常通過(guò)一定量代表性點(diǎn)源觀測(cè)單元獲取，再通過(guò)這些點(diǎn)源觀測(cè)數(shù)據(jù)實(shí)現(xiàn)觀測(cè)數(shù)據(jù)的面源拓展，是一個(gè)用一定量點(diǎn)源觀測(cè)數(shù)據(jù)估算面源上未觀測(cè)單元參數(shù)值的過(guò)程。當(dāng)前，空間插值方法、GIS技術(shù)、估值預(yù)測(cè)模型等常用于解決該問(wèn)題。毛洋洋等[3]利用不同日太陽(yáng)總輻射預(yù)測(cè)模型對(duì)華北地區(qū)6個(gè)站點(diǎn)的逐日太陽(yáng)總輻射數(shù)據(jù)進(jìn)行估算，其估值效果皆可；郭兆夏等[4]利用GIS技術(shù)對(duì)陜西年降水量數(shù)據(jù)進(jìn)行了較準(zhǔn)確的分析與預(yù)測(cè)；Srebotnjak等[5]利用樣條插值法有效完成了全球尺度上水質(zhì)監(jiān)測(cè)并實(shí)現(xiàn)了水質(zhì)TN、TP、DO等數(shù)據(jù)的填充。其中，Kriging空間插值法應(yīng)用較多，實(shí)現(xiàn)了黃土高原區(qū)多年降雨量[6-7]、西部地區(qū)降雨量[8-9]、黃河流域多年降雨量[10]的空間拓展與分析；最近鄰法、反距離加權(quán)法等空間插值方法能很好地預(yù)測(cè)全國(guó)較大區(qū)域范圍、湖南復(fù)雜地形區(qū)的日平均氣溫[11-12]；基于模型的空間插值技術(shù)對(duì)江蘇、安徽逐日氣溫[13]和結(jié)合線性回歸分析等的空間插值方法對(duì)漢江上游多年平均氣溫[14]的預(yù)測(cè)效果顯著。而上述方法面向的基本為氣象數(shù)據(jù)的空間估值拓展，少有涉及對(duì)點(diǎn)源時(shí)間序列數(shù)據(jù)的估值。

點(diǎn)源尺度上時(shí)間序列缺失值的估值，主要是對(duì)一定觀測(cè)時(shí)間段內(nèi)缺失的觀測(cè)數(shù)據(jù)進(jìn)行有效插補(bǔ)。一些研究直接將缺失數(shù)據(jù)的樣本剔除[15]，或采用均值替換所有缺失值[16]，雖然操作簡(jiǎn)單，但會(huì)導(dǎo)致潛在信息丟失，局限性大。鑒于點(diǎn)源時(shí)間序列實(shí)為二維數(shù)據(jù)集，實(shí)際研究中通常運(yùn)用線性內(nèi)插法、樣條插值法等二維曲線擬合的數(shù)學(xué)方法對(duì)缺失數(shù)據(jù)進(jìn)行插補(bǔ)。Ferrari等[17-18]證實(shí)了線性內(nèi)插法可對(duì)降雨量和溫度數(shù)據(jù)的缺失值進(jìn)行有效估算；結(jié)合地形等因素，鄭小波等[19]發(fā)現(xiàn)薄盤(pán)光滑樣條函數(shù)法對(duì)西南地區(qū)溫度和降水?dāng)?shù)據(jù)的插值效果最優(yōu)。當(dāng)前國(guó)內(nèi)外對(duì)點(diǎn)源時(shí)間序列數(shù)據(jù)缺失值的估值問(wèn)題關(guān)注較少，且常集中于某一種估值方法或技術(shù)對(duì)特定類(lèi)型的數(shù)據(jù)集缺失值的估值分析，缺乏不同估值方法對(duì)一類(lèi)或幾類(lèi)數(shù)據(jù)缺失值估值結(jié)果的性能差異比較，也少有分析估值方法對(duì)不同數(shù)據(jù)集的性能響應(yīng)，大多不易推廣和應(yīng)用到其它類(lèi)型數(shù)據(jù)。

為此，本研究選用LIM、KNNM、SIM、PIM和KDEM 5種估值方法，以湖南金井小流域每日氣象數(shù)據(jù)（最高氣溫、最低氣溫、太陽(yáng)輻射量）、水文數(shù)據(jù)（降雨量、地表徑流量）為應(yīng)用實(shí)例，研究上述5種方法對(duì)不同數(shù)據(jù)集的估值性能差異及其影響因素，以期為氣象和水文等領(lǐng)域點(diǎn)源時(shí)間序列數(shù)據(jù)缺失值的估值方法提供選擇，并為提高相關(guān)模型預(yù)測(cè)精度提供參考依據(jù)。

1 材料與方法

1.1 數(shù)據(jù)來(lái)源

數(shù)據(jù)來(lái)源于湖南省長(zhǎng)沙縣金井河流域，流域總面積134.4 km2，位于27°55N-28°40N、112°56E-113°30' E（圖1），屬亞熱帶濕潤(rùn)性季風(fēng)氣候，是典型的亞熱帶紅壤丘陵地貌，年平均降水量為1200～1500mm[20]。

圖1 金井流域水文和氣象觀測(cè)站分布圖

所用水文數(shù)據(jù)為2010-2012年金井河小流域每日降雨量及2010-2014年流域內(nèi)出水口的每日地表徑流量數(shù)據(jù)，氣象數(shù)據(jù)為2010-2013年流域內(nèi)每日最高氣溫、最低氣溫?cái)?shù)據(jù)和太陽(yáng)輻射量數(shù)據(jù)。地表徑流量數(shù)據(jù)采用Simpson's Parabolic Rule方法，用螺旋杯式流速儀實(shí)測(cè)而得，該系統(tǒng)每10min自動(dòng)采集并記錄流量數(shù)據(jù)，據(jù)此計(jì)算流域研究時(shí)段內(nèi)的日地表徑流量。各氣象因子數(shù)據(jù)，則由小型氣象站（Intelimet Advantage，Dynamax Inc.，美國(guó)產(chǎn)）觀測(cè)獲得。

所選取的數(shù)據(jù)集類(lèi)型皆為流域水文和氣象觀測(cè)的基礎(chǔ)類(lèi)型，且各數(shù)據(jù)集間差異明顯，整體上水文數(shù)據(jù)（包括降雨量、地表徑流量）其CV為130.71%～162.57%，較氣象數(shù)據(jù)（最高氣溫、最低氣溫和太陽(yáng)輻射量）變異性大（CV為42.82%～67.51%）。

1.2 估值方法

1.2.1 線性內(nèi)插法

線性內(nèi)插法（LIM）[21]利用時(shí)間與觀測(cè)值之間的等比關(guān)系近似求解時(shí)間序列的缺失值。給定時(shí)間序列集t，已知ti、tk時(shí)刻對(duì)應(yīng)的觀測(cè)值分別為Y（ti）、Y（tk），tj時(shí)刻數(shù)據(jù)樣點(diǎn)值Y（tj）缺失，其中i

由式（1）可見(jiàn)，若數(shù)據(jù)缺失位點(diǎn)處于時(shí)間序列的兩端點(diǎn)，即j=i或j=k，則LIM方法將無(wú)法實(shí)現(xiàn)預(yù)測(cè)。

1.2.2 K-最近鄰插值法

K-最近鄰插值法（KNNM）[22]的核心思想是，搜索與待估算點(diǎn)最鄰近的k個(gè)觀測(cè)點(diǎn)樣本，用這些樣本點(diǎn)觀測(cè)值的加權(quán)和賦予待估值點(diǎn)。樣點(diǎn)之間的鄰近關(guān)系為

時(shí)間序列數(shù)據(jù)的計(jì)算則首先給定與tj鄰近的k個(gè)鄰近點(diǎn)集，然后估算Y（tj）

1.2.3 樣條插值法和多項(xiàng)式插值法

樣條插值法（SIM）是一種特殊的分段3次多項(xiàng)式插值法。相對(duì)普通多項(xiàng)式插值，通常樣條插值方法對(duì)數(shù)據(jù)集的擬合更平滑，輸出的插值誤差更小。給定n+1個(gè)不同的觀測(cè)時(shí)刻ti，并滿足t0＜t1＜…＜tn-1＜tn以及 n+1個(gè)觀測(cè)值Y（ti），樣條插值實(shí)質(zhì)上就是構(gòu)建一個(gè)n階樣條函數(shù)Y（t）逼近觀測(cè)數(shù)據(jù)集，即

多項(xiàng)式插值法（PIM）[23]是用多項(xiàng)式對(duì)一列數(shù)據(jù)進(jìn)行線性擬合，再對(duì)給定待估值點(diǎn)進(jìn)行估值的過(guò)程。給定時(shí)間序列數(shù)據(jù)集Y= {Y（t1），Y（t2），…, Y（tn）}和待估值點(diǎn)Y（tj），用多項(xiàng)式函數(shù)f（t）=β0+β1t+β2t2+…+βntn對(duì)時(shí)間序列數(shù)據(jù)集Y進(jìn)行線性擬合，以求解最優(yōu)的參數(shù)β=（β0，β1，β2，…，βn）。本研究用最小二乘法求解最優(yōu)參數(shù)β。

1.2.4 核密度估值法

核密度估值法（KDEM）[24]是一種從數(shù)據(jù)樣本本身出發(fā)研究數(shù)據(jù)分布特征的密度函數(shù)近似估值算法，不需要有關(guān)數(shù)據(jù)分布的先驗(yàn)知識(shí)。對(duì)給定缺失值Y（tj），核密度估值方法估算式為

式中，K（t）為核函數(shù)；h為核函數(shù)的帶寬；n為參與估值的觀測(cè)值數(shù)目。本研究中，核函數(shù)K（t）采用高斯核函數(shù)；該核函數(shù)是一個(gè)權(quán)函數(shù)，離缺失點(diǎn)tj越近的點(diǎn)對(duì)函數(shù)值的影響越大，其權(quán)值也越大；核函數(shù)帶寬h統(tǒng)一為缺失點(diǎn)tj到其它觀測(cè)點(diǎn)的距離集的中段值。

1.3 缺失值設(shè)置及模型校驗(yàn)

采用的日時(shí)間序列數(shù)據(jù)集（最高氣溫、最低氣溫、太陽(yáng)輻射量、降雨量和地表徑流量）皆為完整數(shù)據(jù)集（即無(wú)缺失值）。通常，時(shí)間序列數(shù)據(jù)集數(shù)據(jù)缺失位點(diǎn)以及數(shù)據(jù)缺失量是隨機(jī)的。缺失量的多少在一定程度上會(huì)影響估值方法性能評(píng)價(jià)的客觀性，目前主流研究以20%～30%缺失量作為研究對(duì)象用于篩選估值方法[25-26]。為有效評(píng)估LIM、KNNM、SIM、PIM和KDEM 5種方法對(duì)缺失值的估值性能的差異，本研究隨機(jī)抽取每個(gè)實(shí)例數(shù)據(jù)集的25%數(shù)據(jù)樣本點(diǎn)為模擬缺失量。

涉及的LIM、KNNM、SIM、PIM和KDEM的代碼實(shí)現(xiàn)以及模型運(yùn)行均在Matlab2011b軟件平臺(tái)完成。其中，LIM、KNNM、SIM和PIM直接調(diào)用Matlab2011b軟件的內(nèi)置包進(jìn)行運(yùn)行；KDEM則為自主編碼實(shí)現(xiàn)。在運(yùn)行模型對(duì)25%抽樣樣本進(jìn)行預(yù)測(cè)前，用交叉校驗(yàn)方法測(cè)試75%的訓(xùn)練樣本，分別為上述5種方法尋找較優(yōu)的模型輸入?yún)?shù)。多次試驗(yàn)證實(shí)，采用12～18的鄰近樣本數(shù)，LIM、KNNM、SIM、PIM及KDEM的估值性能較優(yōu)。考慮后期需要多次進(jìn)行抽樣測(cè)試，故對(duì)25%測(cè)試樣本的估值試驗(yàn)統(tǒng)一鄰近樣本參數(shù)定為15。為消除單次試驗(yàn)帶來(lái)的隨機(jī)誤差，每次試驗(yàn)重復(fù)100次。將100次試驗(yàn)的均方根誤差（RMSE）、絕對(duì)平均誤差（MAE）和Pearson相關(guān)系數(shù)（r）3個(gè)指標(biāo)的平均值作為驗(yàn)證指標(biāo)用于評(píng)估各方法估值性能的優(yōu)劣。

2 結(jié)果與分析

2.1 小流域水文和氣象時(shí)間序列數(shù)據(jù)集的統(tǒng)計(jì)特征

金井河小流域水文和氣象數(shù)據(jù)類(lèi)型中25%缺失數(shù)據(jù)集和75%訓(xùn)練樣本數(shù)據(jù)集數(shù)據(jù)點(diǎn)分布見(jiàn)圖2。由圖可見(jiàn)，所有數(shù)據(jù)從2010-07-01起始，水文數(shù)據(jù)（小流域出水口的地表徑流量）至2013-10-20，樣本數(shù)據(jù)共1206個(gè)，訓(xùn)練樣本數(shù)據(jù)904個(gè)與缺失數(shù)據(jù)302個(gè)隨機(jī)分布在取樣時(shí)段內(nèi)，數(shù)據(jù)點(diǎn)分布趨勢(shì)吻合；氣象要素集的降雨量數(shù)據(jù)共470個(gè)，截至2011-10-13；最高和最低氣溫?cái)?shù)據(jù)共904個(gè)，截至2012-12-20；太陽(yáng)輻射量數(shù)據(jù)總共632個(gè)，截至2012-03-23；缺失數(shù)據(jù)集均隨機(jī)分布在取樣時(shí)間段內(nèi)，總體上與訓(xùn)練樣本數(shù)據(jù)集的數(shù)據(jù)點(diǎn)分布趨勢(shì)吻合。

各指標(biāo)訓(xùn)練樣本數(shù)據(jù)集的統(tǒng)計(jì)特征見(jiàn)表1。由表可見(jiàn)，所選指標(biāo)的數(shù)據(jù)集差異明顯，降雨量數(shù)據(jù)穩(wěn)定性差，數(shù)值變化范圍大，所選時(shí)段內(nèi)最大降雨量為34.62mm，最小降雨量為0.01mm，變異系數(shù)CV最大為162.57%；地表徑流量的主要來(lái)源是降雨，基流匯聚形成地表徑流，最大徑流量為41.93m3，數(shù)據(jù)集的變異系數(shù)CV也較大，僅次于降雨量數(shù)據(jù)，達(dá)130.71%；該兩指標(biāo)均屬?gòu)?qiáng)變異水平[27]。最高氣溫、最低氣溫和太陽(yáng)輻射量數(shù)據(jù)較穩(wěn)定，最高氣溫?cái)?shù)據(jù)集CV最小，為42.82%；最低氣溫和太陽(yáng)輻射量數(shù)據(jù)集CV居中，分別為66.96%、67.51%，均屬弱變異水平。

圖2 各觀測(cè)數(shù)據(jù)日值集中訓(xùn)練數(shù)據(jù)與缺失數(shù)據(jù)的分布

2.2 五種方法對(duì)時(shí)間序列數(shù)據(jù)集中缺失數(shù)據(jù)估值效果的比較

由表2可見(jiàn)，5種估值方法對(duì)不同數(shù)據(jù)集的估值性能具有較大差異。對(duì)于變異系數(shù)較小的氣象數(shù)據(jù)（最高氣溫、最低氣溫及太陽(yáng)輻射量），LIM、KNNM、SIM、PIM及KDEM 5種方法皆表現(xiàn)較佳，估算值與實(shí)測(cè)值相關(guān)性強(qiáng)（r為0.64～0.98，P＜0.05）；其中，LIM估值準(zhǔn)確性最佳，估值結(jié)果誤差最小，其RMSE、MAE分別為1.81～4.58、1.30～3.43，相關(guān)性最高，r為0.78～0.98（P＜0.05）；KDEM和SIM估值效果居中，KDEM對(duì)最高氣溫估值較好，其RMSE、MAE、r為2.91℃、2.12℃、0.95（P＜0.05），SIM對(duì)最低氣溫的估值效果與LIM相同，同為最佳方法，且對(duì)太陽(yáng)輻射量估值也較好；KNNM和PIM兩種方法表現(xiàn)最差，誤差大，相關(guān)性弱。

對(duì)于變異系數(shù)較大的降雨量數(shù)據(jù)，LIM、KNNM、SIM、PIM及KDEM 5種方法估值效果皆不佳，RMSE和MAE偏大，估算值與實(shí)測(cè)值相關(guān)性不顯著（r為0.07～0.13），其中，KDEM相對(duì)較優(yōu)，其RMSE、MAE、r分別為16.75mm、9.22mm、0.13。而受降雨影響的地表徑流量數(shù)據(jù)，SIM的估值性能最優(yōu)，其RMSE、MAE、r分別為12.54m3、3.40m3、0.72，誤差小且相關(guān)系數(shù)較大；LIM和KDEM的性能居中，其RMSE、MAE、r分別為12.66m3、3.60m3、0.71和13.47m3、 3.86m3、0.69；KNNM和PIM的性能最差。

總體上，上述5種估值方法對(duì)日最高氣溫、日最低氣溫、日太陽(yáng)輻射量以及日地表徑流量數(shù)據(jù)的估值結(jié)果較為可靠，但對(duì)日降雨量的估值精度不高，這可能是因?yàn)槿战涤炅繙y(cè)試數(shù)據(jù)集的變異系數(shù)過(guò)大（CV=162.57%）。另外，LIM、SIM和KDEM 3種估值方法對(duì)這5種缺失數(shù)據(jù)的估值效果較好。

表2 五種方法對(duì)水文和氣象數(shù)據(jù)集缺失值的估值效果比較

注：LIM為線性內(nèi)插法、KNNM為K-最近鄰插值法、SIM為樣條插值法、PIM為多項(xiàng)式插值法、KDEM為核密度估值法；RMSE為均方根誤差、MAE絕對(duì)值平均誤差、r為估算值與實(shí)測(cè)值的Pearson相關(guān)系數(shù)。表中數(shù)據(jù)為每種方法重復(fù)100次估算結(jié)果的平均值±標(biāo)準(zhǔn)誤差。

Note: LIM for linear interpolation method, KNNM for K-nearest neighbor interpolation method, SIM for spline interpolation method, PIM for polynomial interpolation method, KDEM for Kernel density estimation method, while RMSE for root mean square error, MAE for absolute mean error, r for Pearson product-moment correlation coefficient between estimated and measured values. The data in the table were mean±standard error values of estimations of repeated 100 times by each of the interpolation methods.

2.3 原數(shù)據(jù)集中變異系數(shù)對(duì)缺失值估算結(jié)果的影響

將最高氣溫、最低氣溫、太陽(yáng)輻射量、降雨量以及地表徑流量訓(xùn)練數(shù)據(jù)集的變異系數(shù)（CV）與交叉驗(yàn)證指標(biāo)值（RMSE、MAE、r）進(jìn)行線性擬合分析，結(jié)果見(jiàn)圖3。由圖3可見(jiàn)，CV與RMSE、MAE、r之間存在明顯的線性相關(guān)關(guān)系。CV與RMSE呈顯著正相關(guān)（P＜0.05），線性擬合方程的決定系數(shù)R2達(dá)0.89；與MAE也呈線性正相關(guān)（P＜0.05），R2達(dá)0.67；與r呈負(fù)線性相關(guān)（P＜0.05），R2達(dá)0.79。說(shuō)明變異系數(shù)是影響缺失值估值結(jié)果的重要因素，變異系數(shù)越大，均方根誤差和絕對(duì)平均誤差越大，相關(guān)性越小；反之，變異系數(shù)越小，誤差越小，相關(guān)性越大。進(jìn)一步分析不同估值方法的估值性能對(duì)5種水文氣象數(shù)據(jù)集CV變化的響應(yīng)相關(guān)性。圖4表明，CV與RMSE和MAE呈線性正相關(guān)，決定系數(shù)R2分別為0.92～0.95和0.69～0.74；與相關(guān)系數(shù)r呈線性負(fù)相關(guān)，決定系數(shù)R2為0.78～0.80。這表明在上述應(yīng)用實(shí)例中5種估值方法輸出的估值誤差，超過(guò)69%是與數(shù)據(jù)集固有的變異性有關(guān)。因此，在本研究中CV是影響估值方法LIM、KNNM、SIM、PIM及KDEM的估值性能的關(guān)鍵因素。數(shù)據(jù)集的變異系數(shù)越大，LIM、KNNM、SIM、PIM及KDEM 5種方法的估值誤差越大，輸出的預(yù)測(cè)值與實(shí)測(cè)值的擬合度越小，對(duì)估值結(jié)果的準(zhǔn)確性影響越大[28]。

圖3 各數(shù)據(jù)集變異系數(shù)與缺失值估值評(píng)估指標(biāo)RMSE、MAE和r的相關(guān)性

圖4 各觀測(cè)數(shù)據(jù)集變異系數(shù)與五種估值方法輸出的估值評(píng)估指標(biāo)的相關(guān)性分析

3 結(jié)論與討論

3.1 討論

以日最高氣溫、日最低氣溫、日太陽(yáng)輻射量、日降雨量及日地表徑流量數(shù)據(jù)為應(yīng)用實(shí)例，模擬和比較了25%樣本缺失量條件下LIM、KNNM、SIM、PIM和KDEM 5種估值方法的性能差異及其主要影響因素。總體上，LIM、SIM和KDEM 3種方法對(duì)氣象數(shù)據(jù)集缺失值的估算性能優(yōu)于其它兩種方法，對(duì)缺失值估算的誤差小且估算值與實(shí)測(cè)值具有線性相關(guān)關(guān)系，尤其是對(duì)氣溫?cái)?shù)據(jù)，其RMSE、MAE分別低至1.81℃、1.30℃，r高達(dá)0.98。

上述估值方法的性能差異與估值方法本身有一定的關(guān)系。LIM運(yùn)算簡(jiǎn)單，適于所有水文氣象數(shù)據(jù)。不論點(diǎn)源還是面源數(shù)據(jù)，KDEM僅從樣本本身出發(fā)，可以估值任何形狀的缺失值概率密度函數(shù)，且連續(xù)性好[29]。KNNM估算時(shí)難以確定的k值易導(dǎo)致估值變化大，穩(wěn)定性不高[30]。PIM受數(shù)據(jù)量大小和運(yùn)算次數(shù)的限制，誤差較大，SIM較PIM更靈活穩(wěn)定，運(yùn)算結(jié)果精度高，不受數(shù)據(jù)量大小影響，運(yùn)算簡(jiǎn)便[31]。文獻(xiàn)研究也證實(shí)了LIM對(duì)點(diǎn)源時(shí)間序列數(shù)據(jù)的估值性能較優(yōu)。例如，Noor等[32]在估算環(huán)境質(zhì)量PM10數(shù)據(jù)集的缺失值時(shí)所表現(xiàn)的高精度和可靠性佐證了LIM的高性能；Saleem等[18]分析發(fā)現(xiàn)LIM對(duì)空氣溫度數(shù)據(jù)缺失值的估值精度最高，r高達(dá)0.99以上（P＜0.01）；唐云輝等[33]基于鄰域特征對(duì)重慶市日最高、日最低氣溫?cái)?shù)據(jù)進(jìn)行缺失填補(bǔ)的擬合精度高，結(jié)果可靠。這可能是氣溫?cái)?shù)據(jù)時(shí)間尺度上變化小，限制因素少，數(shù)據(jù)集相對(duì)穩(wěn)定的原因。

本研究也發(fā)現(xiàn)，對(duì)變異系數(shù)最大的日降雨量數(shù)據(jù)缺失值的估值，5種方法均表現(xiàn)不佳，其相關(guān)性弱（r在0.02～0.11），預(yù)測(cè)誤差大。但相對(duì)其它4種估值方法，線性插值方法對(duì)日降雨量數(shù)據(jù)的估值相對(duì)較好，其RMSE、MAE值分別為8.25mm、5.30mm，但是估值精度低于巴西巴拉那州氣象站[17]和巴基斯坦[18]日降雨量缺失值的估算精度。這歸因于不同研究區(qū)域日降雨量的地理差異。

本研究表明，數(shù)據(jù)集變異系數(shù)小，離散程度小，則5種估值方法對(duì)數(shù)據(jù)集缺失值的估值效果較優(yōu)；反之，數(shù)據(jù)集變異系數(shù)大，離散程度大，5種估值方法的估值效果皆顯著下降。不同估值方法處理后的估值驗(yàn)證指標(biāo)對(duì)5種水文氣象數(shù)據(jù)集CV變化的響應(yīng)關(guān)系也表明：不同估值方法處理下數(shù)據(jù)集CV與RMSE和MAE線性正相關(guān)（P＜0.05），與r線性負(fù)相關(guān)（P＜0.05）。這充分證實(shí)了數(shù)據(jù)集的變異系數(shù)是影響估值方法的估值結(jié)果的重要因素，該結(jié)論與其它學(xué)者研究結(jié)果相吻合。例如，趙彥鋒等[34]發(fā)現(xiàn)有機(jī)質(zhì)數(shù)據(jù)變異系數(shù)小于10%時(shí)對(duì)數(shù)據(jù)集估值結(jié)果的準(zhǔn)確性最高；Yozgatligil等[35]也證實(shí)土耳其降水、溫度數(shù)據(jù)集CV值越小，對(duì)缺失值估值結(jié)果越可靠。

綜合上述研究結(jié)果，數(shù)據(jù)集的變異系數(shù)顯著影響估值方法的估值性能。依據(jù)數(shù)據(jù)集變異系數(shù)CV與估值驗(yàn)證指標(biāo)（RMSE、MAE以及r）之間的線性關(guān)系，可推斷出：數(shù)據(jù)集變異系數(shù)在不超過(guò)45%的情況下，LIM、SIM和KDEM對(duì)數(shù)據(jù)缺失值的估值結(jié)果更可靠。

3.2 結(jié)論

（1）LIM、KNNM、SIM、PIM和KDEM對(duì)點(diǎn)源時(shí)間序列數(shù)據(jù)缺失值的估值效果存在差異，其中LIM、SIM和KDEM的估值性能優(yōu)于KNNM和PIM。

（2）5種估值方法對(duì)氣象數(shù)據(jù)（最高氣溫、最低氣溫、太陽(yáng)輻射量）缺失值的估值效果整體上優(yōu)于水文數(shù)據(jù)（降雨量、地表徑流量）。

（3）數(shù)據(jù)集的變異系數(shù)CV是影響估值性能的主要因素，且 CV與評(píng)估指標(biāo)RMSE、MAE及r線性相關(guān)（P<0.05）；當(dāng)氣象、水文點(diǎn)源時(shí)間序列數(shù)據(jù)集CV不超過(guò)45%時(shí)，推薦使用LIM、SIM和KDEM估算缺失值。

[1]Kantardzic M.Data mining:concepts,models,methods,and algorithms[M].John Wiley & Sons,2011.

[2]關(guān)宏強(qiáng),蔡福,王陽(yáng),等.短時(shí)間序列氣溫要素空間插值方法精度的比較研究[J].氣象與環(huán)境學(xué)報(bào),2007,23(5): 13-16.

Guan H Q,Cai F,Wang Y,et al.Comparison of different spatial interpolation methods for air temperature data of short-time series[J].Journal of Meteorology and Environment,2007,23(5): 13-16.(in Chinese)

[3]毛洋洋,趙艷霞,張祎,等.五個(gè)常見(jiàn)日太陽(yáng)總輻射模型在華北地區(qū)的有效性驗(yàn)證及分析[J].中國(guó)農(nóng)業(yè)氣象,2016,37(5): 520-530.

Mao Y Y,Zhao Y X,Zhang Y,et al.Validation and analysis of five general daily solar radiation estimation models used in Northern China[J].Chinese Journal of Agrometeorology,2016, 37(5):520-530.(in Chinese)

[4]郭兆夏,李星敏,朱琳,等.基于GIS的陜西省年降水量空間分布特征分析[J].中國(guó)農(nóng)業(yè)氣象,2010,31(S1): 121-123.

Guo Z X,Li X M,Zhu L,et al.Research on spatial distribution of annual precipitation in Shanxi Province based on GIS[J].Chinese Journal of Agrometeorology,2010,31(S1):121- 123.(in Chinese)

[5]Srebotnjak T,Carr G,de Sherbinin A,et al.A global water quality index and hot-deck imputation of missing data[J]. Ecological Indicators,2012,17:108-119.

[6]段建軍,王小利,高照良,等.黃土高原地區(qū)50年降水時(shí)空動(dòng)態(tài)與趨勢(shì)分析[J].水土保持學(xué)報(bào),2009,23(5):143-146.

Duan J J,Wang X L,Gao Z L,et al.Dynamics and trends analysis of annual precipitation in the Loess Plateau Region for 50 years[J].Journal of Soil and Water Conservation, 2009,23(5): 143-146.(in Chinese)

[7]馬晶,陳錫云,劉曉燕.地理因素輔助的黃土高塬典型流域面雨量制圖效果比較與評(píng)價(jià)[J].水土保持學(xué)報(bào),2016,30(6): 174-180.

Ma J,Chen X Y,Liu X Y.Comparison and evaluation of areal precipitation mapping effectiveness with consideration of geographic factors in the Loess Plateau[J].Journal of Soil and Water Conservation,2016,30(6):174-180.(in Chinese)

[8]Zhu Q A,Zhang W C,Zhao D Z.Topography-based spatial daily precipitation interpolation by means of PRISM and thiessen polygon analysis[J].Scientia Geographica Sinica, 2005,25(2):233-238.

[9]Gu Z H, Shi P J,Chen J.Precipitation interpolation research over regions with sparse meteorological stations:a case study in Xilingole League[J].Journal of Beijing Normal University (Natural Science),2006,42(2):204-208.

[10]邵曉梅,嚴(yán)昌榮,魏紅兵.基于Kriging插值的黃河流域降水時(shí)空分布格局[J].中國(guó)農(nóng)業(yè)氣象,2006,27(2):65-69.

Shao X M,Yan C R,Wei H B.Spatial and temporal structure of precipitation in the Yellow River Basin based on Kriging method[J].Chinese Journal of Agrometeorology,2006,27(2): 65-69.(in Chinese)

[11]Liu Y,Chen P Q,Zhang W.A spatial interpolation method for surface air temperature and its error analysis[J]. Chinese Journal of Atmospheric Sciences,2006,30(1):146-152.

[12]杜東升,廖玉芳,趙福華.湖南復(fù)雜地形下日平均氣溫空間插值方法探討[J].中國(guó)農(nóng)業(yè)氣象,2011,32(4):607-614.

Du D S,Liao Y F,Zhao F H.Study on the spatial interpolation method for daily mean air temperature over complex terrain in Hunan province[J].Chinese Journal of Agrometeorology, 2011, 32(4):607-614.(in Chinese)

[13]郭建茂,王錦杰,吳越,等.基于衛(wèi)星遙感與氣象站數(shù)據(jù)的水稻高溫?zé)岷ΡO(jiān)測(cè)和評(píng)估模型研究:以江蘇、安徽為例[J].農(nóng)業(yè)現(xiàn)代化研究,2017,38(2):298-306.

Guo J M,Wang J J,Wu Y,et al.Research on monitoring and modeling of rice heat injury based on satellite and meteorological station data:case study of Jiangsu and Anhui[J]. Research of Agricultural Modernization,2017,38 (2): 298- 306. (in Chinese)

[14]任利利,殷淑燕.漢江上游近50多年來(lái)氣溫變化特征與區(qū)域差異[J].農(nóng)業(yè)現(xiàn)代化研究,2013,34(3):348-352.

Ren L L,Yin S Y.Air temperature variation of the upper reaches of Hanjiang River in recent 50 years and its regional differences[J].Research of Agricultural Modernization,2013, 34(3):348-352.(in Chinese)

[15]鮑曉蕾,高輝,胡良平.多種填補(bǔ)方法在縱向缺失數(shù)據(jù)中的比較研究[J].中國(guó)衛(wèi)生統(tǒng)計(jì),2016,33(1):45-48.

Bao X L,Gao H,Hu L P.Comparative study of various imputation methods in dealing with longitudinal missing data[J].Chinese Health Statistics,2016,33(1):45-48.(in Chinese)

[16]楊軍,趙宇,丁文興.抽樣調(diào)查中缺失數(shù)據(jù)的插補(bǔ)方法[J].數(shù)理統(tǒng)計(jì)與管理,2008,27(5):821-832.

Yang J,Zhao Y,Ding W X.On imputation methods of missing data in survey sampling[J].Application of Statistics and Management,2008,27(5):821-832.(in Chinese)

[17]Ferrari G T,Ozaki V.Missing data imputation of climate datasets:implications to modeling extreme drought events[J]. Revista Brasileira de Meteorologia,2014,29(1):21-28.

[18]Saleem M U,Ahmed S R.Missing data imputations for upper air temperature at 24 standard pressure levels over pakistan collected from Aqua satellite[J].Journal of Data Analysis and Information Processing,2016,4(3):132.

[19]鄭小波,羅宇翔,于飛,等.西南復(fù)雜山地農(nóng)業(yè)氣候要素空間插值方法比較[J].中國(guó)農(nóng)業(yè)氣象,2008,29(4):458-462.

Zheng X B,Luo Y X,Yu F,et al.Comparisons of spatial interpolation methods for agro-climate factors in complex mountain areas of southwest China[J].Chinese Journal of Agrometeorology,2008,29(4):458-462.(in Chinese)

[20]孟岑,李裕元,吳金水,等.亞熱帶典型小流域總氮最大日負(fù)荷(TMDL)及影響因子研究:以金井河流域?yàn)槔齕J].環(huán)境科學(xué)學(xué)報(bào),2016,36(2):700-709.

Meng C,Li Y Y,Wu J S,et al.Study on total nitrogen TMDL and its contributing factors in typical subtropical watersheds: a case study of Jinjinghe watershed[J].Acta Scientiae Circumstantiae,2016,36(2):700-709.(in Chinese)

[21]李新,程國(guó)棟,盧玲.空間內(nèi)插方法比較[J].地球科學(xué)進(jìn)展,2000,15(3):260-265.

Li X,Cheng G D,Lu L.Comparison of spatial interpolation methods[J].Advance Earth Sciences,2000,15(3):260-265.(in Chinese)

[22]張曉琴,王敏.基于主成分分析的成分?jǐn)?shù)據(jù)缺失值插補(bǔ)法[J].應(yīng)用概率統(tǒng)計(jì),2016,32(1):101-110.

Zhang X Q,Wang M.Imputation of missing values for compositional data based on principal component analysis[J]. Chinese Journal of Applied Probability and Statistics,2016,32(1): 101-110.(in Chinese)

[23]陳林.基于GIS的流域水文數(shù)據(jù)的時(shí)空分析:以格蘭德河流域徑流數(shù)據(jù)為例[D].青島:山東科技大學(xué),2010.

Chen L.GIS-based spatial-temporal analysis of watershed hydrological data[D].Qingdao:Shandong University of Science and Technology,2010.(in Chinese)

[24]王國(guó)榮,俞耀明,徐兆亮,等.數(shù)值分析(第三版)[M].北京:機(jī)械工業(yè)出版社,2005.

Wang G R,Yu Y M,Xu Z L,et al.Numerical analysis(Third Edition)[M].Beijing:Mechanical Industry Press,2005.(in Chinese)

[25]殷杰,石銳.SAS中處理數(shù)據(jù)集缺失值方法的對(duì)比研究[J].計(jì)算機(jī)應(yīng)用,2007,27(b6):438-439.

Yin J,Shi R.A comparative study on the method of missing value of data set in SAS[J].Computer Applications,2007, 27(b6):438-439.(in Chinese)

[26]花琳琳,施念,楊永利,等.不同缺失值處理方法對(duì)隨機(jī)缺失數(shù)據(jù)處理效果的比較[J].鄭州大學(xué)學(xué)報(bào)(醫(yī)學(xué)版), 2012,47(3):315-318.

Hua L L,Shi N,Yang Y L,et al.Comparison of different methods in dealing with missing values of missing at random[J].Journal of Zhengzhou University(Medical Sciences), 2012,47(3):315-318.(in Chinese)

[27]蔡浩.地質(zhì)統(tǒng)計(jì)學(xué)在地層巖土參數(shù)分布規(guī)律研究中的應(yīng)用[D].蘇州:蘇州科技學(xué)院,2015.

Cai H.Applications of geostatistics to research on the distribution of the geotechnical parameters[D].Suzhou: Suzhou University of Science and Technology,2015.(in Chinese)

[28]Hong T,Kim C J,Jeong J,et al.Framework for approaching the minimum CV(RMSE) using energy simulation and optimization tool[J].Energy Procedia,2016,88:265-270.

[29]張桂銘,朱阿興,楊勝天,等.基于核密度估計(jì)的動(dòng)物生境適宜度制圖方法[J].生態(tài)學(xué)報(bào),2013,33(23):7590-7600.

Zhang G M,Zhu A X,Yang S T,et al.Mapping wildlife habitat suitability using kernel density estimation[J].Acta Ecologica Sinica,2013,33(23):7590-7600.(in Chinese)

[30]于力超,金勇進(jìn),王俊.缺失數(shù)據(jù)插補(bǔ)方法探討:基于最近鄰插補(bǔ)法和關(guān)聯(lián)規(guī)則法[J].統(tǒng)計(jì)與信息論壇,2015, 30(1): 35-40.

Yu L C,Jin Y J,Wang J.The research of missing data imputation method:based on nearest neighbor imputation and association rules[J].Statistic & Information Forum,2015, 30(1):35-40.(in Chinese)

[31]閻洪.薄板光順樣條插值與中國(guó)氣候空間模擬[J].地理科學(xué),2004,24(2):163-169.

Yan H.Modeling spatial distribution of climate in China using thin plate smoothing spline interpolation[J].Scientia Geographica Sinica,2004,24(2):163-169.

[32]Noor N M,Abdullah M M A B,Yahaya A S,et al.Comparison of linear interpolation method and mean method to replace the missing values in environmental data set[J]. Materials Science Forum,2015,(5):10.

[33]唐云輝,高陽(yáng)華.基于鄰域特征的溫度缺失值的填補(bǔ)方法[J].中國(guó)農(nóng)業(yè)氣象,2008,29(4):454-457.

Tang Y H,Gao Y H.Imputation method of missing temperature data based on neighborhood features[J].Chinese Journal of Agrometeorology,2008,29(4):454-457.(in Chinese)

[34]趙彥鋒,陳杰,齊力,等.不同采樣尺度下土壤圖和Kriging法的空間估值精度比較:以砂姜黑土典型地區(qū)的研究為例[J].土壤通報(bào),2011,(4):872-878.

Zhao Y F,Chen J,Qi L,et al.The comparison of soil map and Kriging methods for spatially prediction precision of soil properties with different sample spacings:a case of Shajiang black soil area[J].Chinese Journal of Soil Science,2011,(4): 872-878.(in Chinese)

[35]Yozgatligil C,Aslan S,Iyigun C,et al.Comparison of missing value imputation methods in time series:the case of Turkish meteorological data[J].Theoretical and Applied Climatology, 2013,112(1-2):143-167.

Performance Comparison of Different Interpolation Methods on Missing Values for Time Series Data——A Case Study of Meteorological and Hydrological Data in Subtropical Small Watershed

GAN Lei1, 2, ZHOU Jiao-gen2, SHI Jin2, 3, LI Xi2, SHEN Jian-lin2, LV Dian-qing1, LI Yu-yuan2,WU Jin-shui2

(1. College of Resources and Environmental Sciences, Hunan Normal University, Changsha 410081, China; 2. Key Laboratory of Agro- ecological Processes in Subtropical Region, Institute of Subtropical Agriculture, Chinese Academy of Sciences, Changsha 410125; 3. College of Engineering, Hunan Agricultural University, Changsha 410128)

The effective estimation of the missing values of time series data at the scale of point process could improve its data quality. The meteorological and hydrological data sets (daily maximum air temperature, daily minimum air temperature, daily solar radiation, daily rainfall and daily stream flow) were collected through the long-term field experiments in a typically small subtropical watershed in subtropical zone. The performance differences within five interpolation methods of linear interpolation method(LIM), K-Nearest neighbor interpolation method(KNNM), spline interpolation method(SIM), polynomial interpolation method(PIM) and kernel density estimation method(KDEM) were analyzed on the above-mentioned five data sets. The root mean square error(RMSE), absolute mean error(MAE) and Pearson correlation coefficient(r) were selected to evaluate the advantages and disadvantages of the five methods. The results showed that: (1) The estimation performance of LIM, SIM and KDEM was generally superior to the other two methods. (2) The estimation of the missing values of meteorological data (maximum temperature, minimum temperature and solar radiation) produced the varying values of the three evaluation indices with RMSE values of 1.81-6.35, MAE values of 1.30-4.20 and r values of 0.70-0.98 (P＜0.05), respectively. In contrast, the estimation of missing values of hydrological data (rainfall and stream flow) had relatively high values of RMSE and MAE which were 12.51-26.28 and 3.60-14.21, respectively, and low values of r (0.07-0.72). So the above-mentioned interpolation methods generally produced better estimation of missing values of meteorological data sets than those of hydrological data. (3) Additionally, the coefficient of variation (CV) of the above data sets linearly correlated with the evaluation indices (RMSE, MAE and r) (P＜0.05), and played an important role in affecting the valuation performance of the above-mentioned interpolation methods.

Missing values；Interpolation methods；Coefficient of variance；Time series

10.3969/j.issn.1000-6362.2018.03.007

甘蕾,周腳根,石錦,等.點(diǎn)源時(shí)間序列數(shù)據(jù)缺失值的估值方法比較:以小流域氣象和水文數(shù)據(jù)為例[J].中國(guó)農(nóng)業(yè)氣象,2018,39(3):195?204

收稿日期：2017-07-13

通訊作者。E-mail: zhoujg@isa.ac.cn

國(guó)家科技支撐計(jì)劃項(xiàng)目（2014BAD14B02）；水利部公益性行業(yè)科研專(zhuān)項(xiàng)經(jīng)費(fèi)項(xiàng)目（201501055）；湖南省地理學(xué)重點(diǎn)學(xué)科建設(shè)項(xiàng)目（20110101）

甘蕾（1992-），女，碩士生，主要從事水文生態(tài)與環(huán)境研究。E-mail:805150477@qq.com