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

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

?

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

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

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

對點(diǎn)源時間序列數(shù)據(jù)缺失值進(jìn)行有效估值能提升其數(shù)據(jù)質(zhì)量。為探究不同估值方法對點(diǎn)源時間序列數(shù)據(jù)缺失值的估值效果及其影響因素，以亞熱帶典型小流域長期定位觀測的每日氣象和水文數(shù)據(jù)（最高氣溫、最低氣溫、太陽輻射量、降雨量及地表徑流量）為例，以均方根誤差（RMSE）、絕對平均誤差（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種估值方法對氣象數(shù)據(jù)（最高氣溫、最低氣溫和太陽輻射量）缺失值估值的RMSE為1.81～6.35，MAE為1.30～4.20，r為0.70～0.98（P＜0.05），而對水文數(shù)據(jù)（降雨量和地表徑流量）缺失值估值的RMSE為12.54～26.28，MAE為3.60～14.21，r為0.07～0.72。可見，各估值方法對氣象數(shù)據(jù)的估值性能強(qiáng)于對水文數(shù)據(jù)；（3）上述數(shù)據(jù)集的變異系數(shù)（CV）與估值評估指標(biāo)（RMSE、MAE及r）線性相關(guān)（P＜0.05），是影響估值性能的重要因素。

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

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

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

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

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

1 材料與方法

1.1 數(shù)據(jù)來源

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

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

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

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

1.2 估值方法

1.2.1 線性內(nèi)插法

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

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

1.2.2 K-最近鄰插值法

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

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

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

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

多項(xiàng)式插值法（PIM）[23]是用多項(xiàng)式對一列數(shù)據(jù)進(jìn)行線性擬合，再對給定待估值點(diǎn)進(jìn)行估值的過程。給定時間序列數(shù)據(jù)集Y= {Y（t1），Y（t2），…, Y（tn）}和待估值點(diǎn)Y（tj），用多項(xiàng)式函數(shù)f（t）=β0+β1t+β2t2+…+βntn對時間序列數(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)知識。對給定缺失值Y（tj），核密度估值方法估算式為

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

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

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

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

2 結(jié)果與分析

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

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

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

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

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

由表2可見，5種估值方法對不同數(shù)據(jù)集的估值性能具有較大差異。對于變異系數(shù)較小的氣象數(shù)據(jù)（最高氣溫、最低氣溫及太陽輻射量），LIM、KNNM、SIM、PIM及KDEM 5種方法皆表現(xiàn)較佳，估算值與實(shí)測值相關(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對最高氣溫估值較好，其RMSE、MAE、r為2.91℃、2.12℃、0.95（P＜0.05），SIM對最低氣溫的估值效果與LIM相同，同為最佳方法，且對太陽輻射量估值也較好；KNNM和PIM兩種方法表現(xiàn)最差，誤差大，相關(guān)性弱。

對于變異系數(shù)較大的降雨量數(shù)據(jù)，LIM、KNNM、SIM、PIM及KDEM 5種方法估值效果皆不佳，RMSE和MAE偏大，估算值與實(shí)測值相關(guān)性不顯著（r為0.07～0.13），其中，KDEM相對較優(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種估值方法對日最高氣溫、日最低氣溫、日太陽輻射量以及日地表徑流量數(shù)據(jù)的估值結(jié)果較為可靠，但對日降雨量的估值精度不高，這可能是因?yàn)槿战涤炅繙y試數(shù)據(jù)集的變異系數(shù)過大（CV=162.57%）。另外，LIM、SIM和KDEM 3種估值方法對這5種缺失數(shù)據(jù)的估值效果較好。

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

注：LIM為線性內(nèi)插法、KNNM為K-最近鄰插值法、SIM為樣條插值法、PIM為多項(xiàng)式插值法、KDEM為核密度估值法；RMSE為均方根誤差、MAE絕對值平均誤差、r為估算值與實(shí)測值的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ù)對缺失值估算結(jié)果的影響

將最高氣溫、最低氣溫、太陽輻射量、降雨量以及地表徑流量訓(xùn)練數(shù)據(jù)集的變異系數(shù)（CV）與交叉驗(yàn)證指標(biāo)值（RMSE、MAE、r）進(jìn)行線性擬合分析，結(jié)果見圖3。由圖3可見，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。說明變異系數(shù)是影響缺失值估值結(jié)果的重要因素，變異系數(shù)越大，均方根誤差和絕對平均誤差越大，相關(guān)性越??；反之，變異系數(shù)越小，誤差越小，相關(guān)性越大。進(jìn)一步分析不同估值方法的估值性能對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種估值方法輸出的估值誤差，超過69%是與數(shù)據(jù)集固有的變異性有關(guān)。因此，在本研究中CV是影響估值方法LIM、KNNM、SIM、PIM及KDEM的估值性能的關(guān)鍵因素。數(shù)據(jù)集的變異系數(shù)越大，LIM、KNNM、SIM、PIM及KDEM 5種方法的估值誤差越大，輸出的預(yù)測值與實(shí)測值的擬合度越小，對估值結(jié)果的準(zhǔn)確性影響越大[28]。

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

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

3 結(jié)論與討論

3.1 討論

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

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

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

本研究表明，數(shù)據(jù)集變異系數(shù)小，離散程度小，則5種估值方法對數(shù)據(jù)集缺失值的估值效果較優(yōu)；反之，數(shù)據(jù)集變異系數(shù)大，離散程度大，5種估值方法的估值效果皆顯著下降。不同估值方法處理后的估值驗(yàn)證指標(biāo)對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ù)據(jù)集估值結(jié)果的準(zhǔn)確性最高；Yozgatligil等[35]也證實(shí)土耳其降水、溫度數(shù)據(jù)集CV值越小，對缺失值估值結(jié)果越可靠。

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

3.2 結(jié)論

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

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

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

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

[2]關(guān)宏強(qiáng),蔡福,王陽,等.短時間序列氣溫要素空間插值方法精度的比較研究[J].氣象與環(huán)境學(xué)報,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]毛洋洋,趙艷霞,張祎,等.五個常見日太陽總輻射模型在華北地區(qū)的有效性驗(yàn)證及分析[J].中國農(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].中國農(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年降水時空動態(tài)與趨勢分析[J].水土保持學(xué)報,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]馬晶,陳錫云,劉曉燕.地理因素輔助的黃土高塬典型流域面雨量制圖效果比較與評價[J].水土保持學(xué)報,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插值的黃河流域降水時空分布格局[J].中國農(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].中國農(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)測和評估模型研究:以江蘇、安徽為例[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多年來氣溫變化特征與區(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].中國衛(wèi)生統(tǒng)計,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)計與管理,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].中國農(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é)報,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]李新,程國棟,盧玲.空間內(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)計,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ù)據(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]王國榮,俞耀明,徐兆亮,等.數(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ù)集缺失值方法的對比研究[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]花琳琳,施念,楊永利,等.不同缺失值處理方法對隨機(jī)缺失數(shù)據(jù)處理效果的比較[J].鄭州大學(xué)學(xué)報(醫(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)計學(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].生態(tài)學(xué)報,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)計與信息論壇,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]閻洪.薄板光順樣條插值與中國氣候空間模擬[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]唐云輝,高陽華.基于鄰域特征的溫度缺失值的填補(bǔ)方法[J].中國農(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].土壤通報,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ù)據(jù)缺失值的估值方法比較:以小流域氣象和水文數(shù)據(jù)為例[J].中國農(nóng)業(yè)氣象,2018,39(3):195?204

收稿日期：2017-07-13

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

國家科技支撐計劃項(xiàng)目（2014BAD14B02）；水利部公益性行業(yè)科研專項(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