盧月明，王 亮，仇阿根，張用川，2，趙陽陽

1.中國測繪科學(xué)研究院，北京 100830

2.武漢大學(xué) 資源與環(huán)境科學(xué)學(xué)院，武漢 430079

1 引言

克里金插值法是一種空間最優(yōu)線性無偏估計方法，是基于空間屬性在空間位置上的分布情況，利用半變異函數(shù)確定周圍待插值點的權(quán)重以實現(xiàn)待插值點屬性的估計[1]。該方法綜合考慮了變量的空間結(jié)構(gòu)性與隨機性，通過模擬地理現(xiàn)象空間分布的相關(guān)性和變異性進(jìn)行統(tǒng)計分析，因此，克里金插值方法被廣泛應(yīng)用于氣象[2]、土壤[3]等領(lǐng)域。然而，李杰等的研究結(jié)果表明，無論哪種插值方法，都需要足夠數(shù)量的樣本數(shù)據(jù)才能保證插值結(jié)果準(zhǔn)確可靠[4]。

半監(jiān)督學(xué)習(xí)是介于監(jiān)督學(xué)習(xí)與非監(jiān)督學(xué)習(xí)之間的一種學(xué)習(xí)方式，其學(xué)習(xí)樣本既包括標(biāo)記樣本，又包括未標(biāo)記樣本，既可以利用大量容易獲得的未標(biāo)記樣本，減輕標(biāo)記樣本的工作量，又可以利用標(biāo)記樣本獲得更高效的學(xué)習(xí)模型[5]。協(xié)同訓(xùn)練是一種半監(jiān)督學(xué)習(xí)方法，它利用雙視圖訓(xùn)練兩個分類器來互相標(biāo)記樣本以擴大訓(xùn)練集，以此借助未標(biāo)記樣本提升學(xué)習(xí)性能[6]。Yang等的實驗結(jié)果說明半監(jiān)督協(xié)同訓(xùn)練可利用未標(biāo)記樣本輔助訓(xùn)練，提升只有少量標(biāo)記樣本時模型的學(xué)習(xí)性能[7]。協(xié)同訓(xùn)練法在多視圖數(shù)據(jù)上實驗效果很好，已在理論上得到證明：當(dāng)兩個充分冗余視圖滿足條件獨立時，通過協(xié)同訓(xùn)練可以利用未標(biāo)記樣本把弱分類器的精度提升到任意高[8]。Wang和Zhou證明了協(xié)同訓(xùn)練法的充分必要性定理，結(jié)果表明，協(xié)同訓(xùn)練只關(guān)心權(quán)值矩陣的性質(zhì)，而并不在意權(quán)值矩陣是否通過多視圖得到，這確認(rèn)了基于分歧的學(xué)習(xí)方法并不需要多視圖，僅要求分類器之間存在適當(dāng)?shù)姆制?，其必要性條件是每個未標(biāo)記樣本在聯(lián)合圖中都與標(biāo)記樣本連通[9]。

近年來，對半監(jiān)督學(xué)習(xí)方法的研究，主要聚焦于解決半監(jiān)督學(xué)習(xí)中的分類問題，而對半監(jiān)督學(xué)習(xí)中的回歸問題的研究相對較少，一個主要原因是半監(jiān)督學(xué)習(xí)中的聚類假設(shè)在回歸問題上不成立，且在回歸分析中標(biāo)記置信度的計算也比較困難。對此，Zhou等提出一種協(xié)同回歸計算方法（Co-training Regression，COREG），該方法基于不同的距離度量或不同的k值產(chǎn)生不同的k近鄰回歸模型，然后基于預(yù)測一致性來選擇置信度高的未標(biāo)記樣本進(jìn)行標(biāo)記[10]；馬蕾等利用SVM（Support Vector Machine，SVM）來建立回歸器，實現(xiàn)了基于SVM的半監(jiān)督回歸訓(xùn)練方法[5]；趙陽陽等提出了一種基于半監(jiān)督學(xué)習(xí)的地理加權(quán)回歸方法（Semi-supervised Learning Geographic Weighted Regression，SSLGWR），并分別使用模擬數(shù)據(jù)與真實數(shù)據(jù)說明了SSLGWR的預(yù)測結(jié)果顯著優(yōu)于單純的地理加權(quán)回歸方法[11]；趙陽陽等基于協(xié)同訓(xùn)練，提出了協(xié)同GTWR方法，結(jié)果表明協(xié)同GTWR的性能相對于使用不同核函數(shù)的GTWR均有所提升[12]；馬蕾等提出基于SVM協(xié)同訓(xùn)練的回歸模型，該模型適用于處理大量有輸出的輸入情況，緩解了使用單一回歸模型所造成的錯誤累加問題，提升了回歸模型的泛化能力[5]。綜合上述研究發(fā)現(xiàn)：半監(jiān)督學(xué)習(xí)理論可有效提升模型的精度，廣泛應(yīng)用于樣類分類、語音識別等領(lǐng)域，但隨著訓(xùn)練的進(jìn)行，自動標(biāo)記中的噪音會不斷地累積，其負(fù)作用不斷增大。

在克里金方法建模過程中，標(biāo)記樣本數(shù)據(jù)量的多少直接關(guān)系到模型的精度，當(dāng)標(biāo)記樣本較少時，通常難以構(gòu)建可靠的模型。而在實際應(yīng)用中，常常難以獲取足夠數(shù)量的標(biāo)記樣本，如PM2.5濃度觀測數(shù)據(jù)。本文針對克里金模型在樣本較少時模型精度低這一問題，提出基于半監(jiān)督學(xué)習(xí)的克里金插值方法，即自訓(xùn)練克里金插值模型（Self-Training Kriging，STK）和協(xié)同訓(xùn)練克里金插值模型（Co-Training Kriging，CTK），并使用北京地區(qū)2017年4月和5月的PM2.5濃度數(shù)據(jù)進(jìn)行實驗。通過與普通克里金插值方法進(jìn)行對比實驗，以平均絕對誤差（Mean Absolute Error，MAE）、均方根誤差（Root Mean Square Error，RMSE）作為評價指標(biāo)來說明模型的準(zhǔn)確性。通過以上方法，可獲得PM2.5在不同地點的濃度，并實現(xiàn)其屬性值從點到面的轉(zhuǎn)變，從而為PM2.5空間分布的預(yù)測及可視化提供一種手段。

2 研究方法

2.1 克里金插值

克里金插值作為地統(tǒng)計學(xué)的核心，用來估算未采樣位置的屬性值，其研究對象是區(qū)域化變量，是一種最優(yōu)無偏估計方法。它通過變差函數(shù)來量化觀測數(shù)據(jù)的空間相關(guān)性，建立函數(shù)關(guān)系，將標(biāo)記樣本代入函數(shù)關(guān)系，計算權(quán)重系數(shù)，從而建立插值模型，進(jìn)行分析或預(yù)測。這里的標(biāo)記樣本是指含有自變量和因變量的樣本數(shù)據(jù)，未標(biāo)記樣本指只含有自變量，不含有因變量的樣本數(shù)據(jù)。

克里金插值法可表示為：

其中，λi為權(quán)重系數(shù)，表示各空間樣本點xi處的觀測值Z(xi)對其估計值Z?(x0)的貢獻(xiàn)程度?？死锝鸩逯捣椒ǖ年P(guān)鍵在于求權(quán)重系數(shù)，權(quán)重系數(shù)的計算需要滿足兩個假設(shè)條件：

即

通過求解上述克里金方程組，求出權(quán)重系數(shù)和拉格朗日因子，代入式（1）、（2），即可求出估計值與估計方差。

2.2 半監(jiān)督學(xué)習(xí)

2.2.1 自訓(xùn)練克里金模型

自訓(xùn)練方法最早由Fralick等[14]提出，自訓(xùn)練方法是一種半監(jiān)督學(xué)習(xí)方法，在每一輪的訓(xùn)練過程中反復(fù)運用監(jiān)督學(xué)習(xí)方法，將上一輪標(biāo)記結(jié)果最優(yōu)的樣例和它的類標(biāo)簽一起加入到當(dāng)前訓(xùn)練樣本集中，用自己產(chǎn)生的結(jié)果不斷訓(xùn)練自己[15]。本文將自訓(xùn)練理論應(yīng)用于克里金插值模型，得到自訓(xùn)練克里金插值模型（STK）。

STK模型的算法流程圖如圖1所示。

步驟1確定標(biāo)記樣本集、未標(biāo)記樣本集，初始化克里金插值模型參數(shù)，該克里金插值模型采用高斯核函數(shù)。

步驟2 Kriging模型對未標(biāo)記樣本集進(jìn)行插值估算。

步驟3從插值結(jié)果中選擇置信度最高的未標(biāo)記樣本及其預(yù)測結(jié)果加入到插值模型的標(biāo)記樣本集中，并從未標(biāo)記樣本集中去除該樣本。

步驟4重新訓(xùn)練克里金插值模型，直至訓(xùn)練一定數(shù)量的未標(biāo)記樣本為止。

2.2.2 協(xié)同訓(xùn)練克里金模型

協(xié)同訓(xùn)練是一種半監(jiān)督學(xué)習(xí)方法，可以在少量有標(biāo)記樣本和大量未標(biāo)記樣本的基礎(chǔ)上，通過不斷迭代，使得不同學(xué)習(xí)器互相學(xué)習(xí)[10]。其原理是建立兩個學(xué)習(xí)器，分別在這兩個學(xué)習(xí)器上使用標(biāo)記數(shù)據(jù)與未標(biāo)記數(shù)據(jù)，利用學(xué)習(xí)器和標(biāo)記數(shù)據(jù)來標(biāo)記未標(biāo)記數(shù)據(jù)，不斷更新另一個學(xué)習(xí)器的標(biāo)記數(shù)據(jù)，通過這樣不斷互相學(xué)習(xí)，得到未標(biāo)記數(shù)據(jù)的標(biāo)記，從而擴充標(biāo)記數(shù)據(jù)的樣本量，提升模型性能[12]。

本文基于COREG算法的理論基礎(chǔ)[10]，將克里金插值模型與半監(jiān)督學(xué)習(xí)協(xié)同訓(xùn)練理論相結(jié)合，得到基于半監(jiān)督學(xué)習(xí)的協(xié)同訓(xùn)練克里金插值模型，即協(xié)同訓(xùn)練克里金插值模型（CTK）。協(xié)同訓(xùn)練克里金插值模型不僅集成了協(xié)同訓(xùn)練在小樣本中的獨特優(yōu)勢，彌補了樣本少的不足，也集成了克里金插值模型在地理應(yīng)用中的特點。

CTK模型的算法流程圖如圖2所示。

步驟1確定標(biāo)記樣本集、未標(biāo)記樣本集，初始化兩個Kriging模型，兩個模型分別為基于高斯核函數(shù)的kriging模型和基于指數(shù)核函數(shù)的Kriging模型。

步驟2每個插值模型對其未標(biāo)記樣本集進(jìn)行插值估算，從插值結(jié)果中選擇置信度最高的未標(biāo)記樣本及其插值結(jié)果加入到另一插值模型的標(biāo)記樣本集中，并從未標(biāo)記樣本集中去除該樣本。

步驟3重復(fù)進(jìn)行步驟2，直至訓(xùn)練一定數(shù)量的未標(biāo)記樣本為止。

步驟4最終插值結(jié)果為兩個插值模型插值結(jié)果的平均值。

2.2.3 置信度計算方法

置信度用于從若干未標(biāo)記樣本中選取最優(yōu)的訓(xùn)練結(jié)果，滿足預(yù)測一致性原則，即具有真實標(biāo)記的樣本應(yīng)能夠體現(xiàn)出插值的內(nèi)在規(guī)律。在模型學(xué)習(xí)過程中，每一輪選取的未標(biāo)記樣本都會對新插值模型的精度產(chǎn)生影響。基于預(yù)測一致性原則選取置信度高的結(jié)果，即置信度越高說明越接近真實值。因此，插值模型通過高置信度選擇的樣本應(yīng)該是使插值模型與標(biāo)記樣本最一致的樣本[16]。本文采用均方誤差（Mean Square Error，MSE）作為置信度評判的指標(biāo)，即如果在未標(biāo)記樣本中存在一條數(shù)據(jù)，當(dāng)其加入標(biāo)記樣本集后，使得插值模型的均方誤差變小且變小的幅度最大，則這條數(shù)據(jù)就為置信度最高的未標(biāo)記樣本[10]。置信度計算方法如下：

式中，yL為標(biāo)記樣本的真實值；y?L為標(biāo)記樣本在原插值模型上的估計值；y?′L為標(biāo)記樣本在新插值模型上的估計值，新插值模型是指加入未標(biāo)記樣本后重新構(gòu)建的插值模型。

當(dāng) ξXx∈μ＞0 時 ，令 N(x,u,ν)=arc max(ξXx∈μ) 。N(x,u,ν)即為置信度最高的未標(biāo)記樣本。 ξXx∈μ＞0說明未標(biāo)記樣本的加入使得插值模型性能有所提升。置信度最大說明插值模型性能提升幅度最大，即所選數(shù)據(jù)是參與訓(xùn)練的未標(biāo)記樣本中置信度最高的數(shù)據(jù)。

圖1 自訓(xùn)練克里金插值模型的算法流程圖

圖2 協(xié)同訓(xùn)練克里金插值模型的算法流程圖

3 插值實驗

3.1 研究區(qū)概況

北京位于東經(jīng) 115.7°～117.4°，北緯 39.4°～41.6°，中心位于北緯 39°54′20″，東經(jīng) 116°25′29″，總面積達(dá)16 410.54 km2，全市常住人口達(dá)2 100多萬人，是中國的首都、政治中心、文化中心、科技創(chuàng)新中心。近年來，以PM2.5和PM10為主的大氣顆粒物濃度急劇升高，導(dǎo)致北京的霧霾天氣頻發(fā)，以致頻頻啟動“重霧霾橙色預(yù)警”。

因此，開展大氣污染的相關(guān)研究，對掌控空氣質(zhì)量分布狀況，采取相關(guān)防控舉措具有一定意義。

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

本文選取2017年4月和5月北京地區(qū)35個監(jiān)測站點每小時監(jiān)測數(shù)據(jù)，包括 PM2.5、PM10、NO2、CO、SO2、O3等大氣污染物的濃度數(shù)據(jù)，以及每個監(jiān)測站點的經(jīng)度與緯度，上述監(jiān)測數(shù)據(jù)抓取自網(wǎng)站www.pm25.in。35個監(jiān)測站點遍布北京城區(qū)及其郊縣，其覆蓋范圍基本可以反映整個北京地區(qū)的空氣質(zhì)量狀況。本文將監(jiān)測站點按照1∶1的比例隨機分成標(biāo)記樣本集（標(biāo)記監(jiān)測站點18個）與未標(biāo)記樣本集（未標(biāo)記監(jiān)測站點17個）。標(biāo)記監(jiān)測站點與未標(biāo)記監(jiān)測站點分布，如圖3所示。

圖3 北京地區(qū)空氣質(zhì)量監(jiān)測站點分布圖

3.3 數(shù)據(jù)處理與校驗

首先對每個監(jiān)測站點每天的24組PM2.5數(shù)據(jù)求平均得到PM2.5日均值，再由日均值計算出每個監(jiān)測站點的月均值。然后對監(jiān)測站點的PM2.5濃度值進(jìn)行空間數(shù)據(jù)探索分析，來尋找數(shù)據(jù)內(nèi)在的規(guī)律性，確定數(shù)據(jù)是否適合使用克里金插值法。對于不適宜的數(shù)據(jù)將通過數(shù)據(jù)變換，使原來不適合于插值的數(shù)據(jù)可以進(jìn)行空間插值。由直方圖分析得出4月的峰度為3.21，偏態(tài)為0.26；5月的峰度為3.04，偏態(tài)為-0.57。其中，峰度（Kurtosis）用來描述數(shù)據(jù)分布的高度，標(biāo)準(zhǔn)正態(tài)分布的峰度應(yīng)為3，峰度值越接近3表示該數(shù)據(jù)越接近正態(tài)分布；偏態(tài)（Skewness）用來描述數(shù)據(jù)左右的對稱性，標(biāo)準(zhǔn)正態(tài)分布的偏態(tài)值應(yīng)為0，如果偏態(tài)值大于0，稱正偏態(tài)或右偏態(tài)，此時大部分?jǐn)?shù)據(jù)集中于左邊，如果偏態(tài)值小于0，稱負(fù)偏態(tài)或左偏態(tài)，此時大部分?jǐn)?shù)據(jù)集中于右邊。北京地區(qū)2017年4月、5月的PM2.5濃度Normal QQPlot圖如圖4、圖5所示。

圖4 4月PM 2.5的QQPlot圖

圖5 5月PM 2.5的QQPlot圖

通過對監(jiān)測站點的PM2.5濃度值進(jìn)行直方圖、正態(tài)QQPlot圖以及半變異函數(shù)分析得出，數(shù)據(jù)很接近正態(tài)分布且數(shù)據(jù)的相關(guān)性較高，說明該數(shù)據(jù)進(jìn)行空間插值有意義。

3.4 結(jié)果與分析

本文將35個監(jiān)測站點隨機等分為2組（標(biāo)記監(jiān)測站點集與未標(biāo)記監(jiān)測站點集），分別使用普通克里金插值法（Kriging）、自訓(xùn)練克里金插值法（STK）和協(xié)同訓(xùn)練克里金插值法（CTK）對未標(biāo)記站點集進(jìn)行空間插值分析。為驗證模型的精度，通過平均絕對誤差（MAE）、均方根誤差（RMSE）[17]來對插值結(jié)果進(jìn)行評估。其中平均絕對誤差反映估計值可能的誤差范圍，均方根誤差反映插值函數(shù)的反演靈敏度和極值效應(yīng)，這兩個指標(biāo)均是越小代表模型精度越高。經(jīng)計算，得到3種插值方法的插值精度與精度提升百分比如表1、表2所示。

表1 插值結(jié)果精度分析 μg/m3

表2 插值精度提升百分比 %

由表1、表2中的各項指標(biāo)可知，基于自訓(xùn)練的克里金插值法（STK）相對于普通Kriging法性能沒有提升反而有所下降，這種現(xiàn)象是由于早期加入標(biāo)記樣本集中未標(biāo)記樣本的誤差在后期訓(xùn)練過程中不斷累積放大的結(jié)果?；趨f(xié)同訓(xùn)練的克里金插值方法（CTK）相對于普通Kriging法在4月、5月均有較大幅度的提升，CTK法相對于Kriging插值法（CTK-Kriging）平均絕對誤差（MAE）提升程度在10%左右，均方根誤差（RMSE）提升程度在11%左右，說明采用兩個模型進(jìn)行協(xié)同訓(xùn)練可以有效地削弱早期誤差對后期訓(xùn)練的影響，通過合理地利用未標(biāo)記樣本提升了模型的性能。

克里金插值法、自訓(xùn)練克里金插值法、協(xié)同訓(xùn)練克里金插值法對北京2017年4月、5月PM2.5濃度的插值結(jié)果圖如下所示。其中圖6、圖7和圖8為分別使用克里金插值法、自訓(xùn)練克里金插值法和協(xié)同訓(xùn)練克里金插值法對北京市4月PM2.5濃度的插值結(jié)果圖；圖9、圖10和圖11為分別使用克里金插值法、自訓(xùn)練克里金插值法和協(xié)同訓(xùn)練克里金插值法對北京市5月PM2.5濃度的插值結(jié)果圖。

圖6 克里金模型4月PM 2.5插值圖

圖7 自訓(xùn)練克里金模型4月PM 2.5插值圖

圖8 協(xié)同訓(xùn)練克里金模型4月PM 2.5插值圖

圖9 克里金模型5月PM 2.5插值圖

圖10 自訓(xùn)練克里金模型5月PM 2.5插值圖

圖11 協(xié)同訓(xùn)練克里金模型5月PM 2.5插值圖

由插值結(jié)果圖可看出5月有較明顯的“牛眼”現(xiàn)象。一方面是由于原始測量值中存在奇異值，即孤立點數(shù)據(jù)明顯高于或低于周圍監(jiān)測點數(shù)據(jù)；另一方面是由于插值區(qū)域中監(jiān)測站點分布不均勻，且插值時將距離作為權(quán)重，忽略了方位等其他因素的影響，導(dǎo)致最終結(jié)果圖中形成以插值點為圓心的圈狀現(xiàn)象。此外，5月誤差相較于4月略大，也與該數(shù)據(jù)中存在奇異值有著直接關(guān)系。分別對比4月三種方法的插值結(jié)果圖與5月三種方法的插值結(jié)果圖可發(fā)現(xiàn)，使用自訓(xùn)練克里金方法得到的插值圖的顏色均較深，即屬性值（PM2.5濃度）均較高；克里金方法和協(xié)同訓(xùn)練克里金方法得到的插值圖的顏色存在由淺到深的過渡，更符合實際情況，從這一方面也可看出自訓(xùn)練克里金方法的誤差相較于另外兩種方法偏大。

4 結(jié)束語

本文針對數(shù)據(jù)量較小時，克里金方法插值精度低這一問題，將克里金插值模型與半監(jiān)督學(xué)習(xí)理論相結(jié)合，利用半監(jiān)督學(xué)習(xí)使用未標(biāo)記樣本參與訓(xùn)練來提升回歸模型性能的優(yōu)勢，提出了基于半監(jiān)督學(xué)習(xí)的克里金插值模型，即STK和CTK。這兩個模型既具有半監(jiān)督學(xué)習(xí)的優(yōu)點，適用于解決只有少量標(biāo)記樣本的情況，又可以將離散點的測量數(shù)據(jù)轉(zhuǎn)換為連續(xù)的數(shù)據(jù)曲面，以便與其他空間現(xiàn)象的分布模式進(jìn)行比較。本文采用2017年4月和5月北京地區(qū)的PM2.5濃度數(shù)據(jù)進(jìn)行對比實驗，結(jié)果表明CTK插值法采用兩個協(xié)同訓(xùn)練的回歸模型，削弱了僅使用單一模型的STK法中錯誤累積放大的缺點，提高了插值模型的泛化能力。此外，本文未考慮影響PM2.5濃度的因素，如風(fēng)力、濕度、高程等，未來不僅要在該插值方法上深入研究，還要引入更多的影響因素，以進(jìn)一步提高插值精度。