深圳网站建设公司服务平台,网站设计与编辑,做网站用python好还是PHP好,网站开发建设专业机器学习(课后作业–hw1)
本篇文章全文参考这篇blog 网上找了很多教程#xff0c;这个是相对来说清楚的#xff0c;代码可能是一模一样#xff0c;只是进行了一些微调#xff0c;但是一定要理解这个模型具体的处理方法#xff0c;这个模型我认为最巧妙的它对于数据的处理…机器学习(课后作业–hw1)
本篇文章全文参考这篇blog 网上找了很多教程这个是相对来说清楚的代码可能是一模一样只是进行了一些微调但是一定要理解这个模型具体的处理方法这个模型我认为最巧妙的它对于数据的处理直接把前9天所有的参数参数当作变量最简单粗暴的方法然后再是里面一些矩阵的运算也很重要可能直接看无法理解那就手动创造几个数据进行手动模拟一定会有不小的收获的 作业介绍
题目要求利用每天的前9个小时的18个特征包括PM2.5来预测第10个小时的PM2.5的值。
特征包括 AMB_TEMP, CH4, CO, NHMC, NO, NO2, NOx, O3, PM10, PM2.5, RAINFALL, RH, SO2, THC, WD_HR, WIND_DIREC, WIND_SPEED, WS_HR总计18种。每一种都可以作为一个变量加入你的模型里面。
数据下载Github链接 下载结束以后把文件放入对应的文件夹里面就好了 CSV文件包含台湾丰原地区240天的气象观测资料(取每个月前20天的数据做训练集12月X20天240天每月后10天数据用于测试对学生不可见) 作业任务根据每月的观测数据对未来进行预测
题目分析
Dataset里面我们要用的两个csv文件数据集train.csv和测试集(test.csv)。 数据集中记录了12个月的前20天每天24个小时的18个物质浓度的资料然后测试集就是剩余的数据中再取的。从上述我们可以得到以下结论
模型的输入前9个小时的观测数据每小时18个数据所以输入参数有9*18个参数值模型的输出预测第十个小时的PM2.5浓度模型回归模型Regression Model要求根据前面连续9个小时的观测数据预测接下来一个小的PM2.5浓度 每一天的信息182418个指标24个小时训练数据 训练数据是连续9个小时的观测数据接下来一个小时的PM2.5作为这个训练数据集的label 注意每小时的数据可以作为训练数据也可以作为label一共12个月每个月20天一天24小时每天18个观测数据所以一共有(2024)x18个数据而训练数据集有2024-9471个因为最开始的9个无法作为数据集的label
代码解答
1. Data Preprocessing 我们需要做的第一件事就是数据的读入和预处理 按照文件路径读入去除无用数据预处理数据转化 # 引入必要的包
import pandas as pd
import numpy as np# 读取数据保存到data中路径根据你保存的train.csv位置而有变化
# 数据编码方式这里选择big5台湾是繁体字使用utf-8可能会出现问题
data pd.read_csv(Dataset/train.csv, encodingbig5)
# print(data)# 行保留所有列从第三列开始往后才保留这样去除了数据中的时间、地点、参数等信息
data data.iloc[:, 3:]
# print(data)# 将所有NR的值全部置为0方便之后处理
data[data NR] 0
# print(data)# 将data的所有数据转换为二维数据并用raw_data来保存
raw_data data.to_numpy()
# print(raw_data)
# 可以每一步都打印出结果看到数据的变化2. Feature Engineering 按照train.csv给出的数据肯定是无法满足我们的使用我们在进行提取数据集的时候把一个月的数据全部放在一起只取连续的9个小时不管是那一天的使用我们原先的按天数竖着排列修改为每个月的数据按天数横着排放这个代码进行实现其实不困就是不好理解。 处理以后才更方便我们提取连续的9小时数据作为数据集和接下来1小时PM2.5作为对于数据集的label month_data {}# month 从0-11 共12个月
# 返回一个18行480列的数组保存一个月的data一月20天一天24小时
# day从0-19 共20天
for month in range(12):sample np.empty([18, 480])for day in range(20):# raw的行每次取18行列取全部列。送到sample中sample是18行480列# 行给全部行列只给24列然后列往后增加sample[:, day * 24: (day 1) * 24] raw_data[18 * (20 * month day): 18 * (20 * month day 1), :]month_data[month] sample创建字典month为keysample为value循环遍历sample和raw_data都是两个参数行和列行和列都在变化但是我们都是优先是行然后是列这个建议自己手动演示一下这个代码真的很奇妙 month 0,day 0 遍历raw_data的第0行(0到18-1)的每一列(一共24列)依次赋值给sample的第0行(一共18行)的第0到24-1列。month 0,day 1 遍历raw_data的第18行(18到182-1)的每一列(一共24列)依次赋值给sample的第0行(一共18行)的第24到242-1列。 真的发现自己和zz一样python最大的优点就是神奇python可以直接实现矩阵对应块赋值所以我们之间进行整块转移就好了那就很好理解了 raw_data其实每遍历一次就是去除一天的数据量sample其实也就是取出前一天后面那一块。 经过上面的处理数据已经横向转化为方便我们提取训练集的数据集了
# 每月共480个小时每9个小时一个数据480列最后一列不可以计入如果取到最后一列那么最后一个数据
# 便没有了结果{需要9个小时的输入和第10个小时的第10行作为结果}480-1-91471。
# 12个月共12*471个数据按行排列每一行一个数据数据每小时有18个特征而每个数据9个小时共18*9列
x np.empty([12 * 471, 18 * 9], dtypefloat)# 12个月共12*471个数据每个数据对应一个结果即第10小时的PM2.5浓度也就是提取对应label
y np.empty([12 * 471, 1], dtypefloat)for month in range(12): # month 0-11for day in range(20): # day 0-19for hour in range(24): # hour 0-23if day 19 and hour 14: # 第20天的第23时continue# 取对应month行都要取列取9个依次进行最后将整个数据reshape成一行数据# (列数无所谓)。然后赋给xx内的坐标只是为了保证其从0-471*12# vector dim:18*9x[month * 471 day * 24 hour, :] month_data[month][:, day * 24 hour: day * 24 hour 9].reshape(1, -1)# value,结果对应的行数一直是9行即第10行PM2.5然后列数随着取得数据依次往后进行y[month * 471 day * 24 hour, 0] month_data[month][9, day * 24 hour 9]3. Data Normalization 数据的标准化我在另外一篇博客里面提到过 数据的标准化normalization是将数据按比例缩放使之落入一个小的特定区间。去除数据的单位限制将其转化为无量纲的纯数值便于不同单位或量级的指标能够进行比较和加权。 这里我们采用的是最简单的标准化公式 x* (x - μ ) / σ 其中μ为所有样本数据的均值σ为所有样本数据的标准差。 标准化其实是为了让模型更快收敛 mean_x np.mean(x, axis0) # 18 * 9 求均值axis 0表示对各列求值返回 1* 列数 的矩阵
std_x np.std(x, axis0) # 18 * 9 求标准差axis 0表示对各列求值返回 1* 列数 的矩阵
for i in range(len(x)): # 12 * 471for j in range(len(x[0])): # 18 * 9if std_x[j] ! 0:x[i][j] (x[i][j] - mean_x[j]) / std_x[j]4. Splitting Data
将训练数据按8:2拆成训练数据和验证数据。训练数据用于模型的训练测试数据用于模型的检测这样可以让我们了解自己建立的模型面对新的数据集的时候的表现。
x_train_set x[: math.floor(len(x) * 0.8), :]
y_train_set y[: math.floor(len(y) * 0.8), :]
x_validation x[math.floor(len(x) * 0.8):, :]
y_validation y[math.floor(len(y) * 0.8):, :]5. Model Training:
5.1 Regression Model 采用最普通的线性回归模型如下列公式所展示我们在前面将每一组训练数据集进训练压缩成为一行所以这一行一共有9*18个变量每个变量都有其对应的权重参数同时还设置了一个偏置bias y ∑ i 0 9 ∗ 18 ω i ⋅ x i b y \sum_{i0}^{9*18}{{\omega}_i} \cdot{x_i} b yi0∑9∗18ωi⋅xib xi为对应训练数据集中第i个变量wi为其对应的权重值b为偏置项y为该数据帧样本的预测结果。 5.2 模型训练 数据在进行了前置的处理以后我们就要开始训练自己的模型了我们这里使用的是损失函数Loss和Adagrad算法对参数进行优化。 均方根误差Root Mean Square Error R M S E 1 n ∑ i 1 n ( y i − y ^ i ) 2 RMSE \sqrt{\frac{1}{n}\sum_{i1}^{n}(y_i - \hat{y}_i)^2} RMSEn1i1∑n(yi−y^i)2 Adagrad算法 ω t 1 ω t − η ∑ i 0 t ( g i ) 2 ϵ g t \omega^{t1} \omega^t - \frac{\eta}{\sqrt{\sum_{i0}^t(g^i)^2 \epsilon}}g^t ωt1ωt−∑i0t(gi)2ϵ ηgt 其中 η \eta η是学习率 g i g^i gi是第i次训练得到的梯度 ϵ \epsilon ϵ是一个极小正值防止分母为0 这部分建议先对着代码看一遍然后按照我写的过程理解一下矩阵部分求导有点复杂有兴趣的可以了解一下 # 用来做参数vector的维数加1是为了对bias好处理还有个误差项。即h(x)w1x1w2x2wnxnb
dim 18 * 9 1
# 生成一个dim行1列的数组用来保存参数值
w np.ones([dim, 1])
# np.ones来生成12*471行1列的全1数组np.concatenateaxis1
# 表示按列将两个数组拼接起来即在x最前面新加一列内容之前x是12*471行
# 18*9列的数组新加一列之后变为12*471行18*91列的数组
x np.concatenate((np.ones([12 * 471, 1]), x), axis1).astype(float)
learning_rate 100 # 学习率
iter_time 10000 # 迭代次数# 生成dim行即163行1列的数组用来使用adagrad算法更新学习率
adagrad np.zeros([dim, 1])# 因为新的学习率是learning_rate/sqrt(sum_of_pre_grads**2),
# 而adagradsum_of_grads**2,所以处在分母上而迭代时adagrad可能为0
# 所以加上一个极小数使其不除0
eps 0.0000000001
for t in range(iter_time):# rmse loss函数是从0-n的(X*W-Y)**2之和/(471*12)再开根号# 即使用均方根误差(root mean square error),具体可看公式# /471/12即/N(次数)loss np.sqrt(np.sum(np.power(np.dot(x, w) - y, 2)) / 471 / 12)if (t % 100 0): # 每一百次迭代就输出其损失print(str(t) : str(loss))# dim*1 x.transpose即x的转置后面是X*W-Y,即2*(x的转置*(X*W-Y))是梯度# 具体可由h(x)求偏微分获得.最后生成1行18*91列的数组。转置后的X其每一行# 是一个参数与h(x)-y的值相乘之后是参数W0的修正值同理可得W0-Wn的修正值# 保存到1行18*91列的数组中即gradient# 计算梯度。梯度是损失函数对参数向量w的偏导数gradient 2 * np.dot(x.transpose(), np.dot(x, w) - y)# adagrad用于保存前面使用到的所有gradient的平方进而在更新时用于调整学习率adagrad gradient ** 2w w - learning_rate * gradient / np.sqrt(adagrad eps) # 更新权重
np.save(weight.npy, w) # 将参数保存下来6. Model Evaluation 在前面进行数据集划分的时候我们将数据集划分为了两部分按照8:2的比例分别进行训练和测试。 但是我这里采用的是全部数据进行训练然后使用我们分出来的那部分数据进行测试这样效果相较于其他处理方式效果较好。 一般来说数据集的划分比较随意我们也习惯多次划分尝试然后对模型得到的参数进行评估然后选择最好的数据集划分方式。 w np.load(weight.npy)
# 使用x_validation和y_validation来计算模型的准确率因为X已经normalize了
# 所以不需要再来一遍只需在x_validation上添加新的一列作为bias的乘数即可
x_validation np.concatenate((np.ones([1131, 1]), x_validation), axis1).astype(float)
ans_y np.dot(x_validation, w)
loss np.sqrt(np.sum(np.power(ans_y - y_validation, 2)) / 1131)
print(loss)7. Test Data Prediction
7.1 测试数据预处理 测试数据的预处理和前面的训练数据预处理很类似直接按部就班就好了 testdata pd.read_csv(Dataset/test.csv, headerNone, encodingbig5)
test_data testdata.iloc[:, 2:] # 取csv文件中的全行数即第3列到结束的列数所包含的数据
test_data.replace(NR, 0, inplaceTrue) # 将testdata中的NR替换为0
test_data test_data.to_numpy() # 将其转换为数组
# 创建一个240行18*9列的空数列用于保存textdata的输入
test_x np.empty([240, 18 * 9], dtypefloat)
for i in range(240): # 共240个测试输入数据test_x[i, :] test_data[18 * i: 18 * (i 1), :].reshape(1, -1)
# 下面是Normalize,且必须跟training data是同一种方法进行Normalize
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前面拼接一列全1数组构成240行163列数据
test_x np.concatenate((np.ones([240, 1]), test_x), axis1).astype(float)7.2 数据预测 单纯的矩阵相乘原理很简单得到一个列向量分别对应每一组数据集未来一小时的PM2.5预测值 # 进行预测
w np.load(weight.npy)
ans_y np.dot(test_x, w) # test data的预测值ans_ytest_x与W的积8. Saving Predictions 最后一步当然是把数据存入文件里 # 将预测结果填入文件当中
import csvwith open(submit.csv, modew, newline) as submit_file:csv_writer csv.writer(submit_file)header [id, value]print(header)csv_writer.writerow(header)for i in range(240):row [id_ str(i), ans_y[i][0]]csv_writer.writerow(row)print(row)完整代码 这是我文件的一些对应存放位置调整好了应该就可以直接运行 # 引入必要的包
import pandas as pd
import numpy as np# 读取数据保存到data中路径根据你保存的train.csv位置而有变化
# 数据编码方式这里选择big5台湾是繁体字使用utf-8可能会出现问题
data pd.read_csv(Dataset/train.csv, encodingbig5)
# print(data)# 行保留所有列从第三列开始往后才保留这样去除了数据中的时间、地点、参数等信息
data data.iloc[:, 3:]
# print(data)# 将所有NR的值全部置为0方便之后处理
data[data NR] 0
# print(data)# 将data的所有数据转换为二维数据并用raw_data来保存
raw_data data.to_numpy()
# print(raw_data)
# 可以每一步都打印出结果看到数据的变化month_data {}# month 从0-11 共12个月
# 返回一个18行480列的数组保存一个月的data一月20天一天24小时
# day从0-19 共20天
for month in range(12):sample np.empty([18, 480])for day in range(20):# raw的行每次取18行列取全部列。送到sample中sample是18行480列# 行给全部行列只给24列然后列往后增加sample[:, day * 24: (day 1) * 24] raw_data[18 * (20 * month day): 18 * (20 * month day 1), :]month_data[month] sample# 每月共480个小时每9个小时一个数据480列最后一列不可以计入如果取到最后一列那么最后一个数据
# 便没有了结果{需要9个小时的输入和第10个小时的第10行作为结果}480-1-91471。
# 12个月共12*471个数据按行排列每一行一个数据数据每小时有18个特征而每个数据9个小时共18*9列
x np.empty([12 * 471, 18 * 9], dtypefloat)# 12个月共12*471个数据每个数据对应一个结果即第10小时的PM2.5浓度也就是提取对应label
y np.empty([12 * 471, 1], dtypefloat)for month in range(12): # month 0-11for day in range(20): # day 0-19for hour in range(24): # hour 0-23if day 19 and hour 14: # 第20天的第23时continue# 取对应month行都要取列取9个依次进行最后将整个数据reshape成一行数据# (列数无所谓)。然后赋给xx内的坐标只是为了保证其从0-471*12# vector dim:18*9x[month * 471 day * 24 hour, :] month_data[month][:, day * 24 hour: day * 24 hour 9].reshape(1, -1)# value,结果对应的行数一直是9行即第10行PM2.5然后列数随着取得数据依次往后进行y[month * 471 day * 24 hour, 0] month_data[month][9, day * 24 hour 9]mean_x np.mean(x, axis0) # 18 * 9 求均值axis 0表示对各列求值返回 1* 列数 的矩阵
std_x np.std(x, axis0) # 18 * 9 求标准差axis 0表示对各列求值返回 1* 列数 的矩阵
for i in range(len(x)): # 12 * 471for j in range(len(x[0])): # 18 * 9if std_x[j] ! 0:x[i][j] (x[i][j] - mean_x[j]) / std_x[j]# 将训练数据拆成训练数据验证数据8:2这样用来验证
import mathx_train_set x[: math.floor(len(x) * 0.8), :]
y_train_set y[: math.floor(len(y) * 0.8), :]
x_validation x[math.floor(len(x) * 0.8):, :]
y_validation y[math.floor(len(y) * 0.8):, :]# 用来做参数vector的维数加1是为了对bias好处理还有个误差项。即h(x)w1x1w2x2wnxnb
dim 18 * 9 1
# 生成一个dim行1列的数组用来保存参数值
w np.ones([dim, 1])
# np.ones来生成12*471行1列的全1数组np.concatenateaxis1
# 表示按列将两个数组拼接起来即在x最前面新加一列内容之前x是12*471行
# 18*9列的数组新加一列之后变为12*471行18*91列的数组
x np.concatenate((np.ones([12 * 471, 1]), x), axis1).astype(float)
learning_rate 100 # 学习率
iter_time 10000 # 迭代次数# 生成dim行即163行1列的数组用来使用adagrad算法更新学习率
adagrad np.zeros([dim, 1])# 因为新的学习率是learning_rate/sqrt(sum_of_pre_grads**2),
# 而adagradsum_of_grads**2,所以处在分母上而迭代时adagrad可能为0
# 所以加上一个极小数使其不除0
eps 0.0000000001
for t in range(iter_time):# rmse loss函数是从0-n的(X*W-Y)**2之和/(471*12)再开根号# 即使用均方根误差(root mean square error),具体可看公式# /471/12即/N(次数)loss np.sqrt(np.sum(np.power(np.dot(x, w) - y, 2)) / 471 / 12)if (t % 100 0): # 每一百次迭代就输出其损失print(str(t) : str(loss))# dim*1 x.transpose即x的转置后面是X*W-Y,即2*(x的转置*(X*W-Y))是梯度# 具体可由h(x)求偏微分获得.最后生成1行18*91列的数组。转置后的X其每一行# 是一个参数与h(x)-y的值相乘之后是参数W0的修正值同理可得W0-Wn的修正值# 保存到1行18*91列的数组中即gradient# 计算梯度。梯度是损失函数对参数向量w的偏导数gradient 2 * np.dot(x.transpose(), np.dot(x, w) - y)# adagrad用于保存前面使用到的所有gradient的平方进而在更新时用于调整学习率adagrad gradient ** 2w w - learning_rate * gradient / np.sqrt(adagrad eps) # 更新权重
np.save(weight.npy, w) # 将参数保存下来w np.load(weight.npy)
# 使用x_validation和y_validation来计算模型的准确率因为X已经normalize了
# 所以不需要再来一遍只需在x_validation上添加新的一列作为bias的乘数即可
x_validation np.concatenate((np.ones([1131, 1]), x_validation), axis1).astype(float)
ans_y np.dot(x_validation, w)
loss np.sqrt(np.sum(np.power(ans_y - y_validation, 2)) / 1131)
print(loss)testdata pd.read_csv(Dataset/test.csv, headerNone, encodingbig5)
test_data testdata.iloc[:, 2:] # 取csv文件中的全行数即第3列到结束的列数所包含的数据
test_data.replace(NR, 0, inplaceTrue) # 将testdata中的NR替换为0
test_data test_data.to_numpy() # 将其转换为数组
# 创建一个240行18*9列的空数列用于保存textdata的输入
test_x np.empty([240, 18 * 9], dtypefloat)
for i in range(240): # 共240个测试输入数据test_x[i, :] test_data[18 * i: 18 * (i 1), :].reshape(1, -1)
# 下面是Normalize,且必须跟training data是同一种方法进行Normalize
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前面拼接一列全1数组构成240行163列数据
test_x np.concatenate((np.ones([240, 1]), test_x), axis1).astype(float)# 进行预测
w np.load(weight.npy)
ans_y np.dot(test_x, w) # test data的预测值ans_ytest_x与W的积# 将预测结果填入文件当中
import csvwith open(submit.csv, modew, newline) as submit_file:csv_writer csv.writer(submit_file)header [id, value]print(header)csv_writer.writerow(header)for i in range(240):row [id_ str(i), ans_y[i][0]]csv_writer.writerow(row)print(row)
