2020-08-09

李宏毅机器学习hm1:利用regression预测丰原站pm2.5值

作业概述

输入:9个小时的数据,共18项特征(AMB_TEMP, CH4, CO, NHMC, NO, NO2, NOx, O3, PM10, PM2.5, RAINFALL, RH, SO2, THC, WD_HR, WIND_DIREC, WIND_SPEED, WS_HR)

输出:第10小时的PM2.5数值

模型:线性回归

import sys 
import numpy as np
import pandas as pd

data = pd.read_csv(r"E:\RaoLing\DataSet\李宏毅机器学习数据\hw1\train.csv", encoding="big5") # 读入train.csv,繁体字以big5编码 
data.head(18)

![6CHHJB{K@HE]62L)H{_0C.png](https://upload-images.jianshu.io/upload_images/20148371-91cd7892e51120c1.png?imageMogr2/auto-orient/strip%7CimageView2/2/w/1240)
第1列是日期,第2列是观测站所在地,第3列是观测指标,第4列-第27列是0-23共24小时。

data.shape
输入.png

预处理

可以看到降雨(rainfall)都是字符“NR”,将它变成数值0;从第3列开始是数值数据,提取出这些数值。

data.replace('NR', 0, inplace=True)      # 把降雨的NR字符变成数值0
raw_data = data.to_numpy()               # 把DataFrame转换成numpy array
raw_data
raw_data.png

提取特征(Extract Features)

ExtractFeatures1.png

ExtractFeatures2.png

参考上图,我们需要把数据shape做一个转换,把240天内相同的测量放在同一行,将 数组raw_data:(18 * 20 * 12, 24) => 转字典:{12:[18, 2420]} 实现分组*。

month_data = {}
for month in range(12):
    sample = np.empty([18, 480])
    for day in range(20):
        sample[:, day * 24 : (day + 1) * 24] = raw_data[18 * (20 * month + day ): 18 * (20 * month + day + 1), :]  # 传数据
    month_data[month] = sample
ExtractFeatures3.png

ExtractFeatures4.png

然后我们把10个小时为一个data,把前9个小时的data做训练,第10个小时的值为target,即0-9点的data来预测10点的PM2.5,第1-10点的data来预测11点的PM2.5

最后一笔是471-479来预测480的PM2.5,所以一共是471笔的data,shape从(18 * 5760)变成(12月 * 471笔),(18项 * 9小时)

y就是第10个小时的PM2.5的值,shape为(12 * 471,1)

x = np.empty([12 * 471, 18 * 9], dtype = float)
y = np.empty([12 * 471, 1], dtype = float)
for month in range(12):
    for day in range(20):
        for hour in range(24):
            if day == 19 and hour > 14:
                continue
            x[month * 471 + day * 24 + hour, :] = month_data[month][:, day * 24 + hour: day * 24 + hour + 9].reshape(1,-1) #vector dim:18*9 (9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9)
            y[month * 471 + day * 24 + hour, 0] = month_data[month][9, day * 24 + hour + 9] # value
x
shape.png

标准化(Normalization)

mean_x = np.mean(x, axis = 0) #18 * 9 
std_x = np.std(x, axis = 0) #18 * 9 
for i in range(len(x)): #12 * 471
    for j in range(len(x[0])): #18 * 9 
        if std_x[j] != 0:
            x[i][j] = (x[i][j] - mean_x[j]) / std_x[j]
x   

把训练数据分成训练集train_set和验证集validation,其中train_set用于训练,而validation不会参与训练,仅用于验证。

import math
x_train_set = x[: math.floor(len(x) * 0.8), :]
y_train_set = y[: math.floor(len(y) * 0.8), :]
x_validation = x[math.floor(len(x) * 0.8):, :]
y_validation = y[math.floor(len(y) * 0.8):, :]
validation.png

梯度下降算法:Adagrad

和上图不同处: 下面Loss的代码用到的是 Root Mean Square Error

因为存在常数项b,所以维度(dim)需要多加一列;eps项是极小值,避免adagrad的分母为0.

每一个维度(dim)会对应到各自的gradient和权重w,通过一次次的迭代(iter_time)学习。最终,将训练得到的模型(权重w)存储为.npy格式的文件。


8P8URK0(5KZJ2`KUAAM6UC6.png

adagrad1.png
adagrad2.png
adagrad3.png
dim = 18 * 9 + 1
w = np.zeros([dim, 1])
x = np.concatenate((np.ones([12 * 471, 1]), x), axis = 1).astype(float)
learning_rate = 100
iter_time = 1000
adagrad = np.zeros([dim, 1])
eps = 0.0000000001
for t in range(iter_time):
    loss = np.sqrt(np.sum(np.power(np.dot(x, w) - y, 2))/471/12)#rmse
    if(t%100==0):
        print(str(t) + ":" + str(loss))
    gradient = 2 * np.dot(x.transpose(), np.dot(x, w) - y) #dim*1
    adagrad += gradient ** 2
    w = w - learning_rate * gradient / np.sqrt(adagrad + eps)
np.save('weight.npy', w)
w
87U{999GOO%X8FK0}J83~UX.png

测试

testdata = pd.read_csv(r"E:\RaoLing\DataSet\李宏毅机器学习数据\hw1\test.csv",
                       header=None, encoding="big5")
test_data = testdata.iloc[:, 2:]
test_data.replace("NR", 0, inplace=True)
test_data = test_data.to_numpy()
test_data
test_x = np.empty([240, 18*9], dtype = float)
for i in range(240):
    test_x[i, :] = test_data[18 * i: 18* (i + 1), :].reshape(1, -1)
    
for i in range(len(test_x)):
    for j in range(len(test_x[0])):
        if std_x[j] != 0:
            test_x[i][j] = (test_x[i][j] - mean_x[j]) / std_x[j]          # 用train训练的模型来分类test的数据
test_x = np.concatenate((np.ones([240, 1]), test_x), axis = 1).astype(float)
test_x

载入模型即可对test数据进行预测,得到预测值ans_y

w = np.load('weight.npy')
ans_y = np.dot(test_x, w)
ans_y
res = pd.DataFrame(ans_y)
res['id'] = ['id_' + str(x) for x in res.index]
res.rename(columns={0:'value'}, inplace=True)
res=res[['id', 'value']]  # 调换顺序
res.to_csv('submit.csv', index=False)

output.png

最后,可以调整 learning rate、iter_time (iteration 次数)、
取用 features 的多少(取几个小时,取哪些特征)
甚至是不同的 model 来超越 baseline

参考链接:https://www.jianshu.com/p/2c3b54cf8dbc

最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 213,047评论 6 492
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 90,807评论 3 386
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 158,501评论 0 348
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 56,839评论 1 285
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 65,951评论 6 386
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 50,117评论 1 291
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 39,188评论 3 412
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 37,929评论 0 268
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 44,372评论 1 303
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 36,679评论 2 327
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 38,837评论 1 341
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 34,536评论 4 335
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 40,168评论 3 317
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 30,886评论 0 21
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 32,129评论 1 267
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 46,665评论 2 362
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 43,739评论 2 351