Goal

Today we will be working to build a spam filter using logistic regression. Basically, a filter scans an email and looks for a few key text features. For instance, if the email says “cruise”, are you more or less likely to assume it is spam? Filters take these ideas and convert them into statistical models that use these features to predict the probability that an email is spam. If this probability is high, the email is labeled spam and tucked into the spam folder. If the probability is low, the email appears in your inbox.

Let’s load the libraries we will need for today.

library(dplyr)
library(tidytext)

The Data

To load the data sets, put the following into a chunk in Markdown and press play.

# Load the test data set
test <- read.csv("https://www.dropbox.com/scl/fi/2od0arr2ebclzzkx503pr/newtest.csv?rlkey=gos6y71nc3jdnaqf3nua77y8r&st=sfu5wit6&dl=1")[,1:2]

# Load the training data set 
train <- read.csv("https://www.dropbox.com/scl/fi/ilhmlo9ic8k89bq2on2hp/newtrain.csv?rlkey=els7a72hg490h8t4buvx3nsju&st=3zd96fw6&dl=1")[,1:2]

# Convert to a categorical variable
train$label <- as.factor(train$label)
test$label <- as.factor(test$label)

There are 2 variables in each data set. We have the response variable

• label: Indicator for whether the email is spam or not (“ham”)

and the text:

• text: The text of the email

Our response variable is binary, so before we try to model the data, our first step is always to visualize or summarize our data in some way. With binary data, a good first step is determine how many 0s (ham) and 1s (spam) we actually have in this data set. Are they all spam? Do we have a good mix? The table command is useful for this:

table(train$label)

It is important to check this before we begin an analysis to make sure we have enough emails that are spam to actually find relationships. If there are very few rows with \(y_i=1\) (in this case spam), modeling becomes more difficult.

Feature Engineering

Now that we have \(Y\), we need a feature. There are a lot of different features we could use, but for now, let’s stick with features created by words. In other words, we are going to choose a word that we think might be helpful in distinguishing spam from not spam.

train |>
  filter(label=="spam") |>
  unnest_tokens(word , text) |>
  anti_join(stop_words) |>
  count(word) |>
  slice_max(n=1, n)

train |>
  unnest_tokens(word , text) |>
  anti_join(stop_words) |>
  group_by(label) |>
  count(word) |>
  bind_tf_idf(word, label,n) |>
  slice_max(tf_idf,n=1)

Now that we know which word to use, let’s create a feature that indicates whether or not this word is present in each email. As a reminder, the code we need to do this is below, but change word inside \\bword\\b to be the word you chose in Question 5.

train$topword <- grepl("\\bword\\b", tolower(train$text))
test$topword <- grepl("\\bword\\b", tolower(test$text))

Model 1: Using the Top Word

Now that we have a feature, let’s see how effective a filter we can build using \(X\) = the word we chose in Question 5.

Using a categorical variable as a feature requires converting to an indicator variable. This is a variable that takes on 0 or 1 depending on the value of X. The level (or value) of X that is coded as 0 is called the baseline of X.

When we make a table, the baseline is the value of X that shows up first.

Luckily, R creates indicator variables for us when we put a categorical \(X\) in a model. To see this, let’s fit a logistic regression model using the training data and \(X\) = the indicator variable we chose:

Model1 <- glm( label ~ topword, data = train, family = "binomial")

Running the code above doesn’t actually print anything on the screen. Instead, it stores everything we need, ready to print it out when we ask it to.

To see the coefficients (\(\hat{\beta}\) terms), for instance, we use

summary(Model1)$coefficients

Running this shows us our estimates of \(\hat{\beta}_0\) and \(\hat{\beta}_1\), along with standard errors and other things we could use for inference. It is helpful output, but it is what we call raw R-output. For any kind of formal reports, we want formatted output.

To format the output and make it more professional, we can use

knitr::kable( summary(Model1)$coefficients )

The knitr::kable part of the code is a wrapper than formats any table output (including data sets!) nicely for knitting.

Once we have our coefficients, we can write down the trained (fitted) model.

One of the cool things about logistic regression is that when there is an increase in the log odds of \(y_i = 1\), there is also an increase in the probability that \(y_i = 1\). The same holds for decreases. The increase/decrease is not the same, but the direction is the same.

The answer to Question 10 is what we call an informal interpretation. It gives a general idea of the relationships shown in the model, but it does not actually use any numbers. To incorporate more information from the trained model, we move to a more formal interpretation.

When we interpret a slope, we are used to saying things like “As X goes up by 1…”. However, interpretation is slightly different with categorical features. We are no longer able to start with “If X increases by 1…” because that doesn’t make any sense when an email either has the top word or not. Instead, we start with “If the email contains the word…”

Another way we sometimes discuss categorical features is in terms of the odds ratio. The odds ratio is

\[\frac{Odds(y_i=1|x_i = 1)}{Odds(y_i=1|x_i = 0)}\]

Suppose the odds ratio is 4. We then interpret the result as: “The odds of being spam are predicted to be 4 times higher if X=1”. We replace “X=1” in our interpretation with what that means for our data.

Prediction

Once we have trained our model, we can check to see how well that model performs at prediction. We used our training data train to build our model, and now we can see how well that model can predict for spam status on test, our test data set compromised of 1170 new emails that were not in the training data.

Making predictions for logistic regression requires two steps. First, we need to find \(\hat{\pi}_i\) for all rows in the test data. This means we need to find the probability of being spam for all emails in the test data.

Doing this by hand for all the emails would be tedious. Luckily, we can do this quickly in R using:

pihat <- predict( Model1 , type = "resp", newdata = test)

Here,

• Model1 is the model we want to use to make predictions.
• type = "resp" tells R to predict the probabilities; without this statement, R predicts the log odds.
• newdata = test tells R which data set we want to make predictions on. If you do not include this, R makes predictions using the training data.

Once we have predicted probabilities, we use a threshold to decide whether or not a probability is high enough for us to declare \(\hat{y}_i = 1\). A common choice for threshold is .5, so let’s start with that.

yhat <- as.factor(ifelse( pihat > .5, "spam" , "ham"))

The code ifelse is a handy one. If \(\hat{\pi}_i > .5\), then we set \(\hat{y}_i = 1=spam\). If not, we set \(\hat{y}_i = 0=ham\).

The next step is to see how accurate our predictions are! We usually use a confusion matrix for this.

knitr::kable( table(yhat, test$label) , col.names = c("Truth: Ham", "Truth: Spam"))

Model 2: Numbers

This feature didn’t work particularly well, so let’s try something else. The spam emails have a lot of numbers in them, and this does not seem to be true for the ham emails. What if we create a feature for whether or not an email has numbers?

train$hasNumbers <- ifelse(grepl("[0-9]", train$text), "numbers", "no numbers")
test$hasNumbers <- ifelse(grepl("[0-9]", test$text), "numbers", "no numbers")

library(caret)
library(forcats)

conf_matrix <- confusionMatrix(data = yhat, reference = test$label)

conf_matrix$table

as_tibble(conf_matrix$table) |>
  mutate(freq = n/sum(n)) |>
  ungroup() |>
  ggplot(aes(y = Prediction, x = Reference, fill = freq)) +
  geom_tile(color = "white") +
  geom_text(aes(label = paste0(round(freq * 100, 2),"%\n(",n,")")),
            color = "white", size = 3) +
  scale_fill_gradient(low  = "#9ecae1",high = "#08306b",
                      labels = scales::percent, limits = c(0, 1)) +
  labs(fill = "Percent\n", x = "True Spam Status",
       y = "Predicted Spam Status") +
  coord_equal() +
  theme_linedraw() +
  theme(panel.grid.major = element_blank(),
        axis.ticks = element_line(color = "black"),
        axis.text.x = element_text(angle = 45, hjust = 1),
        axis.text = element_text(color = "black"),
        axis.title = element_text(color = "black", face = "bold"),
        legend.title = element_text(color = "black", face = "bold"),
        legend.text = element_text(color = "black"))

pihat_train <- predict( Model2 , type = "resp", newdata = train)
yhat_train <- as.factor(ifelse( pihat_train > .5, "spam" , "ham"))