Chris Serrano - Assignment 5

Question 2

a

Statement iii is correct: Lasso regression, relative to least squares, is less flexible and will give improved prediction accuracy when its increase in bias is less than its decrease in variance.

b

Statement iii is correct: Ridge regression, relative to least squares, is less flexible and will give improved prediction accuracy when its increase in bias is less than its decrease in variance.

c

Statement ii is correct: A non-linear method, relative to least squares, is more flexible and will give improved prediction accuracy when its increase in variance is less than its decrease in bias.

---

Question 9

a

library(ISLR2)
library(glmnet)
library(pls)

set.seed(1)

# Split dataset 50/50
train_idx <- sample(nrow(College), nrow(College) / 2)
train_set <- College[train_idx, ]
test_set <- College[-train_idx, ]

# Prepare model matrices for glmnet
x_train <- model.matrix(Apps ~ ., data = train_set)[, -1]
y_train <- train_set$Apps
x_test <- model.matrix(Apps ~ ., data = test_set)[, -1]
y_test <- test_set$Apps

b

lm.fit <- lm(Apps ~ ., data = train_set)
lm.pred <- predict(lm.fit, test_set)
lm.err <- mean((test_set$Apps - lm.pred)^2)
cat("OLS Linear Model Test MSE:", lm.err, "\n")

OLS Linear Model Test MSE: 1135758 

c

set.seed(1)
cv.ridge <- cv.glmnet(x_train, y_train, alpha = 0)
bestlambda.ridge <- cv.ridge$lambda.min
cat("Ridge Regression Best Lambda:", bestlambda.ridge, "\n")

Ridge Regression Best Lambda: 405.8404 

ridge.pred <- predict(cv.ridge, s = bestlambda.ridge, newx = x_test)
ridge.err <- mean((y_test - ridge.pred)^2)
cat("Ridge Regression Test MSE:", ridge.err, "\n")

Ridge Regression Test MSE: 976261.5 

d

set.seed(1)
cv.lasso <- cv.glmnet(x_train, y_train, alpha = 1)
bestlambda.lasso <- cv.lasso$lambda.min
cat("Lasso Best Lambda:", bestlambda.lasso, "\n")

Lasso Best Lambda: 1.97344 

lasso.pred <- predict(cv.lasso, s = bestlambda.lasso, newx = x_test)
lasso.err <- mean((y_test - lasso.pred)^2)
cat("Lasso Test MSE:", lasso.err, "\n")

Lasso Test MSE: 1115901 

# Extract coefficients for best lambda
lasso.coef <- predict(cv.lasso, type = "coefficients", s = bestlambda.lasso)
# Count non-zero coefficients excluding the intercept
non_zero_coefs <- sum(lasso.coef[-1, 1] != 0)
cat("Number of Non-zero Coefficients (excluding intercept):", non_zero_coefs, "\n")

Number of Non-zero Coefficients (excluding intercept): 17 

print(lasso.coef)

18 x 1 sparse Matrix of class "dgCMatrix"
                s=1.97344
(Intercept) -7.688896e+02
PrivateYes  -3.127034e+02
Accept       1.762718e+00
Enroll      -1.318195e+00
Top10perc    6.482356e+01
Top25perc   -2.081406e+01
F.Undergrad  7.119149e-02
P.Undergrad  1.246161e-02
Outstate    -1.049091e-01
Room.Board   2.088305e-01
Books        2.926466e-01
Personal     3.955068e-03
PhD         -1.455463e+01
Terminal     5.395858e+00
S.F.Ratio    2.171398e+01
perc.alumni  5.088260e-01
Expend       4.824455e-02
Grad.Rate    7.036148e+00

e

set.seed(1)
pcr.fit <- pcr(Apps ~ ., data = train_set, scale = TRUE, validation = "CV")
validationplot(pcr.fit, val.type = "MSEP")

# Find component size M that minimizes CV error
pcr.cv.err <- RMSEP(pcr.fit)$val[1,, -1] # Exclude 0 components
best_M_pcr <- which.min(pcr.cv.err)
cat("PCR Best Number of Components (M):", best_M_pcr, "\n")

PCR Best Number of Components (M): 17 

pcr.pred <- predict(pcr.fit, test_set, ncomp = best_M_pcr)
pcr.err <- mean((test_set$Apps - pcr.pred)^2)
cat("PCR Test MSE:", pcr.err, "\n")

PCR Test MSE: 1135758 

f

set.seed(1)
pls.fit <- plsr(Apps ~ ., data = train_set, scale = TRUE, validation = "CV")
validationplot(pls.fit, val.type = "MSEP")

# Find component size M that minimizes CV error
pls.cv.err <- RMSEP(pls.fit)$val[1,, -1] # Exclude 0 components
best_M_pls <- which.min(pls.cv.err)
cat("PLS Best Number of Components (M):", best_M_pls, "\n")

PLS Best Number of Components (M): 17 

pls.pred <- predict(pls.fit, test_set, ncomp = best_M_pls)
pls.err <- mean((test_set$Apps - pls.pred)^2)
cat("PLS Test MSE:", pls.err, "\n")

PLS Test MSE: 1135758 

g

test_avg <- mean(test_set$Apps)
test_var <- sum((test_set$Apps - test_avg)^2)

# Compute R-squared values
lm.r2 <- 1 - sum((test_set$Apps - lm.pred)^2) / test_var
ridge.r2 <- 1 - sum((y_test - ridge.pred)^2) / test_var
lasso.r2 <- 1 - sum((y_test - lasso.pred)^2) / test_var
pcr.r2 <- 1 - sum((test_set$Apps - pcr.pred)^2) / test_var
pls.r2 <- 1 - sum((test_set$Apps - pls.pred)^2) / test_var

# Construct comparison table
results <- data.frame(
  Method = c("OLS Linear Model", "Ridge Regression", "Lasso Regression", "PCR", "PLS"),
  Test_MSE = c(lm.err, ridge.err, lasso.err, pcr.err, pls.err),
  Test_R2 = c(lm.r2, ridge.r2, lasso.r2, pcr.r2, pls.r2)
)
print(results)

Interpretation

All models perform exceptionally well, explaining approximately 90% to 93% of the variance in the number of college applications in the test set.

• OLS performs well showing that a linear model fits the data nicely.
• Ridge Regression achieves a test MSE of 976261.5, performing slightly worse than OLS.
• Lasso Regression improves the test MSE showing the benefits of feature selection.
• PLS performs almost identically to Lasso.
• PCR performs the same as PLS.

---

Question 11

a

We will compare:
1. Best Subset Selection
2. Ridge Regression
3. Lasso Regression
4. PCR 

library(leaps)

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] %*% coefi
}

run the 10-fold cross-validation loop.

set.seed(1)
k <- 10

# Partition the Boston data into 10 folds
folds <- sample(1:k, nrow(Boston), replace = TRUE)

err.best <- numeric(k)
err.ridge <- numeric(k)
err.lasso <- numeric(k)
err.pcr <- numeric(k)

# Model matrices for glmnet
x_boston <- model.matrix(crim ~ ., data = Boston)[, -1]
y_boston <- Boston$crim

for (j in 1:k) {
  train_fold <- (folds != j)
  test_fold <- (folds == j)
  
  # 1. Best Subset Selection
  fit.best <- regsubsets(crim ~ ., data = Boston[train_fold, ], nvmax = 13)
  best_size_bic <- which.min(summary(fit.best)$bic)
  pred.best <- predict(fit.best, Boston[test_fold, ], id = best_size_bic)
  err.best[j] <- mean((y_boston[test_fold] - pred.best)^2)
  
  # 2. Ridge Regression
  cv.r <- cv.glmnet(x_boston[train_fold, ], y_boston[train_fold], alpha = 0)
  pred.r <- predict(cv.r, s = cv.r$lambda.min, newx = x_boston[test_fold, ])
  err.ridge[j] <- mean((y_boston[test_fold] - pred.r)^2)
  
  # 3. Lasso Regression
  cv.l <- cv.glmnet(x_boston[train_fold, ], y_boston[train_fold], alpha = 1)
  pred.l <- predict(cv.l, s = cv.l$lambda.min, newx = x_boston[test_fold, ])
  err.lasso[j] <- mean((y_boston[test_fold] - pred.l)^2)
  
  # 4. PCR
  fit.pcr <- pcr(crim ~ ., data = Boston[train_fold, ], scale = TRUE, validation = "CV")
  pcr.cv.err <- RMSEP(fit.pcr)$val[1,, -1]
  best_M <- which.min(pcr.cv.err)
  pred.pcr <- predict(fit.pcr, Boston[test_fold, ], ncomp = best_M)
  err.pcr[j] <- mean((y_boston[test_fold] - pred.pcr)^2)
}

# Print average CV error rate for each method
cat("Best Subset CV MSE:", mean(err.best), "\n")

Best Subset CV MSE: 44.34247 

cat("Ridge CV MSE:", mean(err.ridge), "\n")

Ridge CV MSE: 42.61162 

cat("Lasso CV MSE:", mean(err.lasso), "\n")

Lasso CV MSE: 42.29416 

cat("PCR CV MSE:", mean(err.pcr), "\n")

PCR CV MSE: 42.50125 

b

We propose the Lasso Regression model as the best model for this dataset. Lasso achieves the lowest cross-validation test MSE indicating it has the best predictive accuracy.

c

lasso.full <- cv.glmnet(x_boston, y_boston, alpha = 1)
best_lambda <- lasso.full$lambda.min
coef(lasso.full, s = best_lambda)

14 x 1 sparse Matrix of class "dgCMatrix"
             s=0.04259589
(Intercept)  1.325949e+01
zn           3.710115e-02
indus       -7.429552e-02
chas        -5.986112e-01
nox         -7.539411e+00
rm           2.651568e-01
age          .           
dis         -8.211252e-01
rad          5.178546e-01
tax         -2.017843e-06
ptratio     -2.022884e-01
black       -7.547405e-03
lstat        1.265768e-01
medv        -1.636525e-01

No, the proposed Lasso model does not involve all features.

The Lasso regression model performs automatic variable selection, simplifying the model complexity and making it easier to interpret. It drops variables that are non-informative or redundant, helping to lower the variance of our model predictions without adding significant bias.