張 波,宋國君

大規(guī)?？諝赓|(zhì)量監(jiān)測數(shù)據(jù)缺失處理方法實證研究

張 波,宋國君*

(中國人民大學(xué)環(huán)境學(xué)院,北京 100872)

基于2016年1月至2021年7月的全國1654個國控監(jiān)測點小時級的6種污染物空氣質(zhì)量監(jiān)測數(shù)據(jù),研究缺失值處理方法、效果及其影響.模擬實驗表明交替最小二乘下的低秩矩陣插補算法相比于其他缺失值處理方法擁有更小的均方根誤差、平均百分比誤差,更高的相關(guān)系數(shù)和更快的運算速度,在大規(guī)模數(shù)據(jù)集上性能更優(yōu).實證分析表明應(yīng)用文本方法得到的插補值是有效且合理的,缺失值插補前后污染物濃度評估值會有±10%以內(nèi)的變化,插補后的數(shù)據(jù)集更加準(zhǔn)確和完備.本文建議在基于空氣質(zhì)量監(jiān)測數(shù)據(jù)研究時應(yīng)先采用本文中的缺失數(shù)據(jù)處理方法,對監(jiān)測數(shù)據(jù)中存在的缺失數(shù)據(jù)進行插補,提高研究所使用監(jiān)測數(shù)據(jù)的完整性,保證相關(guān)計算結(jié)果的準(zhǔn)確性和有效性.

監(jiān)測數(shù)據(jù)；大數(shù)據(jù)處理；缺失值處理；矩陣填充；實證研究

我國2016年實施的《環(huán)境空氣質(zhì)量標(biāo)準(zhǔn) (GB3095)》(簡稱“標(biāo)準(zhǔn)”)中對污染物監(jiān)測數(shù)據(jù)的有效性做出明確規(guī)定,一天24h中有至少20h以上數(shù)據(jù)時,計算SO2、NO2、CO、PM2.5和PM10的日均濃度值才是有效日均數(shù)據(jù),O3監(jiān)測數(shù)據(jù)則要求每連續(xù)8h中至少有6h數(shù)據(jù)時才為有效數(shù)據(jù).但在監(jiān)測設(shè)備的日常運行中,由于某些不可預(yù)期因素(如服務(wù)器宕機、網(wǎng)絡(luò)中斷、設(shè)備故障等)導(dǎo)致監(jiān)測數(shù)據(jù)缺失.目前國家公開發(fā)布的空氣質(zhì)量監(jiān)測數(shù)據(jù)以小時數(shù)據(jù)為最小顆粒度,小時數(shù)據(jù)缺失會對后續(xù)日均值、月均值、年度均值、超標(biāo)率、二級標(biāo)準(zhǔn)天數(shù)等指標(biāo)的計算造成直接影響.2018年生態(tài)環(huán)境部印發(fā)的《城市環(huán)境空氣質(zhì)量排名技術(shù)規(guī)定》[1](簡稱“規(guī)定”)中指出SO2、NO2、PM2.5、PM10的評價濃度為日均濃度,O3的評價濃度為日最大8h平均值的90%分位數(shù),CO的評價濃度為日均濃度的95%分位數(shù),不符合有效性規(guī)定的數(shù)據(jù)不能參與計算,為空氣質(zhì)量監(jiān)測數(shù)據(jù)的計算提供了標(biāo)準(zhǔn).大規(guī)?？諝赓|(zhì)量監(jiān)測數(shù)據(jù)是進行污染物時空分布規(guī)律研究的基礎(chǔ)數(shù)據(jù),目前大部分研究均基于中國環(huán)境監(jiān)測總站發(fā)布的連續(xù)在線監(jiān)測數(shù)據(jù)[2-8],數(shù)據(jù)的缺失和不同的缺失處理方法會對研究結(jié)果產(chǎn)生一定影響[9],但空氣質(zhì)量監(jiān)測數(shù)據(jù)中數(shù)據(jù)缺失的影響究竟有多大,使用何種缺失處理方法更加可靠,這方面的研究還非常有限,因而研究如何更加客觀、科學(xué)的處理大規(guī)模環(huán)境監(jiān)測數(shù)據(jù)中缺失問題具有重要的現(xiàn)實應(yīng)用價值.

目前國內(nèi)有關(guān)環(huán)境監(jiān)測數(shù)據(jù)缺失問題的研究主要集中在制度建設(shè)和管理規(guī)范方面,如建立外場監(jiān)測和實驗室分析的兩級質(zhì)量保證體系[10],通過開展專項質(zhì)控工作提升監(jiān)測數(shù)據(jù)質(zhì)量[11],加強對生態(tài)環(huán)境監(jiān)測機構(gòu)的監(jiān)管提升數(shù)據(jù)質(zhì)量[12],通過建立從業(yè)人員監(jiān)管制度確保數(shù)據(jù)質(zhì)量[13],以及對國外經(jīng)驗的總結(jié)和借鑒[14],這些研究從體制機制的角度研究如何提升監(jiān)測數(shù)據(jù)質(zhì)量,但缺少對具體現(xiàn)實具體問題的探討,尤其是存在大規(guī)模數(shù)據(jù)缺失情況下該如何處理,用什么方法處理,處理后的影響如何等問題還鮮有研究.國際上環(huán)境監(jiān)測數(shù)據(jù)的缺失問題一直以來都是研究者非常關(guān)注的問題之一,很多缺失值處理方法也被應(yīng)用于空氣質(zhì)量監(jiān)測數(shù)據(jù)處理和研究中,如最近鄰算法[15]、EM算法[16]、線性回歸[17]、簡單插補[18-20]、多重插補[18,21]、最小核范數(shù)[22]和低秩矩陣插補[23-24]等,多方法對比研究[21,24]發(fā)現(xiàn)多重插補方法和低秩矩陣方法在數(shù)據(jù)缺失較大的情況下更有優(yōu)勢,能夠產(chǎn)生偏差更小的插補值來替代缺失值.以往研究對象面向的是少數(shù)監(jiān)測站點和較短數(shù)據(jù)序列,現(xiàn)實中由于監(jiān)測手段快速改進,空氣質(zhì)量監(jiān)測數(shù)據(jù)已經(jīng)出現(xiàn)爆發(fā)式增加,監(jiān)測點數(shù)量和數(shù)量量都達到很大規(guī)模,如何處理海量環(huán)境數(shù)據(jù)中的缺失問題,已經(jīng)成為越來越需要密切關(guān)注的問題.本文基于全國空氣質(zhì)量國控點的小時級監(jiān)測數(shù)據(jù),利用低秩矩陣插補方法處理海量監(jiān)測數(shù)據(jù)中大量缺失的問題,該方法在保證插補精度的同時能夠大幅度降低大規(guī)模數(shù)據(jù)處理所需時間,得到的插補值服具有較高有效性.

1 模型方法

在缺失值處理方法中,簡單插補和最近鄰算法方法最簡單,但處理效果較差,EM算法只適用于分布為多元正態(tài)分布的數(shù)據(jù),線性回歸、多重插補雖然效果較好但計算量巨大,無法適用于大規(guī)模缺失處理.本文使用交替最小二乘下的低秩矩陣插補(Low Rank SVD via Alternating Least Square, softImpute- ALS)[25]算法,其是在低秩矩陣插補方法[26-27]基礎(chǔ)上為解決大規(guī)模數(shù)據(jù)集處理而提出的.低秩矩陣插補方法最初在Netflix競賽中被提出,用于插補電影評分矩陣中的缺失值,實現(xiàn)電影推薦的目的,因為其在大規(guī)模數(shù)據(jù)集上的優(yōu)異性能而備受關(guān)注.

基本思想是在原始存在缺失值的矩陣基礎(chǔ)上進行重建,即對矩陣中的缺失值進行填充,并保證填充后的完整矩陣與原始矩陣的秩相同.給定數(shù)據(jù)集,且中存在缺失值,由已觀測值的下標(biāo)構(gòu)成集合W,那么softImpute-ALS算法的目標(biāo)是求解以下優(yōu)化問題:

步驟1,初始化矩陣=,其中′r為隨機生成且列正交,′r為單位矩陣.初始化矩陣=,其中′r=0.

步驟2,固定矩陣,通過優(yōu)化式(2)來更新矩陣,

步驟3,固定矩陣,在第二步中交換矩陣和,使用同樣的方法更新矩陣和矩陣.

步驟4,重復(fù)步驟2和步驟3,直到達到收斂條件.

步驟5,計算=′,并對矩陣進行SVD分解,有=R.最終輸出結(jié)果、?和,=diag[(1-)+,…,(-)+],以及填充后的矩陣*.

2 數(shù)據(jù)介紹與分析

2.1 數(shù)據(jù)介紹

本文研究數(shù)據(jù)來自中國環(huán)境監(jiān)測總站實時發(fā)布的全國367個城市,1654個國控監(jiān)測點自2016年1月1日到2021年7月31日的6種空氣污染物小時濃度監(jiān)測數(shù)據(jù),共計48912h,理論上如果每個監(jiān)測點的每種污染物在每個小時都有監(jiān)測數(shù)據(jù)的話,那么監(jiān)測數(shù)據(jù)總量為4.2億條,但由于各種原因?qū)е聰?shù)據(jù)出現(xiàn)缺失(即空值),缺失數(shù)據(jù)總計814萬條.本文將缺失數(shù)據(jù)總量除以理論數(shù)據(jù)總量定義為數(shù)據(jù)缺失率,那么國控點監(jiān)測數(shù)據(jù)的總體缺失率為19.3‰.數(shù)據(jù)中不同污染物缺失率并不相同,SO2缺失率為9.6‰,NO2為10.9‰,CO為13.1‰,O3為14.4‰,PM2.5為15.3‰,PM10為52.4‰,PM10缺失率最高,而SO2缺失率最低.

2.2 監(jiān)測點數(shù)據(jù)缺失情況

不同監(jiān)測點的缺失情況存在很大差異(表1),以PM2.5為例,缺失最嚴(yán)重的監(jiān)測點的缺失率達到331.7‰,即三分之一的數(shù)據(jù)因各種原因缺失,致使監(jiān)測點數(shù)據(jù)有效性大打折扣,在使用該數(shù)據(jù)進行分析時,可能會導(dǎo)致較大誤差.因而有關(guān)部門在使用空氣質(zhì)量監(jiān)測數(shù)據(jù)時,會對監(jiān)測數(shù)據(jù)進行校驗,但當(dāng)前大部分研究中所使用的公開數(shù)據(jù)集往往是未經(jīng)過校驗處理的,其中就存在大量缺失值,缺失值對相關(guān)計算結(jié)果的影信息響究竟有多大,還沒有被深入討論.

表1 不同污染物的監(jiān)測點數(shù)據(jù)缺失統(tǒng)計(‰)

2.3 缺失機制分析

數(shù)據(jù)缺失的內(nèi)在機理可以劃分3種類型[28],分別是完全隨機缺失(MCAR)、隨機缺失(MAR)和非隨機缺失(MNAR).其中MAR相比于MCAR更加常見和符合現(xiàn)實,也是缺失值處理方法最主要的研究對象.MAR指僅在某個特定的組內(nèi)缺失值是隨機產(chǎn)生的,而不同組之間不一定是隨機的.對應(yīng)到空氣質(zhì)量數(shù)據(jù)的缺失機制,由于不同污染物缺失率不同,不同監(jiān)測點的缺失率也有很大差異,顯然污染物缺失并不是完全隨機的,而僅在特定監(jiān)測點和污染物條件下是可以視作是隨機缺失的.此外,空氣質(zhì)量是連續(xù)在線監(jiān)測模式,如果因設(shè)備故障而導(dǎo)致數(shù)據(jù)缺失的話,數(shù)據(jù)往往會連續(xù)缺失一定時段直到設(shè)備故障排除.表2顯示,三分之二左右的缺失間隔僅為1,說明大部分情況下缺失值是偶然出現(xiàn).但也有三分之一左右的缺失間隔大于1,缺失在10h以上的情況占比達到3%以上,這很可能是出現(xiàn)明顯的設(shè)備故障,需要一定的維修時間.其中PM2.5最大連續(xù)缺失200h,PM10最大連續(xù)缺失196h,SO2最大連續(xù)缺失192h,NO2最大連續(xù)缺失197h,O3最大連續(xù)缺失189h,CO最大連續(xù)缺失191h.研究表明隨著缺失間隔的增大,使用簡單缺失值處理方法(如均值填充等)的有效性不斷下降[29-30],更加需要有效方法來處理大規(guī)模缺失的情況,而softImpute-ALS方法在嚴(yán)重缺失的情況也能有很好的插補性能[25],因而適用于處理空氣質(zhì)量監(jiān)測中的缺失值.

表2 不同污染物缺失間隔占比(%)

3 模擬實驗

3.1 實驗設(shè)計

為考察softImpute-ALS方法在大規(guī)模缺失數(shù)據(jù)集上的優(yōu)劣,需要構(gòu)造相應(yīng)數(shù)據(jù)集,并計算各種缺失值處理方法在不同評估指標(biāo)上的表現(xiàn).因此,本文通過以下模擬方案生成算法評估數(shù)據(jù)集.方案包括四個步驟:(1)抽取可能產(chǎn)生缺失的污染物,抽取概率為不同污染物的缺失值數(shù)量占全部缺失值數(shù)量的比例;(2)抽取可能產(chǎn)生缺失的監(jiān)測點,抽取概率為污染物下,監(jiān)測點缺失值數(shù)量占全部該污染物全部缺失值數(shù)量的比例;(3)隨機抽取可能產(chǎn)生缺失的時間點,即該時刻開始出現(xiàn)缺失值;(4)抽取產(chǎn)生缺失的間隔長度,抽取概率為表2中缺失間隔對應(yīng)的占比;(5)記錄缺失值的下標(biāo)(,,,),并將原始數(shù)據(jù)中相同下標(biāo)的數(shù)據(jù)標(biāo)記為缺失值.通過以上四個步驟,在模擬真實缺失機制前提下,本文生成得到兩份數(shù)據(jù)集,一份為原始數(shù)據(jù)集,一份是標(biāo)記了缺失值的數(shù)據(jù)集.接下來,本文對數(shù)據(jù)集運用不同缺失值方法,通過考察其還原數(shù)據(jù)集的程度來評估方法的優(yōu)劣.

3.2 評估指標(biāo)

使用的評估指標(biāo)有均方根誤差(RMSE)、平均絕對百分比誤差(MAPE)和皮爾遜相關(guān)系數(shù)().RMSE指標(biāo)用于測度缺失值和插補值之間絕對偏差的大小,RMSE越小說明插補值約接近原始值.MAPE指標(biāo)用于測度缺失值和插補值之間相對誤差的大小,MAPE可以比較消除量綱后偏差的大小.指標(biāo)用于測度變量間線性相關(guān)關(guān)系,指標(biāo)值介于-1到1之間,越大說明缺失值與插補值之間相關(guān)性越高,二者擬合效果越好.不妨令模擬數(shù)據(jù)中觀測值下標(biāo)集合為W,則缺失值下標(biāo)構(gòu)成的集合是W的補集記為Wc,為原始的模擬數(shù)據(jù),為缺失值填充后的數(shù)據(jù),為缺失值數(shù)量,則評估指標(biāo)的具體計算公式為:

3.3 實驗結(jié)果

選擇簡單填充、線性插值、向前填充、最近鄰插補、EM插補、多重插補方法作為對比,比較這些方法和softImpute-ALS在相同數(shù)據(jù)集上,模擬真實缺失機制下在RMSE、MAPE、和運算時間(time)上的表現(xiàn).重復(fù)實驗10次,最終結(jié)果為10次實驗的平均值,實驗使用平臺配置為2.2GHz的CPU和16G DDR3內(nèi)存,實驗結(jié)果如表3所示(1)softImpute-ALS插補后的值在RMSE、MAPE兩個指標(biāo)上均表現(xiàn)最小,說明softImpute-ALS能夠更加準(zhǔn)確的插補缺失值,同時softImpute-ALS的插補值與真實值有更高的相關(guān)性,其他方法中只有最近鄰填充能夠與softImpute-ALS方法在準(zhǔn)確度方面更為接近.(2)時間消耗方面,softImpute-ALS對不同污染染處理時間介于160～180s之間,雖然相比于簡單填充、線性插值和向前填充耗時顯著增加,但相比最近鄰方法無疑有巨大的優(yōu)勢,最近鄰方法處理每種污染的用時在3h左右,是softImpute-ALS的70倍以上.實驗表明softImpute-ALS方法相比于其他方法有更好的準(zhǔn)確度,并且處理時間可控,因而softImpute-ALS方法在處理大規(guī)模缺失數(shù)據(jù)方面具有顯著的優(yōu)勢.

表3 不同缺失值插補方法比較

4 實證分析

4.1 數(shù)據(jù)有效性分析

“規(guī)定”中指出“SO2、NO2、PM10、PM2.5的評價濃度為評價時段內(nèi)日均濃度的平均值,O3的評價濃度為評價時段內(nèi)日最大8h平均值的第90百分?jǐn)?shù),CO的評價濃度為評價時段內(nèi)日均濃度的第95百分位數(shù)”的計算方法,同時指出城市的日均濃度首先要“計算各監(jiān)測點位日均濃度,然后計算算數(shù)平均值得到城市日均濃度,再由此計算你評價時段內(nèi)城市均值或特定百分位數(shù)”.各監(jiān)測點位的日均濃度由該日24h監(jiān)測值求平均而得,在“標(biāo)準(zhǔn)”中規(guī)定每日至少有20h的平均濃度值”才能計算24h均值,否則該日數(shù)據(jù)無效.

依據(jù)以上標(biāo)準(zhǔn),統(tǒng)計每個監(jiān)測點每日濃度并判斷其有效性,然后計算每個城市平均濃度及其有效性.研究發(fā)現(xiàn),從2016年1月1日至2021年7月31日共計2038d,監(jiān)測點有效天數(shù)最大是O-3的1892d,意味著該監(jiān)測點有多達146d未滿足有效性要求.同時,整體上有超過50%左右的監(jiān)測點數(shù)據(jù)有效時間介于1655～1840d之間,無效時間在198～383d之間(表4),可見由于存在大量缺失值,導(dǎo)致數(shù)據(jù)有效時間顯著減少,每日污染濃度的分析結(jié)果可能存在較大偏差.“規(guī)定”中并未指明在城市某個或某幾個監(jiān)測點當(dāng)日數(shù)據(jù)無效情況下,如何計算該城市當(dāng)日污染物平均濃度,而只是籠統(tǒng)的給出“當(dāng)任何一項污染物不滿足上述有效性規(guī)定且任何一項污染物濃度超過二級標(biāo)準(zhǔn)限值時,以城市當(dāng)日污染物濃度最高點位的數(shù)據(jù)”的規(guī)定,這樣以“最大值”替代的處理方式對于研究污染物時空分布顯然并不適用,更加合理的方式是利用模型算法對缺失值進行準(zhǔn)確的估計.

表4 全部國控監(jiān)測點不同污染物監(jiān)測數(shù)據(jù)有效時間統(tǒng)計

4.2 缺失值影響分析

由于“規(guī)定”中沒有明確給出監(jiān)測點當(dāng)日數(shù)據(jù)無效的處理方法,當(dāng)監(jiān)測點當(dāng)日數(shù)據(jù)無效時,分2種情況分別討論.第1種情況,如果一個城市當(dāng)日有一個或多個監(jiān)測點數(shù)據(jù)無效,則當(dāng)日該城市數(shù)據(jù)也無效;第2種情況,如果監(jiān)測點當(dāng)日數(shù)據(jù)無效,則取該監(jiān)測點當(dāng)日濃度最高值作為當(dāng)日有效數(shù)據(jù).在這2種情況下,缺失值均對后續(xù)當(dāng)日城市污染物濃度的計算均會產(chǎn)生重大影響.

在第1種情況下,以生態(tài)環(huán)境部發(fā)布的《2021年8月全國城市空氣質(zhì)量報告》[31]中空氣質(zhì)量相對較差的20個城市為例,分析這些城市由于數(shù)據(jù)缺失導(dǎo)致的無效數(shù)據(jù)問題.表5顯示,在過去5年多中,這20個城市各種污染物平均每年的無效時間約為60d,即有將近兩個月的時間按照“規(guī)定”是無效數(shù)據(jù),以濟南市PM10為例,無效時間更是達到881d,年均160d,無效時間的比例遠(yuǎn)超按小時統(tǒng)計的1%～5%的數(shù)據(jù)缺失率,其主要原因是每個監(jiān)測點每日須滿足大于等于20個有效監(jiān)測數(shù)據(jù).即使按照規(guī)定中日最大濃度替代日均值濃度的方法,也意味著這些城市每年有將近2個月的時間使用日最高濃度進行排名計算.

表5 空氣質(zhì)量排名后20城市的數(shù)據(jù)無效時間統(tǒng)計

在第2種情況下,污染物濃度在一天內(nèi)存在很大的波動,造成最大濃度與平均濃度之間有很大差異,最大濃度無法反映真實情況.以濟南市機床二廠監(jiān)測點2021年2月3日為例,當(dāng)日有效監(jiān)測數(shù)據(jù)恰為19個,小于20個有效監(jiān)測數(shù)據(jù)的要求,缺失14:00、16:00、17:00、18:00和19:00這5個小時監(jiān)測數(shù)據(jù).圖1顯示一天內(nèi)不同時點污染物濃度存在很大差異,該監(jiān)測點PM2.5最大濃度123mg/m3,19個小時平均值為81mg/m3,相差52%,PM10最大濃度214mg/m3,19小時平均值153mg/m3,相差40%,SO2最大濃度55mg/m3,19小時平均值20mg/m3,相差172%,NO2最大濃度69mg/m3,19小時平均值40mg/m3,相差75%, O3最大濃度87mg/m3,19小時平均值59mg/m3,相差47%,CO最大濃度1.3mg/m3,19小時平均濃度0.9mg/m3,相差41%.“規(guī)定”中的“最大濃度替代”與實際情況存在很大差異,不能反應(yīng)污染物的真實濃度,僅僅是在計算層面對缺失數(shù)據(jù)進行懲罰性處理.當(dāng)存在某些不可抗拒因素導(dǎo)致數(shù)據(jù)缺失時,這樣的處理方法并不合理,而應(yīng)當(dāng)采取更科學(xué)和可靠的插補方法,最大程度近似或還原真實值.

4.3 插補值有效性分析

應(yīng)用softImpute-ALS方法對全國國控點空氣質(zhì)量數(shù)據(jù)進行缺失值的插補處理,未驗證插補值的有效性,為與上文保持一致,選取濟南市機床二廠和科干所兩個監(jiān)測點來驗證插補值的有效性.具體方法是分別計算兩個監(jiān)測點每種污染物在剔除缺失值后的數(shù)據(jù)、插補值數(shù)據(jù)、插補缺失值后的數(shù)據(jù)這3個數(shù)據(jù)集上的相關(guān)系數(shù),并且繪制散點圖(圖2)來分析插補值的有效性與合理性,即如果2個監(jiān)測點在3個數(shù)據(jù)集的相關(guān)系數(shù)比較接近,則說明插補值也符合原有數(shù)據(jù)規(guī)律,插補方法是有效的.

由圖2可見,插補值與剔除缺失的數(shù)據(jù)在污染物相關(guān)關(guān)系上是一致的,插補前兩個監(jiān)測點PM2.5的相關(guān)系數(shù)是0.93,插補后的相關(guān)系數(shù)是0.94,缺失值處理方法得到的插補值也服從監(jiān)測點間相關(guān)的規(guī)律.此外計算濟南市其他任意2個監(jiān)測點之間在3個數(shù)據(jù)集上的相關(guān)系數(shù),均得到大于0.9的相關(guān)系數(shù).

4.4 實證結(jié)果分析

應(yīng)用softImpute-ALS方法對全國國控點空氣質(zhì)量數(shù)據(jù)進行缺失值的插補處理,并依據(jù)“規(guī)定”中城市濃度計算標(biāo)準(zhǔn),分析缺失值插補后的城市濃度值和排名的變化情況.以2021年7月為例,分別對168個城市分別計算PM2.5,PM10,SO2,NO2,CO-95%分位數(shù),O38h-90%分位數(shù)的濃度值及其排名.生態(tài)環(huán)境部在發(fā)布城市空氣質(zhì)量排名時會對監(jiān)測數(shù)據(jù)進行校驗,校驗后數(shù)據(jù)質(zhì)量要優(yōu)于實時監(jiān)測數(shù)據(jù),但校驗后數(shù)據(jù)并未公開發(fā)布,本文實證分析部分基于中國環(huán)境監(jiān)測總站實時在線監(jiān)測數(shù)據(jù),僅用于分析缺失值處理方法對計算結(jié)果的影響,說明在使用公開空氣質(zhì)量監(jiān)測數(shù)據(jù)時對缺失值處理的必要性.

4.4.1 全國城市PM2.5濃度缺失值插補前后計算結(jié)果 20個城市中有18個相同的城市,說明缺失值插補前后排名總體上一致,但局部排名有明顯的變化.主要原因是2021年7月所有城市PM2.5月度均值總體較低,并且不同城市的監(jiān)測值相差非常小,即計算值的微小變化會導(dǎo)致排名的較大變化.缺失高濃度監(jiān)測值會顯著提升城市排名,而缺失低濃度監(jiān)測值會顯著降低城市排名,此時缺失值的處理方法會對城市排名產(chǎn)生關(guān)鍵影響.

表6顯示,不同城市缺失率顯著不同,最高缺失率5.2%,最低缺失率0.8%,而“規(guī)定”中并沒有明確監(jiān)測點小時監(jiān)測數(shù)據(jù)缺失處理方法,只是給出“以城市污染物濃度最高點位的數(shù)據(jù),統(tǒng)計當(dāng)日污染物濃度排名”,這樣的處理方式是不合理的,也沒有考慮污染物在一天的周期性變化特征.使用合理插補方法后,PM2.5月度濃度均值變化在-7.8%～ 10.1%之間,基于公開在線監(jiān)測數(shù)據(jù)進行研究時,應(yīng)當(dāng)采取合理、可靠的缺失值處理方法對缺失值進行插補,減少缺失值對最終計算結(jié)果的影響.

表6 缺失值插補前后城市PM2.5排名比較

注:(1)原排序中相同序號按照城市出現(xiàn)前后的倒序排列;(2)插補前后變動計算公式為:(插補后-插補前)/插補前.

4.4.2 全國城市PM10濃度缺失值插補前后計算結(jié)果 排名后20城市中有16個相同城市,PM10缺失率從0.9%～6.0%之間,排名變化較大的是石家莊市,PM10的月度均值從47mg/m3下降至44mg/m3,從原排名163名上升至152名,上升11名,而運城市的PM10的月度均值從42mg/m3上升至45mg/m3,排名也從150名下降至157名.

表7 缺失值插補前后城市PM10排名比較

注: (1)原排序中相同序號按照城市出現(xiàn)前后的倒序排名;(2)插補前后變動計算公式為:(插補后-插補前)/插補前

4.4.3 全國城市SO2濃度缺失值插補前后計算結(jié)果 排名的后20城市全部相同,SO2缺失率從1.4%～ 5.5%之間.由于SO2月度均值差異較大,插補后排名與原排名的最后7名次序均相同,但在之后排名次序顯著不同,即當(dāng)污染物月度均值有較大差異時,缺失值插補對排名影響較小,但如果污染物月度均值比較接近時,缺失值插補會顯著影響排名情況.

4.4.4 全國城市NO2濃度缺失值插補前后計算結(jié)果 排名后20城市中有17個相同城市,NO2缺失率為1.2%～5.4%之間.與SO2排名類似,NO2月度均值差異較大,插補后排名與原排名的最后5名次序均相同,其他城市的排名情況均有所變化,但排名上升和下降幅度較小.

表8 缺失值插補前后城市SO2排名比較

注:(1)原排序中相同序號按照城市出現(xiàn)前后的倒序排名;(2)插補前后變動計算公式為:(插補后-插補前)/插補前

表9 缺失值插補前后城市NO2排名比較

注:(1)原排序中相同序號按照城市出現(xiàn)前后的倒序排名;(2)插補前后變動計算公式為:(插補后-插補前)/插補前.

表10 缺失值插補前后城市CO日均值95%分位數(shù)排名比較

注:(1)原排序中相同序號按照城市出現(xiàn)前后的倒序排名;(2)插補前后變動計算公式為:(插補后-插補前)/插補前.

表11 缺失值插補前后城市O3日最大8h90%分位數(shù)排名比較

注:(1)原排序中相同序號按照城市出現(xiàn)前后的倒序排名;(2)插補前后變動計算公式為:(插補后-插補前)/插補前.

4.4.5 全國城市CO日均值95%分位數(shù)和O3日最大8h 90%分位數(shù) 通過CO日均值95%分位數(shù)(表10)與O3日最大8h 90%分位數(shù)(表11)的排名比較,在數(shù)值較為接近時,城市排名會受到缺失值的顯著影響,缺失值如何處理將直接影響城市排名前后次序.新余、承德、唐山、運城和忻州這5個城市CO日均值在插補前后變動均超過10%,其中運城和忻州的分別從插值前的157名和151名,在插值后變?yōu)?38名和121名,提升20和30名,可見缺失值對城市評價濃度具有顯著影響.

5 結(jié)論

研究表明數(shù)據(jù)缺失會對計算結(jié)果造成顯著影響程度,尤其是進行日尺度和城市尺度的統(tǒng)計分析,本研究對于提升監(jiān)測數(shù)據(jù)完整率和分析結(jié)果準(zhǔn)確率方面具有一定的現(xiàn)實意義.基于softImpute-ALS方法,對2016年以來的全國國控點小時級的大規(guī)模監(jiān)測數(shù)據(jù)進行缺失值插補處理,模擬實驗表明本文所使用方法能夠得到較好的估計值,RMSE、MAPE和相關(guān)系數(shù)均優(yōu)于其他缺失值常規(guī)處理方法,同時缺失值估計過程具有較快的處理速度.實證分析結(jié)果表明,缺失值插補前后6種污染物評價濃度會有±10%左右的變化,同時監(jiān)測點的缺失值數(shù)量越多,缺失值插補前后濃度的計算值變化也越大,越有必要進行缺失值處理.建議在基于公開的大規(guī)模空氣質(zhì)量監(jiān)測數(shù)據(jù)進行空氣質(zhì)量研究時,首先應(yīng)采用softImpute-ALS方法對監(jiān)測數(shù)據(jù)進行缺失值處理,得到更加準(zhǔn)確的估計值,提升監(jiān)測數(shù)據(jù)完整率,最大程度減少缺失值可能導(dǎo)致的結(jié)果有偏問題,提升相關(guān)研究結(jié)果的準(zhǔn)確性和可靠性.

[1] 中華人民共和國生態(tài)環(huán)境部.城市環(huán)境空氣質(zhì)量排名技術(shù)規(guī)定[R]. 2018.

Ministry of Ecology and Environmental of People's Republic of China. Technical regulations for air quality ranking of cities[R]. 2018.

[2] Deng Q, Yang K, Luo Y. Spatiotemporal patterns of PM2.5in the Beijing–Tianjin–Hebei region during 2013～2016 [J]. Geology, Ecology, and Landscapes, 2017,1(2):95-103.

[3] Li L, Wu A H, Cheng I, et al. Spatiotemporal estimation of historical PM2.5concentrations using PM10, meteorological variables, and spatial effect [J]. Atmospheric Environment, 2017,166:182-191.

[4] Hu M, Wang Y, Wang S, et al. Spatial-temporal heterogeneity of air pollution and its relationship with meteorological factors in the Pearl River Delta, China [J]. Atmospheric Environment, 2021,254:118415.

[5] Li L, Zhang J, Meng X, et al. Estimation of PM2.5concentrations at a high spatiotemporal resolution using constrained mixed-effect bagging models with MAIAC aerosol optical depth [J]. Remote Sensing of Environment, 2018,217:573-586.

[6] Shen Y, Zhang l, Fang X, et al. Spatiotemporal patterns of recent PM2.5concentrations over typical urban agglomerations in China [J]. Science of the Total Environment, 2019,655:13-26.

[7] Zhao S, Yin D, Yu Y, et al. PM2.5and O3pollution during 2015～2019 over 367 Chinese cities: Spatiotemporal variations, meteorological and topographical impacts [J]. Environmental Pollution, 2020,264:114694.

[8] Li K, Jacob D J, Liao H, et al. A two-pollutant strategy for improving ozone and particulate air quality in China [J]. Nature Geoscience, 2019, 12(11):906-910.

[9] Liu J, Li W, Wu J. A framework for delineating the regional boundaries of PM2.5pollution: A case study of China [J]. Environmental Pollution, 2018,235:642-651.

[10] 張 烴,董樹屏,滕 曼,等.區(qū)域大型環(huán)境空氣綜合觀測中外場觀測與實驗室分析數(shù)據(jù)質(zhì)量控制研究[J]. 環(huán)境科學(xué)研究, 2019,32(10): 1664-1671.

Zhang T, Dong S P, Teng M, et al. Quality assurance of field observation and laboratory analysis in regional large scale ambient air joint monitoring campaigns [J]. Research of Environmental Sciences, 2019,32(10):1664-1671.

[11] 師耀龍,呂怡兵,肖建軍.夏季重大活動期間O3監(jiān)測數(shù)據(jù)質(zhì)量提升方法研究[J]. 中國環(huán)境監(jiān)測, 2020,36(2):10-14.

Shi Y L, Lyu Y B, Xiao J J. Data quality control method of ozone monitoring during the guarantee for major events in summer [J]. Environmental Monitoring in China, 2020,36(2):10-14.

[12] 師耀龍,陳傳忠,魏俊山,等.加強生態(tài)環(huán)境監(jiān)測機構(gòu)監(jiān)督管理的思考與分析[J]. 環(huán)境保護, 2018,46(23):56-60.

Shi Yao-long, Chen Chuan-zhong, Wei Jun-shan, et al. The current situation and problem analysis of environmental monitoring organizations' supervision and administration [J]. Environmental Protection, 2018,46(23):56-60.

[13] 劉 媛,彭 溶,張 馳,等.環(huán)境監(jiān)測從業(yè)人員監(jiān)管制度研究[J]. 環(huán)境保護, 2018,46(18):33-35.

Liu Y, Peng R, Zhang C, et al. Research on supervision system of environmental monitoring practitioners [J]. Environmental Protection, 2018,46(18):33-35.

[14] 師耀龍,楊 婧,柴文軒,等.美國環(huán)境空氣監(jiān)測數(shù)據(jù)質(zhì)量核查工作的經(jīng)驗與啟示[J]. 中國環(huán)境監(jiān)測, 2017,33(3):8-14.

Shi Y L, Yang J, Chai W X, et al. Experience and illumination of data quality assessment system for ambient air monitoring in the United States [J]. Environmental Monitoring in China, 2017,33(3):8-14.

[15] Rumaling M I, Chee F pien, Dayou J, et al. Missing value imputation for PM10concentration in Sabah using nearest neighbour method (NNM) and expectation-maximization (EM) algorithm [J]. Asian Journal of Atmospheric Environment, 2020,14:62-72.

[16] Junger W L, Ponce D E Leon A. Imputation of missing data in time series for air pollutants [J]. Atmospheric Environment, 2015,102:96- 104.

[17] Larsen L C, Shah M. A context-intensive approach to imputation of missing values in data sets from networks of environmental monitors [J]. Journal of the Air & Waste Management Association (1995), 2016,66(1):38-52.

[18] Junninen H, Niska H, Tuppurainen K, et al. Methods for imputation of missing values in air quality data sets [J]. Atmospheric Environment, 2004,38(18):2895-2907.

[19] Hadeed S J, O’rourke M K, Burgess J L, et al. Imputation methods for addressing missing data in short-term monitoring of air pollutants [J]. Science of The Total Environment, 2020,730:139140.

[20] Real C, ángel Fernández J, Aboal J R, et al. Substituting missing data in compositional analysis [J]. Environmental Pollution, 2011,159(10): 2797-2800.

[21] Gómez-Carracedo M P, Andrade J M, López-Mahía P, et al. A practical comparison of single and multiple imputation methods to handle complex missing data in air quality datasets [J]. Chemometrics and Intelligent Laboratory Systems, 2014,134:23-33.

[22] Chen X, Xiao Y. A novel method for air quality data imputation by nuclear norm minimization [J]. Journal of Sensors, 2018,2018: e7465026.

[23] Moshenberg S, Lerner U, Fishbain B. Spectral methods for imputation of missing air quality data [J]. Environmental Systems Research, 2015,4(1):26.

[24] Liu X, Wang X, Zou L, et al. Spatial imputation for air pollutants data sets via low rank matrix completion algorithm [J]. Environment International, 2020,139:105713.

[25] Hastie T, Mazumder R, Lee J D, et al. Matrix completion and low-rank SVD via fast alternating least squares [J]. 36.

[26] Candès E J, Recht B. Exact matrix completion via convex optimization [J]. Foundations of Computational Mathematics, 2009, 9(6):717.

[27] Candès E J, Tao T. The power of convex relaxation: near-optimal matrix completion [J]. IEEE Transactions on Information Theory, 2010,56(5):2053-2080.

[28] Buuren S van. Flexible imputation of missing data [M]. 2nd edition. Boca Raton: Chapman and Hall/CRC, 2018.

[29] Liu Y, Dillon T, Yu W, et al. Missing value imputation for industrial IoT sensor data with large gaps [J]. IEEE Internet of Things Journal, 2020,7(8):6855-6867.

[30] Velasco-Gallego C, Lazakis I. A novel framework for imputing large gaps of missing values from time series sensor data of marine machinery systems [J]. Ships and Offshore Structures, 2021,10.1080/ 17445302.2021.1943850.

[31] 中華人民共和國生態(tài)環(huán)境部.2021年8月全國城市空氣質(zhì)量報告[R]. 2021.

Ministry of Ecology and Environmental of People's Republic of China. Air quality reports for cities in China[R]. 2021,8

致謝：本研究受中國人民大學(xué)“雙一流”跨學(xué)科重大創(chuàng)新規(guī)劃平臺—生態(tài)文明跨學(xué)科交叉平臺支持.

Research on the missing value methods for large-scale online air quality monitoring data.

ZHANG Bo, SONG Guo-jun*

(School of Environment and Natural Resources, Renmin University of China, Beijing 100872, China)., 2022,42(5)：2078～2087

Large scale online air quality monitoring data is the basis for air quality research, but there were lots of missing data in large scale online data. In this study, we compared several methods that dealing with the missing values and its impact on the city’s ranking of air quality base on the hourly monitoring data of 1654monitoring sites in China from 1Jan, 2016 to 31July, 2021 of 6types of air pollutants. The simulation results showed that Low Rank SVD via Alternating Least Square had smaller mean squared error, mean absolute percentage error and higher correlation coefficient compared with other traditional methods. The empirical results showed there would be 10% difference before imputation and after imputation for the missing value. The ranking would not change due to the imputation when the air quality assessed value vary greatly, and would change a lot when the assessed value was very close. The study suggested to impute missing value by using the method in this study when analysis the large-scale online air quality monitoring data.

monitoring data；big data；missing value；low rank matrix；empirical research

X323

A

1000-6923(2022)05-2078-10

張 波(1986-),男,內(nèi)蒙古包頭人,講師,博士,研究方向為環(huán)境統(tǒng)計與建模.發(fā)表論文10余篇.

2021-10-12

中國人民大學(xué)科學(xué)研究基金(中央高?；究蒲袠I(yè)務(wù)費專項資金資助)項目成果(22XNF016)

* 責(zé)任作者, 教授, songguojun@ruc.edu.cn