AQI分析与预测

最新推荐文章于 2023-06-23 19:29:10 发布

Danah.F

最新推荐文章于 2023-06-23 19:29:10 发布

阅读量2.4k

点赞数 1

文章标签：空气质量数据分析预测模型线性回归临海城市

本文链接：https://blog.csdn.net/qq_39223444/article/details/119204758

版权

1. 背景

AQl(Ar Qualty Index)，指空气质量指数，用来衡量空气清洁或污染的程度。值越小，表示空气质量越好。近年来，因为环境问题，空气质量越来越受到人们的重视。

2. 任务说明与知识要点

我们期望能够运用数据分析的相关技术，对全国域市空气质量进行研究与分析，希望能够解决如下疑问：

哪些城市的空气质量较好/较差?【描述性统计分析】
空气质量在地理位置分布上，是否具有一定的规律性?【描述性统计分析】。
临海域市的空气质量是否有别于内陆城市?【[推断统计分析】
空气质量主要受哪些因素影响?【相关系数分析】
全国城市空气质量普遍处于何种水平?【区间估计】
怎样预测一个域市的空气质量?【统计建模】

3. 数据集描述

在这里插入图片描述

4. 导入相关库

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import warnings
warnings.filterwarnings("ignore")
import seaborn as sns
# sns.set(style="darkgrid",font="SimHei",rc={"axes.unicode_minus":False})
from pylab import mpl
mpl.rcParams['font.sans-serif'] = ["SimHei"] # 中文字体
mpl.rcParams['axes.unicode_minus'] = False # 坐标轴负号

5. 加载数据集

加载相关的数据集。
可以使用head / tail / sample查看数据的大致情况。

data = pd.read_csv("F:\数据集\AQI\Dataset.csv")
print(data.shape)
print(data.head())# head取前5条，tail取后5条，sample取随机5条
data.sample(5)

6. 数据清洗

6.1 缺失值处理

我们可以使用如下方法查看缺失值:

info
isnull

# data.info()
data.isnull().sum(axis=0)
print(data.isnull().sum(axis=0))
print(data["Precipitation"].skew()) #skew()查找索引轴上数据的偏度
sns.distplot(data["Precipitation"].dropna()) #删除缺失值再显示

# 用中值替换
index = data[data["Precipitation"].isnull()].index
data["Precipitation"][index] = data["Precipitation"].median()
data.fillna({"Precipitation":data["Precipitation"].median()},inplace = True)

6.2 异常值处理

通过describe直看数值信息
可配合箱线图辅助。

print(data.describe())
# 箱线图
plt.figure(figsize=(15,10))
plt.xticks(rotation = 45,fontsize = 15)
sns.boxplot(data=data)

有异常值时的处理

t=data.copy()
for k in t:
	if pd.api.types.is_nunmeric_dtype(t[k]):
		o=t[k].describe()
		IQR=o["75%"]- o["25%"]
		lower = o["25%"]-1.5*IQR
		upper =o["75%"]+1.5*IQR
		t[k][t[k] < lover] = lower
		t[k][t[k] > upper] = upper
plt.figure(figsize=(15,4))
plt.xticks(rotation=45,fontsize=15)
sns.boxplot(data=t)

6.3 重复值处理

使用duplicate检查重复值。可配合keep参数进行调整
使用drop_duplicate删除重复值。

# 发现重复值
data.duplicated().sum()
print("重复值：",data.duplicated().sum())
# 查看重复值
data[data.duplicated()]
# 删除重复值
data.drop_duplicates(inplace=True)

7. 数据分析

7.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])

在这里插入图片描述

plt.figure()
display(t.iloc[-5:])
plt.xticks(rotation=45)
sns.barplot(x="City",y="AQI",data=t.iloc[-5:])

在这里插入图片描述

7.2 全国城市的空气质量

7.2.1 城市空气质量等级统计

国家对空气质量进行等级划分，划分标准如下表所示:

在这里插入图片描述
根据该标准，我们来统计下，全国空气质量每个等级的数量。

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 "五级"
    elif AQI>=300:
        return "六级"
level = data["AQI"].apply(value_to_level)
display(level.value_counts())
plt.figure()
sns.countplot(x=level,order=["一级","二级","三级","四级","五级","六级"])

在这里插入图片描述
可见，我们城市的空气质量主要以一级(优）与二级(良)为主，三级(轻度污染)占一部分，更高污染的城市占少数。

7.2.2 空气质量指数分布

plt.figure()
sns.scatterplot(x="Longitude",y="Latitude",hue="AQI",palette=plt.cm.RdYlGn_r,data = data)

在这里插入图片描述
从结果我们可以发现，从大致的地理位置上看，西部城市好于东部城市，南部城市好于北部城市。

7.2.3 关于空气质量的验证

江湖传闻，全国所有城市的空气质量指数均值在72左右，请问，这个消息可靠吗?
城市平均空气质量指数，我们可以很容易的进行计算。

print(data["AQI"].mean())

75.3343653250774

7.2.3.1 中心极限定理

在这里插入图片描述

7.2.3.2 置信区间

接下来，我们现在根据正态分布的特性，进行概率上的统计，如下图:
在这里插入图片描述

根据中心极限定理，如果多次抽样，则样本均值构成的正态分布。如果我们对总体进行一次抽样，则本次抽样个体的均值有95%的概率会在二倍标准差内，仅有5%的概率会在二倍标准差外。根据小概率事件(很小的概率在一次抽样中基本不会发生），如果抽样的个体均值落在二倍标准差之外，我们就可以认为，本次抽样来自的总体，该总体的均值并非是我们所期望的均值。

通常，我们以二倍标准差作为判定依据，则以均值为中心，正负二倍标准差构成的区间，就是置信区间。而二倍标准差包舍了95%的数据，因此，此时的置信度为95%。换言之，我们有信心认为，总体的均值有95%的概率会在置信区间之内。

7.2.3.3 假设检验——t检验

假设检验，其目的是通过收集到的数据，来验证某个假设是否成立。在假设检验中，我们会建立两个完全对立的假设，分别为原假设（零假设） $H_o$ 与备则假设（对立假设) $H_1$ 。然后根据样本信息进行分析判断，得出P值（概率值)。

假设检验基于小概率反证法，即我们认为小概率事件在一次试验中是不会发生的。如果小概率事件发生，则我们就拒绝原假设，而接受备择假设。否则，我们就没有充分的理由推翻原假设，此时，我们选择去接受原假设。

t检验，就是假设检验的一种，可以用来检验一次抽样中样本均值与总体均值的比较(二者差异是否显著)。其计算方式如下:
在这里插入图片描述

from scipy import stats
r = stats.ttest_1samp(data["AQI"],72)
print("t值：",r.statistic)
print("p值：",r.pvalue)

t值： 1.393763441074581
p值： 0.16435019471704654

我们可以看到，P值大于0.05，故在显著度水平为0.05检验下，我们无法拒绝原假设，因此接受原假设。同样，我们现在可以来计算下，全国所有城市平均空气质量指数的置信区间。

n = len(data)
df = n-1
left = stats.t.ppf(0.025,df=df)
right = stats.t.ppf(0.975,df=df)
print(left,right)
mean = data["AQI"].mean()
std = data["AQI"].std()
print(mean+left*(std/np.sqrt(n)),mean+right*(std/np.sqrt(n)))

70.6277615675309 80.0409690826239
得到全国所有城市平均空气质量指数所在的置信区间大致在70~80之间，置信度为95%

7.2.4 临海域市的空气质量是否有别于内陆城市?

# 统计临海城市和内陆城市的数量
display(data["Coastal"].value_counts())
plt.figure()
sns.countplot(x="Coastal",data=data)

在这里插入图片描述

# 观察临海城市和内陆城市的散点分布
plt.figure()
sns.swarmplot(x="Coastal",y="AQI",data = data)

在这里插入图片描述

# 分组计算均值
display(data.groupby("Coastal")["AQI"].mean())
plt.figure()
sns.barplot(x="Coastal",y="AQI",data = data)

在这里插入图片描述
在柱形图中，仅显示了内陆城市与临海城市空气质量指数（AQl )的均值对比，我们可以使用箱线图来显示更多的信息.

# 箱线图
plt.figure()
sns.boxplot(x="Coastal",y="AQI",data = data)

在这里插入图片描述

# 小提琴图，除了能够展示箱线图的信息外，还可以显示分布密度
plt.figure()
sns.violinplot(x="Coastal",y="AQI",data = data)

在这里插入图片描述

# 将散点图与箱线图或小提琴图结合在一起
plt.figure()
sns.violinplot(x="Coastal",y="AQI",data = data,inner=None)
sns.swarmplot(x="Coastal",y="AQI",color = "g",data = data)

在这里插入图片描述

# 对两样本进行t检验，查看沿海城市与内陆城市的均值差异是否显著
coastal = data[data["Coastal"]==1]["AQI"]
inland = data[data["Coastal"]==0]["AQI"]
# 进行方差齐性检验，为后续的两样本t检验服务
print(stats.levene(coastal,inland))
#result:LeveneResult(statistic=0.08825036641952543, pvalue=0.7666054880248168)
# 如果返回结果的p值远大于0.05，那么我们认为两总体具有方差齐性。
# 方差相等,equal_var=True
# 进行两样本t检验
print(stats.ttest_ind(coastal,inland,equal_var=True))
# result:Ttest_indResult(statistic=-2.7303827520948905, pvalue=0.006675422541012958)
# 两者不相等，说明内陆好于沿海

7.2.5 空气质量主要受哪些因素影响

我们很可能会关注某些问题，例如，我们可能会产生类似如下的疑问:

人口密度大，是否会对空气质量造成负面影响?
绿化本高，是否会提高空气质量?

7.2.5.1 绘制散点图矩阵

通过散点图矩阵，可以显示任意两个变量之间的散点图，我们可以通过散点图规察两个变量之间的关系。

sns.pairplot(data[["AQI","PopulationDensity","GreenCoverageRate"]])

在这里插入图片描述

7.2.5.2 相关系数

相关系数，可以用来体现两个连续变量之间的相关性，最为常用的为皮尔逊相关系数。其定义公式为：
在这里插入图片描述
我们以空气质量(AQI)与阵雨量( Precipitation )为例，计算二者的相关系数,

x = data["AQI"]
y = data["PopulationDensity"]
# 计算协方差
a = (x-x.mean())*(y-y.mean())
cov = np.sum(a)/(len(a)-1)
print(cov)
# 计算相关系数
corr = cov/np.sqrt(x.var()*y.var())
print(corr)

# 自带公式
print(x.cov(y))
print(x.corr(y))

-3319.45192585043
-0.026495735026362967

# 为了更清晰的呈现相关系数数值，使用热图来展示相关系数
print(data.corr())
plt.figure(figsize=(15,12))
sns.heatmap(data.corr(),cmap=plt.cm.RdYlGn,annot=True,fmt=".2f")

在这里插入图片描述

7.2.5.3 结果统计

从结果中可知，空气质星指数主要受阵雨量(-0.40 )与纬度(0.55)影响。

降雨量越多，空气质量越好。
纬度越低，空气质量越好。

此外，我们还能够发现其他一些明显的细节:

GDP(城市生产总值与Incineration(焚烧量）正相关(0.90 )
Temperature(温度)与Precipitation(阵雨量)正相关(0.69 )
Temperature(温度)与Latitude(纬度)负相关( -0.81 ).
Longitude(经度）与Aitude (海拔)负相关( -0.74 )
Latitude(纬度)与Precipitation (阵雨量)负相关(-0.66 )
Temperature(温度）与Atitude (海拔)负相关(-0.46 )
Altitude(海拔)与Precipitation (降雨量)负相关( -0.32 )

7.2.6 对空气质量指数进行预测

对于某域市，如果我们已知降雨量，温度，经纬度等指标，我们是否能够预测该城市的空气质量指数呢?

7.2.6.1 线性回归

答案是肯定的。我们可以通过对以往的数据，去建立一种模式，然后将这种模式去应用于未知的数据，进而预测结果。

from sklearn.linear_model import LinearRegression
from sklearn.model_selection import train_test_split
X = 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)
y_hat = lr.predict(X_test)
print(lr.score(X_train,y_train))
print(lr.score(X_test,y_test))

0.45369066602176233
0.4039626944128524

plt.figure(figsize=(15,10))
plt.plot(y_test.values,"-r",label="真实值",marker="o")
plt.plot(y_hat,"-g",label="预测值",marker="D")
plt.legend(loc="upper left")
plt.title("线性回归预测结果",fontsize=20)

在这里插入图片描述

7.2.6.2 残差分析

y_hat_train = lr.predict(X_train)
residual = y_hat_train-y_train.values
plt.figure()
plt.xlabel("预测值")
plt.ylabel("残差")
plt.axhline(y=0,color="red")
sns.scatterplot(x=y_hat_train,y=residual)

在这里插入图片描述

r = (residual-residual.mean())/residual.std()
print("异常值：",len(y_hat_train[np.abs(r)>2]))
plt.figure()
plt.xlabel("预测值")
plt.ylabel("残差")
plt.axhline(y=0,color="red")
sns.scatterplot(x=y_hat_train[np.abs(r)<=2],y=residual[np.abs(r)<=2],color="b",label="正常值")
sns.scatterplot(x=y_hat_train[np.abs(r)>2],y=residual[np.abs(r)>2],color="r",label="异常值")

在这里插入图片描述

X_train=X_train[np.abs(r)<=2]
y_train=y_train[np.abs(r)<=2]
lr.fit(X_train,y_train)
y_hat = lr.predict(X_test)
print(lr.score(X_train,y_train))
print(lr.score(X_test,y_test))

0.49663662284090554
0.45116945455183965

plt.figure(figsize=(15,10))
plt.plot(y_test.values,"-r",label="真实值",marker="o")
plt.plot(y_hat,"-g",label="预测值",marker="D")
plt.legend(loc="upper left")
plt.title("线性回归预测结果",fontsize=20)