Perform Logistic Regression with R

R for My Class

# Load some packages: 
library(tidyverse)
library(broom)

# Load data: 
default <- as_tibble(ISLR::Default)

# Logistic Curve: 
theme_set(theme_bw())

default %>% 
  mutate(prob = case_when(default == "No" ~ 0, TRUE ~ 1)) %>% 
  ggplot(aes(balance, prob)) + 
  geom_point(aes(color = factor(prob)), show.legend = FALSE, alpha = 0.1) + 
  geom_smooth(method = "glm", method.args = list(family = "binomial")) + 
  labs(title = "An Example of Logistic Function", 
       x = "Balance", 
       y = "Probability of Default")

default %>% 
  mutate(prob = case_when(default == "No" ~ 0, TRUE ~ 1)) %>% 
  ggplot(aes(balance, prob, color = student)) + 
  geom_point(alpha = 0.1) + 
  geom_smooth(method = "glm", method.args = list(family = "binomial"), alpha = 0.15) + 
  labs(title = "An Example of Logistic Function", 
       subtitle = "The S-curve shows that students have a lower\nprobability of default than non-students given the balance", 
       x = "Balance", 
       y = "Probability of Default")

# Splot data: 

train <- default %>% 
  group_by(default) %>% 
  sample_frac(0.5, replace = FALSE)

test <- setdiff(default, train)

# Simple Logistic Regression and results: 
model1 <- glm(default ~ balance, family = "binomial", data = train)
model1 %>% tidy()

## # A tibble: 2 x 5
##   term         estimate std.error statistic  p.value
##   <chr>           <dbl>     <dbl>     <dbl>    <dbl>
## 1 (Intercept) -10.3      0.493        -21.0 1.64e-97
## 2 balance       0.00532  0.000304      17.5 1.40e-68

# ROC curve: 

library(pROC)
my_auc <- roc(test$default, predict(model1, test, type = "response"))
sen_spec_df <- data.frame(TPR = my_auc$sensitivities, FPR = 1 - my_auc$specificities)

sen_spec_df %>% 
  ggplot(aes(x = FPR, ymin = 0, ymax = TPR))+
  geom_polygon(aes(y = TPR), fill = "red", alpha = 0.3)+
  geom_path(aes(y = TPR), col = "firebrick", size = 1.2) +
  geom_abline(intercept = 0, slope = 1, color = "gray37", size = 1, linetype = "dashed") + 
  coord_equal() +
  labs(x = "FPR (1 - Specificity)", 
       y = "TPR (Sensitivity)", 
       title = "Model Performance for Logistic Classifier based on Test Data", 
       subtitle = paste0("AUC Value: ", my_auc$auc %>% round(2)))