第二部分 利用回归预测数值型数据

第 8 章 预测数值型数据：回归

[TOC]

本章内容

• 线性回归
• 局部加权线性回归
• 岭回归和逐步线性回归
• 预测鲍鱼年龄和玩具售价

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

线性回归：

• 优点：结果易于理解，计算熵不复杂。
• 缺点：对非线性的数据拟合不好。
• 适用数类型：数值型和标称型数据。

回归的目的是预测数值型的目标值。最直接的办法是依据输入写成一个目标值的计算公式。

y=a_1*x_1+a_2*x_2

这就是回归方程，其中的 a1 和 a2 称作回归系数（regression weights），求这些回归系数的过程就是回归。一旦有了这些回归系数，再给定输入，做预测就非常容易了，具体的做法是用回归系数乘以输入值，再将结果全部加在一起，就得到了预测值。（这里的回归都是指线性回归（linear regression）。

回归的一般方法：

1. 收集数据：采用任意方法收集数据。
2. 准备数据：回归需要数值型数据，标称型数据将被转成二值型数据。
3. 分析数据：绘制出数据的可视化二维图将有助于对数据做出理解和分析，在采用缩减法求得新回归系数之后，可以将新拟合线绘制在图上作为对比。
4. 训练算法：找到回归系数。
5. 测试算法：使用 R2 或者预测值和数据的拟合度，来分析模型的效果。
6. 使用算法：使用回归，可以在给定输入的时候预测出一个数值，这是对分类方法的提升，因为这样可以预测连续型数值而不仅仅是离散的类别标签。

应当怎样从一大堆数据里求出回归方程呢？假定输入数据存放在矩阵 X 中，而回归系数存放在向量 W 中，那么对于给定的数据 x1，预测结果将会通过

Y_1=X^T_1 W

给出。现在的问题是，有一些 x 个对应的 Y，怎样才能找到 w 呢？一个常用的方法就是找出误差最小的 w。这里的误差是指预测 y 值和真实 y 值之间的差值，使用该差值的简单累加将使得正差值和发差值相互抵消，所以使用平方误差公式：

\sum^m_{i=1}(y_i-x^T_i w)^2

用矩阵表示为：

(y-Xw)^T(y-Xw)

如果对 w 求导，得到

X^T(y-Xw)

令其等于零，解出 w 如下：

\hat{w}=(X^T X)^{-1}X^T y

w 上方的小标记表示，这是当前可以估计出的 w 的最优解。

值得注意的是，上述公式中包含

X^T X^{-1}

，也就是需要对矩阵求逆，因此这个方程只在逆矩阵存在的时候适用（代码中可以用伪逆矩阵）。上述的最佳 w 求解的统计学中的常见问题，除了矩阵方法外还有其他方法可以解决。该方法也称为 OLS，即普通最小二乘法（ordinary least squares）。

对于图8-1的散点图，下面来介绍如何给出该数据的最佳拟合直线。

图8-1

添加 loadDataSet() 函数：

def loadDataSet(fileName):
    """
    读取文本文件
    :param fileName:
    :return:
    """
    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

添加 standRegres 函数：

def standRegres(xArr, yArr):
    """计算最佳拟合直线"""
    xMat = np.mat(xArr)
    yMat = np.transpose(np.mat(yArr))
    xTx = np.transpose(xMat) * xMat
    if np.linalg.det(xTx) == 0.0: # 计算 xTx 的行列式，为零则不能计算逆矩阵
        print("this matrix is singular,cannot do inverse")
        return
    ws = xTx.I * (xMat.T * yMat)
    # ws = np.linalg.solve(xTx, xMat.T * yMat) # 伪逆矩阵，可以不判断行列式是否为零
    return ws

测试：

xArr,yArr = loadDataSet('ex0.txt')
ws = standRegres(xArr,yArr) # 回归系数
# print(ws) # [[ 3.00774324] [ 1.69532264]]
xMat = np.mat(xArr)
yMat = np.mat(yArr)
yHat = xMat * ws # 预测的 y 值： y=ws[0]+ws[1]*X1
print(np.cov(yHat.T, yMat)) # 协方差 [[ 0.24664409  0.24664409] [ 0.24664409  0.25345439]]
print(np.corrcoef(yHat.T, yMat)) # 相关系数 [[ 1.  0.98647356] [ 0.98647356  1.  ]]

# 绘制数据集散点图和最佳拟合直线图
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()

图8-2

几乎任一数据集都可以用上述方法建立模型，那么，如何判断这些模型的好坏呢？比较一下图8-3 的两个子图，如果在两个数据集上分别作线性回归，将得到完全一样的模型（拟合直线）。显然两个数据是不一样的，那么模型分别在二者上的效果如何？我们当如何比较这些效果的好坏？有两种方法可以计算预测值 yHat 序列和真实值 y 序列的匹配程度，那就是计算这两个序列的相关系数。

图8-3

在 numpy 中提供了 corrcoef 来计算预测值与真实值的相关性.

协方差：

Cov(X,Y)=\frac{1}{m}\sum^m_{i=1}(x_i-\bar{x})(y_i-\bar{y})

相关系数：

r=\frac{Cov(X,Y)}{\sqrt{\delta_x}\sqrt{\delta_y}}

其中，

\sqrt{\delta_x}$$、$$\sqrt{\delta_y}

表示 X、Y 的方差，并且 |r|<=1。

该矩阵包含所有两两组合的相关系数。可以看到，对角线上的数据是 1.0，因为 yMat 和自己的匹配是最完美的，而 yHat 和 yMat 的相关系数为 0.98.

2. 局部加权线性回归

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

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

\hat{w}=(X^T WX)^{-1}X^TWy

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

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

w(i,i)=exp(\frac{|x^{(i)}-x|}{-2k^2})

这样就构建了一个只含对角元素的权重矩阵 w，并且点 x 与 x(i) 越近，w(i,i) 将会越大。上述公式包含一个需要用户指定的参数 k，它决定了对附近的点赋予多大的权重，这也就是使用 LWLR 时唯一需要考虑的参数，在图8-4中可以看到参数 k 与权重的关系：

图8-4

添加 lwlr() 函数：

def lwlr(testPoint, xArr, yArr, k=1.0):
    """
    局部加权线性回函数
    :param testPoint:
    :param xArr:
    :param yArr:
    :param k:用于控制衰减速度
    :return:
    """
    xMat = np.mat(xArr)
    yMat = np.mat(yArr).T
    m = np.shape(xMat)[0]
    weights = np.mat(np.eye((m))) # 创建对角矩阵，阶数等于样本点个数
    for j in range(m):
        # 随着样本点与待预测点距离的递增，权重值大小以指数级衰减
        diffMat = testPoint - xMat[j,:]
        weights[j,j] = np.exp(diffMat*diffMat.T/(-2.0*k**2)) # 参数 k，控制衰减的速度
    xTx = xMat.T * (weights * xMat)
    if np.linalg.det(xTx) == 0.0:
        print('this matrix is singular, cannot do inverse')
    ws = xTx.I * (xMat.T * (weights * yMat))
    return testPoint * ws

添加 lwlrTest() 函数：

def lwlrTest(testArr, xArr, yArr, k=1.0):
    """
    用于为数据集中每个点调用 lwlr()，这有助于求解 k 的大小
    :param testArr:
    :param xArr:
    :param yArr:
    :param k:
    :return:
    """
    m = np.shape(testArr)[0]
    yHat = np.zeros(m)
    for i in range(m):
        yHat[i] = lwlr(testArr[i], xArr, yArr, k)
    return yHat

测试：

xArr, yArr = loadDataSet('ex0.txt')
yHat = lwlrTest(xArr, xArr, yArr, 0.01) # k=1, 0.01, 0.003
xMat = np.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], np.mat(yArr).T.flatten().A[0], s=2, c='red')
plt.show()

可以观察到如图8-5所示的效果。图8-5 给出了 k 在不同取值下的结果图。

图8-5

局部加权线性回归也存在一个问题，即增加了计算量，因为它对每个点做预测时都必须使用整个数据集。从图8-5可以看出，k=0.01 时可以得到很好的估计，但是同时看到涂8-4中 k=0.01 的情况，就会发现大多数据点的权重都接近零。如果避免这些计算将可以减少程序运行时间，从而缓解因计算量增加带来的问题。

3. 示例：预测鲍鱼的年龄

在data目录下有一份来自 UCI 数据集合的数据，记录了鲍鱼的年龄。鲍鱼年龄可以从鲍鱼壳的层数推算得到。

预测代码：

abX, abY = loadDataSet('abalone.txt')
yHat01 = lwlrTest(abX[0:99], abX[0:99], abY[0:99], 0.1)
yHat1 = lwlrTest(abX[0:99], abX[0:99], abY[0:99], 1)
yHat10 = lwlrTest(abX[0:99], abX[0:99], abY[0:99], 10)
print(rssError(abY[0:99], yHat01.T)) # 56.7886874305
print(rssError(abY[0:99], yHat1.T)) # 429.89056187
print(rssError(abY[0:99], yHat10.T)) # 549.118170883

可以看到，使用较小的核将得到较低的误差，那么，为什么不在所有数据集上都使用最小的核呢？这是因为 使用最小的核将造成过拟合，对新数据不一定能达到最好的预测效果。下面就来看看它们在新数据上的表现：

yHat01 = lwlrTest(abX[100:199], abX[0:99], abY[0:99], 0.1)
yHat1 = lwlrTest(abX[100:199], abX[0:99], abY[0:99], 1)
yHat10 = lwlrTest(abX[100:199], abX[0:99], abY[0:99], 10)
print(rssError(abY[100:199], yHat01.T)) # 57913.5155016
print(rssError(abY[100:199], yHat1.T)) # 573.52614419
print(rssError(abY[100:199], yHat10.T)) # 517.571190538

从上述结果可以看到，核大小等于 10 时的测试误差最小，但它在训练集上的误差却最大的。接下来和简单的线性回归做个比较：

ws = standRegres(abX[0:99], abY[0:99])
yHat = np.mat(abX[100:199]) * ws
print(rssError(abY[100:199], yHat.T.A)) # 518.636315325

简单线性回归达到了与局部加权线性回归类似的结果，这也表明一点，必须在未知数据上比较效果才能选取到最佳模型。

4. 缩减系数来“理解”数据

如果数据的特征比样本点还多应该怎么办？是否还可以使用线性回归和之前的额方法来做预测？答案是否定的，即不能再使用前面介绍的方法。这是因为在计算

(X^TX)^{-1}

的时候会出错。

如果特征比样本点还多 (n>m) 也就是说输入数据的矩阵 X 不是满秩矩阵。非满秩矩阵在求逆时会出问题。

这节介绍的三种方法就是来解决这个问题。

4.1 岭回归

岭回归（ridge regression）就是在矩阵 $$X^TX$$ 上加一个 $$\lambda I$$ 从而使得矩阵非奇异，进而能对 $$X^T X+\lambda I$$ 求逆。其中矩阵 I 是一个 m*m 的单位矩阵。而 lambda 是一个用户定义数值，后面会做介绍。在这种情况下， 回归系数的计算公式将变成：

\hat{w}=(X^T X+\lambda I)^{-1}X^T y

岭回归最先用来处理特征多余样本数的情况，现在也用于在估计中加入偏差，从而得到更好的估计。这里通过引入 lambda 来限制了所有 w 之和，通过引入该惩罚项，能够减少不重要的参数，这个技术在统计学中也叫做缩减（shrinkage）。

添加 ridgeRegres() 函数：

def ridgeRegres(xMat, yMat, lam=0.2):
    """
    计算回归系数
    :param xMat: 
    :param yMat: 
    :param lam: 
    :return: 
    """
    xTx = xMat.T * xMat
    denom = xTx + np.eye(np.shape(xMat)[1]) * lam # 
    if np.linalg.det(denom) == 0.0:
        print("this matrix is singular,cannot do inverse")
        return
    ws = denom.I * (xMat.T * yMat)
    # ws = np.linalg.solve(xTx + np.eye(np.shape(xMat)[1]) * lam, xMat.T * yMat) # 用伪逆矩阵求解
    return ws

添加 ridgeTest() 函数：

def ridgeTest(xArr, yArr):
    """
    在一组lam上测试结果
    :param xArr:
    :param yArr:
    :return:
    """
    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 = (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))
        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()

这样就得到了 30 个不同 lambda 所对应回归系数。图8-6是所见的额效果图：

图8-6

该图绘制出了回归系数与 log(lam) 的关系。在最左边，即 lambda 最小时，可以得到所有系数的原始值（与线性回归一致）；而在右边，系数全部缩减为 0；在中间部分的某值将可以取得最好的预测效果。为了定量地找到最佳参数值，还需要进行交叉验证。另外，要判断哪些变量对结果预测最具有影响力，在图8-6 中观察它们对应的系数大小就可以。

4.2 lasso

不难证明，在增加如下约束时，普通的最小二乘法回归会得到与岭回归的一样的公式：

\sum^n_{k=1}w^2_k \leq \lambda

上式限定了所有回归系数的平方和不能大于 lambda。使用普通的最小二乘法回归在当两个或更多的特征相关时，可能会得出一个很大的正系数和一个很大的负系数。正式因为上述限制条件的存在，使用岭回归可以避免这个问题。

与岭回归类似，另一个缩减方法 lasso 也对回归系数做了限定，对应的约束条件如下：

\sum^n_{k=1}|w_k| \leq \lambda

唯一的不同点在于，这个约束条件使用绝对值取代了平方和。在 lambda 足够小时，一些系数会因此被迫缩减到 0，这个特性可以帮助我们更好的理解数据，这两个约束条件在公式上看起来相差无几，但细微的变化却极大的增加了计算复杂度（为了在这个新的约束条件下解出回归系数，需要使用二次规划算法）。下面将介绍一个更为简单的方法来得到结果，该方法叫做前向逐步回归。

4.3 前向逐步回归

前向逐步回归算法可以得到与 lasso 差不多的效果，但更加简单。它属于一种贪心算法，即每一步都尽可能减少误差。一开始，所有的权重都设为 1，然后每一步所做的决策时对某个权重增加或减少一个很小的值。

该算法的伪代码如下：

数据标准化，使其分布满足 0 均值和单位方差
在每轮迭代过程中：
    设置当前最小误差啊 lowestError 为正无穷
    对每个特征：
        增大或缩小：
            改变一个系数得到一个新的 W
            计算新 W 下的误差
            如果误差 Error 小于当前最小误差 lowestError：设置 Wbest 等于当前的 W
        将 W 设置为新的 WBest

添加 regularize() 函数：

def regularize(xMat):
    """
    把特征按照均值为0，方差为1进行标准化
    :param xMat:
    :return:
    """
    inMat = xMat.copy()
    inMeans = np.mean(inMat, 0)
    inVar = np.var(inMat, 0)
    inMat = (inMat - inMeans)/inVar
    return inMat

添加 stageWise() 函数：

def stageWise(xArr, yArr, eps=0.01, numIt=100):
    """
    逐步线性回归算法
    :param xArr:
    :param yArr:
    :param eps: 步长
    :param numIt: 迭代次数
    :return:
    """
    xMat =np.mat(xArr)
    yMat = np.mat(yArr).T
    yMean = np.mean(yMat, 0)
    yMat = yMat - yMean
    xMat = regularize(xMat)
    m,n = np.shape(xMat)
    returnMat = np.zeros((numIt,n))
    ws = np.zeros((n,1))
    wsTest = ws.copy() # 为了实现贪心算法建立了ws的两个副本
    wsMax = ws.copy()
    for i in range(numIt): # 迭代
        print('ws.T:',ws.T)
        lowestError = float('inf')
        for j in range(n):
            for sign in [-1,1]: # 两次循环，分别计算增加或减少该特征对误差的影响
                wsTest = ws.copy()
                wsTest[j] += eps*sign
                yTest = xMat*wsTest
                rssE = rssError(yMat.A, yTest.A) # 得到平方误差
                if rssE < lowestError: # 经过与所有的误差比较后取最小的误差
                    lowestError = rssE
                    wsMax = wsTest
        ws = wsMax.copy()
        returnMat[i,:] = ws.T
    return returnMat

测试一下实际效果：

xArr, yArr = loadDataSet('abalone.txt')
# print(stageWise(xArr,yArr,0.01,200))
print(stageWise(xArr,yArr,0.001,5000))

与最小二乘法进行比较，发现在5000次迭代后，逐步线性回归算法与常规的最小二乘法效果类似。使用 0.005 的 epsilon 值并经过 1000 次迭代后的结果参见图8-7：

图8-7

逐步线性回归算法的实际好处并不在于能绘出图8-7这样漂亮的图，主要的优点在于它可以帮助人们理解现有的模型并作出改进。当构建一个模型后，可以运行该算法找出重要的特征，这样就有可能及时停止对那些不重要特征的收集。最后，如果用于测试，该算法每 100 次迭代后就可以构建出一个模型，可以使用类似于 10 折交叉验证比较这些模型，最终选择使误差最小的模型。

当应用缩减方法（如逐步线性回归或岭回归）时，模型也就增加了**偏差（bias），与此同时却减小了模型的方差。

5. 权衡偏差与方差

图8-8给出了训练误差和测试误差的曲线图，上面的曲线图就是测试误差，下面的曲线是训练误差。根据8.3节的实验我们知道：如果降低核的大小，那么训练误差将变小。从图8-8来看，从左到右就表示了核逐渐减小的过程。

图8-8

一般认为，上述两种误差由三个部分组成：偏差、测量误差和随机噪声。在8.2节和8.3节，我们通过引入了三个越来越小的核来不断增大模型方差。

8.4节介绍了缩减法，可以将一些系数缩减成很小的值或直接缩减为0，这是一个增大模型偏差的例子。通过把一些特征的回归系数缩减到 0 ，同时也就减少了模型的复杂度。例子中有8个特征，消除了其中两个后不仅使模型更易理解，同时还降低了预测误差。图8-8左侧是参数缩减过于严厉的结果，右侧是无缩减的效果。

6. 示例：预测乐高玩具套装的价格

一种乐高套装基本上在几年后就会停产，但乐高的收藏着之间仍会在停产后彼此交易，为了给乐高套装估价，下面将用本章的回归技术来建立一个预测模型。

示例：用回归法预测乐高套装的价格：

1. 收集数据：用 google shopping 的 api 收集数据。
2. 准备数据：从返回的 json 数据中抽取价格。
3. 分析数据：可视化并观察数据。
4. 训练算法：构建不同的模型，采用逐步线性回归和直接的线性回归模型。
5. 测试算法：使用交叉验证来测试不同的模型，分析哪个效果更好。
6. 使用算法：这次练习的目的就是生成数据模型。

6.1 收集数据：使用 google 购物的 api

详细 api 介绍可参见：http://code.google.com/apis/shopping/search/v1/gettting_started.html

添加 searchForSet() 函数：

from time import sleep
import json
import urllib

def searchForSet(retX, retY, setNum, yr, numPce, origPrc):
    """
    调用google 购物 api 并保证数据抽取的正确性
    :param retX:
    :param retY:
    :param setNum:
    :param yr:
    :param numPce:
    :param origPrc:
    :return:
    """
    sleep(10) # 防止短时间内有过多的 api 调用
    myAPIstr = 'get from code.google.com'
    searchURL = 'http://www.googleapis.com/shopping/search/v1/public/products?key=%s&country=US&q=lego+%d&alt=json' % (myAPIstr,setNum)
    pg = urllib.request.urlopen(searchURL)
    retDict = json.loads(pg.read())
    for i in range(len(retDict['items'])):
        try:
            currItem = retDict['items'][i]
            if currItem['product']['condition'] == 'new': # 判断是否新产品
                newFlag = 1
            else: newFlag = 0
            listOfInv = currItem['product']['inventories']
            for item in listOfInv:
                sellingPrice = item['price']
                if sellingPrice > origPrc * 0.5: # 如果一个套装的价格比原始价格低一半以上就认为是不完整，过滤掉
                    print("%d\t%d\t%d\t%f\t%f" % (yr, numPce, newFlag, origPrc, sellingPrice))
                    retX.append([yr, numPce, newFlag, origPrc])
                    retY.append(sellingPrice)
        except:
            print('problem with item %d' % i)

测试（无法访问google api）：

def setDataCollect(retX, retY):
    searchForSet(retX, retY, 8288, 2006, 800, 49.99)
    searchForSet(retX, retY, 10030, 2002, 3096, 269.99)
    searchForSet(retX, retY, 10179, 2007, 5195, 499.99)
    searchForSet(retX, retY, 10181, 2007, 3428, 199.99)
    searchForSet(retX, retY, 10189, 2008, 5922, 299.99)
    searchForSet(retX, retY, 10196, 2009, 3263, 249.99)

lgX = []
lgY = []
setDataCollect(lgX, lgY)

6.2 训练算法：建立模型

上一节从网上收集到了一些真实的数据，下面将为这些数据构建一个模型。构建的模型可以对售价作出预测，并帮助我们理解现有数据。

使用常规的线性回归：

# print(np.shape(lgX)) # (58, 4)
lgX1 = np.mat(np.ones((58,5))) # 创建一个全 1 的矩阵
lgX1[:,1:5] = np.mat(lgX) # 将原数据矩阵 lgX 复制到新数据矩阵 lgX1 的第 1 到第 5 列
ws = standRegres(lgX1, lgY) # 回归算法
print(lgX1[0]*ws)
print(lgX1[-1]*ws)
print(lgX1[43]*ws)

使用岭回归确定最佳回归系数：

def crossValidation(xArr, yArr, numVal=10):
    """
    交叉验证测试岭回归
    :param xArr:
    :param yArr:
    :param numVal: 交叉验证的次数
    :return:
    """
    m = len(yArr)
    indexList = range(m)
    errorMat = np.zeros((numVal, 30))
    for i in range(numVal): # 数据分为训练集和测试集
        trainX = []
        trainY = []
        testX = []
        testY = []
        np.random.shuffle(indexList) # 对数据混洗
        for j in range(m):
            if j < m*0.9:
                trainX.append(xArr[indexList[j]])
                trainY.append(yArr[indexList[j]])
            else:
                testX.append(xArr[indexList[j]])
                testY.append(yArr[indexList[j]])
    wMat = ridgeTest(trainX, trainY) # wMat 保存岭回归中的所有回归系数
    for k in range(30):
        # 用训练时的参数将测试数据标准化
        matTestX = np.mat(testX)
        matTrainX = np.mat(trainX)
        meanTrain = np.mean(matTrainX, 0)
        varTrain = np.var(matTrainX, 0)
        matTestX = (matTestX - meanTrain)/varTrain
        yEst = matTestX * np.mat(wMat[k,:]).T + np.mean(trainY)
        errorMat[i,k] = rssError(yEst.T.A, np.array(testY))
    meanErrors = np.mean(errorMat, 0)
    minMean = float(min(meanErrors))
    bestWeights = wMat[np.nonzero(meanErrors==minMean)]
    xMat = np.mat(xArr);
    yMat = np.mat(yArr).T
    meanX = np.mean(xMat, 0)
    varX = np.var(xMat, 0)
    unReg = bestWeights/varX
    print("the best model from Ridge regression is :\n", unReg) # 标准化后的数据还原
    print("with constant term:",-1 * sum(np.multiply(meanX, unReg)) + np.mean(yMat))

运行：

crossValidation(lgX, lgY, 10)

我们本期望找到一个更易于理解的模型，显然没有达到预期结果。为了达到这一点，我们来看一下在缩减过程中回归系数是如何变化的，输入一下的命令：

print(ridgeTest(lgX, lgY)) 

这些系数是经过不同程度的缩减得到的。

这种分析方法似的我们可以挖掘大量数据的内在规律。在特征比较多的时候，可以指出哪些特征时关键的，而哪些特征时不重要的。

7. 本章小结

在回归方程中，求得特征对应的最佳回归系数的方法是最小化误差的平方差。给定输入矩阵 X，如果 $X^TX$ 的逆存在并可以求得的话，回归法对可以直接使用。数据集上计算出的回归方程并不一定意味着它是最佳，可以使用预测值 yHat 和原始值 y 的相关性来度量回归方程的好坏。

• 机器学习实战（笔记）：第 1 章 机器学习基础
• 机器学习实战（笔记）：第 2 章 k-近邻算法
• 机器学习实战（笔记）：第 3 章 决策树
• 机器学习实战（笔记）：第 4 章 基于概率论的分类方法：朴素贝叶斯
• 机器学习实战（笔记）：第 5 章 Logistic 回归
• 机器学习实战（笔记）：第 6 章 支持向量机（未完）
• 机器学习实战（笔记）：第 7 章 利用 AdaBoost 元算法提高分类性能（未完）
• 机器学习实战（笔记）：第 8 章 预测数值型数据：回归