城市空气质量分析与预测
城市空氣質量分析與預測
- 一、AQI分析與預測
- 1、背景信息
- 2、任務說明
- 3、數據集描述
- 二、數據分析流程
- 基本流程
- 三、讀取數據
- 1、導入相關的庫
- 2、加載數據集
- 四、數據清洗
- 1、缺失值
- 1.1、缺失值探索
- 1.2、缺失值處理
- 1、數據分布
- 2、填充數據
- 2、異常值
- 2.1、異常值探索
- 1、describe方法
- 2、3σ的方式
- 3、箱線圖
- 2.2、異常值處理
- 1、對數轉換
- 2、使用邊界值替換
- 3、分箱離散化
- 3、重復值
- 3.1、重復值的探索
- 3.2、重復值處理
- 五、數據分析
- 1、空氣質量最好/最差的5個城市
- 1.1、空氣質量最好的5個城市
- 1.2、空氣質量最差的5個城市
- 2、全國城市的空氣質量
- 2.1、城市空氣質量等級統計
- 2.2、空氣質量分布
- 3、臨海城市是否空氣質量優于內陸城市
- 3.1、數量統計
- 3.2、分布統計
- 3.3、差異檢驗
- 4、空氣質量主要受哪些因素影響
- 4.1、散點圖矩陣
- 4.2、相關系數
- 1、協方差定義
- 2、相關系數定義
- 4.3、結果統計
- 5、關于空氣質量的消息真假驗證
- 6、對空氣質量的預測
- 6.1、數據轉換
- 6.2、基模型
- 6.3、特征選擇
- 1、RFECV
- 特征選擇實例
- 7、異常值的處理
- 7.1、使用臨界值替換
- 7.2、效果對比
- 7.3、分箱操作
- 8、殘差圖分析
- 8.1、異方差性
- 8.2、離群點
- 六、總結
一、AQI分析與預測
1、背景信息
AQI指的是空氣質量指數，用來衡量一個城市的空氣清潔或污染的程度，數值越小則空氣質量越好。近年來，空氣污染問題備受關注，現收集不同城市的數據，運用數據分析的方法來對不同城市的空氣質量進行分析與預測。
2、任務說明
我們期望能夠運用數據分析的相關技術，對全國城市的空氣質量進行分析與預測，同時會根據分析結果來解決一些常見問題，例如：
- 哪些城市的空氣質量比較好？哪些比較差？
- 空氣質量在地理上的分布，是否存在著一定的規律？
- 沿海城市與內陸城市的空氣質量是否存在不同？
- 空氣質量主要受到哪些因素的影響？
- 全國空氣質量的整體情況是怎樣的？
- 怎么來預測一個城市的空氣質量？
3、數據集描述
我們現在獲取了2015年的空氣質量數據集，該數據集包含的是全國主要城市的相關數據及空氣質量指數，數據情況如下：
| City | 城市名 |
| AQI | 空氣質量指數 |
| Precipitation | 降雨量 |
| GDP | 城市生產總值 |
| Temperature | 溫度 |
| Longitude | 緯度 |
| Latitude | 經度 |
| Altitude | 海拔高度 |
| PopulationDensity | 人口密度 |
| Coastal | 是否沿海 |
| GreenCoverageRate | 綠化覆蓋率 |
| Incineration(10000ton) | 焚燒量(10000噸) |
二、數據分析流程
基本流程
在進行數據分析之前，我們清楚整個分析流程，每個流程所對應的步驟與使用的相關分析技術：
- 明確任務需求與分析目標
- 數據收集：
- 內部數據、購買數據、爬取數據、調查問卷、其他收集
- 數據預處理：
- 數據整合（橫向整合、縱向整合）
- 數據清洗（缺失值、異常值、重復值）
- 數據轉換
- 數據分析：
- 描述性分析
- 推斷分析
- 數據建模（特征方程、超參數調整）
- 數據可視化
- 編寫報告
三、讀取數據
1、導入相關的庫
導入相關的庫，同時進行一些初始化的設置
import numpy as np import pandas as pd import matplotlib.pyplot as plt import seaborn as sns import warningssns.set(style='darkgrid') plt.rcParams['font.family']='SimHei' plt.rcParams['axes.unicode_minus']=False warnings.filterwarnings('ignore')2、加載數據集
data=pd.read_csv('data.csv') print(data.shape) （325，12）加載后，可以使用head/tail/sample等方法查看數據的大致情況：
data.head()四、數據清洗
1、缺失值
1.1、缺失值探索
我們可以用以下方法來查看缺失值：
- info
- is null
1.2、缺失值處理
對于缺失值，我們可以使用以下處理方式：
- 刪除缺失值：僅適合于缺失值很少的情況
- 填充缺失值：
- 數值變量：
- 均值填充
- 中值填充
- 類別變量：
- 眾數填充
- 單獨作為一個類別
- 額外處理說明：
- 缺失值小于20%，直接填充
- 缺失值在20%-80%，填充變量后，同時增加一列，標記該列是否缺失，參與后續建模
- 缺失值大于80%，不用原始列，而是增加一列，標記該列是否缺失，參與后續建模
- 數值變量：
1、數據分布
print(data['Precipitation'].skew()) sns.distplot(data['Precipitation'].dropna())2、填充數據
data.fillna({'Precipitation':data['Precipitation'].median()},inplace=True) data.isnull().sum()2、異常值
2.1、異常值探索
我們可以如下方法來發現缺失值：
- describe查看數值信息
- 正態分布3個標準差外(3σ)的方式
- 使用箱線圖輔助
- 相關異常檢測算法
1、describe方法
調用DataFrame對象的describe方法，可以顯示數據的統計信息，但僅為一種簡單的異常探索方式，不夠直觀，局限性強。
2、3σ的方式
根據正態分布的特性，我們可以將3個標準差外（3σ）的數據，視為異常值。現在我們以GDP為例，首先先繪制GDP的分布圖
sns.distplot(data['GDP']) print(data['GDP'].skew())
由圖看出數據呈現嚴重右偏分布，也就是存在很多極大異常值，我們先獲取這些異常值。
3、箱線圖
箱線圖也是一種常見的異常值檢測方式。
sns.boxplot(data=data['GDP'])2.2、異常值處理
對于異常值，我們可以用以下的方式來處理：
- 刪除異常值
- 視為缺失值
- 對數轉換
- 使用臨界值填充
- 使用分箱法離散化處理
1、對數轉換
如果數據中存在較大的異常值，可以通過采用對數來進行轉化，能得到一定的緩解。例如GDP變量呈現右偏分布
import matplotlib.pyplot as plt fig,ax = plt.subplots(1,2) fig.set_size_inches(15,5) sns.distplot(data['GDP'],ax=ax[0]) sns.distplot(np.log(data['GDP']),ax=ax[1])
取對數的方式比較簡單，不過也存在著局限性，比如：
-
對數轉換值針對正數操作，不過可以采用以下方式進行轉換：
*np.sign(X)np.log(np.abs(X)+1)
-
適合右偏分布，不適合左偏分布，不過現實中大多數的數據一般的呈右偏分布較多。
2、使用邊界值替換
對異常值進行‘截斷’處理，即使用臨界值替換異常值。例如在3σ與箱線圖可以這樣處理
3、分箱離散化
有時候特征對目標值存在一定的影響，可能未必是線性的增加，這時可以使用分箱方式，對特征進行離散化處理。
分箱離散化分為兩個階段：
- 分箱：比如將不同年齡段的人進行分化區間，例如0-18歲、18歲-30歲、30歲-40歲、40歲以上。
- 離散化：可以將不同年齡段依次對應為1、2、3、4，作為離散化處理。最后輸入模型中進行分析
3、重復值
3.1、重復值的探索
使用duplicate檢查重復值，可配合keep參數進行調整
#發現重復值 print(data.duplicated().sum()) #查看哪些記錄出現重復值 data[data.duplicated(keep=False)]3.2、重復值處理
重復值對數據分析通常沒有作用，直接刪除即可
data.drop_duplicates(inplace=True) data.duplicated().sum()五、數據分析
1、空氣質量最好/最差的5個城市
1.1、空氣質量最好的5個城市
t=data[['City','AQI']].sort_values('AQI') display(t.iloc[:5]) plt.xticks(rotation=45) sns.barplot(x='City',y='AQI',data=t.iloc[:5])1.2、空氣質量最差的5個城市
display(t.iloc[:-5]) plt.xticks(rotation=45) sns.barplot(x='City',y='AQI',data=t.iloc[-5:])2、全國城市的空氣質量
2.1、城市空氣質量等級統計
對于AQI，可以對空氣質量進行等級劃分，劃分標準如下：
| 0-50 | 一級 | 優 |
| 51-100 | 二級 | 良 |
| 101-150 | 三級 | 輕度污染 |
| 151-200 | 四級 | 中度污染 |
| 201-300 | 五級 | 重度污染 |
| >300 | 六級 | 嚴重污染 |
根據該標準，我們來統計下，全國空氣質量每個等級的數量。
# 編寫函數，將AQI轉換為對應的等級。 def value_to_level(AQI):if AQI>=0 and AQI<=50:return '一級'elif AQI>=51 and AQI<=100:return '二級'elif AQI>=101 and AQI<=150:return '三級'elif AQI>=151 and AQI<=200:return '四級'elif AQI>=201 and AQI<=300:return '五級'else:return '六級' level = data['AQI'].apply(value_to_level) display(level.value_counts()) sns.countplot(x=level,order=['一級','二級','三級','四級','五級'])2.2、空氣質量分布
我們來繪制下全國各城市的空氣質量分布圖。
sns.scatterplot(x='Longitude',y='Latitude',hue='AQI',palette=plt.cm.RdYlGn_r,data=data)
從結果上可以發現，從大致的地理位置，西部城市好于東部城市，南部城市好于北部城市。
3、臨海城市是否空氣質量優于內陸城市
3.1、數量統計
我們先來統計下臨海城市與內陸城市的數量
display(data['Coastal'].value_counts()) sns.countplot(x='Coastal',data=data)3.2、分布統計
然后我們來觀察臨海城市與內陸城市的散點分布
sns.swarmplot(x='Coastal',y='AQI',data=data)
然后我們來分組計算空氣質量的均值:
在柱狀圖中，僅顯示了內陸城市與沿海城市的空氣質量指數（AQI）的均值對比，我們也可以用箱線圖來顯示更多的信息。(上圖均值的豎線是均值的置信區間)
3.3、差異檢驗
這里，我們可以進行兩樣本的t檢驗，來查看臨海城市與內陸城市的均值差異是否顯著
from scipy import stats Coastal = data[data['Coastal']=='是']['AQI'] inland = data[data['Coastal']=='否']['AQI'] #進行方差齊性檢驗，為后續的兩樣本的t檢驗服務。 stats.levene(Coastal,inland)LeveneResult(statistic=0.08825036641952543, pvalue=0.7666054880248168)
這里原假設方差是相等的也就是齊性的，檢驗結果看P值，P值就是對原假設的支持概率，這里pvalue=0.76，則說明兩獨立樣本是齊性的。
Ttest_indResult(statistic=-2.7303827520948905, pvalue=0.006675422541012958)
0.9966622887294936
這里的P值很小則說明，原假設不成立即沿海城市與內陸城市均值不相等。同時怎么判斷兩個樣本均值的大小，看statistic，這里是負值說明，沿海城市均值小于內陸城市均值。
4、空氣質量主要受哪些因素影響
比如我們可能會關注這些問題：
- 人口密度大，是否會對空氣質量造成負面影響？
- 綠化率高，是否會提高空氣質量？
4.1、散點圖矩陣
sns.pairplot(data[['AQI','PopulationDensity','GreenCoverageRate']]) #sns.pairplot(data[['AQI','PopulationDensity','GreenCoverageRate']],kind='reg')4.2、相關系數
1、協方差定義
2、相關系數定義
我們以空氣質量（AQI）與降雨量（Precipitation）為例，計算二者的相關系數：
協方差： -10098.209013903044
相關系數： -0.40184407003013883
不過DataFrame對象提供了計算相關系數的方法。我們可以直接使用
data.corr()為了能夠更清晰的呈現相關系數值，我們可以用熱圖來顯示相關系數
plt.figure(figsize=(15,10)) ax = sns.heatmap(data.corr(),cmap=plt.cm.RdYlGn,annot=True,fmt='.2f') #注意：Matplotlib3.1.1版本的bug，heatmap的首行與末行會顯示不全。 #可手動調整y軸的范圍來進行修復。 a,b=ax.get_ylim() ax.set_ylim(a+0.5,b-0.5)4.3、結果統計
從結果中可知，空氣質量主要受降雨量（-0.4）與緯度（0.55）影響
- 降雨量越多，空氣質量越好。
- 緯度越低，空氣質越好
此外，我們還能發現其他的明顯細節：
- GDP (城市生產總值) 與Incineration (焚燒量)正相關(0.90)。
- Temperature (溫度)與Precipitation (降雨量)正相關(0.69)。
- Temperature (溫度)與Latitude (緯度)負相關(-0.81)。
- Longitude (經度)與Altitude (海拔)負相關(-0.74)。
- Latitude (緯度)與Precipitation (降雨量)負相關(-0.66)。
- Temperature (溫度)與Altitude (海拔)負相關(-0.46)。
- Altitude (海拔)與Precipitation (降雨量)負相關(-0.32)。
5、關于空氣質量的消息真假驗證
網上有傳聞說，全國所有城市的空氣質量指數均值在71左右，請問，這個消息可靠嘛？
我們可以采用假設檢驗來驗證。
首先，城市的平均空氣質量指數，我們肯容易就能進行計算。
75.3343653250774
該需求是要驗證樣本均值是否等于總體均值，根據場景，我們可以使用單樣本t檢驗，置信度為95%
r = stats.ttest_1samp(data['AQI'],71) print('t值：',r.statistic) print('p值：',r.pvalue)t值： 1.8117630617496872
p值： 0.07095431526986647
我們可以看到，P值大于0.05，故我們無法拒絕原假設，因此接受原假設。
我們要清楚，+ -1.96倍的標準差，是正態分布在置信區間為95%下的臨界值，嚴格來說，對 t 分布不是如此，
只不過，當樣本容量較大時，t分布近似正態分布。但當樣本容量較小時，二者會有較大的差異。
我們可以獲取更準確的置信區間：
由此，我們就計算出全國所有城市平均空氣質量指數，95%的可能大致在70.63-80.04之間。
6、對空氣質量的預測
對于某城市，如果我們已知降雨量，溫度，經緯度等指標。我們是否能夠預測該城市的空氣質量指數呢？答案是肯定的。我們可以通過對以往的數據，去建立一種模式，然后將這種模式去應用于未知的數據，進而預測結果。
6.1、數據轉換
對于模型來說，內部進行的都是數學上的運算，在進行建模之前，我們需要首先進行轉換，將類別變量轉換為離散變量。
data['Coastal'] = data['Coastal'].map({'是':1,'否':0}) data['Coastal'].value_counts()6.2、基模型
首先，我們不進行任何處理，建立一個模型，后續的操作，可以在此基礎上進行改進。
from sklearn.linear_model import LinearRegression from sklearn.model_selection import train_test_splitX = data.drop(['City','AQI'],axis=1) y = data['AQI'] X_train,X_test,y_train,y_test = train_test_split(X,y,test_size=0.3,random_state=0)lr=LinearRegression() lr.fit(X_train,y_train) print(lr.score(X_train,y_train)) print(lr.score(X_test,y_test))0.4538897765064036
0.40407705623835266
6.3、特征選擇
剛才，我們使用所有可能的原始數據作為特征，建立模型，然而，特征并非越多越好，有些特征可能對模型質量并沒有什么改善，能夠提高模型訓練速度。
1、RFECV
特征選擇方式有很多，這里我們用RFECV方法來實現特征選擇。RFECV分為兩個部分，分別是：
- RFE：遞歸特征消除，用來對特征進行重要性評級。
- CV：交叉檢驗，在特征評級后，通過交叉檢驗，選擇最佳數量的特征。
具體過程如下： - RFE階段：
1、初始的特征集為所有可用的特征。
2、使用當前特征集進行建模，然后計算每個特征的重要性。
3、刪除最不重要的一個（或多個）特征，更新特征集。
4、跳轉到步驟2，直到完成所有特征的重要性評級。 - CV階段：
1、根據RFE階段確定的特征重要性，依次選擇不用數量的特征。
2、對選定的特征集進行交叉檢驗。
3、確定平均分最高的特征集數量，完成特征選擇。
特征選擇實例
from sklearn.feature_selection import RFECV #estimator:要操作的模型 #step:每次刪除的變量數 #cv:使用的交叉檢驗折數 #n_jobs:并發的數量 #scoring:評估的方式 rfecv = RFECV(estimator=lr,step=1,cv=5,n_jobs=-1,scoring='r2') rfecv.fit(X_train,y_train) #返回經過選擇之后，剩余的特征數量 print(rfecv.n_features_) #返回經過特征選擇后，使用縮減特征后的模型 print(rfecv.estimator) #返回每個特征的等級，數值越小，特征越重要 print(rfecv.ranking_) #返回布爾數組，用來表示特征是否被選擇。 print(rfecv.support_) #返回對應數量特征時，模型交叉檢驗的評分 print(rfecv.grid_scores_)8
LinearRegression(copy_X=True, fit_intercept=True, n_jobs=None, normalize=False)
[1 2 1 1 1 1 3 1 1 1]
[ True False True True True True False True True True]
[-0.06091362 0.1397744 0.2302962 0.22814855 0.22644355 0.21342713
0.24573222 0.26368987 0.25744818 0.25389817]
通過結果可知，我們成功刪除了兩個特征。可以繪制下，在特征選擇過程中，使用交叉檢驗獲取的R平方值。
plt.plot(range(1,len(rfecv.grid_scores_)+1),rfecv.grid_scores_,marker='o') plt.xlabel('特征數量') plt.ylabel('交叉檢驗$R^2$值')
然后，我們對測試集應用這種特征選擇（變換），進行測試，獲取測試結果。
剔除的變量： Index([‘GDP’, ‘PopulationDensity’], dtype=‘object’)
0.45322255556356406
0.39822556329575864
結論：我們發現，經過特征選擇后，消除了GDP與PopulationDensity兩個特征，而使用剩余8個特征訓練的模型，與之前未消除特征訓練的模型（使用全部10個特征訓練的模型），無論在訓練集還是測試集的表現上，都幾乎相同，這就可以證明，我們清楚的這兩個特征，確實對擬合目標（y值）沒有什么幫助，可以去掉。
7、異常值的處理
如果數據中存在異常值，可能會影響模型的效果，因此，我們在建模之前，有必要對異常值進行處理
7.1、使用臨界值替換
我們可以依據箱線圖判斷離群點的原則去探索異常值，然后使用臨界值，替換掉異常值。
#Coastal是類別變量，映射為離散變量，不會有異常值 for col in X.columns.drop('Coastal'):if pd.api.types.is_numeric_dtype(X_train[col]):quartile = np.quantile(X_train[col],[0.25,0.75])IQR = quartile[1]-quartile[0]lower = quartile[0]-1.5*IQRupper = quartile[1]+1.5*IQRX_train[col][X_train[col]<lower]=lowerX_train[col][X_train[col]>upper]=upperX_test[col][X_test[col]<lower]=lowerX_test[col][X_test[col]>upper]=upper7.2、效果對比
去除異常值后，我們使用新的訓練集與測試集來評估模型的效果
lr.fit(X_train,y_train) print(lr.score(X_train,y_train)) print(lr.score(X_test,y_test))0.4631142291492417
0.446142026583966
效果相對之前，似乎有著輕微的改進，不過并不明顯，我們可以使用RFECV在去除異常值的數據上，再次嘗試
rfecv = RFECV(estimator=lr,step=1,cv=5,n_jobs=-1,scoring='r2') rfecv.fit(X_train,y_train) print(rfecv.n_features_) print(rfecv.estimator_) print(rfecv.ranking_) print(rfecv.support_) print(rfecv.grid_scores_)plt.plot(range(1,len(rfecv.grid_scores_)+1),rfecv.grid_scores_,marker='o') plt.xlabel('特征數量') plt.ylabel('交叉檢驗$R^2$值') print('剔除的變量：',X_train.columns[~rfecv.support_]) # X_train_eli = rfecv.transform(X_train) # X_test_eli = rfecv.transform(X_test) X_train_eli = X_train[X_train.columns[rfecv.support_]] X_test_eli = X_test[X_test.columns[rfecv.support_]] print(rfecv.estimator_.score(X_train_eli,y_train)) print(rfecv.estimator_.score(X_test_eli,y_test))剔除的變量： Index([‘PopulationDensity’], dtype=‘object’)
0.46306656191488615
0.4450225589408254
7.3、分箱操作
注意：分箱后，我們不能將每個區間映射為離散數值，而是應當使用One-Hot編碼。
from sklearn.preprocessing import KBinsDiscretizer #KBinsDiscretizer K個分箱離散器。用于將數值（通常是連續變量）變量進行區間離散化操作。 #n_bins:分箱（區間）的個數 #encode:離散化編碼方式。分為：onehot,onehot-dense與ordinal # onehot:使用獨熱編碼，返回稀疏矩陣 # onehot-dense:使用獨熱編碼，返回稠密矩陣 # ordinal:使用序數編碼（0,1,2...） #strategy:分箱的方式。分為：uniform,quantile,kmeans。 # uniform:每個區間長度范圍大致相同 # quantile:每個區間包含的元素大致相同 # kmeans:使用一維kmeans方式進行分箱。 k = KBinsDiscretizer(n_bins=[4,5,14,6],encode='onehot-dense',strategy='uniform') #定義離散化的特征。 discretize = ['Longitude','Temperature','Precipitation','Latitude']r = k.fit_transform(X_train_eli[discretize]) r = pd.DataFrame(r,index=X_train_eli.index) #獲取離散化后特征之外的其他特征 X_train_dis = X_train_eli.drop(discretize,axis =1) #獲取離散化后的特征與其他特征進行重新組合 X_train_dis = pd.concat([X_train_dis,r],axis =1) #對測試集進行同樣的離散化操作。 r = pd.DataFrame(k.transform(X_test_eli[discretize]),index=X_test_eli.index) X_test_dis = X_test_eli.drop(discretize,axis =1) X_test_dis = pd.concat([X_test_dis,r],axis =1) #查看轉換后的格式 display(X_train_dis.head())
這樣，我們就可以對轉換后的數據進行訓練了。
0.6892388692774563
0.6546062348355675
8、殘差圖分析
殘差就是模型預測值與真實值之間的差異，我們可以繪制殘差圖，來對回歸模型進行評估。殘差圖的橫坐標為預測值，縱坐標為殘差值。
8.1、異方差性
對于一個好的回歸模型，誤差應該是隨機分布的，因此殘差也應隨機分布于中心線附近。如果我們從殘差圖中找到規律，這意味著模型遺漏了某些能夠影響殘差的解釋信息。
異方差性，是指殘差具有明顯的方差不一致性，這里我們異常值處理前后的兩組數據，分別訓練模型，然后觀察殘差的效果。
import matplotlib.pyplot as plt fig,ax = plt.subplots(1,2) fig.set_size_inches(15,5) data = [X_train,X_train_dis] title = ['原始數據','處理后數據'] for d,a,t in zip(data,ax,title):model = LinearRegression()model.fit(d,y_train)y_hat_train = model.predict(d)residual = y_hat_train - y_train.valuesa.set_xlabel('預測值')a.set_ylabel('殘差')a.axhline(y=0,color='red')a.set_title(t)sns.scatterplot(x=y_hat_train,y=residual,ax=a)
在右圖中，我們發現隨著預測值的增大，模型的誤差也在增大，對于此情況，我們可以使用對目標y值取對數的方式處理。
此時，異方差性得到解決。同時，模型的效果也可能會得到一定的提升。
8.2、離群點
如果是簡單線性回歸，我們通過繪制回歸線，就可以輕松的看出是否存在一些離散點。然而，對于多元線回歸，其回歸線
已經擴展成為超平面，我們無法通過可視化來進行預測。
然而，我們可以通過繪制殘差圖，通過預測值與實際值之間的關系，來檢測離群點。
0.7354141753913532
0.6302724058812208
最后訓練模型整體效果比之前好很多。
六、總結
1.空氣質量總體分布上來說，南部城市優于北部城市，西部城市優于東部城市。
2.臨海城市的空氣質量整體上好于內陸城市。
3.是否臨海，降雨量與緯度對空氣質量指數的影響較大。
4.我國城市平均空氣質量指數大致在(70.63 - 80.04)這個區間內，在該區間的可能性概率為95%。
5.通過歷史數據，我們可以對空氣質量指數進行預測。
總結
以上是生活随笔為你收集整理的城市空气质量分析与预测的全部內容，希望文章能夠幫你解決所遇到的問題。
- 上一篇： Win11连接投影仪没反应怎么解决？
- 下一篇： 微信支付-免充值立减与折扣(上篇)