中國各城市空氣質(zhì)量的聚類分析

1、中國各城市空氣質(zhì)量聚類和判別分析摘要 中國經(jīng)濟的快速增長導(dǎo)致環(huán)境污染不斷加重，其中空氣污染與人類的健康密切相關(guān)。結(jié)合全國74個城市的空氣污染物濃度數(shù)據(jù)對各城市進行聚類分析，分類方法包括等價關(guān)系法和Kmeans分析兩種方法。結(jié)果表明，?？谑侨珖諝赓|(zhì)量最好的城市，石家莊和邢臺是全國空氣質(zhì)量最差的城市，而武漢、成都和烏魯木齊的空氣質(zhì)量與北京最為接近。關(guān)鍵詞 聚類分析 空氣質(zhì)量 等價關(guān)系 Kmeans1. 介紹隨著中國經(jīng)濟的高速發(fā)展和工業(yè)化、城市化進程的加快，能源的消耗速度也不斷提高。中國的工業(yè)發(fā)展大量依賴煤炭、石油等化石燃料，其燃燒產(chǎn)生的廢氣嚴重污染空氣，導(dǎo)致中國各地區(qū)空氣質(zhì)量不斷下降。90年代初

2、期，中國的500個城市當中，達到國家空氣質(zhì)量I級標準的僅占1%；此外，近年來的數(shù)據(jù)顯示，暴露于未達標空氣中的城市人口占統(tǒng)計城市人口的三分之二1。城市的空氣污染對人體健康構(gòu)成極大威脅，研究表明，即便暴露于污染物密度較低的空氣中也會提高慢性呼吸系統(tǒng)的發(fā)病率以及多種癌癥的患病概率2。因此，有必要對全國各大城市的空氣質(zhì)量進行數(shù)據(jù)收集和分析，確定不同城市的污染程度及相互之間的關(guān)系，為相關(guān)部門制定政策提供有力的數(shù)據(jù)支撐。城市的空氣污染程度主要受經(jīng)濟發(fā)展水平影響，但二者不是呈簡單的倒U型曲線關(guān)系，不同的污染物與經(jīng)濟水平之間有不同的關(guān)系3，因此需要對各種污染物進行綜合分析和評價。而目前對環(huán)境進行綜合評價的方法

3、包括模糊數(shù)學(xué)法、距離判別法和物元分析法4，本文即采取其中的模糊數(shù)學(xué)法對全國74個主要城市的空氣污染數(shù)據(jù)進行聚類和判別分析，以研究目前中國各大城市的空氣污染水平和特點。2. 原始數(shù)據(jù)及聚類分析方法本文所用到的城市空氣污染數(shù)據(jù)來自環(huán)保部相關(guān)統(tǒng)計數(shù)據(jù)5，參見表1。為了便于分析，選取空氣污染指標中量綱相同的三個指標進行考察，分別是SO2濃度、NO2濃度和PM2.5濃度。采用兩種聚類分析方法對這74個城市進行分類，分別是等價關(guān)系法和Kmeans分類方法。表1 全國主要城市空氣污染物濃度城市SO2濃度NO2濃度PM2.5濃度北 京265689天 津595496石 家 莊10568154唐 山1146911

4、5秦 皇 島604765邯 鄲9558139邢 臺11369160保 定6956135張 家 口513240承 德373549滄 州5431102廊 坊4648110衡 水6846122太 原804381呼和浩特564057沈 陽904378大 連343152長 春444473哈 爾 濱445681上 海244862南 京375578無 錫404775徐 州524777常 州414872蘇 州315370南 通283672連 云 港343567淮 安343279鹽 城282965揚 州354370鎮(zhèn) 江304272泰 州262677宿 遷333574杭 州285370寧 波224454溫 州2

5、35158嘉 興304768湖 州295274紹 興384971金 華344170衢 州363768舟 山102233臺 州173453麗 水193249合 肥223988福 州114336廈 門204436南 昌404069濟 南9561110青 島584367鄭 州5952108武 漢336094長 沙334683廣 州205253深 圳114040珠 海133738佛 山325352江 門273351肇 慶283854惠 州162938東 莞234547中 山194248南 寧193857海 口71727重 慶323870成 都316396貴 陽313353昆 明284042拉 薩922

6、26西 安4657105蘭 州333567西 寧484170銀 川774351烏魯木齊2961882.1 等價關(guān)系法將每個城市的三種污染物濃度值作為其數(shù)據(jù)表征，根據(jù)最小最大法6建立74個城市之間的相似關(guān)系矩陣，得74×74矩陣，然后利用平方法求得相似矩陣的傳遞閉包，最后取不同的值將矩陣中的元素分為不同的類。相關(guān)程序代碼見附錄I，部分分類結(jié)果如表2所示。表2 等價關(guān)系法部分分類結(jié)果=0.0423所有城市為一類=0.4362海口北京，··· ，烏魯木齊=0.87石家莊 唐山邢臺=0.89石家莊 邢臺=0.87北京武漢成都烏魯木齊=0.89北京武漢成都烏魯木齊

7、=0.9北京武漢成都烏魯木齊=0.91北京武漢成都烏魯木齊上表中，位于同一單元格內(nèi)的城市為一類，位于同一行不同單元格內(nèi)的城市屬不同的類。從表中可看出，隨著值從0開始增加，?？跒榈谝粋€脫離大類的城市，這意味著?？诘奈廴境潭扰c其它城市的平均污染程度相差最遠。由于污染物濃度值越大意味著污染程度越嚴重，這里定義城市的污染指數(shù)為3種污染物濃度值的代數(shù)平均，則得海口市的污染指數(shù)為17，而其它城市的平均污染指數(shù)為52.5，因此可認為海口是全國空氣質(zhì)量最好的城市。另一方面，當值增大至0.9附近時，石家莊和邢臺依然分在同一類中，說明石家莊和邢臺的空氣污染程度極為接近，而邢臺的污染指數(shù)為114，位于所有城市之首，

8、因此石家莊和邢臺是全國空氣污染最嚴重的兩個城市。當在0.870.91之間時，北京、武漢、成都、烏魯木齊四個城市始終屬于同一類，說明武漢、成都、烏魯木齊的空氣污染狀況與北京比較接近。圖1給出了幾個典型城市的污染物濃度值與全國平均水平的比較情況，從圖中依然可以看出，石家莊、唐山、邢臺和邯鄲的污染程度明顯高于全國城市的平均水平，而北京的PM2.5指標值也高于全國平均水平，北京近年來的霧霾情況反映了這一點。圖1 幾個典型城市的污染物濃度值2.2. Kmeans分類方法Kmeans算法7是一種基于質(zhì)心簇的分類方法，首先在元素當中任意挑選k個元素作為初始質(zhì)心，根據(jù)歐氏距離確定與每個元素相距最近的質(zhì)心，并將

9、每個元素與最近的質(zhì)心歸為一類，接下來再次計算每一簇的質(zhì)心，再次重新歸類，以此循環(huán)直至質(zhì)心不再變化，算法可表示如下：1) 初始化-確定聚類數(shù)和聚類中心2) 計算所有樣本與聚類中心的距離3) 按照最小距離原則對各個樣本進行聚類4) 重新計算聚類中心5) 重復(fù)1-4，直到聚類中心不再改變由于初始質(zhì)心的數(shù)量由用戶自定義，所以最終分成的組數(shù)也是隨意的，表3給出了劃分為9組和5組的部分結(jié)果，相關(guān)程序代碼參見附錄。表3 Kmeans法部分分類結(jié)果承德大連上海寧波溫州臺州麗水舟山福州廈門深圳珠?；葜莺？谑仪f唐山邯鄲邢臺天津滄州廊坊鄭州西安北京哈爾濱南京無錫武漢成都烏魯木齊秦皇島張家口呼和浩特長春徐州青島西寧

10、保定衡水濟南南通連云港淮安鹽城揚州鎮(zhèn)江泰州太原沈陽銀川北京天津滄州廊坊哈爾濱南京成都武漢烏魯木齊石家莊唐山邯鄲邢臺保定衡水濟南秦皇島張家口太原呼和浩特沈陽青島銀川舟山臺州麗水福州廈門深圳珠海江門惠州承德大連長春上海無錫常州蘇州南通連云港從表中可看出，Kmeans方法分類結(jié)果中，北京、武漢、成都、烏魯木齊同樣屬于同一類，類似地，石家莊和邢臺也屬同類，這一結(jié)果與等價關(guān)系法形成互相印證。3. 結(jié)論1) 等價關(guān)系法和kmeans算法能夠有效地將樣本分類2) 聚類分析顯示海口是空氣最好的城市3) 石家莊、邯鄲、邢臺和唐山空氣質(zhì)量較差的城市4) 北京的空氣質(zhì)量同武漢、成都、烏魯木齊相近4. 人員分工數(shù)據(jù)的

11、收集和整理張文平、于鵬程序的編制任文希，彭熾，于鵬結(jié)果分析彭可文、張文平PPT制作彭熾，李冉，任文希報告編寫李冉、彭可文附錄I 等價關(guān)系法程序代碼function r_sol= Cluster(lamda,r_equ)% Written by Wensi Jen, 2016.03.22% 求截矩陣m n=size(r_equ);r_tru=0*ones(m);for i3=1:m for j3=1:m if r_equ(i3,j3)>=lamda r_tru(i3,j3)=1; else r_tru(i3,j3)=0; end endend% 聚類m1,n1=size(r_tru); f

12、or i3=1:m1 var4=0; for j3=1:n1 if(r_tru(i3,j3)=1)%找到每一行的元素1 var4=var4+1; r_find(i3,var4)=j3;%把每一行找到的元素1的列下標順次儲存在矩陣的每一行中 end end end m2,n2=size(r_find); for i4=1:m2 for j4=2:n2 if (r_find(i4,j4)=0) r_find(r_find(i4,j4),:)=0;%將相同行只保留一個其他變?yōu)榱阈?end endendvar5=0;for i5=1:m2 if (sum(r_find(i5,:)=0 var5=var

13、5+1; r_sol(var5,:)=r_find(i5,:); endendendm n=size(u);% 數(shù)據(jù)標準化極差變換for i=1:m for j=1:n u_nor(i,j)=(u(i,j)-min(u(:,j)/(max(u(:,j)-min(u(:,j); endend% 建立相似關(guān)系最大最小法for i1=1:m for j1=1:m var1=0; var2=0; for index=1:3 var1=var1+min(u_nor(i1,index),u_nor(j1,index); var2=var2+max(u_nor(i1,index),u_nor(j1,inde

14、x); end r(i1,j1)=var1/var2; endend% 改造相似關(guān)系為等價關(guān)系平方法求傳遞閉包num=floor(log2(m)+1;r_equ=r;for k2=1:num r=r_equ; for i2=1:m for j2=1:m for index2=1:m var3(index2,1)=min(r(i2,index2),r(index2,j2); end r_equ(i2,j2)=max(var3); end endend% 方法1利用等價關(guān)系聚類% 求lamda% r_equ=1 0.8 0.7;0.8 1 0.3;0.7 0.3 1m_1,n_1=size(r_e

15、qu);tot=m_1*n_1;%矩陣r_equ中總的元素個數(shù)r_row=r_equ(1:tot);%將矩陣r_equ轉(zhuǎn)化為一個行向量r_rows=sort(r_row);%將向量r_row中的元素按照從小到大排序for i_1=1:tot for j_1=i_1+1:tot if (r_rows(i_1)=r_rows(j_1) r_rows(j_1)=0; end endendr_rows1=sort(r_rows);%將r_rows中的元素按照從小到大排序var_1=0;for i_2=1:tot if (r_rows1(i_2)=0) var_1=var_1+1;%統(tǒng)計r_rows1中

16、0的個數(shù) endendr_rows2=r_rows1(tot:-1:var_1+1);%從r_rows2中挑出非零元素并按從大到小的順序排列s=length(r_rows2);for i_3=1:s lamda_row(i_3,1)=r_rows2(i_3); %lamda的列向量end% 直接聚類法% 根據(jù)lamda_row中的lamda值和等價矩陣r_equ，調(diào)用函數(shù)Cluster% 調(diào)用格式 Cluster(lamda,r_equ)Result=Cluster(0.7,r_equ);% 以城市名稱表示結(jié)果% Result=1 3 0;2 0 0;7 0 0m_2 n_2=size(Res

17、ult);for i_4=1:m_2 for j_4=1:n_2 if Result(i_4,j_4)=0 Result1i_4,j_4=cnResult(i_4,j_4); end endend附錄 Kmeans算法程序代碼K=9;%指定分為多少類Idx,C=kmeans(u,K);%Idx N*1的向量，存儲的是每個點的聚類標號m n=size(u);Result=0*ones(K,m);totr torc=size(Idx);for i=1:m for j=1:m if Result(Idx(i),j)=0 Result(Idx(i),j)=i; break end endend% 以城市名稱表示結(jié)果% Result=1 3 0;2 0 0;7 0 0m_2 n_2=size(Result);for i_4=1:m_2 for j_4=1:n_2 if Result(i_4,j_4)=0 Result1i_4,j_4=cnResult(i_4,j_4); end