基于偏最小二乘回歸的空氣質(zhì)量數(shù)據(jù)校準(zhǔn)研究

丁學(xué)利，任 鵬

（阜陽職業(yè)技術(shù)學(xué)院，安徽 阜陽 236031）

0 引言

空氣污染對(duì)生態(tài)環(huán)境和人類健康危害巨大，所以對(duì)空氣質(zhì)量進(jìn)行監(jiān)測(cè)必不可少。目前國家通過監(jiān)測(cè)控制站點(diǎn)（國控點(diǎn)）對(duì)PM2.5、PM10、CO、NO2、SO2、O（3“兩塵四氣”）的濃度進(jìn)行實(shí)時(shí)監(jiān)測(cè)［1］，雖然數(shù)據(jù)較為準(zhǔn)確，但因布控點(diǎn)較少，且數(shù)據(jù)發(fā)布延遲較長、投資較大，無法對(duì)空氣質(zhì)量進(jìn)行實(shí)時(shí)監(jiān)測(cè)和預(yù)報(bào)。通常做法是在國控點(diǎn)近鄰自建一些監(jiān)測(cè)控制站點(diǎn)（自建點(diǎn)），利用微型空氣質(zhì)量檢測(cè)儀對(duì)某一地區(qū)空氣質(zhì)量進(jìn)行實(shí)時(shí)網(wǎng)格化監(jiān)控，并同時(shí)監(jiān)測(cè)溫度、風(fēng)速、氣壓、濕度、降水量等氣象參數(shù)［1，2］。由于所使用的電化學(xué)氣體傳感器在長時(shí)間使用后會(huì)產(chǎn)生一定的零點(diǎn)漂移和量程漂移，非常規(guī)氣態(tài)污染物（氣）濃度變化對(duì)傳感器也存在交叉干擾，以及天氣因素對(duì)傳感器的影響，在自建點(diǎn)上，同一時(shí)間微型空氣質(zhì)量檢測(cè)儀所采集的數(shù)據(jù)與該國控點(diǎn)的數(shù)據(jù)值存在一定的差異。因此，需要利用國控點(diǎn)每小時(shí)的數(shù)據(jù)對(duì)國控點(diǎn)近鄰的自建點(diǎn)數(shù)據(jù)進(jìn)行校準(zhǔn)。

在空氣質(zhì)量數(shù)據(jù)的校準(zhǔn)方面，可通過建立最小二乘擬合或多元線性回歸模型等［3-9］方法進(jìn)行校準(zhǔn)。但由于空氣質(zhì)量數(shù)據(jù)的校準(zhǔn)影響因素較多，且各因素具有較強(qiáng)的非線性，各因素之間存在較強(qiáng)的相關(guān)性，普通的線性回歸模型很難對(duì)數(shù)據(jù)進(jìn)行校準(zhǔn)。因此，本文考慮建立偏最小二乘回歸（PLSR）模型對(duì)空氣質(zhì)量數(shù)據(jù)進(jìn)行校準(zhǔn)。

1 PLSR建模步驟［10-13］

設(shè)有n個(gè)樣本點(diǎn)，矩陣表示含有p個(gè)自變量，矩陣表示含有q個(gè)因變量，PLSR具體步驟如下：

Step 1：將X和Y分別進(jìn)行標(biāo)準(zhǔn)化處理。

Step 2：計(jì)算第1 對(duì)成分t1和u1。記E0的第1 個(gè)成分t1，F(xiàn)0的第1 個(gè)成分u1。為使t1和u1的相關(guān)程度達(dá)到最大，需使如下的內(nèi)積θ1達(dá)到最大。

w1可由矩陣計(jì)算其最大特征值對(duì)應(yīng)的特征向量得到，v1可通過計(jì)算得到。計(jì)算出w1和v1，即可得到第1對(duì)成分：

Step 3：分別建立E0和F0對(duì)t1的回歸方程。

其中E1和F1為殘差矩陣，回歸系數(shù)向量α1和β1如下：

Step 4：用E1和F1分別替代E0和F0，重復(fù)上述步驟。若F1中元素近似為0，則表明第1 個(gè)成分得到的回歸模型精度已達(dá)到要求，可終止成分的抽取；否則，用E1和F1分別替代E0和F0，重復(fù)上述步驟，即可得第2對(duì)成分t2=E1w2，u2=F1v2。則有

Step 5：若E0的 秩 是r，則 存 在r個(gè) 成 分t1,t2,…,tr，于是有

Step 6：交叉有效性檢驗(yàn)［11］。

PLSR 方程一般不需要使用全部的成分t1,t2,…,tr進(jìn)行回歸建模，而是通過截取前h個(gè)成分(h＜r)，即可得到一個(gè)回歸效果較為理想的模型。提取h個(gè)成分的交叉有效性檢驗(yàn)定義如下：

每一次提取成分結(jié)束前，都利用（9）式進(jìn)行檢驗(yàn)。當(dāng)?shù)趆步時(shí)，Q2h＜0.0975，則停止提取成分；否則繼續(xù)提取成分直到達(dá)到精度要求為止。

（9）式中p(h)為預(yù)測(cè)誤差的平方和如下：

其中

（9）式中Ss(h)的表達(dá)式如下：

2 空氣質(zhì)量校準(zhǔn)效果的輔助計(jì)算

2.1 樣本異常值的檢驗(yàn)［10，14］

樣本中異常點(diǎn)的存在會(huì)導(dǎo)致統(tǒng)計(jì)規(guī)律產(chǎn)生較大波動(dòng)，從而使回歸線發(fā)生較大偏離。第i個(gè)樣本點(diǎn)對(duì)第h個(gè)成分th的貢獻(xiàn)率定義為：

2.2 校準(zhǔn)誤差分析［12，13］

為檢驗(yàn)校準(zhǔn)結(jié)果的好壞，需用多個(gè)評(píng)價(jià)指標(biāo)對(duì)校準(zhǔn)效果進(jìn)行整體性的綜合評(píng)價(jià)和衡量。本文應(yīng)用如下評(píng)價(jià)指標(biāo)對(duì)校準(zhǔn)效果進(jìn)行評(píng)價(jià)。

（1）校準(zhǔn)前平方和誤差

n

（2）校準(zhǔn)前均方誤差

（3）校準(zhǔn)前平均絕對(duì)誤差

校準(zhǔn)后平均絕對(duì)誤差

（4）校準(zhǔn)結(jié)果改善百分比物濃度，y?ij為污染物濃度校準(zhǔn)值（j=1,2,…,6）。

其中xij為自建點(diǎn)污染物濃度，yij為國控點(diǎn)污染

3 基于PLSR空氣質(zhì)量數(shù)據(jù)的校準(zhǔn)

3.1 數(shù)據(jù)來源與提取

本文采用2019 年全國大學(xué)生數(shù)學(xué)建模競(jìng)賽D題［1］空氣質(zhì)量數(shù)據(jù)。利用Access 數(shù)據(jù)庫軟件提取自建點(diǎn)符合條件的數(shù)據(jù)，然后將整點(diǎn)附近間隔在5分鐘內(nèi)的自建點(diǎn)數(shù)據(jù)進(jìn)行平均得到與國控點(diǎn)整點(diǎn)匹配的數(shù)據(jù)。刪除不能匹配的數(shù)據(jù)之后，共提取4048 條數(shù)據(jù)，如表1（展示部分?jǐn)?shù)據(jù)）所示。表1 中，x1、x2、x3、x4、x5、x6、x7、x8、x9、x10和x11分別表示自建點(diǎn)的PM2.5 濃度、PM10 濃度、CO 濃度、NO2濃度、SO2濃度、O3濃度、風(fēng)速、壓強(qiáng)、降水量、溫度、濕度；y1、y2、y3、y4、y5和y6分 別 表 示 國 控 點(diǎn) 的PM2.5 濃 度、PM10濃度、CO濃度、NO2濃度、SO2濃度、O3濃度。

3.2 相關(guān)性分析

表2是對(duì)自建點(diǎn)的11個(gè)自變量和國控點(diǎn)的6個(gè)因變量進(jìn)行相關(guān)性分析的結(jié)果。從表2可知，x1與x2、y1、y2；x2與y1；y1與y2之間具有較強(qiáng)的相關(guān)性（大于0.7），而x8與x10的相關(guān)系數(shù)是-0.85，說明具有較強(qiáng)的負(fù)相關(guān)，其他因素之間的相關(guān)系數(shù)均小于0.7。根據(jù)各因素之間的相關(guān)性，可考慮以自建點(diǎn)的11個(gè)變量為自變量，國控點(diǎn)的6個(gè)變量為因變量建立PLSR模型。

3.3 異常數(shù)據(jù)處理

以x1-x11為自變量，y1-y6為因變量作PLSR，提取兩個(gè)主成分t1和t2。根據(jù)（15）式，可在t1-t2平面上畫出散點(diǎn)圖和橢圓圖，如圖1。圖1中的橢圓之外有318個(gè)點(diǎn)，則認(rèn)為這些點(diǎn)是異常點(diǎn)，將其提取并進(jìn)行剔除。

3.4 校準(zhǔn)結(jié)果分析

利用清洗后的數(shù)據(jù)重新進(jìn)行PLSR 建模。首先做交叉有效性檢驗(yàn)（見表3）。當(dāng)h=4 時(shí)，因此，只要抽取前4 個(gè)主成分即可建立合理的PLSR 方程，其回歸方程的系數(shù)如表4所示。此回歸方程即可作為自建點(diǎn)空氣質(zhì)量的校準(zhǔn)方程。

為了對(duì)表4 構(gòu)成的6 個(gè)回歸方程進(jìn)行精度分析，分別作了如下分析。

圖1 t1-t2成分橢圓圖

表1 自建點(diǎn)與國控點(diǎn)匹配數(shù)據(jù)

表2 相關(guān)系數(shù)矩陣

（1）自建點(diǎn)6 種空氣質(zhì)量數(shù)據(jù)的預(yù)測(cè)。圖2 是以自建點(diǎn)的校準(zhǔn)值為橫坐標(biāo)，國控點(diǎn)的實(shí)際值為縱坐標(biāo)，對(duì)6 種空氣質(zhì)量數(shù)據(jù)的樣本點(diǎn)作預(yù)測(cè)圖。在預(yù)測(cè)圖2上，若所有的點(diǎn)都能在直線y=x左右均勻分布，則回歸方程的擬合效果是滿意的，說明校準(zhǔn)值與實(shí)際值差異很小。在圖2中，除SO2的預(yù)測(cè)效果較差外，其余空氣質(zhì)量數(shù)據(jù)均能在直線y=x左右均勻分布，說明校準(zhǔn)效果較滿意。

表3 交叉檢驗(yàn)值

（2）6 種空氣質(zhì)量數(shù)據(jù)的校準(zhǔn)值與實(shí)際值的時(shí)序圖，如圖3。除SO2的前半部分以及O3的后半部分的校準(zhǔn)值與國控點(diǎn)數(shù)據(jù)吻合的較差外，其他數(shù)據(jù)的校準(zhǔn)值都能與國控點(diǎn)數(shù)據(jù)較好地重合，表明利用PLSR取得了較好的校準(zhǔn)效果。

表4 PLSR方程系數(shù)

（3）對(duì)6 種空氣質(zhì)量數(shù)據(jù)校準(zhǔn)前后進(jìn)行誤差分析。由表5 知，各個(gè)校準(zhǔn)誤差指標(biāo)均低于校準(zhǔn)之前的誤差指標(biāo)。另外從整體校準(zhǔn)改善的百分比看，PM10 校準(zhǔn)改善的最好，達(dá)到了63.89%；SO2校準(zhǔn)改善的效果最差，只有16.02%。從整體上看自建點(diǎn)的空氣質(zhì)量數(shù)據(jù)均得到不同程度的校準(zhǔn)，說明基于PLSR模型的校準(zhǔn)方法取得了較好的效果。

圖2 自建點(diǎn)6種空氣質(zhì)量數(shù)據(jù)的預(yù)測(cè)圖

圖3 6種空氣質(zhì)量數(shù)據(jù)的校準(zhǔn)值與國控點(diǎn)實(shí)際值的時(shí)序圖

4 結(jié)語

從數(shù)據(jù)校準(zhǔn)結(jié)果知，自建點(diǎn)空氣質(zhì)量數(shù)據(jù)校準(zhǔn)后的誤差都有明顯降低。同時(shí)發(fā)現(xiàn)PLSR 模型對(duì)具有多重共線性的自變量（如PM2.5 和PM10）校準(zhǔn)效果明顯。但并非對(duì)自建點(diǎn)所有的污染物濃度都能有效地校準(zhǔn)，如對(duì)SO2校準(zhǔn)改善的效果最差，其次是O3，這可能受污染物濃度變化對(duì)傳感器的交叉干擾，以及氣象因素劇烈變化對(duì)傳感器的影響。今后可進(jìn)一步嘗試尋找SO2數(shù)據(jù)校準(zhǔn)的關(guān)鍵影響因素，建立符合數(shù)據(jù)特點(diǎn)的非線性回歸模型。

表5 6種空氣質(zhì)量數(shù)據(jù)校準(zhǔn)前后的誤差對(duì)比