机器学习之线性回归算法

188 阅读7分钟
参与拿奖:本文已参与「新人创作礼」活动,一起开启掘金创作之路
1.相关概念
回归的一般方法:
准备数据:回归需要数值型数据,标称型数据将会转换为二值型数据
分析数据:绘制相关图表(可视化分析)
训练算法:找到回归系数
测试算法:使用R^2或者预测值和数据的拟合度
使用算法:使用训练好的模型进行预测
线性回归
优点:结果易于理解,计算上不复杂
缺点:对非线性数据拟合不好
适用数据类型:数值型和标称型数据
2.案例

2.1 简单例子

利用普通最小二乘法求解最佳回归系数
from numpy import *

# 加载数据
# 数据预处理
def loadDataSet(fileName):      #general function to parse tab -delimited floats
    numFeat = len(open(fileName).readline().split('\t')) - 1 #get number of fields 
    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 standRegres(xArr,yArr):
    xMat=mat(xArr)
    yMat=mat(yArr).T
    xTx=xMat.T*xMat
    # 如果矩阵不可逆
    if linalg.det(xTx)==0.0:
        print("This matrix is singular, cannot do inverse")
        return
    return xTx.I*(xMat.T*yMat)

xArr,yArr=loadDataSet(r'ex0.txt')
# print(xArr,yArr)
ws=standRegres(xArr,yArr)
# shape(ws)
xMat=mat(xArr)
yMat=mat(yArr)
yHat=xMat*ws

# 绘制图像,观察数据
def plotRegres(xArr,yArr):
    import matplotlib.pyplot as plt
    xMat=mat(xArr)
    yMat=mat(yArr)
    ws=standRegres(xArr,yArr)
    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()
    
plotRegres(xArr,yArr)
# 查看相关系数
### corrcoef(yEstimate,yAcutual)
corrcoef(yHat.T,yMat)
# 从数据可以看出两者相关性为 0.98

运行结果:

image.png

2.2 局部加权线性回归

线性回归的一个主要问题是出现欠拟合现象
可以通过引入偏差,从而降低预测的均方误差
一种方法是加权线性回归,可以通过w矩阵对每一个数据点赋予权重
使用核函数来对附近的点赋予更高的权重
### 高斯核函数
def gKernel(x,y,k):
    return exp(abs(x-y)/((-2)*k**2))
# print(gKernel(1,0.5,0.5))


def lwlr(testPoint,xArr,yArr,k=1.0):
    xMat = mat(xArr); yMat = mat(yArr).T
    m = shape(xMat)[0]
    W=[]
    # 创建对数矩阵
    weights = mat(eye((m)))
    for j in range(m):                      #next 2 lines create weights matrix
        diffMat = testPoint - xMat[j,:]
        # 对角线赋值
        # print(diffMat,shape(diffMat))
        weights[j,j] = exp(diffMat*diffMat.T/(-2.0*k**2))
        W.append(weights[j,j])
        # print(W)
    xTx = xMat.T * (weights * xMat)
    # 特征值
    # print(xTx,shape(xTx),linalg.det(xTx))
    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,W

def lwlrTest(testArr,xArr,yArr,k=1.0):  #loops over all the data points and applies lwlr to each one
    m = shape(testArr)[0]
    yHat = zeros(m)
    # print(m)
    for i in range(m):
        # print(lwlr(testArr[i],xArr,yArr,k))
        yHat[i],W = lwlr(testArr[i],xArr,yArr,k)
    return yHat,W

# 测试
xArr,yArr=loadDataSet(r'ex0.txt')
# print(shape(xArr),shape(yArr))
yHat=lwlrTest(xArr,xArr,yArr,k=0.003)
# print(yHat)


# 绘制图像,观察数据
def plotRegres(xArr,yArr,k):
    import matplotlib.pyplot as plt
    xMat=mat(xArr)
    srtInd=xMat[:,1].argsort(0)
    xSort=xMat[srtInd][:,0,:]
    yHat,W=lwlrTest(xArr,xArr,yArr,k)
    fig=plt.figure()
    ax=fig.add_subplot(111)
    ax.plot(xSort[:,1],yHat[srtInd],c="b")
    # ax.plot(xMat[:,1],W,c="g")
    ax.scatter(xMat[:,1].flatten().A[0],yMat.T[:,0].flatten().A[0],s=2,c="r")
    plt.xlabel('x')
    plt.ylabel('y')
    plt.title("k=%.2f"%(k))
    plt.show()

def plotR():
    for k in [1.0,0.1,0.03]:
        plotRegres(xArr,yArr,k)

2.3 预测鲍鱼年龄

### 高斯核函数
def gKernel(x,y,k):
    return exp(abs(x-y)/((-2)*k**2))
# print(gKernel(1,0.5,0.5))


def lwlr(testPoint,xArr,yArr,k=1.0):
    xMat = mat(xArr); yMat = mat(yArr).T
    m = shape(xMat)[0]
    W=[]
    # 创建对数矩阵
    weights = mat(eye((m)))
    for j in range(m):                      #next 2 lines create weights matrix
        diffMat = testPoint - xMat[j,:]
        # 对角线赋值
        # print(diffMat,shape(diffMat))
        weights[j,j] = exp(diffMat*diffMat.T/(-2.0*k**2))
        W.append(weights[j,j])
        # print(W)
    xTx = xMat.T * (weights * xMat)
    # 特征值
    # print(xTx,shape(xTx),linalg.det(xTx))
    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,W

def lwlrTest(testArr,xArr,yArr,k=1.0):  #loops over all the data points and applies lwlr to each one
    m = shape(testArr)[0]
    yHat = zeros(m)
    # print(m)
    for i in range(m):
        # print(lwlr(testArr[i],xArr,yArr,k))
        yHat[i],W = lwlr(testArr[i],xArr,yArr,k)
    return yHat,W

# 测试
xArr,yArr=loadDataSet(r'ex0.txt')
# print(shape(xArr),shape(yArr))
yHat=lwlrTest(xArr,xArr,yArr,k=0.003)
# print(yHat)


# 绘制图像,观察数据
def plotRegres(xArr,yArr,k):
    import matplotlib.pyplot as plt
    xMat=mat(xArr)
    srtInd=xMat[:,1].argsort(0)
    xSort=xMat[srtInd][:,0,:]
    yHat,W=lwlrTest(xArr,xArr,yArr,k)
    fig=plt.figure()
    ax=fig.add_subplot(111)
    ax.plot(xSort[:,1],yHat[srtInd],c="b")
    # ax.plot(xMat[:,1],W,c="g")
    ax.scatter(xMat[:,1].flatten().A[0],yMat.T[:,0].flatten().A[0],s=2,c="r")
    plt.xlabel('x')
    plt.ylabel('y')
    plt.title("k=%.2f"%(k))
    plt.show()

def plotR():
    for k in [1.0,0.1,0.03]:
        plotRegres(xArr,yArr,k)

#%%

# 案例1: 预测鲍鱼年龄
def rssError(yArr,yHatArr):
    return ((yArr-yHatArr)**2).sum()


def lwlrTestPlot(xArr,yArr,k=1.0):  #same thing as lwlrTest except it sorts X first
    yHat = zeros(shape(yArr))       #easier for plotting
    xCopy = mat(xArr)
    xCopy.sort(0)
    for i in range(shape(xArr)[0]):
        yHat[i],W= lwlr(xCopy[i],xArr,yArr,k)
    return yHat,xCopy


def clcE(filename,k):
    abX,abY=loadDataSet(filename)
    # print(shape(abX[0:99]),shape(abY[0:99]))
    resYhat={}
    resYhatN={}
    resErr={}
    resErrN={}
    for i in k:
        # print(lwlrTestPlot(abX[0:99],abY[0:99],0.1))
        resYhat[i],xCopy=lwlrTestPlot(abX[0:99],abY[0:99],i)
        resErr[i]=rssError(abY[0:99],resYhat[i].T)
        resYhatN[i],xCopy=lwlrTestPlot(abX[100:199],abY[100:199],i)
        resErrN[i]=rssError(abY[100:199],resYhatN[i].T)
        #print( resErr[i])
    return resYhat,resYhatN,resErr,resErrN


# 测试
xArr,yArr=loadDataSet(r'abalone.txt')
# print(shape(xArr[1:99]),shape(yArr[1:99]))
yHat,xCopy=lwlrTestPlot(xArr[0:99],yArr[0:99],0.1)
# print(mat(yArr[0:99]-yHat)*(mat(yArr[0:99]-yHat).T))
# print(rssError(yArr[0:99],yHat.T))
resYhat,resYhatN,resErr,resErrN=clcE(r'abalone.txt',[0.1,0.3,0.8]) 
for i in [0.1,0.3,0.8]:
    print("k=%.2f ,resErr[0:99]=%.4f,resErr[100:199]=%.4f"%(i,resErr[i],resErrN[i]),end="\n")


#使用标准线性回归测试
ws=standRegres(xArr[0:99],yArr[0:99])
yHat=mat(xArr[100:199])*ws
k=rssError(yArr[100:199],yHat.T.A)
print("standRegres:yArr[100:199]=%.4f:"%(k))

运行结果

image.png

3.岭回归
如果数据特征多于样本个数,不能使用之前的回归方法进行预测
岭回归:通过再x.T*x上加一个nI使得矩阵非奇异,进而可以对相加结果求逆,nI是一个单位矩阵,对角线元素全为1
缩减:通过引入n来限制所有权重之和,引入该惩罚项,能够减少不重要的参数
岭回归使用了单位矩阵乘以常量n,该矩阵因为对角线为1,可以将其看成一条岭,将上下隔开
缩减方法可以去掉不重要的参数,可以通过预测误差确定n
# 该函数返回回归系数
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
    
# 测试函数
def ridgeTest(xArr,yArr):
    xMat = mat(xArr); yMat=mat(yArr).T
    # 数据标准化
    # 求解均值
    yMean = mean(yMat,0)
    yMat = yMat - yMean     #to eliminate X0 take mean off of Y
    #regularize X's
    xMeans = mean(xMat,0)   #calc mean then subtract it off
    xVar = var(xMat,0)      #calc variance of Xi then divide by it
    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


# 测试
xMat,yMat=loadDataSet(r'abalone.txt')
ridgeWeights=ridgeTest(xMat,yMat)

# 绘图,查看30个不同lama值下的系数
def plotRidgeRegres(ridgeW):
    import matplotlib.pyplot as plt
    fig=plt.figure()
    ax=fig.add_subplot(111)
    ax.plot(ridgeW)
    plt.ylabel('w')
    plt.xlabel('log(lambda)')
    plt.show()
    
plotRidgeRegres(ridgeWeights)

运行结果

image.png

岭回归类似于普通的最小二乘法回归
其限定了所有回归系数的平方和不大于lambda,但是在使用最小二乘法回归当两个或更多的特征相关时,可能会出现一个很大的正系数
和一个很大的负系数,但是使用岭回归可以避免这个问题

其他缩减方法:
lasso:限定条件为,所有系数的绝对值之和不大于lambda

其他解决办法的代码参看本节完整代码[4]

使用前向逐步回归
一种贪心算法,每一步都尽可能减小误差
伪代码:
将数据标准化,使其分布满足0均值和单位方差
在每轮迭代过程中:
     设置当前最小误差loss为正无穷
     对于每个特征:
         增大或减小
             改变一个系数
             计算新系数下的误差
             根据误差来更新系数
                 修改最好的系数(如果误差小于当前误差就修改系数)
            
4.小结
分类预测离散型变量,回归预测连续型变量
回归方程中,求得特征值对应的最佳回归系数的方法是最小化误差的平方和。
当数据样本数比特征少时,可以考虑岭回归
5.参考资料

[1] 机器学习实战

[2] 书籍源码

[3] jupyter版本

[4] 本节代码