The Goal

In our last class, we started to explore classification. This just means models where \(Y\) is categorical rather than numeric. We saw that both KNN and logistic regression can be used if \(Y\) is binary. Today, we are going to implement both of these approaches in R.

The Data

Start by loading the data. There are two data sets today, train and test. Please note that even though this is earthquake data like we used in class, you have a different sample today!! This means that you will likely get different results than we did in the slides.

train <- read.csv("https://www.dropbox.com/scl/fi/rbwypc2pi3fi9g9fuqyhv/train_earthquake.csv?rlkey=rio3x88gq040n0h8zk8s4hom9&st=h33zi0m2&dl=1")
test <- read.csv("https://www.dropbox.com/scl/fi/rsnyt5r94ql4jeqkc0acq/test_earthquake.csv?rlkey=hdy0z7ph9n4l38os2b1kbjc0g&st=kjj7vto0&dl=1")

The training data set train should have 2426 rows and 28 columns. The test data set test should have 1079 rows and 28 columns. The response variable is Damage, which is recorded as 0 or 1. A 0 indicates no to minimal damage from the earthquake and a 1 indicates moderate to severe damage. Our goal is to predict \(Y\) = damage.

KNN

Now that we have checked to see if prediction makes sense, it’s time to choose our method. We are going to try two different predictive approaches today, but we are going to start with KNN.

When we do KNN with a categorical \(Y\), we will use the caret package, as we did with KNN with a numeric \(Y\).

library(caret)

The code is also similar:

yhat <- knn3Train( traindata , testdata , Y, k = )

where:

• traindata: replace with a data set containing the features in your training data.
• testdata: replace with a data set containing the features in your test data.
• Y: replace with a dataset holding the column containing the features in your training data.
• k =: fill in a number after the = to tell R how many nearest neighbors to use.

holder <- table( yhat, test$Damage )
colnames(holder) <- c("True 0", "True 1")
rownames(holder) <- c("Predicted 0", "Predicted 1")

knitr::kable( holder, caption = "K = 5 Confusion Matrix")

So far, we have (1) used KNN to make predictions and (2) assessed the accuracy of those predictions using a few different metrics. However, we have yet to tune our choice of \(K\). Tuning is really important to help approaches like KNN predict as well as possible.

# Find the F1 score
TP <- sum(yhat == 1 & test$Damage == 1)
FP <- sum(yhat == 1 & test$Damage == 0 )
FN <- sum(yhat == 0 & test$Damage == 1 )
F1 <- TP/(TP + .5*(FP+FN))
F1

# Find the GM
TPR <- sum(yhat == 1 & test$Damage == 1) / sum(test$Damage==1)
TNR <- sum(yhat == 0 & test$Damage == 0) / sum(test$Damage==0)
GM <- sqrt( TPR*TNR)
GM 

library(ggplot2)

# If you chose F1
ggplot(Storage, aes( K, F1) ) + 
      geom_line() + 
      geom_vline(xintercept = c(3,5,7,11,13,15)[which.max(Storage$F1)]) + 
      labs(title = "Figure 1", xlab = "K", y = "F1 Score", 
      caption = paste("Maximum F1 Score achieved at K=",seq(from = 3, to = 15, by = 2)[which.max(Storage$F1)]))

# If you chose GM
ggplot(Storage, aes( K, GM) ) + 
      geom_line() + 
      geom_vline(xintercept = c(3,5,7,11,13,15)[which.max(Storage$GM)]) + 
      + labs(title = "Figure 1", xlab = "K", y = "GM", 
      caption = paste("Maximum Geometric Mean of TPR and TNR achieved at K=",seq(from = 3, to = 15, by = 2)[which.max(Storage$GM)]))

holder <- table( yhat, test$Damage )
colnames(holder) <- c("True 0", "True 1")
rownames(holder) <- c("Predicted 0", "Predicted 1")

knitr::kable( holder, caption = "K = ?? Confusion Matrix")

Logistic Regression

Unlike KNN, logistic regression is a model. This means it has the added benefit of being able to used for association as well as prediction. It is also just a different approach to prediction, which means that it can perform better (or worse!) than KNN, depending on the structure of the data.

To train a logistic regression model in R, we use:

logistic <- glm( Damage ~ ., data = train, family = "binomial")

Once the model is trained, you can obtain the probability that \(y_{i^{*}}=1\) for each row in the test data using:

probabilities <- predict( logistic, newdata = test, type = "response")

The trick to making predictions with logistic is that we have to choose a threshold. Probabilities above the threshold are assigned a prediction of 1, while probabilities equal to or below the threshold are assigned a prediction of 0. Like \(K\) in KNN, the threshold is a tuning parameter…so we have to choose it!!

Tuning the threshold requires very similar coding to what we just used to tune \(K\). For these data, I checked for you and using a threshold of .5, which is standard, works very well here - so happy day, you don’t have to code it!

holderL <- table( yhat_logistic, test$Damage )
colnames(holderL) <- c("True 0", "True 1")
rownames(holderL) <- c("Predicted 0", "Predicted 1")

knitr::kable( holderL, caption = "Logistic Regression Confusion Matrix")

Comparing Results