Data_analytics_lab_7
Yagna Praseeda Atmuri
2024-03-19
R Markdown
Answer:
library(ISLR)## Warning: package 'ISLR' was built under R version 4.3.2library(caret)## Warning: package 'caret' was built under R version 4.3.2## Loading required package: ggplot2## Warning: package 'ggplot2' was built under R version 4.3.2## Loading required package: latticelibrary(glmnet)## Loading required package: Matrix## Loaded glmnet 4.1-8library(pls)## Warning: package 'pls' was built under R version 4.3.3## ## Attaching package: 'pls'## The following object is masked from 'package:caret': ## ## R2## The following object is masked from 'package:stats': ## ## loadingslibrary(lars)## Loaded lars 1.3
data("College")
set.seed(123)- Split the data into training and test sets
train_indices <- createDataPartition(College$Apps, p = 0.7, list = FALSE)
train_data <- College[train_indices, ]
test_data <- College[-train_indices, ]- Fit a linear model using least squares on the training set
lm_model <- lm(Apps ~ ., data = train_data)
lm_pred <- predict(lm_model, newdata = test_data)
lm_test_error <- sqrt(mean((test_data$Apps - lm_pred)^2))
cat("Linear Model Test Error:", lm_test_error, "\n")## Linear Model Test Error: 1371.887- Fit a ridge regression model on the training set, with lambda chosen by cross-validation
ridge_model <- cv.glmnet(as.matrix(train_data[, -1]), train_data$Apps, alpha = 0)
ridge_pred <- predict(ridge_model, newx = as.matrix(test_data[, -1]), s = "lambda.min")
ridge_test_error <- sqrt(mean((test_data$Apps - ridge_pred)^2))
cat("Ridge Regression Test Error:", ridge_test_error, "\n")## Ridge Regression Test Error: 796.2182- Fit a lasso model on the training set, with lambda chosen by cross-validation
lasso_model <- cv.glmnet(as.matrix(train_data[, -1]), train_data$Apps, alpha = 1)
lasso_pred <- predict(lasso_model, newx = as.matrix(test_data[, -1]), s = "lambda.min")
lasso_test_error <- sqrt(mean((test_data$Apps - lasso_pred)^2))
num_nonzero_coeffs <- sum(coef(lasso_model, s = "lambda.min") != 0)
cat("Lasso Test Error:", lasso_test_error, "\n")## Lasso Test Error: 140.5866cat("Number of Non-zero Coefficients in Lasso Model:", num_nonzero_coeffs, "\n")## Number of Non-zero Coefficients in Lasso Model: 2- Fit a PCR model on the training set, with M chosen by cross-validation
pcr_model <- train(Apps ~ ., data = train_data, method = "pcr", tuneLength = 10)
pcr_pred <- predict(pcr_model, newdata = test_data)
pcr_test_error <- sqrt(mean((test_data$Apps - pcr_pred)^2))
cat("PCR Model Test Error:", pcr_test_error, "\n")## PCR Model Test Error: 1465.604cat("Value of M selected by PCR cross-validation:", pcr_model$bestTune$npc, "\n")## Value of M selected by PCR cross-validation:- Fit a PLS model on the training set, with M chosen by cross-validation
pls_model <- train(Apps ~ ., data = train_data, method = "pls", tuneLength = 10)
pls_pred <- predict(pls_model, newdata = test_data)
pls_test_error <- sqrt(mean((test_data$Apps - pls_pred)^2))
cat("PLS Model Test Error:", pls_test_error, "\n")## PLS Model Test Error: 1425.504cat("Value of M selected by PLS cross-validation:", pls_model$bestTune$ncomp, "\n")## Value of M selected by PLS cross-validation: 10- Comment on the results
cat("Linear Model Test Error:", lm_test_error, "\n")## Linear Model Test Error: 1371.887cat("Ridge Regression Test Error:", ridge_test_error, "\n")## Ridge Regression Test Error: 796.2182cat("Lasso Test Error:", lasso_test_error, "\n")## Lasso Test Error: 140.5866cat("PCR Model Test Error:", pcr_test_error, "\n")## PCR Model Test Error: 1465.604cat("PLS Model Test Error:", pls_test_error, "\n")## PLS Model Test Error: 1425.504In conclusion, the Lasso model stood out as the most effective in
predicting the number of college applications received, followed by the
Ridge Regression model. The PCR and PLS models did not provide
significant improvements over the basic Linear Model in this
scenario.
8. Answer:
library(ISLR)
library(leaps)## Warning: package 'leaps' was built under R version 4.3.2library(glmnet)
library(pls)
library(caret)a.)
set.seed(123)
X <- rnorm(100, mean = 0, sd = 1)
epsilon <- rnorm(100, mean = 0, sd = 1)b.)
beta0 <- 2
beta1 <- 3
beta2 <- -1
beta3 <- 0.5
Y <- beta0 + beta1*X + beta2*X^2 + beta3*X^3 + epsilon
Y## [1] -0.79399839 1.50727197 5.89335819 1.85918641 1.42060986
## [6] 6.68110649 2.43435978 -6.07579487 -1.07456912 1.33913862
## [11] 4.51558973 3.58123164 1.45599939 2.26491349 0.45705827
## [16] 7.32169336 3.41306992 -12.21116810 2.93496178 -0.71887667
## [21] -2.83486363 0.34590925 -3.16128573 -1.16767185 1.45597679
## [26] -8.95623001 4.34087655 2.51636081 -4.40876921 5.10361586
## [31] 4.58085268 1.46637690 4.28396965 3.77935766 2.01378231
## [36] 4.88631793 1.97926595 2.55026066 2.88328136 -0.75760249
## [41] -0.03259336 1.06632699 -5.98265513 7.38961593 3.44448991
## [46] -3.86993452 -0.86542365 1.01937237 6.06890288 0.45562220
## [51] 3.49165188 2.68257543 2.20171394 4.50610349 1.14650784
## [56] 5.71303429 -6.33935066 3.13953192 3.33414620 2.23164785
## [61] 4.07486176 -0.87185157 -0.38930209 0.61943413 -3.39657052
## [66] 3.13066587 3.82532788 1.67249706 4.82531704 8.62446773
## [71] 0.01121773 -16.35101994 4.48029562 1.31953205 -1.44155426
## [76] 3.46826747 1.09082678 -3.75135143 2.95054278 1.10433011
## [81] 0.95393329 4.29918109 0.37551838 2.78617503 1.04828702
## [86] 2.70635188 5.85716079 3.24210681 1.65271555 4.38544646
## [91] 4.69822591 3.10222798 2.76058889 -1.29712865 4.17931975
## [96] 1.72798281 9.61084819 4.79763616 0.61963122 -3.85897032c.)
data <- data.frame(X = X, X2 = X^2, X3 = X^3, X4 = X^4, X5 = X^5, X6 = X^6, X7 = X^7, X8 = X^8, X9 = X^9, X10 = X^10, Y = Y)
subset_model <- regsubsets(Y ~ ., data = data, nbest = 1, nvmax = 10)
summary(subset_model, Taylor = TRUE)## Subset selection object
## Call: regsubsets.formula(Y ~ ., data = data, nbest = 1, nvmax = 10)
## 10 Variables (and intercept)
## Forced in Forced out
## X FALSE FALSE
## X2 FALSE FALSE
## X3 FALSE FALSE
## X4 FALSE FALSE
## X5 FALSE FALSE
## X6 FALSE FALSE
## X7 FALSE FALSE
## X8 FALSE FALSE
## X9 FALSE FALSE
## X10 FALSE FALSE
## 1 subsets of each size up to 10
## Selection Algorithm: exhaustive
## X X2 X3 X4 X5 X6 X7 X8 X9 X10
## 1 ( 1 ) "*" " " " " " " " " " " " " " " " " " "
## 2 ( 1 ) "*" "*" " " " " " " " " " " " " " " " "
## 3 ( 1 ) "*" "*" "*" " " " " " " " " " " " " " "
## 4 ( 1 ) "*" "*" " " " " "*" " " " " " " "*" " "
## 5 ( 1 ) "*" " " "*" "*" " " "*" " " "*" " " " "
## 6 ( 1 ) "*" " " " " "*" "*" "*" "*" "*" " " " "
## 7 ( 1 ) "*" " " " " "*" "*" "*" "*" "*" " " "*"
## 8 ( 1 ) "*" "*" " " "*" "*" "*" "*" "*" " " "*"
## 9 ( 1 ) "*" "*" " " "*" "*" "*" "*" "*" "*" "*"
## 10 ( 1 ) "*" "*" "*" "*" "*" "*" "*" "*" "*" "*"best_cp_model <- which.min(summary(subset_model)$cp)
cp_model_summary <- summary(subset_model)$outmat[best_cp_model, ]best_bic_model <- which.min(summary(subset_model)$bic)
bic_model_summary <- summary(subset_model)$outmat[best_bic_model, ]best_adjr2_model <- which.max(summary(subset_model)$adjr2)
adjr2_model_summary <- summary(subset_model)$outmat[best_adjr2_model, ]coef(subset_model, id = best_cp_model)## (Intercept) X X2 X3
## 1.9703939 2.9204462 -1.0915431 0.5204363coef(subset_model, id = best_bic_model)## (Intercept) X X2 X3
## 1.9703939 2.9204462 -1.0915431 0.5204363coef(subset_model, id = best_adjr2_model)## (Intercept) X X4 X5 X6 X7
## 1.90032830 3.11707619 -1.53472963 0.30756427 0.57013129 -0.04513887
## X8
## -0.06125841plot(subset_model, scale = "Cp")
plot(subset_model, scale = "bic")
plot(subset_model, scale = "adjr2")
d.)
forward_model <- regsubsets(Y ~ ., data = data, method = "forward", nvmax = 10)
summary(forward_model, scale = TRUE)## Subset selection object
## Call: regsubsets.formula(Y ~ ., data = data, method = "forward", nvmax = 10)
## 10 Variables (and intercept)
## Forced in Forced out
## X FALSE FALSE
## X2 FALSE FALSE
## X3 FALSE FALSE
## X4 FALSE FALSE
## X5 FALSE FALSE
## X6 FALSE FALSE
## X7 FALSE FALSE
## X8 FALSE FALSE
## X9 FALSE FALSE
## X10 FALSE FALSE
## 1 subsets of each size up to 10
## Selection Algorithm: forward
## X X2 X3 X4 X5 X6 X7 X8 X9 X10
## 1 ( 1 ) "*" " " " " " " " " " " " " " " " " " "
## 2 ( 1 ) "*" "*" " " " " " " " " " " " " " " " "
## 3 ( 1 ) "*" "*" "*" " " " " " " " " " " " " " "
## 4 ( 1 ) "*" "*" "*" " " " " "*" " " " " " " " "
## 5 ( 1 ) "*" "*" "*" " " " " "*" " " " " "*" " "
## 6 ( 1 ) "*" "*" "*" " " "*" "*" " " " " "*" " "
## 7 ( 1 ) "*" "*" "*" " " "*" "*" " " " " "*" "*"
## 8 ( 1 ) "*" "*" "*" "*" "*" "*" " " " " "*" "*"
## 9 ( 1 ) "*" "*" "*" "*" "*" "*" " " "*" "*" "*"
## 10 ( 1 ) "*" "*" "*" "*" "*" "*" "*" "*" "*" "*"best_forward_cp_model <- which.min(summary(forward_model)$cp)
forward_cp_model_summary <- summary(forward_model)$outmat[best_forward_cp_model, ]
best_forward_bic_model <- which.min(summary(forward_model)$bic)
forward_bic_model_summary <- summary(forward_model)$outmat[best_forward_bic_model, ]
best_forward_adjr2_model <- which.max(summary(forward_model)$adjr2)
forward_adjr2_model_summary <- summary(forward_model)$outmat[best_forward_adjr2_model, ]backward_model <- regsubsets(Y ~ ., data = data, method = "backward", nvmax = 10)
summary(backward_model, scale = TRUE)## Subset selection object
## Call: regsubsets.formula(Y ~ ., data = data, method = "backward", nvmax = 10)
## 10 Variables (and intercept)
## Forced in Forced out
## X FALSE FALSE
## X2 FALSE FALSE
## X3 FALSE FALSE
## X4 FALSE FALSE
## X5 FALSE FALSE
## X6 FALSE FALSE
## X7 FALSE FALSE
## X8 FALSE FALSE
## X9 FALSE FALSE
## X10 FALSE FALSE
## 1 subsets of each size up to 10
## Selection Algorithm: backward
## X X2 X3 X4 X5 X6 X7 X8 X9 X10
## 1 ( 1 ) "*" " " " " " " " " " " " " " " " " " "
## 2 ( 1 ) "*" " " " " "*" " " " " " " " " " " " "
## 3 ( 1 ) "*" " " " " "*" "*" " " " " " " " " " "
## 4 ( 1 ) "*" " " " " "*" "*" "*" " " " " " " " "
## 5 ( 1 ) "*" " " " " "*" "*" "*" " " "*" " " " "
## 6 ( 1 ) "*" " " " " "*" "*" "*" "*" "*" " " " "
## 7 ( 1 ) "*" " " " " "*" "*" "*" "*" "*" " " "*"
## 8 ( 1 ) "*" "*" " " "*" "*" "*" "*" "*" " " "*"
## 9 ( 1 ) "*" "*" " " "*" "*" "*" "*" "*" "*" "*"
## 10 ( 1 ) "*" "*" "*" "*" "*" "*" "*" "*" "*" "*"best_backward_cp_model <- which.min(summary(backward_model)$cp)
backward_cp_model_summary <- summary(backward_model)$outmat[best_backward_cp_model, ]
best_backward_bic_model <- which.min(summary(backward_model)$bic)
backward_bic_model_summary <- summary(backward_model)$outmat[best_backward_bic_model, ]
best_backward_adjr2_model <- which.max(summary(backward_model)$adjr2)
backward_adjr2_model_summary <- summary(backward_model)$outmat[best_backward_adjr2_model, ]e.)
X_matrix <- model.matrix(Y ~ . - 1, data = data)
lasso_model <- cv.glmnet(X_matrix, Y, alpha = 1)
plot(lasso_model)
optimal_lambda <- lasso_model$lambda.min
optimal_lambda## [1] 0.05387053lasso_model_optimal <- glmnet(X_matrix, Y, alpha = 1, lambda = optimal_lambda)
coef(lasso_model_optimal)## 11 x 1 sparse Matrix of class "dgCMatrix"
## s0
## (Intercept) 1.9277391
## X 2.8850716
## X2 -1.0313866
## X3 0.5046393
## X4 .
## X5 .
## X6 .
## X7 .
## X8 .
## X9 .
## X10 .f.)
set.seed(123)
X <- rnorm(100, mean = 0, sd = 1)
beta0 <- 2
beta7 <- 3
epsilon <- rnorm(100, mean = 0, sd = 1)
Y <- beta0 + beta7 * X^7 + epsilon
subset_model <- regsubsets(Y ~ poly(X, 7), data = data.frame(X = X, Y = Y), nbest = 1)
summary(subset_model)## Subset selection object
## Call: regsubsets.formula(Y ~ poly(X, 7), data = data.frame(X = X, Y = Y),
## nbest = 1)
## 7 Variables (and intercept)
## Forced in Forced out
## poly(X, 7)1 FALSE FALSE
## poly(X, 7)2 FALSE FALSE
## poly(X, 7)3 FALSE FALSE
## poly(X, 7)4 FALSE FALSE
## poly(X, 7)5 FALSE FALSE
## poly(X, 7)6 FALSE FALSE
## poly(X, 7)7 FALSE FALSE
## 1 subsets of each size up to 7
## Selection Algorithm: exhaustive
## poly(X, 7)1 poly(X, 7)2 poly(X, 7)3 poly(X, 7)4 poly(X, 7)5
## 1 ( 1 ) " " " " "*" " " " "
## 2 ( 1 ) "*" " " "*" " " " "
## 3 ( 1 ) "*" " " "*" " " "*"
## 4 ( 1 ) "*" " " "*" " " "*"
## 5 ( 1 ) "*" " " "*" "*" "*"
## 6 ( 1 ) "*" " " "*" "*" "*"
## 7 ( 1 ) "*" "*" "*" "*" "*"
## poly(X, 7)6 poly(X, 7)7
## 1 ( 1 ) " " " "
## 2 ( 1 ) " " " "
## 3 ( 1 ) " " " "
## 4 ( 1 ) " " "*"
## 5 ( 1 ) " " "*"
## 6 ( 1 ) "*" "*"
## 7 ( 1 ) "*" "*"X_matrix <- model.matrix(Y ~ poly(X, 7) - 1, data = data.frame(X = X, Y = Y))
lasso_model <- glmnet(X_matrix, Y, alpha = 1)
plot(lasso_model)
optimal_lambda <- lasso_model$lambda.min
optimal_lambda## NULLlasso_model_optimal <- glmnet(X_matrix, Y, alpha = 1, lambda = optimal_lambda)
coef(lasso_model_optimal)## 8 x 45 sparse Matrix of class "dgCMatrix"## [[ suppressing 45 column names 's0', 's1', 's2' ... ]]##
## (Intercept) 9.993166 9.993166 9.993166 9.993166 9.993166 9.993166
## poly(X, 7)1 . . . 20.737077 100.581662 173.333074
## poly(X, 7)2 . . . . . .
## poly(X, 7)3 . 105.549843 201.722928 289.352262 369.196847 441.948259
## poly(X, 7)4 . . . . . .
## poly(X, 7)5 . . . . . .
## poly(X, 7)6 . . . . . .
## poly(X, 7)7 . . . . . .
##
## (Intercept) 9.993166 9.993166 9.993166 9.993166 9.993166 9.993166
## poly(X, 7)1 239.621451 300.020952 355.054727 405.199454 450.889461 492.520494
## poly(X, 7)2 . . . . . .
## poly(X, 7)3 508.236636 568.636137 623.669912 673.814639 719.504646 761.135679
## poly(X, 7)4 . . . . . .
## poly(X, 7)5 . . . . . 31.853104
## poly(X, 7)6 . . . . . .
## poly(X, 7)7 . . . . . .
##
## (Intercept) 9.993166 9.993166 9.993166 9.993166 9.993166 9.993166
## poly(X, 7)1 530.453140 565.015955 596.508305 625.202961 651.348463 675.171270
## poly(X, 7)2 . . . . . .
## poly(X, 7)3 799.068325 833.631140 865.123490 893.818146 919.963648 943.786455
## poly(X, 7)4 . . . . . .
## poly(X, 7)5 69.785751 104.348566 135.840916 164.535572 190.681074 214.503881
## poly(X, 7)6 . . . . . .
## poly(X, 7)7 . . . . . .
##
## (Intercept) 9.993166 9.993166 9.993166 9.993166 9.993166
## poly(X, 7)1 696.877725 716.655839 734.676919 751.097056 766.058474
## poly(X, 7)2 . . . . .
## poly(X, 7)3 965.492910 985.271024 1003.292104 1019.712241 1034.673659
## poly(X, 7)4 . . . . .
## poly(X, 7)5 236.210336 255.988450 274.009530 290.429667 305.391085
## poly(X, 7)6 . . . . .
## poly(X, 7)7 . . . . .
##
## (Intercept) 9.993166 9.993166 9.993166 9.993166 9.993166
## poly(X, 7)1 779.690760 792.111992 803.429755 813.742080 823.138286
## poly(X, 7)2 . . . . .
## poly(X, 7)3 1048.305945 1060.727177 1072.044940 1082.357265 1091.753471
## poly(X, 7)4 . . . . .
## poly(X, 7)5 319.023371 331.444603 342.762366 353.074691 362.470897
## poly(X, 7)6 . . . . .
## poly(X, 7)7 . . . . .
##
## (Intercept) 9.993166 9.993166 9.993166 9.993166 9.993166
## poly(X, 7)1 831.699759 839.500655 846.608540 853.084981 858.986072
## poly(X, 7)2 . . . . .
## poly(X, 7)3 1100.314944 1108.115840 1115.223725 1121.700166 1127.601257
## poly(X, 7)4 . . . . .
## poly(X, 7)5 371.032370 378.833266 385.941151 392.417592 398.318683
## poly(X, 7)6 . . . . .
## poly(X, 7)7 . . . 3.465811 9.366903
##
## (Intercept) 9.993166 9.993166 9.993166 9.993166 9.993166
## poly(X, 7)1 864.362927 869.262117 873.726076 877.793469 881.499526
## poly(X, 7)2 . . . . .
## poly(X, 7)3 1132.978112 1137.877302 1142.341261 1146.408654 1150.114711
## poly(X, 7)4 . -2.401029 -6.864989 -10.932382 -14.638439
## poly(X, 7)5 403.695538 408.594727 413.058687 417.126080 420.832137
## poly(X, 7)6 . . . . .
## poly(X, 7)7 14.743757 19.642947 24.106906 28.174300 31.880357
##
## (Intercept) 9.993166 9.993166 9.993166 9.993166 9.993166
## poly(X, 7)1 884.876348 887.953182 890.756678 893.311120 895.638632
## poly(X, 7)2 . . . . .
## poly(X, 7)3 1153.491533 1156.568367 1159.371863 1161.926305 1164.253817
## poly(X, 7)4 -18.015261 -21.092095 -23.895591 -26.450033 -28.777545
## poly(X, 7)5 424.208959 427.285793 430.089289 432.643731 434.971243
## poly(X, 7)6 . . . . .
## poly(X, 7)7 35.257178 38.334012 41.137509 43.691950 46.019462
##
## (Intercept) 9.993166 9.993166
## poly(X, 7)1 897.759374 899.691715
## poly(X, 7)2 . .
## poly(X, 7)3 1166.374559 1168.306900
## poly(X, 7)4 -30.898287 -32.830628
## poly(X, 7)5 437.091985 439.024326
## poly(X, 7)6 . .
## poly(X, 7)7 48.140205 50.072546