Yet Another Reason to NOT Use Numerics as Factors

Introduction

When dealing with linear regression models, it is tempting to convert numeric fields to factors. However, caution must be taken when doing this. One of the problems that can arrise when attempting this type of conversion is that it can cloud the detection of outliers in the dataset. This is not surprising. When we think of the numeric values of (1,2,3,36), it is not hard to detect the outliner. However, if we have a dataset like this: (“red”,“blue”,“green”,“yellow”) outlier detection becomes much more difficult if not impossible. The following example demonstrates how this can be manifested in a simple dataset in R.

The following cell constructs a dataset based on the builtin mtcars dataset in R.

rm(list=ls())
library(dplyr)
library(tidyverse)
mynewcars <- mtcars
mynewcars$ID <- seq.int(nrow(mynewcars))
mynewcars$model <- row.names(mynewcars)
row.names(mynewcars) <- mynewcars$ID
mynewcars <- mynewcars %>% select(model, cyl, wt, mpg)
mynewcars <- as_tibble(mynewcars)

The following code introduces an outliers to this dataset. And, what is more of an outlier The Homer?!

mynewcars <- rbind(mynewcars, c(model="The Homer",cyl=36, wt=7.8, mpg=9.3))

mynewcars$cyl <- as.numeric(mynewcars$cyl)
mynewcars$wt <- as.numeric(mynewcars$wt)
mynewcars$mpg <- as.numeric(mynewcars$mpg)
mynewcars

Next we will construct a model with the predictors in our dataset as numeric (dbl) values. We are attempting to predict the fuel efficiency of the cars based on their number of cylinders and weight.

mdl <- lm(mpg~cyl+wt, data=mynewcars)
summary(mdl)

Call:
lm(formula = mpg ~ cyl + wt, data = mynewcars)

Residuals:
    Min      1Q  Median      3Q     Max 
-5.2918 -2.3340 -0.3548  1.7413  7.5235 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  35.2816     1.8984  18.585  < 2e-16 ***
cyl           0.3590     0.1804   1.991   0.0557 .  
wt           -5.3823     0.7902  -6.811 1.48e-07 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 3.373 on 30 degrees of freedom
Multiple R-squared:  0.7246,    Adjusted R-squared:  0.7062 
F-statistic: 39.47 on 2 and 30 DF,  p-value: 3.977e-09

Based on the above we have a resonably good model with an adjusted R-squared value of 0.7062339.

If we look at the outliers for this data, what we get is no surprise.

cd1 <- cooks.distance(mdl)
outliers <- which(cd1>1)
plot(cd1, main="Numeric Cyc Column")

Based on the “rule of thumb” for Cooks distance, Homer’s car appears to be an outlier.

Now, let’s do the same thing, but convert the number of cylinders from a numeric to a factor…

mynewcars_f <- mynewcars
mynewcars_f$cyl <- as.factor(mynewcars_f$cyl)
mdl2 <- lm(mpg~cyl+wt, data=mynewcars_f)
summary(mdl2)

Call:
lm(formula = mpg ~ cyl + wt, data = mynewcars_f)

Residuals:
    Min      1Q  Median      3Q     Max 
-4.5890 -1.1557 -0.4879  1.3370  5.7915 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  33.9908     1.8878  18.006  < 2e-16 ***
cyl6         -4.2556     1.3861  -3.070 0.004718 ** 
cyl8         -6.0709     1.6523  -3.674 0.000999 ***
cyl36         0.3130     4.9411   0.063 0.949942    
wt           -3.2056     0.7539  -4.252 0.000213 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.557 on 28 degrees of freedom
Multiple R-squared:  0.8522,    Adjusted R-squared:  0.8311 
F-statistic: 40.38 on 4 and 28 DF,  p-value: 3.056e-11

cd2 <- cooks.distance(mdl2)
plot(cd2, main="Factor Cyl Column")

Now, all of a sudden, our adjusted R-Squared value has shot way up at o over 0.8311402. However, we should be cautious. If we look at the plot of the cooks distance for this model, no observations have a cooks distance greater than one. We appear to have gotten rid of the outliers in our model, but really we just maked them. If we look very carefully, we will note that the cooks distance was unable to be computed for The Homer.

Why is that?

The cooks distance formula is:

\[ D_i = \frac{\sum_{j=1}^n(\hat Y_j-\hat Y_{j(i)})^2}{(p+1)\hat\sigma^2} \]

Not to get into too much of the gorey details, but note that p on the bottom. Cooks distance summarizes how much the all the values in the model change with the ith observation is removed. In the above formula, p represents the number of parameters in our model. If we remove the The Homer from our dataset, the number of parameters in our dataset changes. We no longer have a predictor for the dummy variable for a car with 36 cylinders. If you want to, you can add more cars with 36 cylinders to the data (maybe the Batmobile or Wonder Woman’s car.) If you do, you will notice that the Cook’s distances for the models with numeric values are going to be different from those with the factors.

Now, just to wrap things up, lets return to model 1 and remove the outlier…

mynewcars2 <- mynewcars %>% slice(-outliers)
mdl3 <- lm(mpg~cyl+wt, data=mynewcars2)
summary(mdl3)

Call:
lm(formula = mpg ~ cyl + wt, data = mynewcars2)

Residuals:
    Min      1Q  Median      3Q     Max 
-4.2893 -1.5512 -0.4684  1.5743  6.1004 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  39.6863     1.7150  23.141  < 2e-16 ***
cyl          -1.5078     0.4147  -3.636 0.001064 ** 
wt           -3.1910     0.7569  -4.216 0.000222 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.568 on 29 degrees of freedom
Multiple R-squared:  0.8302,    Adjusted R-squared:  0.8185 
F-statistic: 70.91 on 2 and 29 DF,  p-value: 6.809e-12

After removing the outlier in our data, notice that the Adjusted R2 value is now 0.8185189, which is line with the factor version of the model. So, converting to a factor didn’t really gain us anything in accuracy anyway, AND the factor hid the fact that one of our observations was an outlier.

Conclusion

This data is clearly made up. The point that I am trying to make, yet again, is that great care should be taken when converting numeric values in models into factors (categorical variables). As seen in the above example, it can have unexpected statistical ramifications including masking outliers in the data.