题目要求:
在上一节的梯度下降法的试验中,大家已经初步通过使用梯度下降法找出最合适的 θ 值实现线性拟合,现在,利用给出的数据集,同时用你上节课实现的梯度下降法和本次试验中的sklearn模块下的LinearRegression对给出的数据进行线性回归,比较两种方法的效果: 本次实验的数据集为2005年至2015年城镇公交车运营数量(Buses,辆)以及人均国民生产总值(PGDP,元),大家可以以2005-2012年的数据为训练集,2013-2015年的数据作为预测集进行模型的建立、预测与评估:
data数据
Year Bus PGDP
0 1996 281516 11232
1 1997 313296 14368
2 1998 315576 16738
3 1999 347969 20505
4 2000 370640 24121
5 2001 371822 26222
6 2002 383161 30876
7 2003 412590 36403
8 2004 432021 40007
9 2005 460970 43852
10 2006 476255 47203
11 2007 502916 50251
sklearn算法实现
%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
from sklearn.linear_model import LinearRegression
# 导入数据集
data = pd.read_csv('./Input/data.csv')
# buses: 城镇公交车运营数量
buses = data['Bus']
# pdgp: 人均国民生产总值
pgdp = data['PGDP']
#对数据进行预处理 参考上方data 可以知道 数据集一共12组
###训练数据集
buses_train = buses[:9] #前9个为训练集数据
pgdp_train = pgdp[:9]
###测试数据集
buses_test = buses[-3:] #后3个为测试数据集
pgdp_test = pgdp[-3:]
#sklearn 算法
#首先定义一个线性回归模型
model = LinearRegression()
#对训练集进行处理
buses_train_X = np.array(buses_train.values).reshape(-1,1)
'''
生成了如下所示的垂直(-1,1)中的-1决定的 是竖直方向维度是1的向量
[281516]
[313296]
[315576]
[347969]
[370640]
[371822]
[383161]
[412590]
[432021]]
'''
pgdp_train_y = np.array(pgdp_train.values).reshape(-1,1)
#填充线性回归模型
model.fit(buses_train_X,pgdp_train_y)
#用剩下的测试集来进行测试
buses_test_X = np.array(buses_test.values).reshape(-1,1)
for i in range(3):
print("201%d"%(i+3),
"年预测GDP为%.f"%model.predict(buses_test_X)[i],"真实值为",pgdp_test.values[i],
"与实际值差了%d"%abs((model.predict(buses_test_X)[i]-pgdp_test.values[i])))
# 首先画出真实数据,用点来描绘
X = np.array(buses.values).reshape(-1,1)
y = np.array(pgdp.values).reshape(-1,1)
plt.plot(X,y,'r.')
plt.plot(X,model.predict(X),'b-')
效果图
梯度下降算法实现:
#梯度下降算法
#设置超参数
lr = 0.1 #学习率
epochs = 900 #迭代次数
SLIP_NUM = int(0.75*len(buses))
#数据归一化
#归一化处理,就是原数据-数据最小值 / (数据最大值-数据最小值) 就是范围
def Normalization_data(data):
return (data-data.min())/(data.max()-data.min())
#定义误差损失函数
def error_function(θ,X,y):
diff = np.dot(X,θ)-y #矩阵做运算 形如 (θ_0 * x_0 + θ_1 * x_1) - y 得到误差
#误差矩阵的转置乘以误差矩阵,相当于误差的平方,然后每一个都乘以1/(2*(数组的个数))
return (1./(2*len(buses))) * np.dot(np.transpose(diff),diff)
#定义J(θ)函数
def J(θ,X,y):
diff = np.dot(X,θ)-y #矩阵做运算 形如 (θ_0 * x_0 + θ_1 * x_1) - y 得到误差
#矩阵X的转置乘以误差,相当于求导后的∇θ矩阵 * X矩阵 就是公式 ∇θ*x^(i)
return (1./SLIP_NUM) * np.dot(np.transpose(X),diff)
#定义梯度下降算法
def gradient_descending(x1,y1):
θ = np.array([1,1]).reshape(2,1) #初始化θ为竖立的一维矩阵
x_matrix_left = np.ones((SLIP_NUM,1))
x_matrix_right = np.array(x1[:SLIP_NUM]).reshape(SLIP_NUM,1)
x_matrix_merge = np.hstack((x_matrix_left,x_matrix_right)) #将两个矩阵左右拼接起来
'''
x_matrix_merge 矩阵形如:
[1. 0. ]
[1. 0.1435411 ]
[1. 0.15383921]
[1. 0.30014905]
[1. 0.40254743]
[1. 0.40788618]
[1. 0.45910117]
[1. 0.59202349]
[1. 0.67978771]
'''
y = y1[:SLIP_NUM].values.reshape(SLIP_NUM,1)
gradient = J(θ,x_matrix_merge,y) #获得梯度
loss = np.zeros(epochs) #定义损失矩阵
for i in range(epochs):
θ = θ - alpha * gradient #迭代θ
loss[i] = error_function(θ,x_matrix_merge,y) #计算损失值
gradient = J(θ,x_matrix_merge,y) #更新梯度
plt.plot(range(0,epochs),loss) #描绘出损失函数图像
plt.figure()
plt.show()
return θ
X = Normalization_data(buses) #数据归一化
Y = Normalization_data(pgdp)
θ = gradient_descending(X,Y) #得到θ
X_test_left = np.ones((len(buses)-SLIP_NUM,1))
X_test_right = X[SLIP_NUM:].values.reshape(-1,1)#获得后半段的测试数据集
X_test = np.hstack((X_test_left,X_test_right)) #拼接矩阵
print("X_test",X_test)
'''
[1. 0.81054201]
[1. 0.87957995]
[1. 1. ]
'''
print("θ",θ)
'''
[-0.04840497]
[ 1.11410302]
'''
Y_test = np.dot(X_test,θ)
print(Y_test)
Y_test = Y_test*(pgdp.max()-pgdp.min()) + pgdp.min()
print("预测值")
print(Y_test)
print("真实值")
print(pgdp[SLIP_NUM:].values)
plt.plot(Y_test,'g-')
plt.plot(pgdp[SLIP_NUM:].values,'r-')
plt.show()
运行结果: