Python机器学习之multiple_linear_regression（多元线性回归）

Python机器学习之multiple_linear_regression（多元线性回归）

实验介绍

1.实验内容
本实验介绍线性回归算法,并通过小实验简单认识一下线性回归算法

实验1：用线性回归找到最佳拟合直线
实验2：局部加权线性回归找到最佳拟合直线
实验3：使用scikit-learn实现线性回归算法
2.实验目标
通过本实验掌握线性回归算法找到最佳拟合直线的方法。
3.实验知识点
线性回归
4.实验环境
python 3.6.5
5.预备知识
概率论与数理统计
Linux命令基本操作
Python编程基础

实验准备

点击屏幕右上方的下载实验数据模块，选择下载multiple_linear_regression.tgz到指定目录下，然后再依次选择点击上方的File->Open->Upload,上传刚才下载的数据集压缩包，再使用如下命令解压：

!tar -zxvf multiple_linear_regression.tgz

用线性回归找到最佳拟合直线

接下来，我们使用普通最小二乘法线性回归来依据给定数据找到最佳拟合直线。
1.数据样例
数据样例为数据集目录下的ex0.txt文件，我们先来看下我们的数据及其分布：

!cat multiple_linear_regression/ex0.txt

其中第一个值总是等于1.0，即x0。我们假定偏移量就是一个常数，第二个值X1，也就是横坐标值，最后一列为Y值，即目标值。
2.数据加载及可视化

import matplotlib.pyplot as plt
import numpy as np

def loadDataSet(fileName):
    xArr=[];yArr=[]
    fr=open(fileName)
    for line in fr.readlines():
        lineArr=line.strip().split()#将每一行的元素变为list，strip删除的字符，按照split中的符号进行每行元素分割为list中的元素
        xArr.append([1.0,float(lineArr[1])])
        yArr.append(float(lineArr[-1]))
    return xArr, yArr
    
def plotDataSet():
    xArr,yArr=loadDataSet('multiple_linear_regression/ex0.txt')
    dataArr=np.array(xArr) #多维数组
    datayArr=np.array(yArr)
    n=np.shape(xArr)[0]
    xcord=[];ycord=[]
    for i in range(n):
        xcord.append(dataArr[i,1])#x坐标
        ycord.append(datayArr[i])#y坐标
    fig=plt.figure()
    ax=fig.add_subplot(111)
    ax.scatter(xcord,ycord,s=20,c='g',marker='s',alpha=.5)
    plt.title('dataset')
    plt.xlabel('x');plt.ylabel('y')
    plt.show()
    pass

if __name__ == '__main__':
    plotDataSet()

3.代码实现

from numpy import *
import matplotlib.pyplot as plt
#标准线性回归算法
def standRegres(xArr,yArr):
    xMat=mat(xArr)#转化为矩阵
    yMat=mat(yArr).T#求转置
    xTx=xMat.T*xMat
    #判断是否行列式为0
    if linalg.det(xTx) == 0.0:
        print("This matrix is singular")
        return
    #行列式不为0，可逆，根据公式计算回归系数
    ws=xTx.I*(xMat.T*yMat)
    return ws
    
def showLinerRegre():
    xArr ,yArr = loadDataSet('multiple_linear_regression/ex0.txt') #加载数据集
    ws = standRegres(xArr,yArr)         #得到回归系数
    xMat=mat(xArr)
    yMat=mat(yArr)
    yHat=xMat*ws
    fig=plt.figure()
    ax=fig.add_subplot(111)
    ax.scatter(xMat[:,1].flatten().A[0],yMat.T[:,0].flatten().A[0])
    xCopy=xMat.copy()
    xCopy.sort(0)
    yHat=xCopy*ws
    ax.plot(xCopy[:,1],yHat)
    plt.show()
    
if __name__ =='__main__':
    showLinerRegre()

几乎任一数据集都可以用上述方法建立模型，那么，如何判断这些模型的好坏呢？

显然我们看到有些数据是没有正确的拟合的，这被称为欠拟合。最佳拟合直线将数据视为直线进行建模，具有十分不错的表现。

但图中的数据似乎还存在其他的潜在模式，下一节我们介绍如何根据数据来局部调整预测。

局部加权线性回归找到最佳拟合直线

1.局部线性回归
线性回归的一个问题是有可能出现欠拟合现象，因为它求的是具有小均方误差的无偏估 计。显而易见，如果模型欠拟合将不能取得好的预测效果。所以有些方法允许在估计中引入一 些偏差，从而降低预测的均方误差。

其中的一个方法是局部加权线性回归（Locally Weighted Linear Regression，LWLR）。在该方法中，我们给待预测点附近的每个点赋予一定的权重。与kNN一样，这种算法每次预测均需要事先选取出对应的数据子集。该算法解除回归系数W的形式如下： 其中W是一个矩阵，这个公式跟我们上面推导的公式的区别就在于W，它用来给每个点赋予权重。

LWLR使用"核"（与支持向量机中的核类似）来对附近的点赋予更高的权重。核的类型可以自由选择，最常用的核就是高斯核，高斯核对应的权重如下：
这样我们就可以根据上述公式，编写局部加权线性回归，我们通过改变k的值，可以调节回归效果
2.代码实现
在当前工作目录下新建文件partialLinerRegre.py，添加如下代码：

def lwlr(testPoint, xArr, yArr, k=1.0):
    #转为矩阵形式
    xMat = mat(xArr); yMat = mat(yArr).T
    #样本个数
    m = shape(xMat)[0]
    weights = mat(eye((m)))  #创建m*m的单位矩阵，作为初始化对角权重矩阵
    for j in range(m):                             #遍历数据集计算每个样本的权重
        #计算预测点与该样本的偏差
        diffMat =testPoint-xMat[j, :]  
        #根据偏差利用函数赋予该样本相应的权重
        weights[j, j] = exp(diffMat * diffMat.T/(-2.0 * k**2))
    #将权重矩阵应用到公式中
    xTx = xMat.T * (weights * xMat)                   #m*m矩阵                     
    #计算行列式值是否为0，即确定是否可逆
    if linalg.det(xTx) == 0.0:
        return
    ws = xTx.I * (xMat.T * (weights * yMat))          #计算回归系数，m*1矩阵
    #返回测试点的预测值
    return testPoint * ws
    return testPoint*ws

def lwlrTest(testArr, xArr, yArr,k=1.0):
    #测试集样本数
    m=np.shape(testArr)[0]
    yHat=np.zeros(m)#用于存储预测结果
    for i in range(m):
        yHat[i]=lwlr(testArr[i],xArr,yArr,k)
    return yHat


def plotlwlrRegression():
    xArr, yArr = loadDataSet('multiple_linear_regression/ex0.txt')                                    #加载数据集
    yHat_1 = lwlrTest(xArr, xArr, yArr, 1.0)                            #根据局部加权线性回归计算yHat
    yHat_2 = lwlrTest(xArr, xArr, yArr, 0.01)                            #根据局部加权线性回归计算yHat
    yHat_3 = lwlrTest(xArr, xArr, yArr, 0.003)                            #根据局部加权线性回归计算yHat
    xMat = mat(xArr)                                                    #创建xMat矩阵
    yMat = mat(yArr)                                                    #创建yMat矩阵
    srtInd = xMat[:, 1].argsort(0)                                        #排序，返回索引值
    xSort = xMat[srtInd][:,0,:]
    fig, axs = plt.subplots(nrows=3, ncols=1,sharex=False, sharey=False, figsize=(10,8))                                        
    axs[0].plot(xSort[:, 1], yHat_1[srtInd], c = 'red')                        #绘制回归曲线
    axs[1].plot(xSort[:, 1], yHat_2[srtInd], c = 'red')                        #绘制回归曲线
    axs[2].plot(xSort[:, 1], yHat_3[srtInd], c = 'red')                        #绘制回归曲线
    axs[0].scatter(xMat[:,1].flatten().A[0], yMat.flatten().A[0], s = 20, c = 'green', alpha = .5)                #绘制样本点
    axs[1].scatter(xMat[:,1].flatten().A[0], yMat.flatten().A[0], s = 20, c = 'green', alpha = .5)                #绘制样本点
    axs[2].scatter(xMat[:,1].flatten().A[0], yMat.flatten().A[0], s = 20, c = 'green', alpha = .5)                #绘制样本点
    #设置标题,x轴label,y轴label
    axs0_title_text = axs[0].set_title(u'lwlr,k=1.0')
    axs1_title_text = axs[1].set_title(u'lwlr,k=0.01')
    axs2_title_text = axs[2].set_title(u'lwlr,k=0.003')
    plt.setp(axs0_title_text, size=8, weight='bold', color='red')  
    plt.setp(axs1_title_text, size=8, weight='bold', color='red')  
    plt.setp(axs2_title_text, size=8, weight='bold', color='red')  
    plt.xlabel('X')
    plt.show()

if __name__ == '__main__':
    plotlwlrRegression()

可以看到，当k越小，拟合效果越好。但是当k过小，会出现过拟合的情况，例如k等于0.003的时候

在算法实现的过程中，我们可以体会到局部加权线性回归也存在一个问题，它对每个点做预测时都必须使用整个数据集，大大的增加了计算量。

使用scikit-learn实现线性回归算法sklearn.linear_model模块提供了很多集成操作，线性回归、岭回归、贝叶斯回归等。本文使用的是LinearRegression。

让我们先看下LinearRegression这个函数，一共有4个参数：
参数说明如下：

fit_intercept：可选参数，布尔值，默认为True。是否对数据进行中心化。如果该变量为false，则表明数据已经进行了中心化，在下面的过程中不进行中心化处理，否则对输入数据进行中心化处理。

normalize：可选参数，布尔值，默认为False。是否对数据进行标准化处理。当fit_intercept设置为false时，该参数会被忽略。如果该变量为true，则对数据进行标准化处理。需要在使用normalize=False的estimator调用fit之前使用sklearn.preprocessing.StandardScaler。

copy_X:可选参数，布尔值，默认为True。该值为true，复制X；否则，重写X。

n_jobs:可选参数，整型，默认为1。计算时设置的任务个数(number of jobs)。如果选择-1则代表使用所有的CPU。这一参数的对于目标个数>1（n_targets>1）且足够大规模的问题有加速作用。

返回值：

coef_ :数组型变量，shape为（n_feature）或者(n_targets,n_features)。线性回归问题中的回归系数。如果输入为多目标问题，即fit二维数据，则返回一个二维数组，shape为（n_targets,n_features）；如果输入为单目标问题，返回一个一维数组。

intercept_：数组型变量。线性模型中的独立项。
方法：

decision_function(X) 对训练数据X进行预测

fit(X, y[, n_jobs]) 对训练集X, y进行训练。是对scipy.linalg.lstsq的封装

get_params([deep]) 得到该估计器(estimator)的参数。

predict(X) 使用训练得到的估计器对输入为X的集合进行预测（X可以是测试集，也可以是需要预测的数据）。

score(X, y[,]sample_weight) 返回对于以X为samples，以y为target的预测效果评分。

set_params(**params) 设置估计器的参数

了解了这些，接下来，我们使用sklearn库来实现我们第一个实验中的线性回归模型。添加如下代码：

#!/usr/bin/env python
#-*- coding:utf-8 -*-
import matplotlib.pyplot as plt
from numpy import *
from sklearn import datasets, linear_model
from sklearn.linear_model import LinearRegression
from sklearn.metrics import mean_squared_error, r2_score


def loadDataSet(fileName):  
    numFeat=len(open(fileName).readline().split('\t'))-1
    with open(fileName) as fr:
        #~ numFeat=len(fr.readline().split('\t'))-1
        dataMat=[];labelMat=[]
        for line in fr.readlines():
            lineArr=[]
            curLine=line.strip().split('\t')
            for i in range(numFeat):
                lineArr.append(float(curLine[i]))
            dataMat.append(lineArr)
            labelMat.append(float(curLine[-1]))
        return dataMat,labelMat
    return dataMat,labelMat

if __name__ =='__main__':
    dataX, dataY =loadDataSet('multiple_linear_regression/ex0.txt')
    matX=mat(dataX);matY=mat(dataY).T #将数据保存到矩阵中
    regr = LinearRegression()  #生成线性回归模型
    regr.fit(matX,matY)#对训练集X, y进行训练。
    #填充训练数据 matX(n_samples,n_features);matY(n_samples,n_targets)

    xCopy = matX.copy()
    xCopy.sort(0)
    predictY = regr.predict(xCopy) #得到模型预测值 predict使用训练得到的估计器对输入为X的集合进行预测

    plt.scatter(matX[:,1].flatten().A[0],matY[:,0].flatten().A[0],s=20,color='green',alpha=.5) #绘制散点图
    plt.plot(xCopy[:,1],predictY,color='red',linewidth=1) #绘制最佳拟合直线
 
    plt.xticks(())
    plt.yticks(())

    plt.show()