Python实现线性回归算法，只使用python中的pandas和numpy库

提示：写完文章后，目录可以自动生成，如何生成可参考右边的帮助文档

---

前言

这一篇来介绍以下最简单的模型就是线性回归，这是很多机器学习的起步算法，虽然简单，但是通过Python最基础的数据结构去实现，只是用Numpy,Matplotlib 基础的库包去实现，有利于提高我们的编程能力，接下来让我们一起来编写代码。
当然，你想要读懂以下内容，你需要学习过Python基础，并有认识机器学习的一些基础概念，如目标函数、标准化、过拟合、梯度下降、评价指标等。也许你也可以在阅读过程中一步一步去认识它们。

---

一、线性回归介绍

线性回归就是用一条线来解释自变量与因变量之间的关系。目标函数可以采用两种方式进行求解：最小二乘法和梯度下降法

二、最小二乘法实现线性回归

最小二乘法知识链接：最小二乘法

1.引入库

import numpy as np
import pandas as pd

2.读入数据

数据来源：boston.csv

data = pd.read_csv(r"boston.csv")
data.head()

data= data.drop('Unnamed: 0',axis=1)
data.head()

波士顿房价数据集字段书名

• CRIM 房屋所在镇犯罪率
• ZN 面积大于25000平方英尺住宅所占比例
• INDUS 房屋所在镇非零售区所占比例
• CHASS 房屋是否是位于河边的，河边为1，否则为0
• NOX 一氧化氮浓度
• RM 平均房间数量
• AGE 1940年前建成房屋所占的比例
• DIS 房屋距离波士顿五大就业中心的加权距离
• RAD 距离房屋最近的公路
• TAX 税收额度
• PTRATIO 房屋所在镇师生比例
• BLACK 计算公式： 1000*(房屋所在镇非美籍人口所在比例 - 0.63)**2
• LSTAT 弱势群体人口所占比例
• MEDV 房间平均价格

3.编写线性回归实现类

class LinearRegression:
    """
    使用python是实现线性回归(最小二乘法)
    
    """
    def fit(self,X,y):
        
        
        """
        Parameters
        ----
        X:类数组类型。形状：[样本数量，特征数量]
        特征矩阵，用来对模型进行训练。
        y:类数组类型，形状：[样本数量]
        """
        #转化为矩阵
        X = np.asmatrix(X.copy()) #说明，如果X是数组对象的一部分，而不是完整的对象数据（例如，X是由其他对象通过切片传递过来的）
        #注意，我们现在要进行矩阵的运算，因此需要的是二维的结构，通过reshape方法进行转换
        y = np.asmatrix(y).reshape(-1,1)  #y为一维结构（行向量或列向量），可以不进行拷贝
        #通过最小二乘法公式，求解出最佳的权重值
        self.v_ = (X.T*X).I*X.T*y
        
    def predict(self,X):
        """
        根据参数传递的样本X,对样本数据进行预测
        X:类数组类型。形状：[样本数量，特征数量]
          待预测的样本特征（属性）
          
        Returns
        ----
        result:数组类型
        预测的结果
        """
        X = np.asmatrix(X.copy())
        result = X*self.v_
        #将矩阵转换为ndarray数组，进行扁平化处理,使用ravel()
        return np.array(result).ravel()
    
        
    
    

4.割分数据集，训练模型，测试模型

#1.不考虑截距情况
t = data.sample(len(data),random_state = 0)
train_X = t.iloc[:400,:-1]
train_y = t.iloc[:400,-1]
test_X = t.iloc[400:,:-1]
test_y = t.iloc[400:,-1]

lr = LinearRegression()
lr.fit(train_X,train_y)
result = lr.predict(test_X)
result

5.评价指标：均方误差

#均方误差
display(np.mean((result - test_y)**2))

#查看模型权重值
display(lr.v_)

#2.考虑截距，增加一列，该列所有值为1
t = data.sample(len(data),random_state = 0)
#可以这样增加一列
# t["Intercept"] = 1
new_columns = t.columns.insert(0, "Intercept")
t = t.reindex(columns = new_columns,fill_value=1)

train_X = t.iloc[:400,:-1]
train_y = t.iloc[:400,-1]
test_X = t.iloc[400:,:-1]
test_y = t.iloc[400:,-1]

lr = LinearRegression()
lr.fit(train_X,train_y)
result = lr.predict(test_X)
result

6.可视化

import matplotlib as mpl
import matplotlib.pyplot as plt
#可显示中文
mpl.rcParams['font.family'] = "SimHei"   #中文显示
mpl.rcParams['axes.unicode_minus'] = False  #中文字体下，可以显示中文符号

#设置画布大小
plt.figure(figsize=(10,10))
plt.plot(result,"ro-",label="预测值")
plt.plot(test_y.values,"go--",label="真实值")
plt.title("线性回归预测-最小二乘法")
plt.xlabel("样本序号")
plt.ylabel("房价")
plt.legend()
plt.show()

三、梯度下降法实现线性回归

梯度下降知识：梯度下降

数据处理过程与上述大部分步骤一致，只是对于目标函数的优化不一样。
因此，只保留重要步骤，以上代码进行复用即可‘

如果你学过线性代数，你就会知道，

• 需要说明的是，梯度下降法和最小二乘法不同，因为最小二乘法需要求出矩阵的逆，而不是所有矩阵都可以求出逆矩阵。相对最小二乘法，梯度下降法局限性更小一些。

1.创建实现类

class LinearRegression2:
    """
    使用Python语言实现线性回归算法（梯度下降）
    """
    def __init__(self,alpha,times):
        """
        初始化方法
        Parameter
        ----
        alpha:float
             学习率。用来控制步长（权重调整的幅度）
        times:int
             循环迭代的次数。
        
        """
        self.alpha = alpha
        self.times = times
    def fit(self,X,y):
        """
        根据提供的数据，对模型进行训练。
        
        Parameter
        ----
        X:类数组类型。形状：[样本数量，特征数量]
        待训练的样本特征属性（特征矩阵）
        
        y:类数组类型。形状：[样本数量]
        目标值（标签信息）
        """
        X = np.asarray(X)
        y = np.asarray(y)
        #创建权重的向量，初始值为0（或任何其他值）：长度比特征数量多一
        self.v_=np.zeros(1+X.shape[1])
        #创建损失列表，来保存每次迭代后的损失值，损失值计算：（预测值-真实值）的平方和除以2
        self.loss_=[]
        #进行循环，多次迭代。在每次迭代过程中，不断去调整权重值，使得损失值不断减少。
        for i in range(self.times):
            #计算预测值
            y_hat = np.dot(X,self.v_[1:])+self.v_[0]
            #计算真实值与预测值之间的差距
            error = y-y_hat
            #将损失值加入到损失列表当中
            self.loss_.append(np.sum(error**2)/2)
            #根据差距调整权重v_,根据公式，调整为权重（j） = 权重（j）+学习率*sum((y-y_hat)*x(j))
            self.v_[0]+=self.alpha*np.sum(error)
            self.v_[1:]+=self.alpha*np.dot(X.T,error)
    def predict(self,X):
        """
        根据参数传递的样本，对样本数据进行预测。
        Parameter
        ----
        X:类数组类型。形状：[样本数量，特征数量]
        待预测的样本特征属性（特征矩阵）
        Returns
        ----
        result:数组类型
        预测的结果
        """
        X = np.asarray(X)
        result = np.dot(X,self.v_[1:])+self.v_[0]
        return result

2.分割数据，训练模型，测试模型

lr1 = LinearRegression2(alpha=0.001,times=20)
t = data.sample(len(data),random_state=0)
train_X = t.iloc[:400,:-1]
train_y = t.iloc[:400,-1]
test_X = t.iloc[400:,:-1]
test_y = t.iloc[400:,-1]

lr1.fit(train_X,train_y)
result = lr1.predict(test_X)
display(np.mean((result-test_y)**2))

那问题来了？为什么数据如此之大？？？？
我们忽略了数据标准化，每一个字段之间的数据差异太大了
所以告诉我们，在梯度下降中，如果每一列的数据大小差异较大，要进行标准化

3.添加数据标准化过程

#数据标准化类
class StadardScaler:
    """
    该类对数据进行标准化处理
    """
    def fit(self,X):
        
    
        """
        根据传递的样本，计算每一个特征列的均值和标准差
        Parameter
        ----
        X；类数组类型
         训练数据，用来计算均值和标准差
        """
        X= np.asarray(X)
        self.std_ = np.std(X,axis=0)  #按列计算
        self.mean_ = np.mean(X,axis=0)
    def transform(self,X):
        """
        对给定的数据X，进行标准化处理。（将X的每一列都变成标准正态分布）
        
        Parameter
        ----
        X：类数组类型
           待转换的数据
           
        Returns
        ----
        result:类数组类型。
        参数X转换成标准正态分布后的结果
        """
        
        return (X-self.mean_)/self.std_
    def fit_transform(self,X):
        """
        对数据进行训练，并转换，返回转换之后的结果
        
        Parameter
        ----
        X：类数组类型
        待转换的数据
        
        Returns
        ----
        result:类数组类型
        参数转换成标准正态分布
        """
        
        self.fit(X)
        return self.transform(X)

重新进行模型训练和测试

#为了避免每一个特征数量级的不同，从而在梯度下降的过程中带来的影响
#我们现在考虑对每一个特征进行标准化处理
lr2 = LinearRegression2(alpha=0.0005,times=20)
t = data.sample(len(data),random_state=0)
train_X = t.iloc[:400,:-1]
train_y = t.iloc[:400,-1]
test_X = t.iloc[400:,:-1]
test_y = t.iloc[400:,-1]

#对数据进行标准化处理
s = StadardScaler()
train_X = s.fit_transform(train_X)
test_X = s.transform(test_X)

s2 = StadardScaler()
train_y = s2.fit_transform(train_y)
test_y = s2.transform(test_y)
lr2.fit(train_X,train_y)
result = lr2.predict(test_X)
display(np.mean((result-test_y)**2))

这次均方误差数据就正常了。

4.可视化

#查看权重
display(lr2.v_)

#查看损失值
display(lr2.loss_)

#设置画布大小
plt.figure(figsize=(10,10))
plt.plot(result,"ro-",label="预测值")
plt.plot(test_y.values,"go--",label="真实值")
plt.title("线性回归预测-梯度下降法")
plt.xlabel("样本序号")
plt.ylabel("房价")
plt.legend()
plt.show()

#累计的误差值 loss_
plt.plot(range(1,lr2.times+1),lr2.loss_,"o-")

# 因为房价更新涉及多个维度，不方便可视化
# 为了可视化，我们只选择其中一个维度（RM），并画出直线，进行拟合
lr2 = LinearRegression2(alpha=0.0005, times=20)
t = data.sample(len(data), random_state=0)
train_X = t.iloc[:400, 5:6]
train_y = t.iloc[:400, -1]
test_X = t.iloc[400:, 5:6]
test_y = t.iloc[400:, -1]
# 标准化
ss = StadardScaler()
train_X = ss.fit_transform(train_X)
test_X = ss.transform(test_X)
ss2 = StadardScaler()
train_y = ss2.fit_transform(train_y)
test_y = ss2.transform(test_y)
lr2.fit(train_X, train_y)
result = lr2.predict(test_X)
# display(result)
display(np.mean((result - test_y)**2))

意思是我们只选取一个自变量，一个因变量。因为本来自变量是很多的，就是有很多列的指标，只选一列来训练模型。

plt.scatter(train_X["rm"],train_y)
#查看方程
display(lr2.v_)
#构建方程 y = -2.6134650e-16 +6.4744239e-01 *x
x = np.arange(-5,5,0.1)
y = -2.6134650e-16 +6.4744239e-01 *x
# plt.plot(x,y,"r")
#也可以这样绘制
# x.reshape(-1,1) 把一维转为二维 
plt.plot(x,lr2.predict(x.reshape(-1,1)),"r")