机器学习学习笔记（5）线性回归

线性回归

给定数据集D={(x1,y1),(x2,y2),...,(xm,ym)}，其中xi=(xi1;xi2;xi3;...;xid)，yi是实数。线性回归试图学得一个线性模型以尽可能准确地预测实值输出标记。

对离散属性，若属性值之间存在“序”关系，可通过连续化将其转化为连续值；若属性值之间不存在序关系，假定有k个属性值，则通常转化为k维向量。

线性回归试图学得

，使得

要确定w和b，使得均方误差最小化：

用均方误差最小化来进行模型求解的方法称为“最小二乘法”。在线性回归中，最小二乘法就是试图找到一条直线，使所有样本到直线上的欧式距离之和最小。

求解w和b使得

最小化的过程，称为线性回归模型的最小二乘参数估计，可以将E(w,b)分别对w和b求导，得到

，

，令这两个式子都为0，得到w和b的闭式解。

，

其中

为x的均值。

当每个x有若干个属性时，试图学得

。

最小化损失函数为：

，求导得

，最终得到

。

但是当XX^T不是满秩矩阵时，会求出多个w，选择哪一个解作为输出将由学习算法的归纳偏好决定，常见的做法是引入正则化。

代码语言：javascript复制

# 代码和数据集源自于机器学习实战，可见https://github.com/AnnDWang/MachineLearning/blob/master/thirdbook/ch8/regression.py
from numpy import *

def loadDataSet(fileName):
    numFeat=len(open(fileName).readline().split('t'))-1
    dataMat=[]
    labelMat=[]
    fr=open(fileName)
    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

def standRgres(xArr,yArr):
    xMat=mat(xArr)
    yMat=mat(yArr).T
    xTx=xMat.T*xMat
    if linalg.det(xTx)==0.0:
        print('the matrix is singular, cannot do inverse')
        return
    ws=xTx.I*(xMat.T*yMat)
    return ws

xArr,yArr=loadDataSet('ex0.txt')
ws=standRgres(xArr,yArr)
print(ws)
xMat=mat(xArr)
yMat=mat(yArr)
yHat=xMat*ws

# 绘制出数据集散点图和最佳拟合直线图
import  matplotlib.pyplot as plt
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()

实验结果如图：

局部加权线性回归

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

其中的一个方法是局部加权线性回归（Locally Weighted Linear Regression，LWLR），在该计算方法中，给每个点赋予一定的权重，在这个子集上基于最小均方误差来进行普通的回归。这种算法每次预测均需要事先取出对应的数据子集，该算法解出的回归系数w形式如下：

，

是一个矩阵，用来给每个数据点赋予权重。

LWLR使用“核”来对附近的点赋予更高的权重，核的类型可以自由选择，最常用的就是高斯核，高斯核对应的权重如下：

这样就构成了一个只包含对角元素的权重矩阵，点x与x(i)越近，w(i,i)越大。k决定了对附近的点赋予多大的权重，是使用LWLR时唯一需要考虑的参数。

代码如下：

代码语言：javascript复制

# 代码和数据集源自于机器学习实战
# 局部加权线性回归函数
def lwlr(testPoint,xArr,yArr,k=1.0):
    xMat=mat(xArr)
    yMat=mat(yArr).T
    m=shape(xMat)[0]
    weights=mat(eye(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)
    if linalg.det(xTx)==0.0:
        print('this matrix is singular, cannot do inverse')
        return
    ws=xTx.I*(xMat.T*(weights*yMat))
    return testPoint*ws

def lwlrTest(testArr,xArr,yArr,k=1.0):
    m=shape(testArr)[0]
    yHat=zeros(m)
    for i in range(m):
        yHat[i]=lwlr(testArr[i],xArr,yArr,k)
    return yHat

xArr,yArr=loadDataSet('ex0.txt')
print(lwlr(xArr[0],xArr,yArr,1.0))
print(lwlr(xArr[0],xArr,yArr,0.001))
yHat=lwlrTest(xArr,xArr,yArr,0.003)
xMat=mat(xArr)
srtInd=xMat[:,1].argsort(0)
xSort=xMat[srtInd][:,0,:]

fig=plt.figure()
ax=fig.add_subplot(111)
ax.plot(xSort[:,1],yHat[srtInd])
ax.scatter(xMat[:,1].flatten().A[0],mat(yArr).T.flatten().A[0],s=2,c='red')
plt.show()

结果如图：

局部加权线性回归存在的问题，增加了计算量，因为对每个点做预测时都需要用到整个数据集

岭回归

简单来说，岭回归就是在

上增加一个

使得矩阵非奇异，从而能够进行求逆，矩阵I大小为mxm，对角线上元素为1，其他元素全为0。岭回归的计算公式为：

岭回归最先用来处于特征多于样本数的情况，现在也用于在估计中加入偏差，从而得到更高的估计，这里通过引入

来限制了所有的w的和，通过引入该惩罚项，能够减少不重要的参数，这个技术在统计学中也叫做缩减。缩减方法可以去掉不重要的参数，因此能够更好地理解数据。

代码如下：

代码语言：javascript复制

# 代码和数据集主要源自于机器学习实战
# 岭回归

# 用于计算回归系数，
# 实现了给定了lambda下的岭回归求解
# 如果没有指定lambda，默认为0.2
def ridgeRegres(xMat,yMat,lam=0.2):
    xTx=xMat.T*xMat
    denom=xTx eye(shape(xMat)[1])*lam
    if linalg.det(denom)==0.0:
        print('this matrix is singular, cannot do inverse')
        return
    ws=denom.I*(xMat.T*yMat)
    return ws

# 用于在一组lambda上测试结果
# 为了使用岭回归和缩减技术，首先需要对特征进行标准化处理
# 具体的做法是所有特征都减去各自的均值并处以方差
def ridgeTest(xArr,yArr):
    xMat=mat(xArr)
    yMat=mat(yArr).T
    yMean=mean(yMat,0)
    yMat=yMat-yMean
    xMeans=mean(xMat,0)
    xVar=var(xMat,0)
    xMat=(xMat-xMeans)/xVar
    numTestPts=30
    wMat=zeros((numTestPts,shape(xMat)[1]))
    for i in range(numTestPts):
        ws=ridgeRegres(xMat,yMat,exp(i-10))
        wMat[i,:]=ws.T
    return wMat
abX,abY=loadDataSet('abalone.txt')
ridgeWeights=ridgeTest(abX,abY)

fig=plt.figure()
ax=fig.add_subplot(111)
ax.plot(ridgeWeights)
plt.show()

岭回归的系数变化图如下：