-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathridge_regression_and_lasso.R
58 lines (42 loc) · 1.16 KB
/
ridge_regression_and_lasso.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
library(glmnet)
# data
library(ISLR)
summary(Hitters)
# the "-1" is to remove the intercept
x = model.matrix(Salary~.-1,data=Hitters)
y = Hitters$Salary
####################
# Ridge Regression #
####################
fit.ridge = glmnet(x,y,alpha = 0)
plot(fit.ridge,xvar="lambda",label=TRUE)
# do cross validation
cv.ridge = cv.glmnet(x,y,alpha = 0)
plot(cv.ridge)
##################
# Lasso #
##################
fit.lasso = glmnet(x,y) # default alpha=1 lasso
plot(fit.lasso, xvar="lambda",label=TRUE)
plot(fit.lasso, xvar="dev",label=TRUE)
# do cross validation
cv.lasso = cv.glmnet(x,y)
plot(cv.lasso)
# coeff for the best model
coef(cv.lasso)
# use train/validation
set.seed(1)
train = sample(seq(nrow(Hitters)), 180, replace = FALSE)
# train lasso
lasso.tr = glmnet(x[train,],y[train])
# compute predictions on test set
pred = predict(lasso.tr,x[-train,])
dim(pred)
# compute test error
rmse = sqrt(apply((y[-train]-pred)^2,2,mean))
plot(log(lasso.tr$lambda), rmse, type="b", xlab="Log(lambda)")
# get best lambda looking for minimum test error
lam.best = lasso.tr$lambda[order(rmse)[1]]
lam.best
# get coefficients for best lambda
coef(lasso.tr, s=lam.best)