作业一-手写Linear-Regression&Gradient Descent预测PM2.5

前情提要:

image

image
本文刚过strong baseline,入门水平


第一个作业是手写线性回归+梯度下降预测PM2.5,某acm退役人终于不能拿不会sklearn当做不写的理由了(((

当然纯C++选手之前完全没有写过python,也没写过机器学习算法,需要参考课程范例以及一些前人的经验

主要参考课程样例,差不多是跟着这个做的(为什么进不去懂的都懂)

部分参考这个博客


作业要求
image

image

数据集

image
每个月记录20天,每天记录20个小时,总共有18个特征(PM2.5包括在其中)
image
需要用前9个小时预测第十个小时的PM2.5监测值

数据预处理

首先读取数据集

import pandas as pd
import numpy as np
data = pd.read_csv(r"./train.csv")

发现读了个倒着的进来(读了繁体字),但.columns和.values没有问题,那数据也没问题,可以用朴素的代码去解决这个问题(

image

尝试着把有问题的列名改了,然后就好了:

data.rename(columns={data.columns[0]:'日期',data.columns[1]:'测站',data.columns[2]:'测项'},inplace=True)

image


首先把RAINFALL里面的数据看看:

s = set()
for i in range(0,len(data.values)):
    if data.values[i][2]=='RAINFALL':
        for j in range(3,len(data.values[i])):
            s.add(data.values[i][j])

image

发现除了'NR'以外其他均为实数,把'NR'变成0

#data[data == 'NR'] = 0
for i in range(0,len(data.values)):
    if data.values[i][2] == 'RAINFALL':
        for j in range(3,len(data.values[i])):
            if data.values[i][j] == 'NR':
                data.values[i][j] = 0

这样就把RAINFALL里面的非实数消除了
image


然后就是对表格中特征的提取了。

9个小时里总共的特征数为\(9*18=162\),目标值为第10个小时的PM2.5

对于每一天我们有24个小时的观测数据\([0h,23h]\)

一个月观测20天,也就是480小时(\([0h,479h]\)),总共有特征\(18*480 = 23040\)个。

那么只取\([0h,8h]\&[10h,18h]\&[20h,28h]\)。。。作为\(x\),而选取\(9h\)\(19h\)\(29h\)。。。的PM2.5作为\(y\)这样选取样本是不是太少了?每天只有2组样本(最多3组)

解决方案是通过滑动窗口选取:\([0h,8h] + 9h\)\([1h,9h] + 10h\)\([2h,10h] + 11h\)\([3h,11h] + 12h\)。。。。。。

这样一个月480个小时,我们就可以选取\(480-9 = 471\)组样本了(由于只观测20天,与下个月时间上不连续,所以只取到\([470,478] + 479h\))

官方的范例图示就是上面表达意思(一年12个月,也就是\(18*480*12\),数据只给了2014年):
image

接下来就是代码实现了

想到要把数据分成\(12\)组,每组\(18*480\),原表格中格式好像完全规整,那就没时间、特征名啥事了

1、只留下数据列

data = data.iloc[:,3:]
raw_data = data.to_numpy()

image

2、分组(把每个月的数据块(二维数组)放入字典)

开一个字典,保存每个月的观测值(过程如下图)

month_data={}
for month in range(12):
    tmp = np.empty([18,480]) # 创建18*480的二维数组
    for day in range(20):
        tmp[:, 24 * day : 24 * (day + 1)] = raw_data[18 * (month * 20 + day):18 * (month * 20 + day + 1),:]
    month_data[month] = tmp

image

image

3、滑动窗口取特征

首先构造数组:

对于x数组:有18个特征,每组有9个小时,因此有\(18 * 9 = 162\)列;总共有12个月,每个月可以取471组(480-9=471,之前写了为什么),因此有\(12 * 471\)

对于y数组:当然只有1列,即目标值"第十小时的PM2.5";行数与x数组相同

x = np.empty([471 * 12, 18 * 9], dtype = float)
y = np.empty([471 * 12, 1], dtype = float)
col_cnt = 0
for month in range(12):
    for day in range(20):
        for hour in range(24):
            if day==19 and hour > 14: 
                break
            x[col_cnt,:] = month_data[month][:, day * 24 + hour : day * 24 + hour + 9].flatten() # 按列切片(month_data[month]为18*480)
            y[col_cnt,0] = month_data[month][9, day * 24 + hour + 9] # 第10行为PM2.5的数据
            col_cnt += 1

image

这样我们就完成了数据的预处理,得到x数组和y数组

image

Normalize

特征工程中的「归一化」有什么作用?- 知乎
采用的是Z-score(零均值规范化): \(X_{after} = \frac{X - X_{mean}}{\sigma}\)
标准化后的变量值围绕0上下波动,大于0说明高于平均水平,小于0说明低于平均水平。

image

mean_x = np.mean(x, axis=0) # 按列(竖着求均值)
std_x = np.std(x, axis=0) # 标准差同上

image
image

然后进行normalize

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

image

Split training data into "train" and "validation"

有现成的sklean.train_test_split可以用

from math import floor
# 分前80%为train_set, 后20%为validation_set
x_train = x[:floor(len(x)*0.8),:]
y_train = y[:floor(len(y)*0.8),:]
x_vali = x[floor(len(x)*0.8):,:]
y_vali = y[floor(len(y)*0.8):,:]

image

Train

先解释这个图什么意思:

对于每次迭代,\(y' = x \cdot w\),然后得到\(Loss = y' - y\)
image
gradient为\(2x^T \cdot Loss\)(如下),
image
然后让\(w_{i+1} = w_i - {learning\_rate} * gradient\)

不过样例演示使用的是adgrad算法:

该算法的思想是独立地适应模型的每个参数:具有较大偏导的参数相应有一个较大的学习率,而具有小偏导的参数则对应一个较小的学习率 具体来说,每个参数的学习率会缩放各参数反比于其历史梯度平方值总和的平方根

伪代码如下:
image
对于每次迭代,前面的步骤一直到5.1都和之前没有区别(但代码里Loss改为了均方根误差rmse:\(rmse = \sqrt{\frac{1}{m} \sum\limits_{i=1}^{m}(y_i-\hat{y_i})^2}\)

然后令\(prev\_grad\) += \(gradient^2\)\(ada = \sqrt{prev\_gra}\)

最后\(w_{i+1} = w_i - \frac{{learning\_rate} * gradient}{ada}\)

其中\(prev\_grad\)初始化为一个全1的与\(gradient\)维数相同的向量

dim = 18 * 9 + 1 # 维度dimension,加上1个常数项
w = np.zeros([dim,1]) # dim行1列,全部填0
pregrad = np.zeros([dim, 1])
x = np.concatenate((np.ones([12 * 471, 1]), x), axis = 1).astype(float) # 同样每一行补上常数项
learning_rate = 100 # 学习率
iter_time = 1000 # 迭代次数
eps = 1e-8 # 避免分母ada为0而加上的极小项
for t in range(iter_time):
    loss = np.sqrt(np.sum(np.power(np.dot(x, w) - y, 2))/471/12) # rmse = y - y_hat求和的平方取平均开根号
    if(t % 100 == 0):
        print(str(t) + ":" + str(loss)) # 查看每100次迭代当前的loss值
    gradient = 2 * np.dot(x.transpose(), np.dot(x,w) - y) #套公式
    pregrad += gradient ** 2 #套公式
    ada = np.sqrt(pregrad + eps) #套公式
    w = w - learning_rate * gradient / ada #套公式

np.save('weight.npy', w) # 可以将w保存在本地

发现loss(rmse)从27逐渐减小到7

image

同时得到train结果w:

image

Testing

首先预处理test_set,流程同train的预处理:

testdata = pd.read_csv("./test.csv",header=None) # 没有列index
test_data = testdata.iloc[:, 2:] # 把id_0和AMB_TEMP去掉
test_data[test_data=='NR'] = 0
test_data = test_data.to_numpy()
test_x = np.empty([240, 18 * 9], dtype = float) # id_0~id_239, 240组测试数据
for i in range(240):
    test_x[i, :] = test_data[18 * i:18 * (i + 1), :].flatten()
    
# Normalize,用之前的mean_x和std_x
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)

得到test_x:

image

Predict

最后一步就是预测了,这一步直接\(y\_test = x\_test \cdot w\)即可

ansy = np.dot(test_x,w)

image

最后把得到的结果用csv保存即可(不写了直接复制):

import csv
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)

这样直接提交public只能达到simple baseline
image

凭借数学直觉,主观感觉需要增加迭代次数和降低学习率,然后过了public的strong baseline:

image

image

但不知道private会怎么样,感觉还会有很多优化空间,多上几节课学点东西再回来看看。。。。。。

推荐这些文章:

Conjugate Gradient Descent (1) - Python实现

算法特征:①. 将线性方程组等效为最优化问题; ②. 以共轭方向作为搜索方向.
算法推导:Part Ⅰ 算法细节现以如下线性方程组为例进行算法推导,
\begin{equation}Ax = b\label{eq_1}\end{equation}如​上式$\eqref{eq_1}$解存在, 则等效如下最优化问题,
\begin{equation}\min\quad \frac{1}{2}\|Ax - b\|_2^2 \quad \Rightarrow\quad \min\quad \frac{1}{2}x^\mathrm{T}A^\mathrm{T}Ax - b^\mathrm{T}Ax...

手写linear-regression预测PM2.5

作业要求 数据集 每个月记录20天,每天记录20个小时,总共有18个特征(PM2.5包括在其中) 需要用前9个小时预测第十个小时的PM2.5监测值 数据预处理 首先读取数据集 import pandas as pd import numpy as np data = pd.read_csv(r"C:UsersAdministratorDesktop机器学习深度学习day1hw1hw1train.csv") 发现读了个倒着的进来(读了繁体字),但.columns和.values没有问题,那数据也没问题,可以用朴素的代码去解决这个问题( 尝试着把有问题的列名改了,然后就好了: dat...

手写linear-regression预测PM2.5

作业要求 数据集 每个月记录20天,每天记录20个小时,总共有18个特征(PM2.5包括在其中) 需要用前9个小时预测第十个小时的PM2.5监测值 数据预处理 首先读取数据集 import pandas as pd import numpy as np data = pd.read_csv(r"C:UsersAdministratorDesktop机器学习深度学习day1hw1hw1train.csv") 发现读了个倒着的进来(读了繁体字),但.columns和.values没有问题,那数据也没问题,可以用朴素的代码去解决这个问题( 尝试着把有问题的列名改了,然后就好了: dat...

手写linear-regression预测PM2.5

作业要求 数据集 每个月记录20天,每天记录20个小时,总共有18个特征(PM2.5包括在其中) 需要用前9个小时预测第十个小时的PM2.5监测值 数据预处理 首先读取数据集 import pandas as pd import numpy as np data = pd.read_csv(r"C:UsersAdministratorDesktop机器学习深度学习day1hw1hw1train.csv") 发现读了个倒着的进来(读了繁体字),但.columns和.values没有问题,那数据也没问题,可以用朴素的代码去解决这个问题( 尝试着把有问题的列名改了,然后就好了: dat...

深度学习1:Regression-CaseStudy(预测宝可梦进化后CP值)

Step1:设立model   例如Y=b+w*Xcp     A set of  function
 
Step2:Goodness of Function   定义loss function(损失函数L)来评估function的好坏
 
通过training data评估后挑选出最好的function(f*),用testing data来评估误差 ,注意可能会出现overfitting(过度拟合)的情况。
 
 
Step3:Gradient Descent ...

人工智能作业2

import pandas as pdimport numpy as npimport matplotlib.pyplot as pltdef sigmoid(x): # 定义网络激活函数 return 1/(1+np.exp(-x))data_tr = pd.read_csv('3.3 data_tr.txt') # 训练集样本data_te = pd.read_csv('3.3 data_te.txt') # 测试集样本n = len(data_tr)yita = 0.1 # 自己设置学习率out_in = np.array([0.0, 0, 0, 0, -1]) ...

机器学习中的优化 Optimization Chapter 2 Gradient Descent(1)

1. Step of Gradient descent
\[\begin{equation}
x_{t+1} = x_t-\gamma \nabla f(x_t)
\end{equation}
\]2. Vanilla Analysis
\(\text{Let }{\bf g_t} = \nabla f(x_t)\),\(\text{ therefore we can get:}\)
\[\begin{equation}
g_t = (x_t-x_{t+1})/\gamma
\end{equation}
\]\(\text{hence we get:}\)
\[\begin{...

机器学习中的优化 Optimization Chapter 2 Gradient Descent(2)

\(\large \bf{Theorem }\ 2.7:\)
\(f:\mathbb{R^d}\rightarrow\mathbb{R}\text{ be convex and differentiable with a global minimum }x^*;\text{ Suppose }f\text{ is smooth with parameter }L.\text{ Choosing stepsize: }\gamma = \frac{1}{L},\text{ gradients descent yields:}\)
\[\begin{align}
f(x_T)-f(x^*...

第12次作业--你的生日

题目:利用Calendar类计算自己的出生日期距今天多少天,再将自己的出生日期利用SimpleDateFormat类设定的格式输出显示。
 

package com;

import java.text.SimpleDateFormat;
import java.util.*;

public class A {

public static void main(String[] args) {
Scanner input = new Scanner(System.in);
Calendar calendar = Calendar.getInstance();
Sy...

数据分析课后作业--企业所得税分析预测模型(代码)

1.求取企业所得税各特征间的相关系数(1)求取原始数据特征之间的Pearson相关系数。(2)判断各特征之间的相关性。

#求取企业所得税各特征间的相关系数
import numpy as np
import pandas as pd
inputfile = 'income_tax.csv' #读取数据文件
data = pd.read_csv(inputfile) #读取数据
#输出Pearson相关系数,并保留两位小数
print('相关系数矩阵为:','\n',np.round(data.iloc[1:,1:].corr(method = 'pearson'), 2))

由...

文章标题:作业一-手写Linear-Regression&Gradient Descent预测PM2.5
文章链接:https://www.dianjilingqu.com/4171.html
本文章来源于网络,版权归原作者所有,如果本站文章侵犯了您的权益,请联系我们删除,联系邮箱:saisai#email.cn,感谢支持理解。
THE END
< <上一篇
下一篇>>