机器学习实战【9】（岭回归和Lasso）

最新推荐文章于 2024-04-25 17:26:24 发布

0-00-0

最新推荐文章于 2024-04-25 17:26:24 发布

阅读量1.4k

点赞数

CC 4.0 BY-SA版权

分类专栏：机器学习文章标签：机器学习 python

本文链接：https://blog.youkuaiyun.com/qq_29846039/article/details/78467109

机器学习专栏收录该内容

9 篇文章

订阅专栏

本文介绍了机器学习中的岭回归和Lasso两种正则化方法，详细阐述了它们的原理与区别。岭回归通过在最小二乘法基础上引入正则项来解决多重共线性问题，而Lasso不仅解决共线性，还能进行特征选择。内容包括理论解释和Python实现，使用交叉验证确定最佳λ参数。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

本博客记录《机器学习实战》（MachineLearningInAction）的学习过程，包括算法介绍和python实现。

岭回归

岭回归（Ridge Regression），是一种线性回归方法，在最小二乘法的基础上加入一个正则项，以解决样本数少于特征数以及特征存在线性相关（多重共线性）时产生的问题。

最小二乘法

多元线性模型 $y=X\beta+\epsilon$ 中， $\epsilon$ 为残差项，残差项越小，模型也就越拟合数据。最小二乘法就是求解使得残差平方和最小的参数的方法。残差平方和RSS为：

Q(β)=(y−Xβ)T(y−Xβ) $Q(\beta)=(y-X\beta)^T(y-X\beta)$

把这个式子求导计算一波就可以推出 $\beta$ 的最小二乘估计：

β^=(XTX)−1XTy $\hat\beta=(X^TX)^{-1}X^Ty$

如果X存在线性相关关系，或者p小于n，这时逆矩阵就无法求解。

岭回归

在最小二乘估计的基础上，岭回归增加了一项，称为岭回归估计：

β^=(XTX+λI)−1XTy $\hat\beta=(X^TX+\lambda I)^{-1}X^Ty$

此时求解的式子变成了：

β^=argmin{(∑i=1nyi−β0−∑j=1pXijβj)2+λ∑j=1pβ2j}

$\hat\beta=argmin\{(\sum_{i=1}^ny_i-\beta_0-\sum_{j=1}^pX_{ij}\beta_j)^2+\lambda\sum_{j=1}^p\beta_j^2\}$

其中新增的最后一项称为惩罚函数。进一步，这个问题又可以转化成：

β^= a r g m i n (\sum i = 1 n y i - β 0 - \sum j = 1 p X i j β j) 2 s . t . \sum j = 1 p β 2 j < t

$\hat\beta=argmin(\sum_{i=1}^ny_i-\beta_0-\sum_{j=1}^pX_{ij}\beta_j)^2\\ s.t. \sum_{j=1}^p\beta_j^2<t$

相当于在最小二乘法的基础上添加了对参数 $\beta$ 的约束。岭参数 $\lambda$ 接近0时，岭回归变为最小二乘法；岭参数 $\lambda$ 变大时，参数 $\beta$ 趋近于0。

LASSO

LASSO算法与岭回归类似，但是其惩罚函数是一阶的，即：

β^=argmin{(∑i=1nyi−β0−∑j=1pXijβj)2+λ∑j=1p|βj|}

$\hat\beta=argmin\{(\sum_{i=1}^ny_i-\beta_0-\sum_{j=1}^pX_{ij}\beta_j)^2+\lambda\sum_{j=1}^p|\beta_j|\}$

相比岭回归，LASSO更易求解，并且更容易筛选变量。

python实现：

程序中计算了从小到大30个 $\lambda$ 参数的结果，通过交叉验证的方法，取使得最终的平方误差最小的一组参数作为结果。需要注意的是数据集在处理之前要先进行标准化以去除公式中的常数项。

# 数据标准化
def regularize(xMat):  # regularize by columns
    inMat = xMat.copy()
    inMeans = np.mean(inMat, 0)  # calc mean then subtract it off
    inVar = np.var(inMat, 0)  # calc variance of Xi then divide by it
    inMat = (inMat - inMeans) / inVar
    return inMat
# 岭回归
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("This matrix is singular, cannot do inverse")
        return
    ws = denom.I * (xMat.T * yMat)
    return ws
# 岭回归测试，返回一系列（30个）lambda下计算出的回归参数，从e-10到e20
def ridgeTest(xArr, yArr):
    xMat = np.mat(xArr)
    yMat = np.mat(yArr).T

    # 归一化数据
    yMean = np.mean(yMat, 0)
    yMat = yMat - yMean  # to eliminate X0 take mean off of Y
    xMat = regularize(xMat)

    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
# 交叉验证
def crossValidation(xArr, yArr, numVal=10):
    m = len(yArr)
    indexList = range(m)
    errorMat = np.zeros((numVal, 30))  # create error mat 30columns numVal rows
    for i in range(numVal):
        trainX = []
        trainY = []
        testX = []
        testY = []
        # 随机取90%的数据训练，剩下10%的数据测试
        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]])

        # 计算30个岭回归lambda参数下的平方误差
        wMat = ridgeTest(trainX, trainY)  # get 30 weight vectors from ridge
        for k in range(30):  # loop over all of the ridge estimates
            matTestX = np.mat(testX)
            matTrainX = np.mat(trainX)
            meanTrain = np.mean(matTrainX, 0)
            varTrain = np.var(matTrainX, 0)
            matTestX = (matTestX - meanTrain) / varTrain  # regularize test with training params
            yEst = matTestX * np.mat(wMat[k, :]).T + np.mean(trainY)  # test ridge results and store
            errorMat[i, k] = rssError(yEst.T.A, np.array(testY))
            print(errorMat[i, k])

    # 取多次迭代的误差平均，使用平均误差最小的岭回归参数计算出的结果
    meanErrors = np.mean(errorMat, 0)  # calc avg performance of the different ridge weight vectors
    minMean = float(min(meanErrors))
    bestWeights = wMat[np.nonzero(meanErrors == minMean)]
    # can unregularize to get model
    # when we regularized we wrote Xreg = (x-meanX)/var(x)
    # we can now write in terms of x not Xreg:  x*w/var(x) - meanX/var(x) +meanY
    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))