思路主要参考:https://www.cnblogs.com/HL-space/p/10676637.html
一、作业描述![](https://i-blog.csdnimg.cn/blog_migrate/10c31e961d30233cce5cc5e38b513343.png)
数据集描述:
- 取12个月每个月前20天共240天的气象数据
- 每天的数据有18个指标
二、思路分析
2.1 数据预处理
利用pd.read_csv 读取标题行之外的数据
数据中共有很多NR数据,该指标RAINFALL 为是否降雨,将NR替换为0
如图所示,取每天从0时开始到8时的所有数据为训练数据,第9个小时的PM2.5值为label,以此类推,则每天共有15个label,240天共有240X15 = 3600个label
#数据预处理,取出每天24h的数据
def dataProcessing(df):
df = df.iloc[:,3:].replace(['NR'],[0.0]) #将空数据替换为0
array = np.array(df).astype(float) # astype转换array中的元素类型
x_list = []
y_list = []
for i in range(0,4320,18):
for j in range(0,15):
mat = array[i:i+18,j:j+9] #每一天中每10个小时的前9个小时当作参数
label = array[i+9,j+9] #第10列是pm2.5的值,第10个小时的pm2.5作为label
x_list.append(mat) #插入参数和label列表
y_list.append(label)
x = np.array(x_list)
y = np.array(y_list)
return x, y
2.2 模型建立
2.2.1 回归模型
简单起见,仅利用线性函数考虑前9个小时的PM2.5的对第10小时PM2.5的影响
其中
为第i小时的PM2.5数据,
为其对应的权重值,
为偏移量,
为第10小时的预测值
2.2.2 损失函数
利用预测值和每个label的平均欧氏距离来衡量损失
其中,为label的值,
为函数的预测值,
为label的数量,
为label的损失
未加入正则项
2.2.3 梯度下降
要想使loss最小,则需要求loss在b和w上的偏微分
更新参数
为w更新时的学习率,
为b更新时的学习率
2.2.4 更新学习率
为累计梯度平方和的二次平方根,
为更新前的学习率,
为更新后的学习率
#调整参数,训练模型
def train(x_train,y_train,iteration):
# ydata = b + w0 * xdata0 +...+ w8 * xdata8
b = 0 #initial b
w = np.ones(9) #initial w
lr = 5 #learning rate
pre_gra_b = 0 # adagrad history b
pre_gra_w = np.zeros(9) # adagrad history w
ly = len(y_train)
#迭代开始, 用前9个小时的pm2.5数据计算预测第10个小时
for i in range(iteration):
b_grad = 0.0
w_grad = np.zeros(9)
#利用adagrad更新学习率
for n in range(ly):
b_grad = b_grad - 2.0 * (y_train[n] - w.dot(x_train[n,9,:])-b)
for m in range(9):
w_grad[m] = w_grad[m] - 2.0 * (y_train[n] - w.dot(x_train[n,9,:]) - b) * (x_train[n,9,m])
#mean
b_grad = b_grad/ly
w_grad = w_grad/ly
#adagrad
pre_gra_b = pre_gra_b + b_grad ** 2
pre_gra_w = pre_gra_w + w_grad ** 2
#update parameters
b = b - (lr * b_grad) / np.sqrt(pre_gra_b)
for k in w:
k = k - (lr * w_grad ) / np.sqrt(pre_gra_w)
#loss on training data, per 200
if i % 200 == 0:
loss = 0
for j in range(ly):
loss = loss + (y_train[j] - w.dot(x_train[j,9,:]) - b) ** 2
print('iteration:{},loss:'.format(i),loss/ly)
return w,b
2.3 在测试集上测试
#testing on test data
def test(x_test,y_test,w,b):
lt = len(y_test)
loss = 0
for j in range(lt):
loss = loss + (y_test[j] - w.dot(x_test[j,9,:]) - b) ** 2
print('the loss on test data is:', loss/lt)
return loss/lt
三、完整代码
# -*- coding: utf-8 -*-
"""
Created on Fri Nov 1 19:49:52 2019
@author: Lenovo
"""
import pandas as pd
import numpy as np
#数据预处理,取出每天24h的数据
def dataProcessing(df):
df = df.iloc[:,3:].replace(['NR'],[0.0]) #将空数据替换为0
array = np.array(df).astype(float) # astype转换array中的元素类型
x_list = []
y_list = []
for i in range(0,4320,18):
for j in range(0,15):
mat = array[i:i+18,j:j+9] #每一天中每10个小时的前9个小时当作参数
label = array[i+9,j+9] #第10列是pm2.5的值,第10个小时的pm2.5作为label
x_list.append(mat) #插入参数和label列表
y_list.append(label)
x = np.array(x_list)
y = np.array(y_list)
return x, y
#调整参数,训练模型
def train(x_train,y_train,iteration):
# ydata = b + w0 * xdata0 +...+ w8 * xdata8
b = 0 #initial b
w = np.ones(9) #initial w
lr = 5 #learning rate
pre_gra_b = 0 # adagrad history b
pre_gra_w = np.zeros(9) # adagrad history w
ly = len(y_train)
#迭代开始, 用前9个小时的pm2.5数据计算预测第10个小时
for i in range(iteration):
b_grad = 0.0
w_grad = np.zeros(9)
#利用adagrad更新学习率
for n in range(ly):
b_grad = b_grad - 2.0 * (y_train[n] - w.dot(x_train[n,9,:])-b)
for m in range(9):
w_grad[m] = w_grad[m] - 2.0 * (y_train[n] - w.dot(x_train[n,9,:]) - b) * (x_train[n,9,m])
#mean
b_grad = b_grad/ly
w_grad = w_grad/ly
#adagrad
pre_gra_b = pre_gra_b + b_grad ** 2
pre_gra_w = pre_gra_w + w_grad ** 2
#update parameters
b = b - (lr * b_grad) / np.sqrt(pre_gra_b)
for k in w:
k = k - (lr * w_grad ) / np.sqrt(pre_gra_w)
#loss on training data, per 200
if i % 200 == 0:
loss = 0
for j in range(ly):
loss = loss + (y_train[j] - w.dot(x_train[j,9,:]) - b) ** 2
print('iteration:{},loss:'.format(i),loss/ly)
return w,b
#testing on test data
def test(x_test,y_test,w,b):
lt = len(y_test)
loss = 0
for j in range(lt):
loss = loss + (y_test[j] - w.dot(x_test[j,9,:]) - b) ** 2
print('the loss on test data is:', loss/lt)
return loss/lt
def main():
file = 'HW1-train.csv'
data = pd.read_csv(file, encoding = 'gb18030')# open data
iteration = 2000 #迭代次数
x,y= dataProcessing(data)
x_train, y_train = x[:3200], y[:3200]
x_test,y_test = x[3200:3600], y[3200:3600]
weights,bias = train(x_train, y_train,iteration)
test(x_test,y_test,weights,bias)
if __name__ == '__main__':
main()
mark:不知道为什么损失总是特别大,还需要继续调一下
四、待改进
- 只是按天进行每十个小时的切割,应该跨天更合理一点
- 增加正则项,防止过拟合
- 考虑其他天气因素的影响,增加参数
- 随机抽取训练集和测试集,而不是简单粗暴的直接前后划分