Homework 1: Linear Regression
本次目標:由前 9 個小時的 18 個 features (包含 PM2.5)預測的 10 個小時的 PM2.5。
目录
Load ‘train.csv’
train.csv 的資料為 12 個月中,每個月取 20 天,每天 24 小時的資料(每小時資料有 18 個 features)。
import sys
import pandas as pd
import numpy as np
data = pd.read_csv('./data/train.csv', encoding = 'big5')
data.iloc[0:4,0:12] #看看数据
日期 | 測站 | 測項 | 0 | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | |
---|---|---|---|---|---|---|---|---|---|---|---|---|
0 | 2014/1/1 | 豐原 | AMB_TEMP | 14 | 14 | 14 | 13 | 12 | 12 | 12 | 12 | 15 |
1 | 2014/1/1 | 豐原 | CH4 | 1.8 | 1.8 | 1.8 | 1.8 | 1.8 | 1.8 | 1.8 | 1.8 | 1.8 |
2 | 2014/1/1 | 豐原 | CO | 0.51 | 0.41 | 0.39 | 0.37 | 0.35 | 0.3 | 0.37 | 0.47 | 0.78 |
3 | 2014/1/1 | 豐原 | NMHC | 0.2 | 0.15 | 0.13 | 0.12 | 0.11 | 0.06 | 0.1 | 0.13 | 0.26 |
Preprocessing
將 ‘RAINFALL’ 欄位全部補 0。数据RAINFALL中有一些值是NR,把他变成0
data=data.iloc[:,3:] #去掉数据的前面不包含数据的三列
data.head(3)
0 | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 | ... | 14 | 15 | 16 | 17 | 18 | 19 | 20 | 21 | 22 | 23 | |
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
0 | 14 | 14 | 14 | 13 | 12 | 12 | 12 | 12 | 15 | 17 | ... | 22 | 22 | 21 | 19 | 17 | 16 | 15 | 15 | 15 | 15 |
1 | 1.8 | 1.8 | 1.8 | 1.8 | 1.8 | 1.8 | 1.8 | 1.8 | 1.8 | 1.8 | ... | 1.8 | 1.8 | 1.8 | 1.8 | 1.8 | 1.8 | 1.8 | 1.8 | 1.8 | 1.8 |
2 | 0.51 | 0.41 | 0.39 | 0.37 | 0.35 | 0.3 | 0.37 | 0.47 | 0.78 | 0.74 | ... | 0.37 | 0.37 | 0.47 | 0.69 | 0.56 | 0.45 | 0.38 | 0.35 | 0.36 | 0.32 |
3 rows × 24 columns
data[data=='NR']=0 #把值为NR的置为0
raw_data=data.to_numpy()
raw_data[0:3,0:6]
array([['14', '14', '14', '13', '12', '12'],
['1.8', '1.8', '1.8', '1.8', '1.8', '1.8'],
['0.51', '0.41', '0.39', '0.37', '0.35', '0.3']], dtype=object)
提取特征
原始的数据是(4320, 24),现在把它分为12个月。每个月18(features) * 480 (hours) .只有480=20天*24个小时
raw_data.shape
(4320, 24)
month_data={} #声明字典
for month in range(12):
sample=np.empty([18,480]) #暂存一个月的数据
for day in range(20): #一个月只要20天
sample[:,day*24:(day+1)*24]=raw_data[18 * (20 * month + day) : 18 * (20 * month + day + 1), :]
month_data[month]=sample
每個月會有 480hrs,每 9 小時形成一個 data,每個月會有 471 個 data,故總資料數為 471 * 12 筆,而每筆 data 有 9 * 18 的 features (一小時 18 個 features * 9 小時)。
對應的 target 則有 471 * 12 個(第 10 個小時的 PM2.5)
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
print(x.shape)
print(y.shape)
(5652, 162)
(5652, 1)
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]
Normalize
np.set_printoptions(suppress=True) # 使输出不是科学计数法的形式,而是浮点型的数据输出。
print(mean_x[0:5])
np.mean(x,axis=0)[0:5]
[22.53651805 22.54324133 22.54978769 22.55368011 22.55598018]
array([ 0., -0., 0., -0., -0.])
划分数据集 “train_set” and “validation_set”
這部分是針對作業中 report 的第二題、第三題做的簡單示範,以生成比較中用來訓練的 train_set 和不會被放入訓練、只是用來驗證的 validation_set。
import math
x_train_set=x[:math.floor(len(x)*0.8),:] #train set是前80%的数据
y_train_set=y[:math.floor(len(x)*0.8),:]
x_validation=x[math.floor(len(x)*0.8):,:]
y_validation=y[math.floor(len(x)*0.8):,:]
训练
dim=18*9+1
w=np.zeros([dim,1])
x=np.concatenate( (np.ones([12*471,1]) , x), axis=1).astype(float) #增加b
learning_rate=100
iter_time=1000
adagrad=np.zeros([dim,1])
eps = 0.0000000001
loss_all=[]
for i in range(iter_time):
loss= np.sqrt(np.sum((np.dot(x,w)-y)**2)/471/12)
if(i%50==0):
loss_all.append(loss)
gradient=2*np.dot(x.T, np.dot(x,w)-y)
adagrad+=gradient**2
w=w-learning_rate*gradient/np.sqrt(adagrad+eps)
np.save('weight.npy', w)
import matplotlib.pyplot as plt
plt.plot(loss_all)
np.mean(np.dot(x,w)-y)
-3.3389977108613344e-15
dim=18*9+1
w=np.zeros([dim,1])
x=np.concatenate( (np.ones([len(x_train_set),1]) , x_train_set), axis=1).astype(float) #增加b
learning_rate=100
iter_time=1000
adagrad=np.zeros([dim,1])
eps = 0.0000000001
loss_all=[]
y=y_train_set
for i in range(iter_time):
loss= np.sqrt(np.sum((np.dot(x,w)-y)**2)/471/12)
if(i%50==0):
loss_all.append(loss)
gradient=2*np.dot(x.T, np.dot(x,w)-y)
adagrad+=gradient**2
w=w-learning_rate*gradient/np.sqrt(adagrad+eps)
np.save('weight.npy', w)
np.mean(np.dot(x,w)-y)
np.sqrt(np.sum((np.dot(x,w)-y)**2)/471/12)
5.342550298780461
x_validation=np.concatenate( (np.ones([len(x_validation),1]) , x_validation), axis=1).astype(float) #增加b
np.mean(np.dot(x_validation,w)-y_validation)
np.sqrt(np.sum((np.dot(x_validation,w)-y_validation)**2)/471/12)
2.6447202742907496