利用主成分分析法及地理加權回歸模型分析AOD數據

李廣超，李如仁，盧月明，趙陽陽，余 博

(1. 沈陽建筑大學交通學院，遼寧 沈陽 110168； 2. 中國測繪科學研究院，北京 100830)

全球大氣氣溶膠類型和含量變化與氣候變化和大氣環境污染密切相關，是氣象學、環境學和醫學研究關注的熱點問題[1]。氣溶膠光學厚度作為大氣氣溶膠最重要的參數之一，是表征大氣混濁度或大氣中氣溶膠總含量的重要物理量[2]。許多學者對氣溶膠光學厚度影響因素的相關性進行了研究，發現地理數據(高程等)和氣象數據(濕度等)是影響AOD數據的顯著因素[3-6]。因此，可采用氣象數據和地理數據來估算AOD數據。在地學空間分析中，觀測數據是在不同的空間位置上獲取的，全局空間回歸模型就是假定回歸參數與樣本數據的地理位置無關，而在實際問題研究中經常發現回歸參數隨地理位置變化，這時如果仍采用全局空間回歸模型，得到的回歸參數估計將是回歸參數在整個研究區域內的平均值，不能反映回歸參數的真實空間特征[7-8]。當空間數據存在自相關時，地理加權回歸(geographically weighted regression，GWR)模型提供了一個優于傳統OLS模型的估計方法，傳統的OLS方法僅僅提供全局參數的估計，GWR容許分解成局部參數估計，深刻地闡釋了地理空間數據的某類指標和空間影響因子之間的關系，是傳統OLS無法比擬的[9]。GWR模型是一種典型的局部模型，認為回歸系數隨著空間位置的變化而變化，具有空間非平穩性[10]。

近年來，許多研究者對GWR進行了研究。在GWR模型應用方面，武文娟等[11]以四川省縣公立醫院床位數為例，利用空間分析技術研究了區縣床位的時空變化特征，并采用GWR解釋了經濟、人口、交通、地形等因素對其時空變異性的影響。董沖亞等[12]應用地理加權回歸模型探討了氣象因素和大氣污染因素影響我國女性肺癌發病的空間屬性，結果表明東北地區尤其是遼東半島為我國女性肺癌的高發區域，具有明顯的地區集聚性。孫偉偉[13]以長三角經濟區為例，采用MODIS/Terra AOT時序衛星影像產品為數據源，構建AOD-PM2.5季節地理加權回歸模型，為PM2.5區域聯防聯控提供了信息支撐和科學依據。陳輝等[14]采用GWR進行回歸分析，構建了我國區域范圍內近地面的PM2.5遙感反演模型，結果表明利用GWR進行PM2.5遙感估算，既能體現PM2.5時空分布的全局變化特征，又能從局部體現全國PM2.5組分、污染程度及垂直分布結構的空間差異性特性。龐瑞秋等[15]以吉林省各縣域為基本單元，結合GWR和空間相關分析方法，討論了人口城鎮化水平與國有動力、農業動力及外向動力等因素的空間相關關系，并依次解釋人口城鎮化分縣域差異的影響因素。在GWR模型改進方面，趙陽陽等[16]將半監督學理論與GWR相結合，提出基于半監督學習的地理加權回歸方法，并通過試驗驗證了該方法在模型性能上的優勢。覃文忠等[17]以迭代算法為基礎，推導出混合GWR的常系數(全局參數)和變系數(局部參數)的計算方法，并以上海市樓盤價格進行驗證，結果表明混合GWR的計算量略大于GWR，但對樣本數據的擬合效果更好，局部參數估計更穩健。趙陽陽等[18]基于貪心算法，通過引入Akaike信息法則，設計了適用于GWR的特征變量選擇方法。

由上述研究可以得出，AOD數據也可以利用GWR模型進行建模，但如果只是簡單地加入到模型計算中，沒有考慮各個影響因素之間的共線性，且考慮的影響因素較多會嚴重影響模型的計算效率。針對上述問題，劉蓓[19]將非線性主成分分析法(nonlinear principal component analysis，NLPCA)應用于巖體質量等級分類問題中，構建了巖體質量等級分類的非線性主成分綜合評價模型，結果表明，該評價模型有效且預測精度高。高青松等[20]提出的非線性檢驗度量可有效地檢測出給定數據集中各變量之間的線性或非線性關系，為是否對數據采用非線性主成分分析提供了依據。檀菲菲等[21]提出了利用非線性主成分分析法和施密特正交化(NLPCA-GSO)相結合的方法評價區域的可持續發展水平來彌補傳統方法的不足，結果表明基于NLPCA-GSO的可持續發展水平模型很好地彌補了傳統主成分分析及對各子系統結果的綜合評價的不足。周永正等[22]提出了非線性主成分分析法與神經網絡算法的融合模型，并將非線性主成分神經網絡融合模型應用于水泥強度的預測研究，得到的結果表明預測誤差很小。孫康等[23]提出了一種基于非線性主成分分析的高光譜圖像目標探測算法，試驗表明，基于神經網絡的非線性主成分分析法可以將線性不可分的目標與背景分離，使用非線性特征和原始特征的組合可以獲得更好的目標探測效果。

綜上論述，本文提出一種基于主成分分析的地理加權回歸方法(PCA-GWR)，該方法檢驗了AOD影響因素之間的共線性；通過非線性主成分分析法對影響AOD值的若干相關變量進行處理，既消除了相關變量彼此之間的多重共線性，又可以起到降維的作用；再利用得到的較少幾個綜合指標，通過地理加權回歸模型對AOD濃度進行預測；最后通過與常規GWR模型對比，采用MAE、RMSE、AIC、R2作為評價指標，評價了本文方法的有效性。

1 研究方法

1.1 非線性主成分分析

非線性主成分分析是在線性主成分分析基礎上的擴展。傳統的主成分分析方法在一些情況下降維不好，一般只能處理線性問題。綜合評價的實際結果與評價指標間的相關程度高低成正比，評價指標間相關程度越高，主成分分析的結果越好，指標間的相關性越小，每一個主成分承載的信息量就越少，為了滿足累積貢獻度達到一定水平(通常為85%以上)，可能需要選取較多的主成分，此時主成分分析法的降維效果將會不明顯。非線性主成分分析通過對原始數據作中心化對數比變換，將主成分表示為原始數據的非線性組合，可以較好地保留數據本身的非線性特征；分析的出發點是協方差矩陣，而不是之前的相關系數矩陣，這樣會明顯提升降維效果，用更少的主成分反映更多的原始指標信息，并且評價的穩定性與合理性也有所提高。非線性主成分分析原理如下：

設有P維向量x=(x1,x2,…,xp)的樣本資料(xij)n×p

對原始數據作中心化對數比變換

(1)

計算中心化對數比樣本協方差矩陣

S=(Sij)p×p

(2)

從S出發求樣本主成分如下：

設λ1>λ2…>λp是S的P個特征根，(a1,a2,…,ap)是相應的標準化特征向量，則第i個主成分為

(3)

1.2 地理加權回歸

地理加權回歸是英國圣安德魯斯大學的Fortheringham等[24]在空間變異系數回歸的基礎上利用局部光滑思想提出的。地理加權回歸是對普通線性回歸的擴展，即將樣本點的地理位置引入到回歸參數中，其公式為

(4)

式中，(yi,xi1,xi2,…,xid)為因變量y和自變量(x1,x2,…,xd)在數據點(ui,vi)處的n組數據值；βk(ui,vi)(k=1,2,…,d)為第i個觀測點(ui,vi)處的未知參數；εi(i=1,2,…,n)為獨立同分布的誤差項，通常假定其服從N(0,δ2)分布。

地理加權回歸模型中的回歸參數與樣本數據的地理位置有關，其影響程度(空間權重)可用一個距離函數表示，該函數簡稱核函數。常用的核函數有高斯核函數、Bisquare核函數等。GWR模型的關鍵在于選擇核函數，并確定其最優帶寬。研究發現，不同核函數的帶寬敏感度不同，而帶寬的變化對結果影響較大[16]，因此，可通過核函數和帶寬來區分回歸模型。

核函數的帶寬過大會導致回歸參數估計偏大，過小則導致回歸參數估計偏小[7]。為減小帶寬不適造成的誤差，本文采用CV交叉驗證法[25]來計算最優帶寬。CV法計算公式為

(5)

1.3 PCA-GWR原理

本文充分結合了非線性主成分分析法與地理加權回歸模型的優勢，提出了一種基于主成分分析的地理加權回歸方法(PCA-GWR)。該方法在使用GWR預測之前，首先檢測數據之間的相關性，然后利用主成分分析法去除數據之間的多重共線性，并得到幾個綜合指標，在綜合指標中選取累積貢獻度超過85%的前幾個指標作為GWR模型的輸入變量，從而對京津冀地區的AOD進行分析預測。PCA-GWR的主要流程為：①利用Pearson相關系數檢驗主變量相關變量數據之間的相關性；②對存在相關性的數據進行非線性主成分分析預處理，去除數據之間的多重共線性，并得到幾個綜合指標；③從綜合指標中選取累積貢獻度超過85%的前幾個指標作為GWR模型的輸入變量，對主變量(PM2.5濃度)進行分析或預測；④將分析或預測結果與其他模型對比，采用MAE、RMSE、AIC、R2作為評價指標，來驗證本文方法的預測精度。

2 試 驗

2.1 研究區概況

研究區域選用地理范圍為35.5°N—43°N、113°E—120°E的京津冀地區，土地面積約為21.8萬km2，包含北京、天津與河北。北京位于華北平原的北部，背靠燕山，有永定河流經老城西南，毗鄰天津市和河北省；天津位于華北平原海河五大支流匯流處，東臨渤海，北依燕山，海河在城中蜿蜒而過，是天津的母親河；河北東臨渤海、內環京津，轄保定、唐山、石家莊、邢臺、邯鄲、衡水、滄州、廊坊、秦皇島、張家口、承德等11個地級市。京津冀是中國北方經濟規模最大、最具活力的地區，2015年底人口約為11 143萬人，地區生產總值約為69 358.89億元。隨著京津冀經濟的快速發展，生態環境也隨之破壞嚴重，因此優化生態環境成為京津冀地區的重要任務。大氣氣溶膠類型和含量變化與氣候變化和大氣環境污染密切相關，開展大氣氣溶膠的研究對大氣環境污染分析及防控具有重要意義。

2.2 數據來源

本文以地理空間數據(高程、坡度、坡向)和氣象數據(風速、氣溫、濕度、氣壓)為自變量，以AOD數據為因變量，估算AOD數據的值。其中地理空間數據來自地理空間數據云網站(http:∥www.gscloud.cn/)，選用SRTMEDM 90 m分辨率原始高程數據、SRTMSLOPE 90 m分辨率坡度數據、SRTMASPECT 90 m分辨率坡向數據；氣象數據來自中國氣象科學數據共享服務網(http:∥www.escience.gov.cn)，共110個氣象監測站點地理位置信息，頻率為每天；AOD數據來自 (https:∥ladsweb.nascom.nasa.gov/data)網站的Terra MODIS C06二級氣溶膠產品(代號MYD04_3K)，頻率為每天，空間分辨率為3 km，本文采用MODIS Collection 6 MYD04_3K數據集中參數名為“Optical_ Depth_ Land_ And_ Ocean”，波段為550 nm的2級AOD數據。本文選擇2015年5月期間的氣象和AOD數據作為研究對象。京津冀AOD數據采樣點分布如圖1 所示。(注：1 mile=1.61 km)

2.3 數據預處理

為了數據的時空一致性，對地理空間數據、氣象數據、 AOD數據進行預處理， 流程如圖2所示。首先對氣象數據進行處理，采用SQL Server求110個氣象站點處日均值，再利用Java語言程序求氣象數據的月均值；其次對地理空間數據進行處理，將高程、坡度、坡向影像數據進行投影坐標轉換，再利用ArcGIS10.2對高程、坡度、坡向影像數據進行氣象站點空間數據數值的提取；最后對AOD數據處理，針對京津冀地區創建一個覆蓋全區域的5 km×5 km網格，并提取網格中心點的位置坐標，以網格中心點代表該網格的空間位置，利用C#、ArcGIS Engine、Visual Studio 2013程序對MODIS影像數據進行批量處理，通過重采樣獲取網格中心點AOD的值，利用Java語言程序求出網格中心點AOD數據月均值，再利用ArcGIS10.2對網格中心點AOD數據的月均值進行克里金插值，然后提取氣象站點處AOD的值。

圖1 京津冀AOD數據采樣點分布

圖2 數據預處理流程

2.3.1 影響因素相關性檢驗

通過相關系數矩陣檢驗變量間是否存在多重共線性，相關矩陣見表1。由表1可看出，濕度、氣壓和氣溫與AOD成正相關，高程、風速、坡度和坡向與AOD成負相關，也可看出本文考慮的這7項指標(濕度、氣壓、高程、風速、氣溫、坡度、坡向)均與AOD存在一定程度上的相關性，大部分指標的相關系數均在0.30以上，其中高程和氣壓與AOD的相關系數更是在0.57以上。可以看出影響主變量的相關變量之間存在多重共線性，通過主成分分析法來消除相關變量之間的多重共線性[19]。(以下數據均保留5位小數)

表1 Pearson相關系數

2.3.2 非線性主成分分析

針對AOD數據的預測，對7項指標作非線性主成分分析，本文前3個主成分的特征根及其相應的貢獻度與累積貢獻度，見表2。(以下數據均保留5位小數)

表2 非線性主成分分析結果

由表2的數據預測分析結果可知，第一主成分的特征根、貢獻度、累積貢獻度分別為1.553 74、0.646 51、0.646 51；第二主成分的特征根、貢獻度、累積貢獻度分別為1.378 96、0.155 97、0.802 47；第三主成分的特征根、貢獻度、累積貢獻度分別為1.193 89、0.135 04、0.937 52。前3個主成分的累積貢獻度高達93.752%，已接近100%，表明前3個主成分足以代表原始數據的絕大部分信息，因此本文選取前3主成分作為模型的輸入變量。

3 結果與分析

為了評估本文研究方法的預測效果，將其與常規的GWR模型進行對比，由表1可知，氣壓、高程、氣溫與AOD相關系數在7項指標中最高，以AOD為因變量，以氣壓、高程、氣溫為自變量，作GWR回歸分析，如圖3所示。通過非線性主成分分析后，以第一主成分、第二主成分、第三主成分為自變量，AOD為因變量作PCA-GWR回歸分析，如圖4所示。同時，分別計算各個模型的平均絕對誤差(MAE)、均方根誤差(RMSE)、Akaike信息量(AIC)、擬合優度(R2)等4項評價指標，以及模型之間的提升度，計算結果見表3。

表3 預測結果誤差對比

圖3 GWR分析結果

圖4 PCA-GWR分析結果

從表3可看出，在4項評價指標MAE、RMSE、AIC、R2中，PCA-GWR方法所得值分別為0.235 06、0.322 25、116.033 14、0.607 26；GWR方法所得值分別為0.258 84、0.340 95、123.207 71、0.576 29。PCA-GWR模型分析所得結果優于常規GWR模型，充分說明本文方法預測結果的優越性。其中，PCA-GWR模型比常規GWR模型的MAE提升9.19%、RMSE提升5.48%、AIC提升5.82%、R2提升5.37%。證明在使用GWR模型預測之前，使用非線性主成分分析法對數據預處理的可行性與優越性；同時，也可以看出PCA-GWR模型比GWR模型的預測效果具有顯著的提升。

4 結 論

本文方法在變量個數相同的前提下對AOD濃度的預測精度有明顯的提升，說明該方法不僅可以有效地降低計算工作量，減少原始數據信息損失，簡化數據結構，還可以消除各個影響因素之間的多重共線性，提高了預測精度。

(1) 相關分析得出濕度、氣壓和氣溫與AOD呈正相關，高程、風速、坡度和坡向與AOD呈負相關，這7個影響因素與AOD相關性大小依次為氣壓>高程>氣溫>坡度>風速>濕度>坡向。

(2) 通過交叉驗證，說明本文采用非線性主成分分析處理后綜合指標進行GWR模型預測的精度較高，且減輕了計算量，提高了運算效率。經過對比，本文方法所得4項評價指標MAE、RMSE、AIC、R2均優于常規GWR模型，MAE提升9.19%、RMSE提升5.48%、AIC提升5.82%、R2提升5.37%。

在后續的研究中，將進一步加入其他考慮因子，如污染氣體(如CO、SO2)、人類活動、地表覆蓋等因子，進而對AOD數據精度進行修正，得到高精度的AOD數據。

參考文獻：

[1] 李曉靜,高玲,張興贏,等.衛星遙感監測全球大氣氣溶膠光學厚度變化[J].科技導報,2015,33(17):30-40.

[2] 劉浩，高小明，謝志英,等.京津冀晉魯區域氣溶膠光學厚度的時空特征[J].環境科學學報，2015，35(5):1506-1511.

[3] 張磊,江洪,陳誠,等.廣東省MODIS氣溶膠光學厚度時空分布及其影響因素[J].地理空間信息,2017,15(1):46-49.

[4] 王浩洋.遙感反演安徽地區氣溶膠光學厚度及其時空特征分析[D].合肥：安徽大學,2015.

[5] 黎麗莉.廣東氣溶膠光學厚度及穗深空氣污染物時空特征和影響因素研究[D].廣州：中國科學院研究生院(廣州地球化學研究所),2015.

[6] 葉瑜,李秀央,陳坤,等.大氣氣溶膠光學厚度與大氣污染物及氣象因素關系的時間序列研究[J].氣候與環境研究,2011,16(2):169-174.

[7] 覃文忠.地理加權回歸基本理論與應用研究[D].上海：同濟大學,2007.

[8] 楊毅.顧及時空非平穩性的地理加權回歸方法研究[D].武漢：武漢大學,2016.

[9] 湯慶園，徐偉，艾福利.基于地理加權回歸的上海市房價空間分異及其影響因子研究[J].經濟地理，2012，32(2):52-58.

[10] HUANG B，WU B，BARRY M.Geographically and Temporally Weighted Regression for Modeling Spatio-temporal Variation in House Prices[J].International Journal of Geographical Information Science，2010，24(3):383-401.

[11] 武文娟，徐京華，時進,等.基于GWR的四川省醫院床位數時空分布及其影響因素研究[J].測繪通報，2016(4):49-53.

[12] 董沖亞，康曉平.基于地理加權回歸模型的我國女性肺癌發病空間影響因素分析[J].環境與健康雜志，2014，31(9):769-772.

[13] 孫偉偉.基于GWR模型的PM2.5遙感估算模型研究——以長三角為例[C]∥浙江省地理學會.浙江省地理學會2016年學術年會暨浙江省第三屆地理名師名校長聯盟高峰論壇論文摘要集.杭州:浙江省地理學會，2016:1.

[14] 陳輝，厲青，張玉環,等.基于地理加權模型的我國冬季PM2.5遙感估算方法研究[J].環境科學學報，2016，36(6):2142-2151.

[15] 龐瑞秋，騰飛，魏冶.基于地理加權回歸的吉林省人口城鎮化動力機制分析[J].地理科學，2014，34(10):1210-1217.

[16] 趙陽陽，劉紀平，徐勝華,等.一種基于半監督學習的地理加權回歸方法[J].測繪學報，2017,46(1):123-129.

[17] 覃文忠，王建梅，劉妙龍,等.混合地理加權回歸模型算法研究[J].武漢大學學報(信息科學版)，2007，32(2):115-119.

[18] 趙陽陽，劉紀平，張福浩,等.貪心算法的地理加權回歸特征變量選擇方法[J].測繪科學，2016，41(7):41-46.

[19] 劉蓓.非線性主成分分析法在巖體質量等級分類中的應用[J].水電能源科學,2011,29(12)：78-80.

[20] 高青松,薛付忠.非線性主成分分析中數據非線性特征的檢驗方法[J].中國衛生統計,2011,28(5):488-491,494.

[21] 檀菲菲,陸兆華.基于NLPCA-GSO可持續發展評價——以環渤海區域為例[J].生態學報,2016,36(8):2403-2412.

[22] 周永正,袁曉輝,周勇.基于非線性主成分神經網絡水泥強度預測研究[J].數學的實踐與認識,2013,43(3):83-91.

[23] 孫康,耿修瑞,唐海蓉,等.一種基于非線性主成分分析的高光譜圖像目標檢測方法[J].測繪通報,2015(1):105-108.

[24] FOTHERINGHAM A S，CHARLTON M，BRUNSDON C.Measuring Spatial Variations in Relationships with Geogr Aphically Weighted Regression[M].[S.l.]：Springer，1997：60-63.

[25] CLEVELAND W S.Robust Locally Weighted Regression and Smoothing Scatterplots[J].Journal of the American Statistical Association，1979，74(368):829-836.