李宏毅机器学习2020 作业1：Regression

Shayndel-

已于 2022-07-19 09:12:57 修改

阅读量606

点赞数

分类专栏：课程学习文章标签：机器学习 python

于 2022-07-18 14:58:56 首次发布

本文链接：https://blog.csdn.net/Chral317/article/details/125848103

版权

课程学习专栏收录该内容

1 篇文章 0 订阅

订阅专栏

前置课程

学完梯度下降（Gradient Descent1，2，3）

任务描述

本次作业的资料是从行政院环境环保署空气品质监测网所下载的观测资料。
希望大家能在本作业实作 linear regression 预测出 PM2.5 的数值

数据描述

本次作业使用丰原站的观测记录，分成 train set 跟 test set。

train set 是丰原站每个月的前 20 天所有资料。test set 则是从丰原站剩下的资料中取样出来。
train.csv: 每个月前 20 天的完整资料。
test.csv : 从剩下的资料当中取样出连续的 10 小时为一笔，前九小时的所有观测数据当作feature，第十小时的 PM2.5 当作 answer。一共取出 240 笔不重複的 test data，请根据 feature 预测这 240 笔的 PM2.5。
Data 含有 18 项观测数据 AMB_TEMP, CH4, CO, NHMC, NO, NO2, NOx, O3, PM10, PM2.5, RAINFALL, RH, SO2, THC, WD_HR, WIND_DIREC, WIND_SPEED, WS_HR。

未解决问题

Adagrad优化算法（梯度下降的课还没看
调整learning rate，iter_time (iteration 次数)来调整最后的值。
为什么答案就是权重矩阵w和测试数据test_x的乘积

思路及代码

代码非原创，是作业资料中给的。

读入数据文件（Load ’train.csv')

train.csv的资料是12个月每个月取20天每天24小时每小时有18个特征值的资料

采用pandas库里自带的read_csv函数读入

#打开训练文件，编码方式为：big5，根据不同数据可以改为UTF-8或GBK
data = pd.read_csv('./train.csv',encoding='big5')

预处理(preprocessing)

打开train.csv文件我们会发现，有几行的数据内容是NR，需要将降水这里的'NR'变为0

#iloc函数.iloc[a,b]a表示行b表示列 设定范围
#如果是范围则a:b,c:d表示[a,b-1]，从0开始，没有数值表示到边界
#将‘RAINFALL’全部变成0，原本是NR
data = data.iloc[:,3:]
data[data == 'NR'] = 0
raw_data = data.to_numpy() #将数据转换为一个numpy数组

原本train.csv中的文件排列方式为4320（12*20*18)*24的形式，为了方便后续计算，我们将其变为12*18*480（24*20）的形式

month_data = {}
for month in range(12):
    sample = np.empty([18,480]) #sample[18][480] = {0}
    for day in range(20):
        sample[:, day * 24 : (day + 1) * 24] = raw_data[18 * (20 * month + day):18 * (20 * month + day + 1),:]
    month_data[month] = sample

根据要求我们可以知道，令9个小时为feature，预测第十个小时的值。为了不浪费数据，我们连续计算。令9个小时为一组，类似滑动窗口的形式往下走，这样的话每个月480小时，每9个小时形成一个data，每个月就有471个data（以最后9个小时为data的第一位的无法成立，后续数据缺失）这样一共有471*12个data，并且每个data有9*18个feature（每小时有18个feature）

target是第10个小时的PM2.5，因为有471*12个data，就有471*12个target

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]
#print(x)
#print(y)

标准化(Normalize)

Feature Scaling也是在Gradient Descent中提到的方法，有利于提高Gradient Descent的效率

$z= \frac{x-\mu }{\sigma }$

是数据标准化的结果， $\mu$ 是数据的平均值， $\sigma$ 是数据的标准差（这里用无偏即除以N，若是有偏除以N-1）

$\sigma = \sqrt{\frac{1}{N}\sum_{i=1}^{N}(x_{i}-\mu )^{2}}$

axis说明：

axis不设置，对m*n个数求均值，返回一个实数

axis=0 压缩行，对各列取均值，返回1*n的矩阵

axis=1 压缩列，对各行求均值，返回m*1矩阵

#mean函数，取均值
#np.std 无偏标准差（一般用于总体） 
#pd.std有偏标准差（一般用于样本）
mean_x = np.mean(x, axis = 0) 
std_x = np.std(x, axis = 0)

#len(x)如果x是一个矩阵，返回值是多少行
for i in range(len(x)):
    for j in range(len(x[0])):
        if std_x[j] != 0:
            x[i][j] = (x[i][j] - mean_x[j]) / std_x[j]

训练(Training)

Adagrad优化算法（自适应学习率优化算法）比Gradient Descent快

这个算法在Gradient Descent的课中有讲到

本质上是为了得到权重矩阵（不知道原理，根据结果推出来的结论）


#zeros创建全0数组 [x,y]x行 y列 
#ones全1数组
dim = 18 * 9 + 1
w = np.zeros([dim, 1])

#concatenate沿现有的某个轴对一系列数组进行拼接 .astype(float)规定类型会float
x = np.concatenate((np.ones([12 * 471, 1]), x), axis = 1).astype(float)
learning_rate = 100 #学习率 参数(可以通过调整这个值调整模型）
iter_time = 1000 #迭代次数（也可以改）
adagrad = np.zeros([dim, 1]) #adagrad动态调整学习率的方法 初始化定义为0矩阵
eps = 0.0000000001 #精度

#w为权重矩阵 初始值是全0 后续根据学习率和梯度进行迭代
#x是train_x训练矩阵
#y是train_y全1矩阵 每次不变
#y' = x * w 
#L = y' - y
#dot函数 获取两个矩阵的乘积

#Adagrad算法 自适应学习率优化算法 加快深层神经网络训练速度 学完Gradient Descent
#核心：计算参数的梯度和平方

for t in range(iter_time):

    L = np.dot(x,w) - y #Loss矩阵
    loss = np.sqrt(np.sum(np.power(L, 2))/471/12)#rmse均方根误差的公式
    if(t % 100 == 0):#每100次输出
        print(str(t) + ":" + str(loss))
    #梯度 transpose()行列转置
    gradient = 2 * np.dot(x.transpose(), L) #dim*1 gt
    adagrad += gradient ** 2 #**幂运算 求gi的平方
    w = w - learning_rate * gradient / np.sqrt(adagrad + eps)
np.save('weight.npy', w)

测试(Testing)

载入测试集，进行计算

testdata = pd.read_csv('./test.csv', header = None, encoding = 'big5')

test_data = testdata.copy()
test_data[test_data == 'NR'] = 0 #将降雨置为0
test_data = test_data.iloc[:, 2:]
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)
#reshape(1,-1)表示将一个二维数组重组为一个一列二维数组

#此时的mean_x[j]和std_x[j]就是前面training里的
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)

预测（Prediction）

答案就是权重矩阵w和测试数据test_x的乘积（为啥？）

w = np.load('weight.npy')
ans_y = np.dot(test_x, w)

保存结果（Save Prediction to CSV File）

with open('submit.csv', mode='w', 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 sys
import csv
import pandas as pd
import numpy as np
import math as ms
import matplotlib.pyplot as plt

#load
data = pd.read_csv('./train.csv',encoding='big5')

#preprocessing
data = data.iloc[:,3:]
data[data == 'NR'] = 0
raw_data = data.to_numpy()

month_data = {}
for month in range(12):
    sample = np.empty([18,480]) #sample[18][480] = {0}
    for day in range(20):
        sample[:, day * 24 : (day + 1) * 24] = raw_data[18 * (20 * month + day):18 * (20 * month + day + 1),:]
    month_data[month] = sample

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]

#normalize 

mean_x = np.mean(x, axis = 0)
std_x = np.std(x, axis = 0)

for i in range(len(x)):
    for j in range(len(x[0])):
        if std_x[j] != 0:
            x[i][j] = (x[i][j] - mean_x[j]) / std_x[j]

#training
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 #精度

#Adagrad算法 自适应学习率优化算法
for t in range(iter_time):
    L = np.dot(x,w) - y 
    loss = np.sqrt(np.sum(np.power(L, 2))/471/12)#rmse均方根误差
    if(t % 100 == 0):#每100次输出
        print(str(t) + ":" + str(loss))
    gradient = 2 * np.dot(x.transpose(), L) 
    adagrad += gradient ** 2 #**幂运算
    w = w - learning_rate * gradient / np.sqrt(adagrad + eps)
np.save('weight.npy', w)


#testing
testdata = pd.read_csv('./test.csv', header = None, encoding = 'big5')
print(testdata)
#test_data = testdata.iloc[:, 2:]
test_data = testdata.copy()
test_data[test_data == 'NR'] = 0 #将降雨置为0
test_data = test_data.iloc[:, 2:]
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)

#prediction
w = np.load('weight.npy')
ans_y = np.dot(test_x, w)

#save the prediction
with open('submit.csv', mode='w', 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)