当前位置:   article > 正文

R语言实现LASSO回归——自己编写LASSO回归算法

R语言实现LASSO回归——自己编写LASSO回归算法

93ff826e1db77836f1ebff56d32cc96b.png

  1. 本文约500字,建议阅读5分钟
  2. 这篇文章中我们可以编写自己的代码来计算套索(lasso)回归。

相关视频

2921660848f35183252eb26175786923.png

我们必须定义阈值函数

1860d0185e9f612ddc3d72690d64289a.png

R函数是:

  1. thresh = function(x,a){
  2. sign(x) * pmax(abs(x)-a,0)
  3. }

要解决我们的优化问题,设置

a3b00b34c98d1d7681fca0427c25db2b.png

这样就可以等效地写出优化问题
55a94fdab20a89c6d4b5caa042d58e08.png

因此

82c519e1177422b720f28df4833a6021.png

一个得到

2f7bf49adf46aae33fa79e4be16144fb.png

d406f7a33713ab2e017e73a3bbddda61.png

同样,如果有权重ω=(ωi),则按坐标更新将变为

a261d0cc462bebdd94ccb80c90660d6a.png

计算此分量下降的代码是:

  1. lasso = function(X,y,beta,lambda,tol=1e-6,maxiter=1000){
  2. beta0 = sum(y-X%*%beta /(length(y))
  3. beta0list[1] = beta0
  4. for (j in 1:maxiter){
  5. for (k in 1:length beta)){
  6. r = y - X[,-k]%*%beta[-k] - beta0*rep(1,length(y )
  7. beta[k] = (1/sum(omega*X[,k]^2) *
  8. threshog(t(omega*r)%*%X[,k ,length(y *lambda)
  9. }
  10. beta0 = sum(y-X%*%beta)/(length(y))
  11. obj[j] = (1/2)*(1/length(y))*norm(omega*(y - X%*%beta -
  12. beta0*rep(1,length(y))),'F')^2 + lambda*sum(abs(beta))
  13. if (norm(rbind(beta0list[j],betalist[[j]]) -
  14. rbind(beta0,beta),'F') ) { break }

例如,考虑以下(简单)数据集,其中包含三个协变量:

chicago = read.table("data.txt",header=TRUE,sep=";")

我们可以“标准化”

  1. for(j in 1:3) X[,j] = (X[,j]-mean(X[,j]))/sd(X[,j])
  2. y = (y-mean(y))/sd(y)

要初始化算法,使用OLS估算

lm(y~0+.,)$coef

例如

  1. lasso(X,y,beta_init,lambda=.001)
  2. $obj
  3. [1] 0.001014426 0.001008009 0.001009558 0.001011094 0.001011119 0.001011119
  4. $beta
  5. [,1]
  6. X_1 0.0000000
  7. X_2 0.3836087
  8. X_3 -0.5026137
  9. $intercept
  10. [1] 2.060999e-16

我们可以通过循环获得标准的lasso图

6a7d15b16b74428d309f4fda2d51bff1.png


原文链接:http://tecdat.cn/?p=18840

编辑:于腾凯

校对:林亦霖

945cb6c4452b19015ae36b0959bfba34.png

声明:本文内容由网友自发贡献,版权归原作者所有,本站不承担相应法律责任。如您发现有侵权的内容,请联系我们。转载请注明出处:【wpsshop博客】
推荐阅读
  

闽ICP备14008679号