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.0     ✔ stringr   1.5.1
## ✔ ggplot2   4.0.0     ✔ tibble    3.3.0
## ✔ lubridate 1.9.4     ✔ tidyr     1.3.1
## ✔ purrr     1.1.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
logistic <- glm(Outcome ~ Glucose + BMI + Age, data=data, family="binomial")

summary(logistic)
## 
## 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

R²

r_square <- 1 - (logistic$deviance/logistic$null.deviance)

r_square
## [1] 0.25626

The R² value is 0.256. This means 25.6% of the variation in the onset of diabetes (outcome) is explained by glucose, BMI, and age.

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

The intercept, -9.032377, is the log-odds of diabetes when every other predictor is zero.

We can calculate the probability from log odds using the sigmoid curve function: 1/(1+e^-x)

1/(1+e^-(-9.032377)) = 0.00011946388247

So this intercept represents approximately 0.012% probability of diabetes onset when every other predictor is zero. This is not realistic as BMI, Age, and Glucose can’t be zero.

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)?

A one-unit increase raises the odds of diabetes because all of the coefficients are positive. The coefficients represent the change in log odds of diabetes onset, so if it’s positive, then they increase the odds of diabetes onset. If they’re negative, they decrease the odds of diabetes onset.

All variables are significant because they have a p-value less than 0.05.

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.data <- data.frame(
  probability.of.outcome=logistic$fitted.values,
  outcome=data_subset$Outcome)

predicted.data <- predicted.data[
  order(predicted.data$probability.of.outcome, decreasing=FALSE),]
predicted.data$rank <- 1:nrow(predicted.data)

head(predicted.data)
##     probability.of.outcome outcome rank
## 681             0.01422697       0    1
## 63              0.01490161       0    2
## 618             0.01550448       0    3
## 98              0.01718930       0    4
## 91              0.02040067       0    5
## 590             0.02132939       0    6
tail(predicted.data)
##     probability.of.outcome outcome rank
## 580              0.9460555       1  747
## 207              0.9496269       1  748
## 488              0.9505221       0  749
## 547              0.9547556       1  750
## 155              0.9602228       1  751
## 446              0.9681720       1  752
# Predicted classes

predicted.probs <- logistic$fitted.values

predicted.classes <- ifelse(predicted.probs > 0.5, 1, 0)



# Confusion matrix

confusion <- table(
  Predicted = factor(predicted.classes, levels = c(0, 1)),
  Actual = factor(data_subset$Outcome_num, levels = c(0, 1))
)

confusion
##          Actual
## Predicted   0   1
##         0 429 114
##         1  59 150
#Extract Values:
TN <- 429
FP <- 59
FN <- 114
TP <- 150

#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’s accuracy is 77%. For a medical diagnosis, this isn’t really a good model because diabetes is a serious condition and a common threshold is usually way higher for medical fields, like an accuracy greater than 90%.

This model is better at detecting non-diabetes (specificity) than diabetes (sensitivity). This might matter for medical diagnosis because medical researchers and doctors don’t want to inaccurately identify a patient with diabetes as having no diabetes. Even though the the model can detect non-diabetic patients well, a model still needs to be able to detect diabetes at least the same or better because a huge purpose of the model is to detect diabetic patients.

Question 3: ROC Curve, AUC, and Interpretation

#Enter your code here

library(pROC)
## Warning: package 'pROC' was built under R version 4.5.3
## Type 'citation("pROC")' for a citation.
## 
## Attaching package: 'pROC'
## The following objects are masked from 'package:stats':
## 
##     cov, smooth, var
# ROC curve & AUC on full data
roc_obj <- roc(response = data_subset$Outcome_num,
               predictor = logistic$fitted.values,
               levels = c(0, 1),
               direction = "<")  # smaller prob = non-diabetic

# Print AUC value
auc_val <- auc(roc_obj); auc_val
## Area under the curve: 0.828
# Plot ROC with AUC displayed
plot.roc(roc_obj, print.auc = TRUE, legacy.axes = TRUE,
         xlab = "False Positive Rate (1 - Specificity)",
         ylab = "True Positive Rate (Sensitivity)")

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

The area under the curve is 82.8%. This means the model is better at predicting than just chance (50%), however, it’s not perfect.

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

For diabetes diagnosis, medical researchers and professionals should prioritize sensitivity because the purpose of this model is to detect the sensitivity (diabetes onset). For such serious cases like diabetes diagnosis, the threshold should be very high. The sensitivity rate should at least be around 95% because if the model incorrectly detects someone with diabetes when they don’t actually have it, it could ruin their life.