最近开始学习机器学习,听同学介绍,选择了台大李宏毅老师的视频,https://www.bilibili.com/video/BV1JE411g7XF
老师布置了很多作业,这里记录下自己写作业的过程和一些参考资料
本篇是作业一,根据丰源站历史的空气数据,预测下一小时的pm2.5的值。
本次作业说明视频参见 https://www.bilibili.com/video/BV1gE411F7td
import numpy as np
import pandas as pd
本次作业要求只用numpy,pandas手写线性回归,所以只导入pandas和numpy两个包,
学习这两个包推荐一个gitbook, 利用python进行数据分析
另外这本书有本实体版的,是另外一个人翻译的,翻译的非常差,非常不建议购买,看这个gitbook版的就好。
本次作业代码参考 https://zhuanlan.zhihu.com/p/115335196
上面知乎专栏的大佬写的比较简单,我添加了更多的细节,便于新人理解。
data = pd.read_csv('data/train.csv', encoding='big5')
首先使用pandas读取训练文件,并做一些预处理,李宏毅老师是湾湾人,
用的是繁体,所以要指定下编码为big5,否则里面的中文为乱码。
看一下文件的结构,有27列,第一列为日期,第三列为检测项,后面的是0-23点每一点的数值。
然后每天有18行,即18个检测项,其中第10行是pm2.5。
一共12月*20天*18项=4320行
data.head(20)
然后前三项是没用的,我们用切片去掉。shape就变为了4320*24
数据中还有很多是NR,我们替换为0。
接下来就是numpy的处理范围了,使用to_numpy转换为矩阵。
data = data.iloc[:,3:]
data.replace('NR',0,inplace=True)
row_data = data.to_numpy()
row_data
参考上图,我们需要把数据shape做一个转换,把240天内相同的测量项放在同一行。
也就是从(12月*20天*18项),24小时变为(18项,(24小时*20天),12月
下面的两个循环就是在做转换,每个月的数据写在以月份为键名的字典里的。
month_data={}
for month in range(12):
sample = np.empty([18, 480])
for day in range(20):
sample[:,day*24:(day+1)*24] = row_data[18*(month*20+day):18*(month*20+(day+1)),:]
month_data[month] = sample
然后我们以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从变为
y就是第10个小时的pm2.5的值,shape为
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)
y[month * 471 + day * 24 + hour, 0] = month_data[month][9, day * 24 + hour + 9]
x
接下来做Normalize,教学视频里说明了为什么要做Normalize,不知道的请复习视频
按照上面的公式,需要求平均值和标准差,然后x处理对应位置的值
mean_x = np.mean(x, axis = 0)
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
接下来就是最重要的部分了,我们开始做梯度下降
我们有9小时*18项的变量,所以需要9*18个w值,还有一个bias:b
注意下面代码的第三行,这里是整个矩阵前拼接了一列1值,这一列在点乘后就是b
所以这就是为什么dim有18*9+1项
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
然后定义loss函数,这里知乎作者用的是均方根误差来做loss函数
贴一段百度百科的解释:
均方根误差是预测值与真实值偏差的平方与观测次数n比值的平方根,
在实际测量中,观测次数n总是有限的,真值只能用最可信赖(最佳)值来代替。
标准误差对一组测量中的特大或特小误差反映非常敏感,
所以,标准误差能够很好地反映出测量的精密度。
这正是标准误差在工程测量中广泛被采用的原因。
因此,标准差是用来衡量一组数自身的离散程度,
而均方根误差是用来衡量观测值同真值之间的偏差,
它们的研究对象和研究目的不同,但是计算过程类似。
因为常数项的存在,所以 dimension (dim) 需要多加一栏;
eps 项是避免 adagrad 的分母为 0 而加的极小数值。
rmse公式:
然后adagrid的公式:
这里对gt的计算我没有看懂,根据公式应该是loss对w的偏导,以后搞懂了再来补说明吧
adagrid这一行就是上面的累计求和
然后对w做更新,这样就可以实现梯度下降了
最后保存计算出的w矩阵
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和b就计算完了,接下来拿测试集做测试
还是和train.csv一样先做预处理和Normalize,然后用w和测试集做点乘就可以得到预测值了
记得在拼接1列1作为b
testdata = pd.read_csv('./data/test.csv', header = None, encoding = 'big5')
test_data = testdata.iloc[:, 2:].copy()
test_data.replace('NR', 0, inplace=True)
test_data = test_data.to_numpy()
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]
test_x = np.concatenate((np.ones([240, 1]), test_x), axis = 1).astype(float)
w = np.load('weight.npy')
ans_y = np.dot(test_x, w)
ans_y[:10]
ans_y就是预测值了,然后按照作业要求写成csv文件。
这里转回pandas处理,因为pandas的to_csv方法很好用
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)
最后,可以调整 learning rate、iter_time (iteration 次数)、
取用 features 的多少(取几个小时,取哪些特征)
甚至是不同的 model 来超越 baseline。
网友评论