This is an R Markdown Notebook. When you execute code within the notebook, the results appear beneath the code.
Try executing this chunk by clicking the Run button within the chunk or by placing your cursor inside it and pressing Ctrl+Shift+Enter.
library(ISLR2)
## Warning: package 'ISLR2' was built under R version 4.2.3
library(caret)
## Loading required package: ggplot2
## Loading required package: lattice
library(pROC)
## Type 'citation("pROC")' for a citation.
##
## Attaching package: 'pROC'
## The following objects are masked from 'package:stats':
##
## cov, smooth, var
library(MASS)
##
## Attaching package: 'MASS'
## The following object is masked from 'package:ISLR2':
##
## Boston
# Load Default dataset
data(Default)
# Split data into training and testing sets
set.seed(123)
train_index <- createDataPartition(Default$default, p = 0.7, list = FALSE)
train_data <- Default[train_index, ]
test_data <- Default[-train_index, ]
# Fit logistic regression model
glm_model <- glm(default ~ balance + income, data = train_data, family = "binomial")
glm_predictions <- predict(glm_model, newdata = test_data, type = "response")
glm_table <- table(test_data$default, glm_predictions > 0.5)
glm_accuracy <- sum(diag(glm_table))/sum(glm_table)
glm_precision <- glm_table[2,2]/sum(glm_table[,2])
glm_recall <- glm_table[2,2]/sum(glm_table[2,])
glm_f1_score <- 2*glm_precision*glm_recall/(glm_precision + glm_recall)
glm_auc <- roc(test_data$default, glm_predictions)$auc
## Setting levels: control = No, case = Yes
## Setting direction: controls < cases
glm_predicted_classes <- ifelse(glm_predictions > 0.5, "Yes", "No")
glm_roc <- roc(test_data$default, glm_predictions)
## Setting levels: control = No, case = Yes
## Setting direction: controls < cases
# Fit LDA model
lda_model <- lda(default ~ balance + income, data = train_data)
lda_predictions <- predict(lda_model, newdata = test_data)
lda_table <- table(test_data$default, lda_predictions$class)
lda_accuracy <- sum(diag(lda_table))/sum(lda_table)
lda_precision <- lda_table[2,2]/sum(lda_table[,2])
lda_recall <- lda_table[2,2]/sum(lda_table[2,])
lda_f1_score <- 2*lda_precision*lda_recall/(lda_precision + lda_recall)
lda_auc <- roc(test_data$default, lda_predictions$posterior[,2])$auc
## Setting levels: control = No, case = Yes
## Setting direction: controls < cases
lda_predicted_classes <- as.character(lda_predictions$class)
lda_roc <- roc(test_data$default, lda_predictions$posterior[,2])
## Setting levels: control = No, case = Yes
## Setting direction: controls < cases
# Fit QDA model
qda_model <- qda(default ~ balance + income, data = train_data)
qda_predictions <- predict(qda_model, newdata = test_data)
qda_table <- table(test_data$default, qda_predictions$class)
qda_accuracy <- sum(diag(qda_table))/sum(qda_table)
qda_precision <- qda_table[2,2]/sum(qda_table[,2])
qda_recall <- qda_table[2,2]/sum(qda_table[2,])
qda_f1_score <- 2*qda_precision*qda_recall/(qda_precision + qda_recall)
qda_auc <- roc(test_data$default, qda_predictions$posterior[,2])$auc
## Setting levels: control = No, case = Yes
## Setting direction: controls < cases
qda_predicted_classes <- as.character(qda_predictions$class)
qda_roc <- roc(test_data$default, qda_predictions$posterior[,2])
## Setting levels: control = No, case = Yes
## Setting direction: controls < cases
# Print out-of-sample performance measures for logistic regression, LDA, and QDA
cat("Logistic regression performance on testing set:\n")
## Logistic regression performance on testing set:
cat(paste0("Accuracy: ", round(glm_accuracy, 3), "\n"))
## Accuracy: 0.975
cat(paste0("Precision: ", round(glm_precision, 3), "\n"))
## Precision: 0.717
cat(paste0("Recall: ", round(glm_recall, 3), "\n"))
## Recall: 0.384
cat(paste0("F1 score: ", round(glm_f1_score, 3), "\n"))
## F1 score: 0.5
cat(paste0("AUC-ROC: ", round(glm_auc, 3), "\n"))
## AUC-ROC: 0.957
cat("LDA performance on testing set:\n")
## LDA performance on testing set:
cat(paste0("Accuracy: ", round(lda_accuracy, 3), "\n"))
## Accuracy: 0.973
cat(paste0("Precision: ", round(lda_precision, 3), "\n"))
## Precision: 0.757
cat(paste0("Recall: ", round(lda_recall, 3), "\n"))
## Recall: 0.283
cat(paste0("F1 score: ", round(lda_f1_score, 3), "\n"))
## F1 score: 0.412
cat(paste0("AUC-ROC: ", round(lda_auc, 3), "\n"))
## AUC-ROC: 0.957
cat("QDA performance on testing set:\n")
## QDA performance on testing set:
cat(paste0("Accuracy: ", round(qda_accuracy, 3), "\n"))
## Accuracy: 0.973
cat(paste0("Precision: ", round(qda_precision, 3), "\n"))
## Precision: 0.721
cat(paste0("Recall: ", round(qda_recall, 3), "\n"))
## Recall: 0.313
cat(paste0("F1 score: ", round(qda_f1_score, 3), "\n"))
## F1 score: 0.437
cat(paste0("AUC-ROC: ", round(qda_auc, 3), "\n"))
## AUC-ROC: 0.957
# Print confusion matrix and ROC curve for logistic regression
cat("Confusion matrix for logistic regression:\n")
## Confusion matrix for logistic regression:
glm_predictions <- predict(glm_model, newdata = test_data, type = "response")
glm_predicted_classes <- ifelse(glm_predictions > 0.5, "Yes", "No")
glm_conf_mat <- table(test_data$default, glm_predicted_classes)
print(glm_conf_mat)
## glm_predicted_classes
## No Yes
## No 2885 15
## Yes 61 38
cat("ROC curve for logistic regression:\n")
## ROC curve for logistic regression:
plot(glm_roc, main = "ROC Curve - Logistic Regression")
# Print confusion matrix and ROC curve for LDA
cat("Confusion matrix for LDA:\n")
## Confusion matrix for LDA:
lda_predictions <- predict(lda_model, newdata = test_data)
lda_predicted_classes <- as.character(lda_predictions$class)
lda_conf_mat <- table(test_data$default, lda_predicted_classes)
print(lda_conf_mat)
## lda_predicted_classes
## No Yes
## No 2891 9
## Yes 71 28
cat("ROC curve for LDA:\n")
## ROC curve for LDA:
plot(lda_roc, main = "ROC Curve - LDA")
# Print confusion matrix and ROC curve for QDA
cat("Confusion matrix for QDA:\n")
## Confusion matrix for QDA:
qda_predictions <- predict(qda_model, newdata = test_data)
qda_predicted_classes <- as.character(qda_predictions$class)
qda_conf_mat <- table(test_data$default, qda_predicted_classes)
print(qda_conf_mat)
## qda_predicted_classes
## No Yes
## No 2888 12
## Yes 68 31
cat("ROC curve for QDA:\n")
## ROC curve for QDA:
plot(qda_roc, main = "ROC Curve - QDA")