基于INLA-SPDE方法的區(qū)域污染物模擬與預(yù)測

袁 澤，陳 斌

基于INLA-SPDE方法的區(qū)域污染物模擬與預(yù)測

袁 澤，陳 斌

(北京大學(xué)地球與空間科學(xué)學(xué)院，北京 100871)

采用傳統(tǒng)的空間插值方法對區(qū)域污染物進(jìn)行模擬與預(yù)測，針對源數(shù)據(jù)分布不均，效果一般的問題，提出了采用INLA-SPDE模型來模擬與預(yù)測區(qū)域污染物的方法。模型的空間分量使用隨機(jī)偏微分方程表達(dá)，時間分量則采用一階時序自相關(guān)模型，同時還包含氣象參數(shù)等10種協(xié)變量，以2019年度京津冀地區(qū)日均PM2.5濃度為例，逐月建立了時空模擬與預(yù)測模型。實(shí)驗(yàn)結(jié)果表明，與經(jīng)典的克里金插值方法相比，在區(qū)域污染物分布的模擬上具有更好的效果，尤其在高值污染的預(yù)測上精度效果提升明顯，同時可得到區(qū)域污染風(fēng)險等級等多種結(jié)果。進(jìn)一步基于模型的預(yù)測結(jié)果實(shí)現(xiàn)了京津冀地區(qū)日均PM2.5濃度時空可視化和虛擬仿真系統(tǒng)，為普通民眾的出行或政府相關(guān)部門決策提供支持，驗(yàn)證了模型的實(shí)用性和價值。

細(xì)顆粒物PM2.5；貝葉斯時空建模；INLA算法；仿真系統(tǒng)；決策支持

空氣污染一直是一個全球性問題[1]，隨著我國高人口密度的增加和快速工業(yè)化的結(jié)合，與經(jīng)濟(jì)增長、工業(yè)化和城市化相關(guān)的嚴(yán)重空氣污染事件頻繁發(fā)生。當(dāng)下大氣顆粒物(particulate matter，PM)，即懸浮在空氣當(dāng)中的固體顆粒或液滴，是最主要的污染物之一。而細(xì)顆粒物(空氣動力學(xué)直徑小于等于2.5 μm的顆粒物稱為PM2.5)因其粒徑更為微小，且具有最為明顯的污染效應(yīng)。同時許多流行病學(xué)和毒理學(xué)研究表明，PM2.5暴露會導(dǎo)致許多不良健康影響[2]。基于此，政府已經(jīng)將污染監(jiān)測和預(yù)警工作提到了一個新的高度。

目前我國已經(jīng)建立了一定規(guī)模的地面空氣質(zhì)量監(jiān)測站點(diǎn)[3]，其監(jiān)測結(jié)果可以有效地反映監(jiān)測站周邊的空氣質(zhì)量狀況。但受限于建設(shè)成本以及空間占地等問題，地面監(jiān)測站數(shù)量還比較有限，且多位于發(fā)達(dá)大城市周邊，無法覆蓋大陸全域。同時近年來隨著監(jiān)測網(wǎng)絡(luò)的完善，產(chǎn)生了大量的時空數(shù)據(jù)，這些數(shù)據(jù)大多具有較高的時間維度分辨率，一般為小時級甚至分鐘級，而在空間維度上則相對稀疏，空間數(shù)據(jù)大多會隨著時間發(fā)生變化，對聯(lián)合時間與空間的數(shù)據(jù)模型的研究日益增多。

在多種時空模型中，基于貝葉斯框架的時空建模方法近年來引起了越來越多學(xué)者的關(guān)注，其原因是時空數(shù)據(jù)具有內(nèi)在的不確定性，而貝葉斯方法從概率的角度出發(fā)，能夠?qū)r空數(shù)據(jù)的不確定性進(jìn)行建模，最終給出概率上最優(yōu)的結(jié)果[4]。在研究空氣質(zhì)量時空分布模型時，SARTO等[5]使用分層貝葉斯模型對意大利中部佩魯賈鎮(zhèn)的PM2.5濃度分布進(jìn)行了模擬，使用溫度、濕度等氣象因子作為協(xié)變量，同時添加空間效應(yīng)和時間自回歸項(xiàng)，獲得了較好的擬合效果。

目前對于貝葉斯時空模型參數(shù)的后驗(yàn)分布求解，大都采用馬爾可夫鏈蒙特卡洛方法(Markov chain Monto Carlo，MCMC)，但由于到達(dá)馬爾可夫鏈的平穩(wěn)態(tài)往往需要很長時間，在實(shí)際應(yīng)用中，MCMC方法并非作為一種常規(guī)分析工具。

2009年，RUE等[6]針對一類特殊的時空隱高斯隨機(jī)場的貝葉斯推斷，提出了一種稱為集成嵌套拉普拉斯(integrated nested Laplace approximations，INLA)的快速求解算法。該方法可以求解一類常見的結(jié)構(gòu)化加法回歸模型，這類模型符合被一些超參數(shù)約束的隱高斯場，同時具有非高斯響應(yīng)變量。通過使用集成的INLA逼近及其簡化版本，該算法可以直接計算出后驗(yàn)分布非常精確的逼近，同時在MCMC算法需要運(yùn)行數(shù)小時或數(shù)天的情況下，該算法可在數(shù)秒或數(shù)分鐘內(nèi)提供更精確的估計。2011年，LINDGREN等[7]提出了一種隨機(jī)偏微分方程(stochastic partial differential equations，SPDE)算法，其可有效地鏈接連續(xù)高斯隨機(jī)場和離散的高斯馬爾科夫隨機(jī)場(Gaussian Markov random fields，GMRF)，從而實(shí)現(xiàn)對基于Matern相關(guān)系數(shù)的連續(xù)分布的空間效應(yīng)的快速解算，以進(jìn)一步加速INLA算法對貝葉斯時空模型后驗(yàn)分布的逼近。

目前，該方法在國外已經(jīng)開始普及，在與空氣質(zhì)量有關(guān)的研究方面，CAMELETTI等[8]針對意大利皮埃蒙特地區(qū)的健康狀況，采用INLA-SPDE方法獲得了區(qū)域污染物濃度，并開展了對污染物濃度分布和區(qū)域健康風(fēng)險關(guān)系的研究。FIORAVANTI等[9]基于意大利全國410個監(jiān)測站點(diǎn)2015年的日均PM10濃度數(shù)據(jù)，使用INLA算法得到了意大利全域時空貝葉斯插值結(jié)果。但國內(nèi)，該方法使用較少，幾乎只用于流行病時空分布的分析與模擬[10]。

我國目前的空氣質(zhì)量地面監(jiān)測站主要分布在城市中心地帶，尚未形成空間分布均衡的監(jiān)測網(wǎng)絡(luò)，因此難以對大陸區(qū)域空氣質(zhì)量狀況的時空分布進(jìn)行合理地評估。因此本文以京津冀地區(qū)為例，綜合利用研究區(qū)域內(nèi)的氣象數(shù)據(jù)和地理要素數(shù)據(jù)，建立了區(qū)域PM2.5濃度時空分布貝葉斯模型，并采取文獻(xiàn)[6-7]提出的INLA-SPDE加速算法，實(shí)現(xiàn)了對區(qū)域日均PM2.5濃度的快速模擬與預(yù)測。并以該預(yù)測結(jié)果為基礎(chǔ)，構(gòu)建區(qū)域PM2.5濃度虛擬仿真決策支持系統(tǒng)，實(shí)現(xiàn)京津冀區(qū)域PM2.5濃度時空分布圖、污染風(fēng)險等級圖的可視化，支持與行政區(qū)劃圖、高程圖、人口圖等進(jìn)行疊加顯示，從而驗(yàn)證模型的結(jié)果實(shí)用性。該系統(tǒng)還支持交互放置特定的高污染源于地圖上，快速生成新情景下的空氣質(zhì)量空間分布圖，對比觀察污染物分布的變化，為污染物減排及控制等相關(guān)的空間問題提供決策支持。

1 原始數(shù)據(jù)處理與分析

本文使用的原始PM2.5濃度數(shù)據(jù)來自京津冀地區(qū)236個國控空氣質(zhì)量監(jiān)測站，數(shù)據(jù)的原始時間分辨率為1 h，缺失值由站點(diǎn)當(dāng)前時刻前后3 h污染物濃度數(shù)據(jù)取均值代替。同時為了保證數(shù)據(jù)分布的正態(tài)性，去除了大于650的極端污染值數(shù)據(jù)，日均濃度數(shù)據(jù)則由當(dāng)日24 h污染物濃度數(shù)據(jù)求平均值獲得，時間跨度從2019年1月到12月。最終得到月度PM2.5濃度日均值箱線圖(圖1)，圖中的橘色虛線和紅色虛線分別對應(yīng)的PM2.5濃度為75 μg/m3和150 μg/m3，這是WHO標(biāo)準(zhǔn)中的第三級和第一級日均PM2.5濃度標(biāo)準(zhǔn)，也是我國分級標(biāo)準(zhǔn)中的污染和重度污染分界線。

圖1 京津冀地區(qū)2019年月度PM2.5濃度日均值箱線圖

在2019年中，PM2.5濃度日均值在0～626 μg/m3之間，年均值為49.9 μg/m3，中位數(shù)為37 μg/m3，上下四分位數(shù)分別為22.2 μg/m3和60.8 μg/m3。該箱線圖與第2節(jié)研究一致，呈現(xiàn)出季節(jié)性規(guī)律：12月、1月、2月的PM2.5濃度日均值中位數(shù)大于50 μg/m3，夏季則小于30 μg/m3，其他月份在40 μg/m3。標(biāo)準(zhǔn)偏差也具有類似的規(guī)律，1月和2月的標(biāo)準(zhǔn)偏差大于60 μg/m3，夏季小于20 μg/m3，其他月份在35 μg/m3左右。

參與構(gòu)建貝葉斯模型的協(xié)變量在經(jīng)過反復(fù)試驗(yàn)，最終選定了10個時空變量作為模型的預(yù)報因子。具體的數(shù)據(jù)信息見表1。為了避免數(shù)據(jù)自身的量級對研究的影響，所有的數(shù)據(jù)均進(jìn)行了標(biāo)準(zhǔn)化處理，即全部轉(zhuǎn)換為了均值為0，標(biāo)準(zhǔn)差為1的正態(tài)分布。

表1 貝葉斯時空模型使用的協(xié)變量

使用的6種氣象因子全部來自由美國國家環(huán)境中心(National Centers for Environment-Mental Prediction，NCEP)提供的6 h氣象再分析資料。該數(shù)據(jù)提供了每日0點(diǎn)、6點(diǎn)、12點(diǎn)、18點(diǎn)的空間分辨率為0.25度的格網(wǎng)化氣象要素，下載后以日為尺度重新計算了均值。氣溶膠光學(xué)厚度(aerosol optical depth，AOD)數(shù)據(jù)則來自歐空局(Copernicus Atmosphere Monitoring Service，CAMS)氣象再分析資料，是在空間分辨率為10 km的格網(wǎng)上對波長為550 nm的AOD數(shù)值模擬得到的結(jié)果。

京津冀地區(qū)各級道路數(shù)據(jù)獲取來自O(shè)SM (OpenStreetMap Project)，從Geofabrik網(wǎng)站提供的服務(wù)下載得到。本文選取了OSM數(shù)據(jù)集中標(biāo)記為highway的線狀地物中的前四級標(biāo)簽作為基礎(chǔ)道路數(shù)據(jù)，包括了京津冀地區(qū)的國道、省道以及縣道。最終投入模型的路網(wǎng)密度數(shù)據(jù)則是通過將京津冀地區(qū)劃分為10 km×10 km的格網(wǎng)，計算每個格網(wǎng)內(nèi)道路數(shù)據(jù)的總長度(km)來獲得的。海拔數(shù)據(jù)則來自USGS提供的GTOPO30全球數(shù)字高程模型，分辨率為30弧秒(約1 km)。人口數(shù)據(jù)采用WorldPop datasets提供的中國大陸1 km分辨率的人口密度數(shù)據(jù)。

2 京津冀地區(qū)PM2.5日均濃度時空貝葉斯建模

2.1 貝葉斯時空建模方案設(shè)計

對2019年全年的京津冀地區(qū)日均PM2.5濃度進(jìn)行建模，按月份建立12個模型，具體模型如下：

圖2 區(qū)域SPDE空間效應(yīng)三角網(wǎng)示意圖

2.2 預(yù)測結(jié)果精度評估

首先檢驗(yàn)?zāi)Ｐ偷恼军c(diǎn)預(yù)報精度，表2為12個月度模型在訓(xùn)練集和測試集上3次交叉驗(yàn)證后的平均結(jié)果，這些指標(biāo)數(shù)值均在轉(zhuǎn)換為原始數(shù)據(jù)尺度后計算獲得。總體來說，12個模型在訓(xùn)練集和測試集上的預(yù)測均表現(xiàn)良好，而且效果相當(dāng)，說明模型并不存在過擬合。RMSE和MAE的值在冬季要高于夏季，這與冬季PM2.5濃度值整體偏高有關(guān)。

表2 模型交叉驗(yàn)證結(jié)果(μg/m3)

觀測值和預(yù)測值的關(guān)聯(lián)性如圖3所示，為了避免過多的散點(diǎn)圖，本文將12個月度模型按季度展示其真實(shí)值與預(yù)測值的一致性。圖中顏色越淺的部分代表處于此濃度的觀測數(shù)據(jù)的數(shù)量越多，紅色實(shí)線對應(yīng)1﹕1參考線。無論是在訓(xùn)練集還是測試集上，點(diǎn)幾乎都均勻地分布于參考線兩側(cè)，同時淺色區(qū)域幾乎圍繞著參考線兩側(cè)最近處分布，可以認(rèn)為模型在4個季度對于大部分樣本數(shù)據(jù)預(yù)測效果表現(xiàn)良好。不過從圖3中也可以看到，在測試集上，當(dāng)觀測值大于250 μg/m3時，出現(xiàn)了較多偏小的預(yù)測值。

圖3 模型預(yù)測值與真實(shí)值對比圖((a)訓(xùn)練集；(b)測試集)

圖4以污染最嚴(yán)重的1月和相對最輕的8月為例，繪制了區(qū)域內(nèi)236個站點(diǎn)月度的日均觀測值與預(yù)測置信區(qū)間的時序關(guān)聯(lián)圖，圖中黃色區(qū)域代表預(yù)測濃度值的95%置信區(qū)間。該時序圖表明京津冀地區(qū)PM2.5濃度日均值的時間變異性可以很好地被模型捕捉，觀測值濃度越高，預(yù)測值的置信區(qū)間范圍越大。但真實(shí)值總是落在預(yù)測的置信區(qū)間內(nèi)部，證明了模型對于時序特征的復(fù)現(xiàn)能力。

2.3 預(yù)測結(jié)果對比分析與可視化

本文選取具有代表性的日期，將模型的預(yù)測結(jié)果與傳統(tǒng)的空間插值通用方法——克里金插值方法進(jìn)行比較。即在相同的訓(xùn)練集上建立INLA-SPDE模型和進(jìn)行普通克里金空間插值，最終將測試集上的插值結(jié)果與INLA-SPDE模型在測試集上的預(yù)報結(jié)果進(jìn)行對比，其中貝葉斯模型的預(yù)測結(jié)果均采用模型的預(yù)測后驗(yàn)均值數(shù)據(jù)。

圖4 PM2.5日均濃度時序驗(yàn)證圖((a)1月；(b)8月)

具體選取的實(shí)驗(yàn)數(shù)據(jù)從2019年1月7日—18日，這一時段內(nèi)京津冀地區(qū)經(jīng)歷了一次污染從嚴(yán)重到消散的過程，具有一定的代表性。最終在完全相同的測試集上試驗(yàn)得到的RMSE數(shù)值見表3，可以發(fā)現(xiàn)使用的INLA-SPDE模型在空間上的預(yù)測效果要明顯優(yōu)于傳統(tǒng)的克里金插值。根據(jù)表3的數(shù)據(jù)繪制得到圖5，可以發(fā)現(xiàn)2種模型預(yù)測精度的差距在污染濃度越大時越明顯，證明了本文模型在高污染值上的預(yù)測能力遠(yuǎn)強(qiáng)于傳統(tǒng)的克里金插值方法。從PM2.5濃度預(yù)測的實(shí)際應(yīng)用出發(fā)，對捕捉嚴(yán)重污染出現(xiàn)時機(jī)和位置的需求更大，因此本文模型具有較強(qiáng)的實(shí)用價值與現(xiàn)實(shí)意義。

圖6以京津冀地區(qū)5 km×5 km的格網(wǎng)為基礎(chǔ)，展示了克里金插值方法與本文模型的對比結(jié)果。選取了1月13日和15日分別為污染嚴(yán)重和開始消散的代表日期，值得注意的是左右2幅圖顏色分級尺度不同。可以發(fā)現(xiàn)，本文模型對于空間效應(yīng)的模擬與克里金插值相似，大尺度上的空間分布兩者保持著高度一致性。該模型對于污染濃度偏高的值更為敏感，能更為準(zhǔn)確地反映嚴(yán)重污染地帶的分布，2幅圖中極值亮點(diǎn)均出現(xiàn)在左圖，相比與傳統(tǒng)的克里金方法，INLA-SPDE模型對于小尺度上的空間變異性復(fù)現(xiàn)更為精確。

表3 INLA模型與克里金插值預(yù)測RMSE結(jié)果(μg/m3)

圖5 INLA模型與克里金插值預(yù)測RMSE對比圖

圖6 INLA模型與克里金插值預(yù)測結(jié)果對比圖((a)1月13日；(b)1月15日)

最后相較于傳統(tǒng)的克里金空間插值，因?yàn)樵撃Ｐ偷念A(yù)測結(jié)果本質(zhì)上是所有待預(yù)測格網(wǎng)點(diǎn)處污染濃度值的后驗(yàn)概率密度，因此還可以給出區(qū)域污染風(fēng)險等級的空間分布，這對于人們?nèi)粘５某鲂邪才藕驼块T空間決策提供了支持。京津冀地區(qū)典型日期的污染風(fēng)險如圖7所示。

圖7 京津冀地區(qū)1月12日(左)和11月22日(右)污染風(fēng)險等級空間分布圖

3 區(qū)域PM2.5濃度虛擬仿真決策支持原型系統(tǒng)

本節(jié)基于第2節(jié)建立的模型構(gòu)建京津冀地區(qū)PM2.5濃度虛擬仿真決策支持原型系統(tǒng)。主要目標(biāo)是驗(yàn)證該模型的實(shí)用性，實(shí)現(xiàn)任意時間模型結(jié)果的區(qū)域可視化，呈現(xiàn)時空預(yù)測結(jié)果中蘊(yùn)含的有效信息。并通過交互控制實(shí)現(xiàn)人工模擬大氣環(huán)境的改變來觀測區(qū)域效應(yīng)，從而為普通民眾的出行或政府部門相關(guān)人員的決策提供支持。

3.1 系統(tǒng)模塊設(shè)計

原型系統(tǒng)包含3個基礎(chǔ)模塊，分別為模型層、數(shù)據(jù)層和交互層，如圖8所示。

圖8 原型系統(tǒng)工作模塊圖

將數(shù)據(jù)層中監(jiān)測站的實(shí)時數(shù)據(jù)與有關(guān)地理要素數(shù)據(jù)合并成數(shù)據(jù)表，并存儲到Postgres數(shù)據(jù)庫中。氣象數(shù)據(jù)為格網(wǎng)化數(shù)據(jù)，以GeoTIFF格式存儲。

模型層輸入為選定日期的測站歷史數(shù)據(jù)，調(diào)用相關(guān)月度的INLA-SPDE模型，將結(jié)果數(shù)據(jù)存儲到數(shù)據(jù)層后，再對區(qū)域PM2.5濃度預(yù)報的空間結(jié)果可視化。

交互模塊通過Web前端界面加入相應(yīng)交互按鈕，支持人工模擬新污染源的添加與修改，求解得到交互后的模擬空氣質(zhì)量時空分布結(jié)果并可視化。

3.2 模型的原型系統(tǒng)驗(yàn)證

用于驗(yàn)證的原型系統(tǒng)基本界面主要包括，交互和結(jié)果展示2個區(qū)域。交互區(qū)如圖9左側(cè)邊欄所示，包含日期選擇、想要修改污染值的點(diǎn)位經(jīng)緯度信息、污染輻射范圍和濃度值的錄入，其中位置的經(jīng)緯度信息也可以通過地圖區(qū)的繪制功能交互點(diǎn)選實(shí)現(xiàn)。

圖9為主區(qū)域?yàn)榫┙蚣降貐^(qū)PM2.5濃度時空分布結(jié)果地圖，系統(tǒng)可根據(jù)所選可視化區(qū)域PM2.5濃度空間分布模擬結(jié)果。地圖區(qū)通過左側(cè)的功能欄，支持顯示污染風(fēng)險等級圖和空間分布結(jié)果圖與多種驗(yàn)證數(shù)據(jù)圖的疊加。該結(jié)果展示區(qū)還包括若干功能，如上方顯示的是PM2.5濃度超標(biāo)的測站個數(shù)、預(yù)測PM2.5濃度超標(biāo)的區(qū)域占京津冀地區(qū)的比例和重污染風(fēng)險濃度在75%以上的區(qū)域比例，對于整個區(qū)域這一時刻的污染情況可以有直觀且清晰的感知。右側(cè)則顯示了污染具體的行政區(qū)劃分布情況，包括站點(diǎn)PM2.5濃度排序和以市為基本單位計算的高污染風(fēng)險區(qū)域所占的比例，從而為特定城市的相關(guān)人群提供有價值的參考信息。

界面下方為交互模擬結(jié)果展示，根據(jù)交互區(qū)輸入的相關(guān)信息得到新的空間分布結(jié)果并可視化，與主地圖區(qū)對比，效果如圖10所示。圖中更改濃度值的點(diǎn)位為原結(jié)果圖中的藍(lán)色圓圈，將此處的污染值向下修正到了50 μg/m3，并設(shè)定輻射范圍為100 km。在變動結(jié)果圖中，可以清楚地看到藍(lán)色圓圈周圍的監(jiān)測站值全部降低到75 μg/m3以下(圖中的紅色代表污染值在75 μg/m3以上的站點(diǎn)，綠色代表在75 μg/m3以下的站點(diǎn))。從區(qū)域分布上可知，整個北京以北的地區(qū)污染狀況顯著改善，為了量化這種污染改善的程度，系統(tǒng)右下角還展示了不同污染等級的區(qū)域格網(wǎng)數(shù)目變動百分比，可以看到重度污染以上的區(qū)域數(shù)目大幅減少，而轉(zhuǎn)為良好等級的區(qū)域數(shù)目則明顯上升。

圖9 原型系統(tǒng)基礎(chǔ)界面圖

(a) (b)

因此，通過該交互控制功能，可以清楚地展現(xiàn)在特定區(qū)域減排或設(shè)立新的污染源的影響程度，為相關(guān)政策的針對性制定和部分建筑的空間選址問題提供了重要的參考，真正實(shí)現(xiàn)了空間決策與支持功能，進(jìn)一步驗(yàn)證了模型的實(shí)用價值。

4 結(jié)論與反思

本文提出了基于INLA-SPDE模型的區(qū)域污染物模擬與預(yù)測方法，實(shí)現(xiàn)了京津冀地區(qū)PM2.5濃度日均值的時空建模與預(yù)測，驗(yàn)證了模型預(yù)測精度，得到了京津冀地區(qū)5 km×5 km格網(wǎng)的日均PM2.5濃度時空分布預(yù)報圖和污染風(fēng)險等級圖等。結(jié)果表明，INLA-SPDE方法在區(qū)域污染模擬與預(yù)測上相比傳統(tǒng)空間插值方法具有更好的效果，尤其對于高污染值的預(yù)測更為精確，并能給出區(qū)域污染風(fēng)險等級預(yù)測結(jié)果，在區(qū)域污染的預(yù)警和防控上具有重要的現(xiàn)實(shí)意義。基于該模型構(gòu)建了京津冀地區(qū)PM2.5濃度預(yù)報與控制決策支持原型系統(tǒng)，解決了區(qū)域級別的PM2.5濃度時空可視化的問題。通過將人工模擬得到的仿真結(jié)果與原始結(jié)果對比，為人類的空間決策行為提供了有效地支持，驗(yàn)證了該模型的實(shí)用性與價值。

[1] Report from Intergovernmental Panel on Climate Change (IPCC). Climate change 2007: the physical science basis[EB/OL]. [2021-04-01]. http://www.doc88.com/p- 293361412559.html.

[2] APTE J S, BRAUER M, COHEN A J, et al. Ambient PM2.5reduces global and regional life expectancy[J]. Environmental Science & Technology Letters, 2018, 5(9): 546-551.

[3] CAI S Y, MA Q, WANG S X, et al. Impact of air pollution control policies on future PM2.5concentrations and their source contributions in China[J]. Journal of Environmental Management, 2018, 227: 124-133.

[4] GAUDARD M, KARSON M, LINDER E, et al. Bayesian spatial prediction[J]. Environmental and Ecological Statistics, 1999, 6(2): 147-171.

[5] SARTO S D, RANALLI M G, BAKAR K S, et al. Bayesian spatiotemporal modeling of urban air pollution dynamics[M]// Topics on Methodological and Applied Statistical Inference. Cham: Springer International Publishing, 2016: 95-103.

[6] RUE H, MARTINO S, CHOPIN N. Approximate Bayesian inference for latent Gaussian models by using integrated nested Laplace approximations[J]. Journal of the Royal Statistical Society: Series B: Statistical Methodology, 2009, 71(2): 319-392.

[7] LINDGREN F, RUE H, LINDSTR?M J. An explicit link between Gaussian fields and Gaussian Markov random fields: the stochastic partial differential equation approach[J]. Journal of the Royal Statistical Society: Series B: Statistical Methodology, 2011, 73(4): 423-498.

[8] CAMELETTI M, GóMEZ-RUBIO V, BLANGIARDO M. Bayesian modelling for spatially misaligned health and air pollution data through the INLA-SPDE approach[J]. Spatial Statistics, 2019, 31: 100353.

[9] FIORAVANTI G, MARTINO S, CAMELETTI M, et al. Spatio-temporal modelling of PM10daily concentrations in Italy using the SPDE approach[J]. Atmospheric Environment, 2021, 248: 118192.

[10] 丁書珍, 張輝國, 胡錫健, 等. 利用R-INLA方法分析宏觀因素對艾滋病疫情的影響[J]. 中國艾滋病性病, 2018, 24(12): 1192-1196.

DING S Z, ZHANG H G, HU X J, et al. Influence of macroscopic factors on HIV/AIDS epidemic based on R-INLA method[J]. Chinese Journal of AIDS & STD, 2018, 24(12): 1192-1196 (in Chinese).

[11] BLANGIARDO M, CAMELETTI M. Spatial and spatio-temporal Bayesian models with R-INLA[M]. Chichester: John Wiley & Sons, Ltd, 2015: 268-283.

[12] SIMPSON D, RUE H, RIEBLER A, et al. Penalising model component complexity: a principled, practical approach to constructing priors[EB/OL]. [2021-04-04]. https://arxiv.org/ pdf/1403.4630.pdf.

[13] S?RBYE S H, RUE H. Penalised complexity priors for stationary autoregressive processes[J]. Journal of Time Series Analysis, 2017, 38(6): 923-935.

Simulation and prediction of regional pollutants based on INLA-SPDE method

YUAN Ze, CHEN Bin

(School of Earth and Space Sciences, Peking University, Beijing 100871, China)

The simulation and prediction of regional pollutants generally use the traditional spatial interpolation method, which cannot obtain accurate results when the source data is not uniformly distributed. To address these problems, a method for simulation and prediction of regional pollutants based on the INLA-SPDE model was proposed. The interpolation model was based on a Bayesian hierarchical model where the spatial-component was represented through the stochastic partial differential equation (SPDE) approach, with a lag-1 temporal autoregressive component (AR1). In addition, the model included 10 spatial and spatio-temporal predictors such as meteorological variables. By building 12 models for each month with the integrated nested Laplace approximation (INLA), this research realized the spatio-temporal simulation and prediction of PM2.5concentration at daily resolution in the Beijing-Tianjin-Hebei region in 2019. Experiments show that compared with traditional Kriging interpolation methods, the proposed model can yield a better prediction of air pollutants at regional scale. Particularly, the prediction accuracy of high-value pollutants was improved significantly, and air pollutants exceedance probabilities can also be generated. Furthermore, a system for regional PM2.5concentration simulation and decision support was established, the system can provide support for the travel of ordinary people or the decision-making of government officials, and verify the practicability and value of the proposed model.

PM2.5; spatio-temporal Bayesian hierarchical model; integrated nested Laplace approximation; simulation system; decision support

18 May，2021；

TP 391

10.11996/JG.j.2095-302X.2022010125

A

2095-302X(2022)01-0125-08

2021-05-18；

2021-08-21

21 August，2021

袁 澤(1992–)，男，碩士研究生。主要研究方向?yàn)樘摂M地理環(huán)境。E-mail：180121089@pku.edu.cn

YUAN Ze (1992–), master student. His main research interest covers virtual geographic environments. E-mail：180121089@pku.edu.cn

陳 斌(1973–)，男，教授，博士。主要研究方向?yàn)樘摂M地理環(huán)境。E-mail：gischen@pku.edu.cn

CHEN Bin (1973–), professor, Ph.D. His main research interest covers virtual geographic environments.E-mail：gischen@pku.edu.cn