自定义error、北京污染物数据处理、可视化

最新推荐文章于 2024-08-04 17:54:27 发布

Damon Tong

最新推荐文章于 2024-08-04 17:54:27 发布

阅读量937

点赞数 1

分类专栏：现代程设文章标签：数据挖掘数据分析

本文链接：https://blog.csdn.net/qq_46658082/article/details/127394559

版权

现代程设专栏收录该内容

13 篇文章 2 订阅

订阅专栏

在这里插入图片描述
数据资源已上传在资源库

自定义error、北京污染物数据处理、可视化

README

原始数据说明：
各类污染物都是微克/立方米，但是在计算AQI时CO是使用毫克的
原始数据的格式分布：
['No', 'year', 'month', 'day', 'hour', 'PM2.5', 'PM10', 'SO2', 'NO2', 'CO', 'O3', 'TEMP', 'PRES', 'DEWP', 'RAIN', 'wd', 'WSPM', 'station']
分别为各类污染物，温度、气压、露点温度、雨、风向、风力
Region类的数据格式
{地区:[时间,污染物数据]}
Province类的数据格式关于时间的输入
["year-month-day-hour","year-month-day-hour"]
Province类返回的数据
{region:{time:,各类数值，AQI的值}
region包含的地区
['Aotizhongxin','Changping','Dingling','Dongsi','Guanyuan','Gucheng','Huairou','Nongzhanguan','Shunyi','Tiantan','Wanliu','Wanshouxigong']
有几个地方地图上貌似找不到，用区来试一试,区也不太好使，有几个地方一个区orz。那么我就对这个区的AQI取平均。
关于AQI的计算方式连接如下：
https://blog.csdn.net/m0_58032776/article/details/117912575

cor_relation返回的是每个地区的污染属性与天气情况的相关系数的dataframe构成的字典


关于Visualization模块，我设计了一个Visual类，关于visual类有以下函数：
pie 饼图，对每种属性在这个时间段内的均值做一个饼图？有待考量，貌似不是那么的合理
line 线图 对于某地某污染物的数值做线图
river 河流图 对各地的污染物变化水平做河流图,暂时不太想做，不急
correlation 每个地区的相关系数热力图
map 对每个地区的污染物分布情况做一个地图
给污染情况做一个权值的处理，然后在地图上展示？
计算AQI，然后对每个时刻的AQI做一个均值，然后将这个值在北京地图上表示出来。
其实还可以做一个主成分分析关于各污染物的分类，或者做k-means分析，画个碎石图之类的，但是那个用datafram格式的数据就可以直接画，本身数据处理的时候，用了字典来分类。如果重新读入补充异常值以后的csv文件似乎又没什么必要，下面我把方法链接贴在下面：
https://blog.csdn.net/qq_25990967/article/details/121366143
https://blog.csdn.net/wyn1564464568/article/details/125898241 sklearn主成分分析

关于data输出的值分别为某种污染物的所有地区数值字典，北京某时间段的各地污染物数值字典，北京各地的污染物相关系数字典

数据分析类 Analaysis

import csv
import openpyxl
from Error import *
import csv
import datetime
import numpy as np
import pandas as pd

#一开始不如用pandas库里的read_csv，直接对dataframe类型进行处理
class Region:#该类分析的是某区域的污染物数据
    def __init__(self,region,pollutant):
        self.region = region
        self.pollutant = pollutant

    def read(self):
        with open("./PRSA_Data_20130301-20170228/{:}.csv".format(self.region), "r") as f:
            file = csv.DictReader(f)  # 转为字典的格式列表，每一行表示一个字典
            file=list(file)#很多类型都不支持一些调用，先转为列表比较方便
            data=[]
            if self.pollutant in file[0].keys():
                for row in file:
                    time=datetime.datetime(year=int(row["year"]),month=int(row['month']),day=int(row['day']),hour=int(row['hour']))
                    li=[time,row[self.pollutant]]
                    data.append(li)#deposit data
            else:
                print("you enter the wrong pollutant")
        return(data)#返回一个列表包含了时间和污染物的数据


class Province:
    def __init__(self,*args):#输入一个时间或者时间段
        self.args = args

    def __AQI(self,PM2_5, PM10, SO2, NO2, CO, O3):  # 注意在拼盘IQA时，CO单位用的是mg
        """此为1h的AQI计算代码，24h的只需要更换Idata即可"""
        qua = [0, 50, 100, 150, 200, 300, 400, 500]
        data = [PM2_5, PM10, SO2, NO2, CO / 1000, O3]
        IAQI = list(np.zeros(6))
        Idata = [
            [0, 35, 75, 115, 150, 250, 350, 500],  # PM2.5 1小时平均
            [0, 50, 150, 250, 350, 420, 500, 600],  # （PM10）1小时平均
            [0, 150, 500, 650, 800],  # so2
            [0, 100, 200, 700, 1200, 2340, 3090, 3840],  # no2
            [0, 5, 10, 35, 60, 90, 120, 150],  # co
            [0, 160, 200, 300, 400, 800, 1000, 1200]  # o3
        ]
        for i in range(len(Idata)):
            POL = data[i]  # 提取第i种污染物的数值
            # 寻找区间
            for j in range(len(Idata[i])):
                if Idata[i][j] >= POL:
                    break  # 此时找到污染物的上界,且记录上界的位置，此时对应分值表的上界
            iqa = round((qua[j] - qua[j - 1]) / (Idata[i][j] - Idata[i][j - 1]) * (POL - Idata[i][j - 1]) + qua[j - 1])
            IAQI[i] = iqa
        IAQI = np.array(IAQI)
        AQI = max(IAQI)
        return (AQI)

    def Time_Region(self):
        region=['Aotizhongxin','Changping','Dingling','Dongsi','Guanyuan','Gucheng','Huairou','Nongzhanguan','Shunyi','Tiantan','Wanliu','Wanshouxigong']
        region_dic={}
        for place in region:
            with open("./PRSA_Data_20130301-20170228/{:}.csv".format(place),"r") as f:
                file = csv.DictReader(f)  # 转为字典的格式列表，每一行表示一个字典
                file = list(file)  # 很多类型都不支持一些调用，先转为列表比较方便
                shuju=[]#该地区所有数据
                if len(self.args)==1:#表示时间点
                    ar=self.args#注意调用的时候是self.args
                    args=ar[0].split('-')
                    for row in file:
                        content={}#每一个数据存在一个字典里
                        if row['year']==args[0] and row['month']==args[1] and row['day']==args[2] and row['hour']==args[3]:#匹配时间点
                            content['time']=datetime.datetime(year=int(row["year"]),month=int(row['month']),day=int(row['day']),hour=int(row['hour']))
                            content['PM2.5']=row['PM2.5']
                            content['PM10']=row['PM10']
                            content['SO2']=row['SO2']
                            content['NO2']=row['NO2']
                            content['CO']=row['CO']
                            content['O3']=row['O3']
                            content['TEMP']=row['TEMP']
                            content['PRES']=row['PRES']
                            content['DEWP']=row['DEWP']
                            content['WSPM']=row['WSPM']
                            content['RAIN']=row['RAIN']
                            content['wd']=row['wd']
                            AQI=self.__AQI(float(row['PM2.5']),float(row['PM10']),float(row['SO2']),float(row['NO2']),float(row['CO']),float(row['O3']))
                            content['AQI']=AQI
                            shuju.append(content)#插入每条数据
                elif len(self.args) == 2:  # 表示时间段
                    start=self.args[0].split("-")#直接转数字计算
                    end=self.args[1].split("-")
                    start=int(start[0])*1000000+int(start[1])*10000+int(start[2])*100+int(start[3])
                    end=int(end[0])*1000000+int(end[1])*10000+int(end[2])*100+int(end[3])
                    for row in file:
                        content ={}
                        #以2013 03 10 0为例
                        time=int(row['year'])*1000000+int(row['month'])*10000+int(row['day'])*100+int(row['hour'])
                        if time<=end and time>=start:
                            content['time'] = datetime.datetime(year=int(row["year"]), month=int(row['month']), day=int(row['day']),
                                                                hour=int(row['hour']))
                            content['PM2.5'] = row['PM2.5']
                            content['PM10'] = row['PM10']
                            content['SO2'] = row['SO2']
                            content['NO2'] = row['NO2']
                            content['CO'] = row['CO']
                            content['O3'] = row['O3']
                            content['TEMP'] = row['TEMP']
                            content['PRES'] = row['PRES']
                            content['DEWP'] = row['DEWP']
                            content['WSPM'] = row['WSPM']
                            content['RAIN'] = row['RAIN']
                            content['wd'] = row['wd']
                            AQI = self.__AQI(float(row['PM2.5']), float(row['PM10']), float(row['SO2']),
                                             float(row['NO2']), float(row['CO']), float(row['O3']))
                            content['AQI'] = AQI
                            shuju.append(content)#放在里面是为了防止加入空的字典
                else:
                    print("wrong number of time")
            region_dic[place]=shuju
        return(region_dic)#返回各个地区在这个时间段的数据

class Cor:
    def __init__(self):
        self.region=['Aotizhongxin','Changping','Dingling','Dongsi','Guanyuan','Gucheng','Huairou','Nongzhanguan','Shunyi','Tiantan','Wanliu','Wanshouxigong']

    def analysis(self):
        relation_cor={}
        for i in self.region:
            data = pd.read_csv(f"./PRSA_Data_20130301-20170228/{i:}.csv")
            new = data.iloc[:, [i for i in range(5, 17)]]#提取属性所在的列的序号，利用iloc返回部分数据框
            matrix = new.corr("spearman")
            relation_cor[i]=matrix
        return relation_cor#返回各地的构成,字典形式展示

可视化类 Visualization

#-*- coding: utf-8 -*
from Error import *
import matplotlib.pyplot as plt
import seaborn as sns
import numpy as np
import pandas as pd
from seaborn.matrix import heatmap  #绘制热力图
from pyecharts.charts import Map  #pyecharts绘制map,pie,themeriver,treemap
from pyecharts import options as opts
from pyecharts.charts import ThemeRiver
from pyecharts.charts import Geo, Pie, Timeline
from pyecharts.charts import TreeMap
from pyecharts import options as opts
from pyecharts.globals import ThemeType #引入主题
class visual:
    def __init__(self,pollutant):
        self.region = ['Aotizhongxin', 'Changping', 'Dingling', 'Dongsi', 'Guanyuan', 'Gucheng', 'Huairou', 'Nongzhanguan',
                  'Shunyi', 'Tiantan', 'Wanliu', 'Wanshouxigong']
        self.region_CH=['朝阳区','昌平区','昌平区','东城区','西城区','石景山区','怀柔区','朝阳区','顺义区','东城区','海淀区','宣武区']
        self.pollutant = pollutant
        print("let's do it")

    def pie(self,region_dic):
        sum_li=[]
        name=[]
        for reg in self.region:
            data = region_dic[reg]
            da = list(zip(*data))  # 解码，分成时间序列和数值,对于tuple是无法做索引的
            shuzhi = np.array([float(i) for i in da[1]])#转为tuple类型才能求和
            sum_li.append(sum(shuzhi))
            name.append(reg)#以防提取名字的时候出问题，因为字典是不按顺序的
        pie=Pie()
        pie.add("",[list(z) for z in zip(name,sum_li)])
        pie.set_series_opts(label_opts=opts.LabelOpts(formatter="{b}:{c} {d}%"))#设置格式
        pie.set_global_opts(title_opts=opts.TitleOpts(title=f"各地{self.pollutant:}占比"))
        pie.render(f"./{self.pollutant:} pie.html")

    def line(self,region_dic):#传入的是北京各地的某污染物的线图
        plt.figure(figsize=(24,18))
        i=1
        for reg in self.region:
            data=region_dic[reg]
            da=list(zip(*data))#解码，分成时间序列和数值,对于tuple是无法做索引的
            time=da[0]
            shuzhi=[float(i) for i in da[1]]
            plt.subplot(3,4,i)
            plt.plot(time,shuzhi)
            plt.title(f"{reg:}\'s {self.pollutant:}")
            plt.xticks(fontsize=6)
            i+=1
        plt.savefig(f"./beijing regions of {self.pollutant:} line")
        plt.show()#这个放保存后面不然存不进去,不知道为啥


    def river(self):
        pass

    def map(self,map_dic):
        region=self.region
        data={}
        for j in range(len(region)):
            da=map_dic[region[j]]
            med=0
            for i in range(len(da)):
                med+=float(da[i]['AQI'])
            mean=med/len(da)
            if data.get(self.region_CH[j],0)!=0:#有些区有重复就很烦
                data[self.region_CH[j]]=(mean+data.get(self.region_CH[j],0))/2
            else:
                data[self.region_CH[j]]=mean
        data=list(zip(data.keys(),data.values()))
        map=Map(init_opts=opts.InitOpts(width='900px',height='800px',theme=ThemeType.DARK))
        map.add('AQI',data,maptype='北京')
        map.set_global_opts(title_opts=opts.TitleOpts(title='北京各地AQI'),visualmap_opts=opts.VisualMapOpts(max_=150,min_=90,split_number=6,range_text="AQI颜色区间"))
        map.render("./Map.html")


    def correlation(self,C_Matrix):#传入一个dataframe的字典列表
        key=list(C_Matrix.keys())
        plt.figure(figsize=(24,18))
        for i in range(len(C_Matrix)):
            plt.subplot(3,4,i+1)
            heatmap(C_Matrix[key[i]],annot=False)
            plt.title(f"{key[i]:}heatmap of correlation",fontsize=10)
            cax = plt.gcf().axes[-1]
            cax.tick_params(labelsize=5)  # colorbar刻度字体大小
            plt.tick_params(labelsize=5)#设置刻度大小
        plt.savefig("./北京各地热力图.png")
        plt.show()

检查缺失值 examination

def examine():#检查是否有空值
    region = ['Aotizhongxin', 'Changping', 'Dingling', 'Dongsi', 'Guanyuan', 'Gucheng', 'Huairou', 'Nongzhanguan',
              'Shunyi', 'Tiantan', 'Wanliu', 'Wanshouxigong']
    for item in region:
        data = pd.read_csv(f"./PRSA_Data_20130301-20170228/PRSA_Data_{item:}_20130301-20170228.csv")
        name = data.columns.values.tolist()  # 列名称
        x, y = np.where(data.isnull())
        x = list(x)
        y = list(y)
        if len(x) > 0:  # 说明有异常值
            try:
                pollutant = name[y[0]]
                year = data['year'][x[0]]
                region = "Aotizhongxin"
                month = data['month'][x[0]]
                day = data['day'][x[0]]
                hour = data['hour'][x[0]]
                raise NotNumError_new(region, year, month, day, hour,pollutant)  # 之前使用NotNumError的时候，一直报错，说不显示属性，可能存在重名等问题，因此改了个名字就成功了

            except NotNumError_new as e:  # 认为捕获错误，同时对类实例化
                print(e.message)
            # #下列函数块能够把所有的na的位置打印出来，但是如果全部打印出来会很乱，因此我只打印第一个na
            # for i in range(len(x)):
            #     try:
            #         pollutant=name[y[i]]
            #         year=data['year'][x[i]]
            #         region="Aotizhongxin"
            #         month=data['month'][x[i]]
            #         day=data['day'][x[i]]
            #         hour=data['hour'][x[i]]
            #         raise NotNumError_new(region, year, month, day,hour,pollutant)#之前使用NotNumError的时候，一直报错，说不显示属性，可能存在重名等问题，因此改了个名字就成功了
            #
            #     except NotNumError_new as e:#认为捕获错误，同时对类实例化
            #         print(e.message)

        data.fillna(method='pad', inplace=True)#利用前后值补足,对于还没补全的直接变成0，此时变成0的已经很少了
        data.fillna('0',inplace=True)#此时完美解决空值问题
        try:
            data.to_csv(f"./PRSA_Data_20130301-20170228/{item:}.csv", index=False,
                        sep=',')  # 在调试的时候，由于该函数不支持覆写，因此采用try，当之前已经完成了na值补全的文件时，则跳过报错，继续运行下一个
        except:
            continue

数据分析实现 Data

def data(pollutant):#展示数据
    river={}
    place=['Aotizhongxin','Changping','Dingling','Dongsi','Guanyuan','Gucheng','Huairou','Nongzhanguan','Shunyi','Tiantan','Wanliu','Wanshouxigong']
    for reg in place:
        B_R=Region(reg,pollutant)
        river[reg]=B_R.read()#将所有数值保存在字典中
    Ana_P=Province('2013-3-1-0','2015-3-1-0')
    Ana_C=Cor()
    Ana_P_data=Ana_P.Time_Region()
    Ana_C_data=Ana_C.analysis()
    print(river)
    print(Ana_P_data)
    print(Ana_C_data)
    return river,Ana_P_data,Ana_C_data

数据可视化实现 Vis

def Vis(R,P,C,pollutant):#可视化
    Graph=visual(pollutant)
    Graph.correlation(C)#相关系数热力图
    Graph.line(R)
    Graph.pie(R)
    Graph.map(P)#传入map的数据

主函数

def main():
    #examine()
    R,P,C=data("SO2")#以二氧化硫为例，我们希望输出的是单个地区，所有地区的
    Vis(R,P,C,"SO2")

if __name__ == '__main__':
    main()

北京各地的SO2含量关于时间的波动

可以发现SO2的含量波动成周期性波动，每年的一月份左右处于峰值，每年的七月份左右处于峰谷。下面通过饼图来反映各地的SO2含量占比和大小。

在这里插入图片描述

北京各地so2含量占比饼图

右下方的饼图可以看出，定陵、怀柔、顺义地区的SO2含量较少，可以看出，相比于市中心，偏远地区、景区的SO2含量较少，至于这些地区的环境是否更好还有待考量。
在这里插入图片描述

北京各地污染物，气象环境相关系数热力图

通过分析相关系数热力图，可以发现PM10和PM2.5强正相关、各类污染物除了O3之间都成正相关，O3与其他污染物成负相关，其中NO2与O3成强负相关（NO2与O3会出现热化学反应，NO2中的可能会对O3起催化作用）。温度、露点温度强正相关，但与大气压强成负相关，温度越高，气压越小符合物理规律。是否下雨也许可以影响PM值，但似乎影响不大。NO2与风速有比较强的负相关性。

在这里插入图片描述

北京各地AQI数值分布地图

其实对于各地的AQI分析，由于我已经将数值分析出来了，其实还可以按照风向来画图，看看AQI数值是否跟风向有关，如果说是北风，也许南边的污染程度更严重，如果是南风，也许北边的污染程度更严重。

那么根据2013——2017年的AQI均值

我们可以发现海淀等在北京市中心的AQI值比周边地区的AQI高一些，也就是说，通过均值的统计，市中心地区的污染比较严重。北京市长期处于轻度污染的环境中。

那么在污染物治理的时候，就可以着重关注污染情况较为严重的地区，分析其局部气象环境，污染物含量等进行分析、治理。

在这里插入图片描述

某地SO2的数据分析

[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-6ZUKRw5b-1666090992906)(C:\Users\kerrla\AppData\Roaming\Typora\typora-user-images\image-20221018002007020.png)]

某时间段各地的污染物气象数据，包括每小时的AQI

[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-Zqd87rbY-1666090992907)(C:\Users\kerrla\AppData\Roaming\Typora\typora-user-images\image-20221018002242137.png)]

热力图的相关系数数据框

[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-OkTNvxwV-1666090992908)(C:\Users\kerrla\AppData\Roaming\Typora\typora-user-images\image-20221018002333489.png)]

缺失值位置报错

在这里插入图片描述

生成的新文件

[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-kenEruya-1666090992910)(C:\Users\kerrla\AppData\Roaming\Typora\typora-user-images\image-20221018002548851.png)]
那么在污染物治理的时候，就可以着重关注污染情况较为严重的地区，分析其局部气象环境，污染物含量等进行分析、治理。

[外链图片转存中…(img-PRyfeFFs-1666090992906)]