机器学习：线性回归

标准线性回归函数系数推导：

2

代码实现：

import numpy as np

#加载数据并将变量与标签转换为数组
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 standRegress(xArr,yArr):
    xMat = np.mat(xArr)
    yMat = np.mat(yArr).T
    xTx = xMat.T * xMat
    if np.linalg.det(xTx) == 0.0:
        print u'矩阵不可逆'
        return
    ws = xTx.I * (xMat.T * yMat)
    return ws

#数据文件路径及变量与标签值
path = 'Datas/Ch08/ex0.txt'
xArr,yArr = loadDataSet(path)

#计算预测值与真实值的相关系数
ws = standRegress(xArr,yArr)
xMat = np.mat(xArr)
yMat = np.mat(yArr)
yHat = xMat * ws

#求相关系数计算预测值与真实值的匹配程度
np.corrcoef(yHat.T,yMat)

#输出 array([[ 1.        ,  0.98647356],
#      [ 0.98647356,  1.        ]])

#画图
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)
yHa = xCopy * ws
ax.plot(xCopy[:,1],yHa)
plt.show()

标准线性回归拟合直线：

标准线性回归拟合直线

 

#计算R^2评价拟合度
import math
def CalcRScore(y,Py):
    m = len(y)
    avgy = np.average(y)
    rss = 0.0
    tss = 0.0
    for i in range(m):
        rss += math.pow(y[i]-Py[i],2)
        tss += math.pow(y[i]-avgy,2)
    rScore = 1.0 - 1.0*rss/tss
    return rScore

y = yMat.T.getA()
Py = yHat.getA()
print CalcRScore(y,Py)
# 0.973130088986

 

局部加权线性回归(LWLR)

此算法中我们给待预测点附近的每个点赋予一定的权重W，回归系数形式如下：

LWLR使用“核”来对附近的点赋予更高的权重，最常用的是高斯核函数，公式如下;

代码如下：

#利用高斯核函数求某一点局部加权后的预测值
def lwlr(testPoint,xArr,yArr,k=1.0):
    lw_xMat = np.mat(xArr)
    lw_yMat = np.mat(yArr).T
    m = np.shape(xArr)[0]
    weights = np.eye((m))
    for i in range(m):
        diffMat = testPoint - lw_xMat[i,:]
        weights[i,i] = np.exp(diffMat*diffMat.T/(-2.0*k**2)) #高斯核函数
    xTx = lw_xMat.T * (weights * lw_xMat)
    if np.linalg.det(xTx) == 0.0:
        print u"This matrix is singular,cannot do inverse!"
        return
    ws = xTx.I *(lw_xMat.T * (weights*lw_yMat))
    return testPoint * ws

# 应用于数据集
def lwlrTest(testArr,xArr,yArr,k=1.0):
    m = np.shape(testArr)[0]
    lw_yHat = np.zeros(m)
    for i in range(m):
        lw_yHat[i] = lwlr(testArr[i],xArr,yArr,k)
    return lw_yHat

#对单点进行估计
yArr[0]
#3.176513
lwlr(xArr[0],xArr,yArr,1.0)
#matrix([[ 3.12204471]])
lwlr(xArr[0],xArr,yArr,0.001)
#matrix([[ 3.20175729]])

绘图：

#调用lwlrTest函数求得数据集中所有点的估计
lw_yHat = lwlrTest(xArr,xArr,yArr,0.01)

#绘图函数需将数据点按序排列
lw_xMat = np.mat(xArr)
srtInd = lw_xMat[:,1].argsort(0) #argsort 返回数组中元素从小到大的索引值
xSort = lw_xMat[srtInd][:,0,:]


import matplotlib.pyplot as plt
fig = plt.figure()
ax = fig.add_subplot(111)
ax.plot(xSort[:,1],lw_yHat[srtInd])
ax.scatter(lw_xMat[:,1].flatten().A[0],np.mat(yArr).T.flatten().A[0],s=2,c='red')
plt.show()

k=0.01时的拟合效果：

局部线性回归--k=0.01

岭回归：

岭回归就是在矩阵XTX 上加上一个λI从而使得矩阵非奇异，进而能对XTX + λI求逆。其中矩阵I是一个m*m的单位矩阵，对角线上元素全为1，其他元素均为0，λ为用户定义的数值，回归系数的计算公式变为：

代码清单：

#岭回归
def ridgeRegres(xMat,yMat,lam=0.2):
    xTx = xMat.T *xMat
    denom = xTx + np.eye(np.shape(xMat)[1])*lam
    if np.linalg.det(denom) == 0.0:
        print u"This matrix is singular,cannot do inverse"
        return
    ws = denom.I * (xMat.T * yMat)
    return ws

#获取30个不同lambada下求得的回归系数
def ridgeTest(xArr,yArr):
    xMat = np.mat(xArr)
    yMat = np.mat(yArr).T
    yMean = np.mean(yMat,0)
    yMat = yMat - yMean
    xMeans = np.mean(xMat,0)
    xVar = np.var(xMat,0)
    xMat = (xMat - xMeans)/xVar #数据标准化：所有特征值都减去各自的均值并除以方差
    numTestPts = 30
    wMat = np.zeros((numTestPts,np.shape(xMat)[1]))
    for i in range(numTestPts):
        ws = ridgeRegres(xMat,yMat,np.exp(i-10)) #lambda 以指数级变化，可以看出在取非常小与非常大的值时对结果造成的影响
        wMat[i,:]=ws.T
    return wMat