一、问题描述

希望用线性回归解决问题：
给定某个地区连续9小时包括PM2.5在内的18项污染物每小时的数据，预测第10个小时的PM2.5的值。

二、设计简要描述

机器学习的三个基本步骤——

程序设计思路——

三、程序清单

import pandas as pd
import numpy as np# 1 训练集处理
# 1.1 读取训练集数据，big5是针对于文档中存在繁体字的编码
train_data = pd.read_csv('./train.csv', usecols=range(2, 27), encoding='big5')  # 2~27，从测项到最后一小时数据
row = train_data['測項'].unique()  # 得到的row是一个list，存放了各个污染物的名字# 1.2训练数据规整化
# 定义18维空列表，每一个维度处理一种污染物
# DataFrame是一个二维表格,让它有18行，24*240列(240=12个月x每个月前20天)
new_train_data = pd.DataFrame(np.zeros([18, 24 * 240]))
# 计算维度？
n = 0
# 遍历18类污染物，将数据处理为 [18, 24*240]
for i in row:# 依照污染物名取出相同的行train_data1 = train_data[train_data['測項'] == i]  # 得到的train_data1是一个二维表格，每一行是一天24小时某项污染物的值，共240天# 删去多出的“测项”那一列（存放污染物的名字）train_data1.drop(['測項'], axis=1, inplace=True)# 格式处理，赋值到 new_train_datatrain_data1 = np.array(train_data1)  # 刚刚的二维表格依然是二维表格，只不过可以处理了train_data1[train_data1 == 'NR'] = '0'  # 将降雨量NR替换为0train_data1 = train_data1.astype('float')# 24*240=5760# 二维列表原先包含240个含24个元素的一维列表，现在包含5760个含1个元素的一维列表train_data1 = train_data1.reshape(5760, 1)# 进行转置，二维列表现在包含1个含5760个元素的一维列表train_data1 = train_data1.T# 将刚刚在一起的5760个相同污染物的量存放进预先设定好的DataFrame的一行new_train_data.loc[n] = train_data1# 下一类污染物处理n += 1# 可以将现在的new_train_data和train_array看成是包含18个含5760个元素的二维表格
train_array = np.array(new_train_data).astype(float)# 1.3 处理训练集
# 训练样本features集合
X_train = []
# 训练样本目标PM2.5集合
y_train = []
for i in range(new_train_data.shape[1] - 9):# 每次取9个小时的数据作训练集# 每次数据遍历每行前9个数据，全部加入到训练集中，18 X 9X_temp = np.ones(18 * 9)# 记录count = 0  # count的取值是0~161(18*9=162)for j in range(18):x = train_array[j, i:i + 9]for k in range(9):X_temp[count] = x[k]count += 1# 将样本分别存入X_train中X_train.append(X_temp)# 取本次第10个小时的 PM2.5 的值作为训练的真实值y = int(train_array[9, i + 9])  # 第一个9是PM2.5，第二个+9是找到第十个小时# 将样本分别存入X_train、y_train中y_train.append(y)
X_train = np.array(X_train)  # X_train里面是由5751个含有172个元素的一维列表组成的二维列表
y_train = np.array(y_train)  # X_train里面是由5751个含有1个元素的一维列表组成的二维列表# 2.训练模型
# 训练轮数
epoch = 1000(第一次)
epoch = 2000(第二次)
# 开始训练
# 更新参数，训练模型
# 偏置值初始化
bias = 0
# 权重初始化
weights = np.ones(18 * 9)
# 初始学习率
learning_rate = 1
# 用于存放偏置值的梯度平方和
bg2_sum = 0
# 用于存放权重的梯度平方和
wg2_sum = np.zeros(18 * 9)for i in range(epoch):b_g = 0w_g = np.zeros(18 * 9)# 在所有数据上计算Loss_label的梯度for j in range(len(X_train)):# weights.dot(X_train[j])是172个权重和对应的输入x相乘(dot的本意是点乘)b_g += (y_train[j] - weights.dot(X_train[j]) - bias) * (-1)for k in range(18 * 9):w_g[k] += (y_train[j] - weights.dot(X_train[j]) - bias) * (-X_train[j, k])# 求平均b_g /= len(X_train)  # 这个数字是5751w_g /= len(X_train)# adagradbg2_sum += b_g ** 2wg2_sum += w_g ** 2# 更新权重和偏置bias -= learning_rate / bg2_sum ** 0.5 * b_gweights -= learning_rate / wg2_sum ** 0.5 * w_g# 每训练10轮，输出一次在训练集上的损失if i % 10 == 0:loss = 0for j in range(len(X_train)):loss += (y_train[j] - weights.dot(X_train[j]) - bias) ** 2print('经过{}轮训练,训练集上的损失值为:'.format(i), loss / len(X_train))# 测试机处理
# 3.1读取测试数据
test_data = pd.read_csv('./train2.csv', usecols=range(2, 27), encoding='big5')  # 用train的最后三个月作为test 3242-4321
row2 = test_data['測項'].unique()
# 3.2 测试数据规整化
new_test_data = pd.DataFrame(np.zeros([18, 24 * 60]))
# 计算维度
m = 0
# 遍历18类污染物，将数据处理为[18，24*60]
for i in row2:# 依照污染物名字取出相同的行test_data1 = test_data[test_data['測項'] == i]# 删去多出的“测项”那一列test_data1.drop(["測項"], axis=1, inplace=True)# 格式处理，赋值到new_train_datatest_data1 = np.array(test_data1)test_data1[test_data1 == "NR"] = "0"test_data1 = test_data1.astype("float")# 24*60=1440# 二维列表原先包含60个含24个元素的一维列表，现在包含1440个含1个元素的二维列表test_data1 = test_data1.reshape(1440, 1)# 转置test_data1 = test_data1.T# 将刚刚在一起的1440个相同污染物的量存放进预先设定好的new_test_data的一行new_test_data.loc[m] = test_data1# 下一类污染物处理m += 1
# 现在可以将new_test_data和test_array看成史包含18个含有1440个元素的二维表格
test_array = np.array(new_test_data).astype(float)
# 3.3 处理测试集
# 测试样本features集合
X_test = []  # 二维列表
# 测试样本PM2.5集合
y_test = []  # 一维列表
for i in range(new_test_data.shape[1] - 9):# 每次取9个小时的数据做测试集的输入# 每次数据遍历每行前9个数据，全部加入到测试集中，18x9X_temp = np.ones(18 * 9)# 记录count = 0for j in range(18):x = test_array[j, i:i + 9]for k in range(9):X_temp[count] = x[k]count += 1# 将样本分别存入X_test中X_test.append(X_temp)# 取本次第十个小时的PM2.5值作为测试的真实值y = int(test_array[9, i + 9])# 将样本存入y_test中y_test.append(y)
X_test = np.array(X_test)
y_test = np.array(y_test)# 4.测试训练好的模型
loss = 0
for i in range(len(X_test)):loss += (y_test[i] - weights.dot(X_test[i]) - bias) ** 2
print('测试集上的损失值为:', loss / len(X_test))

四、结果分析

1. 训练轮数为1000时，在测试集上的损失函数结果为7.208519807515013
   
2. 训练轮数为2000时，在测试集上的损失函数结果为测试集上的损失值为: 6.622474294012304
   
   分析 ：从训练集和测试集上的损失值相近可知过程无误，但结果不是特别理想，可能原因有：处理训练集时把每一天的数据连在一起，增大了利用率，但未考虑到train.csv的数据是截取的每个月前20天的，也就是连接的过程把上个月的第二十天跟下个月第一天连在一起了，这是不合理的。
   weights是一个18*9=162维的向量，参数太多，可能过拟合了。

五、调试报告

1. 数据集的编码方式一定要注意，如果pycharm一开始显示不出数据集编码方式，可以先用notepad++打开数据集，查看其编码方式。
2. 阅读实验指导书上的代码时难以理解，可以通过确定数字的最大值来解决，例如这一段

count = 0  
for j in range(18):
x = train_array[j, i:i + 9]for k in range(9):X_temp[count] = x[k]count += 1

起初无法理解X_temp是什么形状，但看到count经过两轮循环，数值变化范围是0~18*9。可知X_temp是把18项数据在9个小时里的连续值进行了首尾拼接，便于后续与权重72维向量相乘。
3. 打桩调试
训练轮数太多，想确定程序执行到哪一步了，可以对在一些数值处进行输出。例如

# 每训练10轮，输出一次在训练集上的损失
if i % 10 == 0:loss = 0for j in range(len(X_train)):loss += (y_train[j] - weights.dot(X_train[j]) - bias) ** 2
print('经过{}轮训练,训练集上的损失值为:'.format(i), loss / len(X_train))

六、实验总结

1. 学到了一些数据处理知识
   ① .astype(‘float’)将数据转化为浮点型
   ② 对数据进行转置 train_data1 = train_data1.T
   ③ 根据列名得到表格中的某一项数据
   row = train_data[‘測項’].unique()
   ④ 删除某一列数据train_data1.drop([‘測項’], axis=1, inplace=True)
2. 需要改进的地方
   ① 对测试集和验证集(后三个月)的处理出现大量的重复代码，应该编写函数，但是对于Python掌握还不够，遂放弃
   ② 把每天的数据连在一起增加了数据的利用率，但是每个月的应该分开
   ③ 采用前9个小时18项数据，权重作为一个72维的向量，参数过多，增加了训练成本，并且容易过拟合