作业任务:
本次作业使用丰原站的观测记录,分成train set跟test set,train set是丰原站每个月的前20天所有数据。test set则是从丰原站剩下的数据中取样出来。
数据的样例如下
每18行为一日中18个feature的变化情况,一天共有24小时,因此有24列数据。
首先进行数据预处理
需要将前两列不需要的数据删除并将所有RAINFALL行为NR的数据项替换成0
import pandas as pd
import numpy as np
import math
# 取需要的数值部分,将 'RAINFALL' 栏位全部置0
data = pd.read_csv('train.csv')
data = data.iloc[:, 3:]
data[data == 'NR'] = 0
raw_data = data.to_numpy()
将原始4320(18 * 20 * 12,18个特征,20天,12个月) * 24(24小时)的数据依照每个月分重组成12个18(features)* 480(hours)的数据。
month_data = {}
for month in range(12):
temp_data = np.empty((18, 480))
for day in range(20): # 将一个月内的数据往后面插入
temp_data[:, day * 24:(day + 1)*24] = raw_data[(month * 20 + day) * 18:(month * 20 + day + 1)*18, :]
month_data[month] = temp_data
每个月会有480小时,每9个小时形成一个data(每次窗口向后滑动一个单位),每个月会有471个data,故总数据数为471 * 12笔,而每笔data有9 * 18的features。对应的target则有471 * 12个(第10个小时的PM2.5)
X = np.empty((471 * 12, 9 * 18), dtype=float)
Y = np.empty((471 * 12, 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)
Y[month * 471 + day * 24 + hour, 0] = month_data[month][9, day * 24 + hour + 9]
正则化
其转化函数为:
x* = (x - μ ) / σ
其中μ为所有样本数据的均值,σ为所有样本数据的标准差。
mean_x = np.mean(X, axis=0)
std_x = np.std(X, axis=0)
for i in range(len(X)):
for j in range(len(X[0])):
if std_x[j] != 0:
X[i][j] = (X[i][j] - mean_x[j]) / std_x[j]
训练模型
loss选用的是均方根误差(Root Mean Square Error)
dim = 18 * 9 + 1 # 共有18个feature,多一个是去掉参数中的b,用w0替代
w = np.zeros([dim, 1])
x = np.concatenate((np.ones([471 * 12, 1]), X), axis=1).astype(float) # 将x0的与其他x拼接起来
learning_rate = 100
iter_time = 10000
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)
if t % 100 == 0:
print(str(t) + ":" + str(loss))
gradient = 2 * np.dot(x.transpose(), np.dot(x, w) - Y)
adagrad += gradient ** 2 # adagrad中,学习率要除以历史所有的梯度平方和
w = w - learning_rate * gradient / np.sqrt(adagrad + eps)
np.save('weight.npy', w)
adagrad算法
为什么gradient = 2 * np.dot(x.transpose(), np.dot(x, w) - Y)?
最后验证阶段
w = np.load('weight.npy')
x_validation = np.concatenate((np.ones([1131, 1]), x_validation), axis=1).astype(float)
ans_y = np.dot(x_validation, w)
loss = np.sqrt(np.sum(np.power(ans_y - y_validation, 2)) / 1131)
print(loss)