一.項目背景

AQI（air Quality Index）指空氣質(zhì)量指數(shù)，用來衡量空氣清潔或者污染程度。值
越小，表示空氣質(zhì)量越好。近年來因為環(huán)境問題，空氣質(zhì)量越來越受到人們重視。

上篇文章[AQI分析與預(yù)測（一）](http://mp.weixin.qq.com/s?__biz=MzIzODI4ODM2MA==&mid=2247486525&idx=1&sn=a92d40e0d715f9ab63563b4e8b054a0d&chksm=e93ae0bade4d69ac57407216484cdbb6d7a8a812baef1c0e4ab0d5aeae0c7f0aff9fd274127b&scene=21#wechat_redirect)我們進(jìn)行了初步分析，主要分析了空氣質(zhì)量最好/差城市和臨海城市是否空氣質(zhì)量優(yōu)于內(nèi)陸城市這兩個問題，本篇我們在之前基礎(chǔ)上繼續(xù)研究如下問題。

1.空氣質(zhì)量受那些因素影響
2.關(guān)于空氣質(zhì)量驗證
3.構(gòu)建空氣質(zhì)量預(yù)測模型

二.實現(xiàn)過程

1.空氣質(zhì)量受那些因素影響
  》》指標(biāo)：協(xié)方差和相關(guān)系數(shù)
  》》圖形：熱力圖

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
import warnings


sns.set(style="darkgrid", font_scale=1.2)
plt.rcParams["font.family"] = "SimHei"
plt.rcParams["axes.unicode_minus"] = False
warnings.filterwarnings("ignore")


#讀取文件
data = pd.read_csv("data.csv")
#查看數(shù)據(jù)
data.head()

#kind：繪制圖像的類型。可選值：
#scatter：散點圖（默認(rèn)值）。
#reg：帶有回歸線的散點圖。
#vars：顯示哪些變量之間的兩兩關(guān)系，默認(rèn)為顯示所有變量。
sns.pairplot(data, vars=["AQI", "PopulationDensity", "GreenCoverageRate"])

#計算相關(guān)數(shù)據(jù)
x = data["AQI"]
y = data["Precipitation"]
#計算AQI與Precipitation的協(xié)方差。
a = (x - x.mean()) * (y - y.mean())
#計算協(xié)方差
cov = np.sum(a) / (len(a) - 1)
print("協(xié)方差：", cov)

#計算AQI與Precipitation的相關(guān)系數(shù)。
corr = cov / np.sqrt(x.var() * y.var())
print("相關(guān)系數(shù)：", corr)

#pandas封裝了相關(guān)方法
print("協(xié)方差：", x.cov(y))
print("相關(guān)系數(shù)：", x.corr(y))

#初始化畫布
plt.figure(figsize=(15, 10))
#繪制熱力圖
ax = sns.heatmap(data.corr(),
                 cmap=plt.cm.RdYlGn, 
                 annot=True,
                 fmt=".2f")
  
 結(jié)論：降雨量越多，空氣質(zhì)量越好；維度越低，空氣質(zhì)量越好

2.關(guān)于空氣質(zhì)量驗證
  》》問題：全國所有城市的空氣質(zhì)量指數(shù)均值在71左右，請問此結(jié)論準(zhǔn)確嗎？
  》》方法：假設(shè)檢驗

該需求是驗證樣本均值是否等于總體均值，根據(jù)條件，我們可以使用單樣本t檢驗，
置信度為95%。

#進(jìn)行單樣本t檢驗
r = stats.ttest_1samp(data["AQI"], 71)
#輸出檢驗統(tǒng)計量
print("t值：", r.statistic)
#輸出p值
print("p值：", r.pvalue)


結(jié)論：我們可以看到P值是大于0.05的，y因此我們無法拒絕原假設(shè)，因此接受原假設(shè)

#計算均值
mean = data["AQI"].mean()
#計算標(biāo)準(zhǔn)差
std = data["AQI"].std()
#計算置信區(qū)間
stats.t.interval(0.95, df=len(data) - 1, loc=mean, scale=std / np.sqrt(len(data)))


結(jié)論：我們就計算出全國城市平均空氣質(zhì)量指數(shù)，95%的可能大致在70.63~80.04之間

3.構(gòu)建空氣質(zhì)量預(yù)測模型
  》》對于一些城市，如果能夠已知降雨量，溫度，經(jīng)緯度等指標(biāo)，我們是否能夠預(yù)測該
城市的空氣質(zhì)量指數(shù)呢？因此我們需要構(gòu)建模型，預(yù)測新的數(shù)據(jù)。
  》》過程：基模型構(gòu)建，異常值處理后構(gòu)建模型和特征選擇后構(gòu)建模型，將結(jié)果與基模
型進(jìn)行對比，看看是否進(jìn)行優(yōu)化。

#進(jìn)行類別轉(zhuǎn)換
data["Coastal"] = data["Coastal"].map({"是": 1, "否": 0})
#統(tǒng)計類別數(shù)目
data["Coastal"].value_counts()

from sklearn.linear_model import LinearRegression
from sklearn.model_selection import train_test_split


#City(城市名稱)對預(yù)測毫無用處，刪掉。
X = data.drop(["City","AQI"], axis=1)
#目標(biāo)值
y = data["AQI"]
#分離測試集和訓(xùn)練集
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.3, random_state=0)
#構(gòu)建基礎(chǔ)線性回歸模型
lr = LinearRegression()
#訓(xùn)練模型
lr.fit(X_train, y_train)
#輸出訓(xùn)練集模型評分
print(lr.score(X_train, y_train))
#輸出測試集模型評分
print(lr.score(X_test, y_test))

#繪制圖形觀察


#預(yù)測y值
y_hat = lr.predict(X_test)
#初始化畫布
plt.figure(figsize=(15, 5))
#繪制真實值折線圖
plt.plot(y_test.values,
         "-r", 
         label="真實值",
         marker="o")
#繪制預(yù)測值折線圖
plt.plot(y_hat,
         "-g",
         label="預(yù)測值",
         marker="D")
#設(shè)置圖例
plt.legend(loc="upper left")
#設(shè)置標(biāo)題
plt.title("線性回歸預(yù)測結(jié)果", fontsize=20)

# Coastal是類別變量，映射為離散變量，不會有異常值。


#遍歷列
for col in X.columns.drop("Coastal"):
    #對數(shù)值型數(shù)據(jù)進(jìn)行判斷
    if pd.api.types.is_numeric_dtype(X_train[col]):
        #獲取分位數(shù)
        quartile = np.quantile(X_train[col], [0.25, 0.75])
        #計算IQR
        IQR = quartile[1] - quartile[0]
        #計算正常數(shù)值下限
        lower = quartile[0] - 1.5 * IQR
        #計算正常數(shù)值上限
        upper = quartile[1] + 1.5 * IQR
        #用邊界值進(jìn)行填充異常值
        X_train[col][X_train[col] < lower] = lower
        X_train[col][X_train[col] > upper] = upper
        X_test[col][X_test[col] < lower] = lower
        X_test[col][X_test[col] > upper] = upper

#訓(xùn)練模型
lr.fit(X_train, y_train)
#去除異常值后評估模型效果
print(lr.score(X_train, y_train))
print(lr.score(X_test, y_test))


結(jié)論：去除異常值后所構(gòu)建模型效果相比之前有所改進(jìn)

#對數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化處理，


from sklearn.preprocessing import StandardScaler
#數(shù)據(jù)標(biāo)準(zhǔn)化
s = StandardScaler()
#對訓(xùn)練集進(jìn)行標(biāo)準(zhǔn)化
X_train_scale = s.fit_transform(X_train)
#對測試集進(jìn)行標(biāo)準(zhǔn)化
X_test_scale = s.transform(X_test)

#對數(shù)據(jù)進(jìn)行特征選擇，目的是提高模型準(zhǔn)確率和訓(xùn)練速度


REFCV方法
  》》RFE（REcursive feature elimination）：遞歸特征消除，用來對特征進(jìn)行重要性評級
  》》CV（Cross Validation）：交叉驗證，通過交叉驗證，選擇最佳數(shù)量特征
具體過程如下：
RFE階段：
  1.初始特征集為所有可用特征
  2.使用當(dāng)前特征集進(jìn)行建模，然后計算每個特征的重要性
  3.刪除最不重要的一個或者多個特征，更新特征集
  4.跳到步驟2，直到完成所有特征集重要性評級
CV階段
  1.根據(jù)REF階段確定的特征重要性，依次選擇不同數(shù)量特征
  2.對選定的特征集進(jìn)行交叉驗證
  3.確定平均分最高的特征數(shù)量，完成特征選擇。

from sklearn.feature_selection import RFECV
#estimator：要操作的模型。
#step：每次刪除的變量數(shù)。
#cv：使用的交叉驗證折數(shù)。
#n_jobs：并發(fā)的數(shù)量。
#scoring: 評估的方式。
rfecv = RFECV(estimator=lr, 
              step=1,
              cv=5, 
              n_jobs=-1, 
              scoring="r2")
rfecv.fit(X_train_scale, y_train)
#返回經(jīng)過選擇之后，剩余的特征數(shù)量。
print(rfecv.n_features_)
#返回經(jīng)過特征選擇后，使用縮減特征訓(xùn)練后的模型。
print(rfecv.estimator_)
#返回每個特征的等級，數(shù)值越小，特征越重要。
print(rfecv.ranking_)
#返回布爾數(shù)組，用來表示特征是否被選擇。
print(rfecv.support_)
#返回對應(yīng)數(shù)量特征時，模型交叉驗證的評分。
print(rfecv.grid_scores_)

#繪制圖形
plt.plot(range(1, len(rfecv.grid_scores_) + 1), 
         rfecv.grid_scores_,
         marker="o")
#設(shè)置x軸標(biāo)簽
plt.xlabel("特征數(shù)量")
#設(shè)置y軸標(biāo)簽
plt.ylabel("交叉驗證$R^2$值")

print("剔除的變量：", X.columns.values[~rfecv.support_])
#應(yīng)用到訓(xùn)練集
X_train_eli = rfecv.transform(X_train_scale)
#應(yīng)用到測試集
X_test_eli = rfecv.transform(X_test_scale)
#輸出模型評分
print(rfecv.estimator_.score(X_train_eli, y_train))
print(rfecv.estimator_.score(X_test_eli, y_test))

#獲取列名與對應(yīng)的權(quán)重，構(gòu)成一個元組，作為列表的元素。
li = list(zip(X.columns.values[rfecv.support_], rfecv.estimator_.coef_))
#根據(jù)權(quán)重的絕對值，對列表進(jìn)行降序排列。
li.sort(key=lambda x: abs(x[1]), reverse=True)
#轉(zhuǎn)換為Series
s = pd.Series(dict(li))
#初始化畫布
plt.figure(figsize=(15, 5))
#繪制柱狀圖
ax = sns.barplot(y=s.index, x=s.values)
for y, x in enumerate(s):
    #繪制標(biāo)注
    t = ax.text(x / 2, y, round(x, 3))
    #設(shè)置居中對齊
    t.set_ha("center")
#顯示圖形
plt.show()

from sklearn.preprocessing import KBinsDiscretizer
#KBinsDiscretizer K個分箱的離散器。用于將數(shù)值（通常是連續(xù)變量）變量進(jìn)行區(qū)間離散化操作。
#n_bins：分箱（區(qū)間）的個數(shù)。
#encode：離散化編碼方式。分為：onehot，onehot-dense與ordinal。
#   onehot：使用獨熱編碼，返回稀疏矩陣。
#   onehot-dense：使用獨熱編碼，返回稠密矩陣。
#   ordinal：使用序數(shù)編碼（0,1,2……）。
#strategy：分箱的方式。分為：uniform，quantile，kmeans。
#uniform：每個區(qū)間的長度范圍大致相同。
#quantile：每個區(qū)間包含的元素個數(shù)大致相同。
#kmeans：使用一維kmeans方式進(jìn)行分箱。


#對數(shù)據(jù)進(jìn)行分箱操作
k=KBinsDiscretizer(n_bins=[4, 5, 10, 6],
                     encode="onehot-dense", 
                     strategy="uniform")
#定義離散化的特征。
discretize=["Longitude", "Temperature", "Precipitation", "Latitude"]
#訓(xùn)練集數(shù)據(jù)轉(zhuǎn)換為DataFrame
X_train_eli=pd.DataFrame(data=X_train_eli, columns=X.columns[rfecv.support_])
#測試集數(shù)據(jù)轉(zhuǎn)換為DataFrame
X_test_eli=pd.DataFrame(data=X_test_eli, columns=X.columns[rfecv.support_])
#應(yīng)用到訓(xùn)練集
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)
#將離散化后的特征與其他特征進(jìn)行重新組合。
X_train_dis=pd.concat([X_train_dis, r], axis=1)
#對測試集進(jìn)行同樣的離散化操作。
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)
#查看轉(zhuǎn)換之后的格式。
display(X_train_dis.head()

#訓(xùn)練模型
lr.fit(X_train_dis, y_train)
#去除異常值后評估模型效果
print(lr.score(X_train_dis, y_train))
print(lr.score(X_test_dis, y_test))


結(jié)論：離散化后模型效果進(jìn)一步提升