Chapter 6 question #2

  1. The lasso is less flexible and will improve prediction accuracy when its increase in bias is smaller than its decrease in variance. This because lasso will shrink coefficient estimates toward zero. This will make the lasso less flexible than least squares.

b)Ridge regression is less flexible and will improve prediction accuracy when its increase in bias is smaller than its decrease in variance. Unlike the lasso, ridge will generally doesn’t set coefficients exactly to zero but it will still constrain the model, to when we compare to least squares it makes it less flexible.

c)Non-linear methods are more flexible and will improve prediction accuracy when their increase in variance is smaller than their decrease in bias.Non-linear method can model more complicated relations between the response and predictors which makes it more flexible.

#Chapter 6 question #9
library(ISLR2)
## Warning: package 'ISLR2' was built under R version 4.5.3
library(tidyverse)
## Warning: package 'ggplot2' was built under R version 4.5.2
## Warning: package 'dplyr' was built under R version 4.5.2
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.4     ✔ readr     2.1.5
## ✔ forcats   1.0.0     ✔ stringr   1.5.2
## ✔ ggplot2   4.0.1     ✔ tibble    3.3.0
## ✔ lubridate 1.9.4     ✔ tidyr     1.3.1
## ✔ purrr     1.1.0     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(glmnet)
## Warning: package 'glmnet' was built under R version 4.5.3
## Loading required package: Matrix
## 
## Attaching package: 'Matrix'
## 
## The following objects are masked from 'package:tidyr':
## 
##     expand, pack, unpack
## 
## Loaded glmnet 4.1-10
library(pls)
## Warning: package 'pls' was built under R version 4.5.3
## 
## Attaching package: 'pls'
## 
## The following object is masked from 'package:stats':
## 
##     loadings
data(College)

#a)
set.seed(1)

train <- sample(
  1:nrow(College),
  nrow(College)/2
)

test <- (-train)

x <- model.matrix(Apps ~ ., College)[,-1]
y <- College$Apps

x.train <- x[train,]
x.test  <- x[test,]

y.train <- y[train]
y.test  <- y[test]
#b
lm.fit <- lm(
  Apps ~ .,
  data = College,
  subset = train
)

lm.pred <- predict(
  lm.fit,
  newdata = College[test,]
)

lm.mse <- mean(
  (y.test - lm.pred)^2
)

lm.mse
## [1] 1135758
#c
set.seed(1)

cv.ridge <- cv.glmnet(
  x.train,
  y.train,
  alpha = 0
)

best.lambda.ridge <- cv.ridge$lambda.min

best.lambda.ridge
## [1] 405.8404
ridge.pred <- predict(
  cv.ridge,
  s = "lambda.min",
  newx = x.test
)

ridge.mse <- mean(
  (y.test - ridge.pred)^2
)

ridge.mse
## [1] 976261.5
#d
set.seed(1)

cv.lasso <- cv.glmnet(
  x.train,
  y.train,
  alpha = 1
)
best.lambda.lasso <- cv.lasso$lambda.min
best.lambda.lasso
## [1] 1.97344
#prediction

lasso.pred <- predict(
  cv.lasso,
  s = "lambda.min",
  newx = x.test
)

lasso.mse <- mean(
  (y.test - lasso.pred)^2
)
lasso.mse
## [1] 1115901
#non-zero
lasso.coef <- predict(
  cv.lasso,
  type = "coefficients",
  s = "lambda.min"
)
sum(lasso.coef != 0)
## [1] 18
#e
set.seed(1)

pcr.fit <- pcr(
  Apps ~ .,
  data = College,
  subset = train,
  scale = TRUE,
  validation = "CV"
)

validationplot(
  pcr.fit,
  val.type = "MSEP"
)

pcr.rmsep <- RMSEP(pcr.fit)

best.M.pcr <- which.min(
  pcr.rmsep$val[1,1,-1]
)

best.M.pcr
## 17 comps 
##       17
pcr.pred <- predict(
  pcr.fit,
  College[test,],
  ncomp = best.M.pcr
)

pcr.mse <- mean(
  (y.test - pcr.pred)^2
)

pcr.mse
## [1] 1135758
#f
set.seed(1)

pls.fit <- plsr(
  Apps ~ .,
  data = College,
  subset = train,
  scale = TRUE,
  validation = "CV"
)

validationplot(
  pls.fit,
  val.type = "MSEP"
)

pls.rmsep <- RMSEP(pls.fit)

best.M.pls <- which.min(
  pls.rmsep$val[1,1,-1]
)

best.M.pls
## 17 comps 
##       17
pls.pred <- predict(
  pls.fit,
  College[test,],
  ncomp = best.M.pls
)

pls.mse <- mean(
  (y.test - pls.pred)^2
)

pls.mse
## [1] 1135758
#g
results <- tibble(
  Model = c(
    "Least Squares",
    "Ridge",
    "LASSO",
    "PCR",
    "PLS"
  ),
  Test_MSE = c(
    lm.mse,
    ridge.mse,
    lasso.mse,
    pcr.mse,
    pls.mse
  )
)

results
## # A tibble: 5 × 2
##   Model         Test_MSE
##   <chr>            <dbl>
## 1 Least Squares 1135758.
## 2 Ridge          976261.
## 3 LASSO         1115901.
## 4 PCR           1135758.
## 5 PLS           1135758.
#chapter 6 question 11

library(ISLR2)
library(tidyverse)
library(leaps)
## Warning: package 'leaps' was built under R version 4.5.3
library(glmnet)
library(pls)

data(Boston)

set.seed(101)

train <- sample(
  1:nrow(Boston),
  nrow(Boston)/2
)

test <- (-train)

x <- model.matrix(crim ~ ., Boston)[,-1]
y <- Boston$crim

x.train <- x[train,]
x.test  <- x[test,]

y.train <- y[train]
y.test  <- y[test]
#11

lm.fit <- lm(
  crim ~ .,
  data = Boston,
  subset = train
)

lm.pred <- predict(
  lm.fit,
  newdata = Boston[test,]
)

lm.mse <- mean(
  (y.test - lm.pred)^2
)

lm.mse
## [1] 35.14697
#11 reg

regfit <- regsubsets(
  crim ~ .,
  data = Boston[train,],
  nvmax = 13
)
summary(regfit)$which
##    (Intercept)    zn indus  chas   nox    rm   age   dis  rad   tax ptratio
## 1         TRUE FALSE FALSE FALSE FALSE FALSE FALSE FALSE TRUE FALSE   FALSE
## 2         TRUE FALSE FALSE FALSE FALSE FALSE FALSE FALSE TRUE FALSE   FALSE
## 3         TRUE FALSE FALSE FALSE FALSE  TRUE FALSE FALSE TRUE FALSE   FALSE
## 4         TRUE FALSE FALSE FALSE FALSE  TRUE FALSE  TRUE TRUE FALSE   FALSE
## 5         TRUE  TRUE FALSE FALSE FALSE  TRUE FALSE  TRUE TRUE FALSE   FALSE
## 6         TRUE  TRUE  TRUE FALSE FALSE  TRUE FALSE  TRUE TRUE FALSE   FALSE
## 7         TRUE  TRUE FALSE FALSE  TRUE  TRUE FALSE  TRUE TRUE FALSE    TRUE
## 8         TRUE  TRUE FALSE FALSE  TRUE  TRUE FALSE  TRUE TRUE  TRUE    TRUE
## 9         TRUE  TRUE  TRUE FALSE  TRUE  TRUE FALSE  TRUE TRUE  TRUE    TRUE
## 10        TRUE  TRUE  TRUE  TRUE  TRUE  TRUE FALSE  TRUE TRUE  TRUE    TRUE
## 11        TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE TRUE  TRUE    TRUE
## 12        TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE TRUE  TRUE    TRUE
##    lstat  medv
## 1  FALSE FALSE
## 2  FALSE  TRUE
## 3  FALSE  TRUE
## 4  FALSE  TRUE
## 5  FALSE  TRUE
## 6  FALSE  TRUE
## 7  FALSE  TRUE
## 8  FALSE  TRUE
## 9  FALSE  TRUE
## 10 FALSE  TRUE
## 11 FALSE  TRUE
## 12  TRUE  TRUE
predict.regsubsets <- function(object, newdata, id) {
  
  form <- as.formula(object$call[[2]])
  mat <- model.matrix(form, newdata)
  
  coefi <- coef(object, id = id)
  xvars <- names(coefi)
  
  mat[, xvars, drop = FALSE] %*% coefi
}

max_size <- nrow(summary(regfit)$which)

val.errors <- rep(NA, max_size)

for (i in 1:max_size) {
  
  pred <- predict.regsubsets(
    object = regfit,
    newdata = Boston[test, ],
    id = i
  )
  
  val.errors[i] <- mean(
    (Boston$crim[test] - pred)^2
  )
}

best.size <- which.min(val.errors)

best.size
## [1] 2
val.errors[best.size]
## [1] 33.71521
coef(regfit, id = best.size)
## (Intercept)         rad        medv 
##   2.0413737   0.5639283  -0.1719006
#11 ridge
set.seed(101)

cv.ridge <- cv.glmnet(
  x.train,
  y.train,
  alpha = 0
)

ridge.pred <- predict(
  cv.ridge,
  s = "lambda.min",
  newx = x.test
)

ridge.mse <- mean(
  (y.test-ridge.pred)^2
)

ridge.mse
## [1] 34.00884
#11 lasso

set.seed(101)

cv.lasso <- cv.glmnet(
  x.train,
  y.train,
  alpha = 1
)

lasso.pred <- predict(
  cv.lasso,
  s = "lambda.min",
  newx = x.test
)

lasso.mse <- mean(
  (y.test-lasso.pred)^2
)

lasso.mse
## [1] 34.47705
#11 pcr 
set.seed(101)

pcr.fit <- pcr(
  crim ~ .,
  data = Boston,
  subset = train,
  scale = TRUE,
  validation = "CV"
)

pcr.rmsep <- RMSEP(pcr.fit)

best.M <- which.min(
  pcr.rmsep$val[1,1,-1]
)

best.M
## 12 comps 
##       12
pcr.pred <- predict(
  pcr.fit,
  Boston[test,],
  ncomp = best.M
)

pcr.mse <- mean(
  (y.test-pcr.pred)^2
)

pcr.mse
## [1] 35.14697
results <- tibble(
  Model = c(
    "Least Squares",
    "Best Subset",
    "Ridge",
    "LASSO",
    "PCR"
  ),
  Test_MSE = c(
    lm.mse,
    min(val.errors),
    ridge.mse,
    lasso.mse,
    pcr.mse
  )
)

results
## # A tibble: 5 × 2
##   Model         Test_MSE
##   <chr>            <dbl>
## 1 Least Squares     35.1
## 2 Best Subset       33.7
## 3 Ridge             34.0
## 4 LASSO             34.5
## 5 PCR               35.1