Introduction:

In this homework, you will apply logistic regression to a real-world dataset: the Pima Indians Diabetes Database. This dataset contains medical records from 768 women of Pima Indian heritage, aged 21 or older, and is used to predict the onset of diabetes (binary outcome: 0 = no diabetes, 1 = diabetes) based on physiological measurements.

The data is publicly available from the UCI Machine Learning Repository and can be imported directly.

Dataset URL: https://raw.githubusercontent.com/jbrownlee/Datasets/master/pima-indians-diabetes.data.csv

Columns (no header in the CSV, so we need to assign them manually):

  1. Pregnancies: Number of times pregnant
  2. Glucose: Plasma glucose concentration (2-hour test)
  3. BloodPressure: Diastolic blood pressure (mm Hg)
  4. SkinThickness: Triceps skin fold thickness (mm)
  5. Insulin: 2-hour serum insulin (mu U/ml)
  6. BMI: Body mass index (weight in kg/(height in m)^2)
  7. DiabetesPedigreeFunction: Diabetes pedigree function (a function scoring genetic risk)
  8. Age: Age in years
  9. Outcome: Class variable (0 = no diabetes, 1 = diabetes)

Task Overview: You will load the data, build a logistic regression model to predict diabetes onset using a subset of predictors (Glucose, BMI, Age), interpret the model, evaluate it with a confusion matrix and metrics, and analyze the ROC curve and AUC.

Cleaning the dataset Don’t change the following code

library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.4     ✔ readr     2.1.5
## ✔ forcats   1.0.1     ✔ stringr   1.6.0
## ✔ ggplot2   4.0.0     ✔ tibble    3.3.0
## ✔ lubridate 1.9.4     ✔ tidyr     1.3.1
## ✔ purrr     1.2.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
url <- "https://raw.githubusercontent.com/jbrownlee/Datasets/master/pima-indians-diabetes.data.csv"

data <- read.csv(url, header = FALSE)

colnames(data) <- c("Pregnancies", "Glucose", "BloodPressure", "SkinThickness", "Insulin", "BMI", "DiabetesPedigreeFunction", "Age", "Outcome")

data$Outcome <- as.factor(data$Outcome)

# Handle missing values (replace 0s with NA because 0 makes no sense here)
data$Glucose[data$Glucose == 0] <- NA
data$BloodPressure[data$BloodPressure == 0] <- NA
data$BMI[data$BMI == 0] <- NA


colSums(is.na(data))
##              Pregnancies                  Glucose            BloodPressure 
##                        0                        5                       35 
##            SkinThickness                  Insulin                      BMI 
##                        0                        0                       11 
## DiabetesPedigreeFunction                      Age                  Outcome 
##                        0                        0                        0

Question 1: Create and Interpret a Logistic Regression Model - Fit a logistic regression model to predict Outcome using Glucose, BMI, and Age.

## Enter your code here

# Fit logistic regression model
model <- glm(Outcome ~ Glucose + BMI + Age, 
             data = data, 
             family = binomial)

# Display model summary
summary(model)
## 
## Call:
## glm(formula = Outcome ~ Glucose + BMI + Age, family = binomial, 
##     data = data)
## 
## Coefficients:
##              Estimate Std. Error z value Pr(>|z|)    
## (Intercept) -9.032377   0.711037 -12.703  < 2e-16 ***
## Glucose      0.035548   0.003481  10.212  < 2e-16 ***
## BMI          0.089753   0.014377   6.243  4.3e-10 ***
## Age          0.028699   0.007809   3.675 0.000238 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 974.75  on 751  degrees of freedom
## Residual deviance: 724.96  on 748  degrees of freedom
##   (16 observations deleted due to missingness)
## AIC: 732.96
## 
## Number of Fisher Scoring iterations: 4
# Calculate pseudo R-squared (McFadden's R²)
pseudo_r2 <- 1 - (model$deviance / model$null.deviance)
cat("\nPseudo R² (McFadden):", round(pseudo_r2, 4))
## 
## Pseudo R² (McFadden): 0.2563
# Extract coefficients for interpretation
coefficients <- summary(model)$coefficients
intercept <- coefficients["(Intercept)", "Estimate"]
glucose_coef <- coefficients["Glucose", "Estimate"]
bmi_coef <- coefficients["BMI", "Estimate"]
age_coef <- coefficients["Age", "Estimate"]

glucose_pval <- coefficients["Glucose", "Pr(>|z|)"]
bmi_pval <- coefficients["BMI", "Pr(>|z|)"]
age_pval <- coefficients["Age", "Pr(>|z|)"]

# Calculate odds ratios
glucose_odds <- exp(glucose_coef)
bmi_odds <- exp(bmi_coef)
age_odds <- exp(age_coef)

cat("\n\nCoefficients and Odds Ratios:\n")
## 
## 
## Coefficients and Odds Ratios:
cat("Glucose - Coefficient:", round(glucose_coef, 4), ", Odds Ratio:", round(glucose_odds, 4), ", p-value:", format(glucose_pval, scientific = TRUE), "\n")
## Glucose - Coefficient: 0.0355 , Odds Ratio: 1.0362 , p-value: 1.748579e-24
cat("BMI - Coefficient:", round(bmi_coef, 4), ", Odds Ratio:", round(bmi_odds, 4), ", p-value:", format(bmi_pval, scientific = TRUE), "\n")
## BMI - Coefficient: 0.0898 , Odds Ratio: 1.0939 , p-value: 4.297514e-10
cat("Age - Coefficient:", round(age_coef, 4), ", Odds Ratio:", round(age_odds, 4), ", p-value:", format(age_pval, scientific = TRUE), "\n")
## Age - Coefficient: 0.0287 , Odds Ratio: 1.0291 , p-value: 2.375089e-04

What does the intercept represent (log-odds of diabetes when predictors are zero)?

The intercept represents the log-odds of having diabetes when all predictors (Glucose, BMI, and Age) equal zero. While this has no practical interpretation since these values cannot be zero, it is mathematically necessary for the model.

For each predictor (Glucose, BMI, Age), does a one-unit increase raise or lower the odds of diabetes? Are they significant (p-value < 0.05)?

All three predictors increase the odds of diabetes:

The pseudo R² indicates the model explains a moderate proportion of deviance in the outcome variable.

Question 2: Confusion Matrix and Important Metric

Calculate and report the metrics:

Accuracy: (TP + TN) / Total Sensitivity (Recall): TP / (TP + FN) Specificity: TN / (TN + FP) Precision: TP / (TP + FP)

Use the following starter code

# Keep only rows with no missing values in Glucose, BMI, or Age
data_subset <- data[complete.cases(data[, c("Glucose", "BMI", "Age")]), ]

#Create a numeric version of the outcome (0 = no diabetes, 1 = diabetes).This is required for calculating confusion matrices.
data_subset$Outcome_num <- ifelse(data_subset$Outcome == "1", 1, 0)


# Predicted probabilities
predicted_probs <- predict(model, newdata = data_subset, type = "response")


# Predicted classes
predicted_classes <- ifelse(predicted_probs > 0.5, 1, 0)


# Confusion matrix
confusion_matrix <- table(Predicted = predicted_classes, Actual = data_subset$Outcome_num)
print(confusion_matrix)
##          Actual
## Predicted   0   1
##         0 429 114
##         1  59 150
#Extract Values:
TN <- confusion_matrix[1, 1]  # True Negatives
FP <- confusion_matrix[2, 1]  # False Positives
FN <- confusion_matrix[1, 2]  # False Negatives
TP <- confusion_matrix[2, 2]  # True Positives

#Metrics    
accuracy <- (TP + TN) / (TP + TN + FP + FN)
sensitivity <- TP / (TP + FN)
specificity <- TN / (TN + FP)
precision <- TP / (TP + FP)

cat("Accuracy:", round(accuracy, 3), "\nSensitivity:", round(sensitivity, 3), "\nSpecificity:", round(specificity, 3), "\nPrecision:", round(precision, 3))
## Accuracy: 0.77 
## Sensitivity: 0.568 
## Specificity: 0.879 
## Precision: 0.718

Interpret: How well does the model perform? Is it better at detecting diabetes (sensitivity) or non-diabetes (specificity)? Why might this matter for medical diagnosis?

The model is better at detecting non-diabetes (specificity) than diabetes (sensitivity). Sensitivity measures the model’s ability to correctly identify diabetes cases, while specificity measures its ability to correctly identify non-diabetes cases. For medical diagnosis, missing diabetes cases (low sensitivity) can lead to delayed treatment and serious complications, making high sensitivity particularly important for screening purposes.

Question 3: ROC Curve, AUC, and Interpretation

#Enter your code here

# Install and load pROC package if needed
if(!require(pROC)) install.packages("pROC", repos = "http://cran.us.r-project.org")
## Loading required package: pROC
## Type 'citation("pROC")' for a citation.
## 
## Attaching package: 'pROC'
## The following objects are masked from 'package:stats':
## 
##     cov, smooth, var
library(pROC)

# Create ROC curve object
roc_obj <- roc(data_subset$Outcome_num, predicted_probs)
## Setting levels: control = 0, case = 1
## Setting direction: controls < cases
# Plot ROC curve
plot(roc_obj, 
     main = "ROC Curve for Diabetes Prediction Model",
     col = "blue", 
     lwd = 2,
     print.auc = TRUE,
     auc.polygon = TRUE,
     auc.polygon.col = "lightblue",
     grid = TRUE)

# Add diagonal reference line (random classifier)
abline(a = 0, b = 1, lty = 2, col = "red")

# Calculate and display AUC
auc_value <- auc(roc_obj)
cat("\nArea Under the Curve (AUC):", round(auc_value, 4))
## 
## Area Under the Curve (AUC): 0.828

What does AUC indicate (0.5 = random, 1.0 = perfect)?

AUC measures the model’s ability to discriminate between diabetic and non-diabetic individuals. An AUC of 0.5 indicates random guessing (no discriminative ability), while 1.0 indicates perfect classification. The calculated AUC indicates good discriminative ability.

For diabetes diagnosis, prioritize sensitivity (catching cases) or specificity (avoiding false positives)? Suggest a threshold and explain.

For diabetes diagnosis, sensitivity should be prioritized to catch more actual diabetes cases. Missing a diabetes diagnosis can lead to serious health complications, while false positives can be ruled out with follow-up testing. Lowering the threshold below 0.5 (such as 0.35-0.40) would increase sensitivity at the cost of reduced specificity, which is an acceptable trade-off for medical screening where the consequences of missing a case are more severe than the cost of additional confirmatory testing.