基于python实现梯度下降法（GD）在线性回归中的应用（Boston房价预测数据集）

【环境准备】

用到的包主要有两个：numpy和sklearn，都是机器学习常用的库。

 

【数据集介绍】

波士顿房价数据集（Boston House Price Dataset）

使用sklearn.datasets.load_boston即可加载相关数据。该数据集是一个回归问题。每个类的观察值数量是均等的，共有 506 个观察，13 个输入变量和1个输出变量。

在这里我将数据集视作线性关系，即

所以项目的目标就是求出该数据集十三个参数的适配系数。

 

【方法介绍】

假如我们要优化一个函数 f(x) ，即找到它的最小值, 常用的方法叫做Gradient Descent (GD), 也就是最速下降法. 

注：图中可以认为是含有两个未知数的损失函数即loss=loss(a,b)

 

宏观的想一想，我们现在有十三个待求系数，因此这十三个数是我们构建的函数的未知数，而函数的形式则取决于现有的数据集。

比如待求关系为f(x,y,z)=x+2y-z（但是未知）, 我们此时有（1,1,1，2），（1，0,0,1），（1,0,1,0），（0,0,0,0）等一些已知数据，需要求出函数关系中的系数，就需要构建一个loss函数，如

就可以求出loss对于a/b/c的偏导，即最大梯度的概念（具体梯度下降法的思想可以在其他文章中学习）

因此我在这里就是构建了一个含有十三个未知数的损失函数，使用梯度下降法收敛到最小值。

这里提一句个人理解，梯度下降法的缺陷是有可能陷入局部最优而非全局最优（比如上面图中不同的“坑”，可能进了一个坑就出不来了收敛在这个坑的底部，但并非全局的最优解），不过这种限于较复杂的loss函数，我个人理解是对于这样的简单线性关系，loss采用二阶损失时，并不存在多个梯度为0的点，因此这种情况下求出来的应该就是全局最优（如果有不对的地方请一定指出！）

 

【代码实现】

import numpy as np
import sklearn.datasets as datasets

#全局常量定义
max_iteration=100000#最大循环次数


#数据载入 
data=datasets.load_boston()
#print(data.feature_names)
X=data.data
y=data.target
scale=data.data.shape[0]

print(scale)#数据条数
print(data.feature_names)#数据意义
print(X[0])

def loss(param):
    #param在这里代表预测模型的所有参数
    l=0
    for i in range(0,X.shape[0]):
        l=l+(np.dot(X[i], param)-y[i])**2
    return l

def lossPar(param,paramNum):
        par=0
        i=0 
        for i in range(scale):
            par=par+2*(np.dot(param, X[i])-y[i])*X[i][paramNum]
        return par

if __name__ == "__main__":
    param=[0 for _ in range(X.shape[1])]
    #param=[-0.08773350967026737, 0.08390719227658551, -0.002136708832591626, 0.26581396762237985, 0.12646983121925526, 2.9558609343418913, 0.05770653081727748, -0.4251674577386858, 0.1526937339473097, -0.008546558464844483, 0.20400742795045249, 0.02092363236316196, -0.6701843284790404]
    print(param)
    epoch=0.00000001
    l=loss(param)
    i=0
    while l > 100 and i < max_iteration:#损失函数误差阈值
        j=0
        for j in range(X.shape[1]):
            param[j]=param[j]-epoch*lossPar(param,j)
            
        print("param is:",param)    
        l=loss(param)
        print("count: ",i)
        print("loss: ",l)
        i=i+1
    #print(loss)

【结果分析】

代码中有一段被我注释掉的部分，就是10w次后得到的参数结果

[-0.08773350967026737, 0.08390719227658551, -0.002136708832591626, 0.26581396762237985, 0.12646983121925526, 2.9558609343418913, 0.05770653081727748, -0.4251674577386858, 0.1526937339473097, -0.008546558464844483, 0.20400742795045249, 0.02092363236316196, -0.6701843284790404]

这里关注几个数值比较大的，说明算法认为这几个数据对房价影响较大

序号（从0开始计）	参数意义	运行结果
5	RM：住宅平均房间数	2.9558609343418913
7	DIS：到波士顿五个中心区域的加权距离。	-0.4251674577386858
3	LSTAT：人口中地位低下者的比例。	-0.6701843284790404

 

定性分析判断结果是否可信：

RM（住宅平均房间数）与房价应呈现较大关系的正相关，合理

DIS：到波士顿五个中心区域的加权距离与房价应呈现负相关，合理

LSTAT：人口中地位低下者的比例与房价应呈负相关，合理

 

因此基本认为运行结果具有可信度。

 

【问题反思】

loss跳出循环条件问题：这里的loss值小于100的判断条件出不去，循环10w次后最佳loss接近于16000左右，距离目标还差很多。我的想法是或许应该选择别的损失函数或者添加一定的缩放条件

训练率问题：epoch选择的值为0.00000001，其实按以往的经验是不太合理的，不过可能是数据本身的范围导致，学习率更大时会导致迅速“步子过大”而无法找到最佳点，loss会迅速指数级爆炸（梯度下降法常见问题），如果学习率更小会使收敛速度更慢。