In supervised learning, we want to use input feature values to predict the value of a response variable. In the training dataset, we usually organize the input feature values in a matrix \(\mathbf{X} \in \mathbb{R}^{n \times p}\), where \(n\) is the number of observations and \(p\) is the number of input features. Each row of \(\mathbf{X}\) corresponds to an observation in the training dataset while each column corresponds to an input feature. We usually organize the response values in a vector \(\mathbf{y}\).
The core model-fitting functions in glmnet (the glmnet and cv.glmnet functions) take in \(\mathbf{y}\) via the function argument y. It is typically a vector, but the type of the vector elements depends on the type of model you are fitting. For example, the elements should be real numbers for linear models, and should be censored responses for survival models. (We go into more detail for what is needed for each type of model in Parts 2 and 3 of the book.)
We focus on the input matrix \(\mathbf{X}\) in the rest of this chapter. The glmnet and cv.glmnet functions take in \(\mathbf{X}\) via the function argument x, which must be a matrix of real numbers. However, this is not always the case for the raw input data:
The glmnet package includes a function makeX that makes it easy to address these two issues, creating an input matrix x which glmnet can accept. In the next sections, we show how to use makeX to format factor variables and to impute missing data.
Here is an example of a data frame with both numeric and factor variables:
library(knitr)
set.seed(101)
n_train <- 5
X1 <- matrix(rnorm(n_train), nrow = n_train)
X2 <- sample(letters[1:3], n_train, replace = TRUE)
X3 <- sample(LETTERS[1:3], n_train, replace = TRUE)
df <- data.frame(X1, X2, X3)
kable(head(df))| X1 | X2 | X3 |
|---|---|---|
| -0.3260365 | c | C |
| 0.5524619 | b | C |
| -0.6749438 | c | B |
| 0.2143595 | c | C |
| 0.3107692 | a | C |
makeX takes a dataframe as input, detects which columns are factor variables, and performs one-hot encoding to form indicator matrices these variables:
makeX(df)
#> X1 X2a X2b X2c X3B X3C
#> 1 -0.3260365 0 0 1 0 1
#> 2 0.5524619 0 1 0 0 1
#> 3 -0.6749438 0 0 1 1 0
#> 4 0.2143595 0 0 1 0 1
#> 5 0.3107692 1 0 0 0 1For each category level in a factor variable, makeX creates a column for it. In our example above, variable X2 and X3 had 3 and 2 levels respectively, and so the resulting x matrix has 3 columns corresponding to X2 and 2 columns corresponding to X3.
Finally, notice in the example above that makeX leaves the numeric variables remain unchanged, as we would expect.
Note that the behavior of creating one column for each category level is not a universal default. For example, stats::model.matrix drops one level from each factor variable when creating the one-hot encoding matrices:
model.matrix(~ X2 + X3, data = df)
#> (Intercept) X2b X2c X3C
#> 1 1 0 1 1
#> 2 1 1 0 1
#> 3 1 0 1 0
#> 4 1 0 1 1
#> 5 1 0 0 1
#> attr(,"assign")
#> [1] 0 1 1 2
#> attr(,"contrasts")
#> attr(,"contrasts")$X2
#> [1] "contr.treatment"
#>
#> attr(,"contrasts")$X3
#> [1] "contr.treatment"This is to avoid multicollinearity, which can cause issues in unpenalized linear regression. However, multicollinearity is not an issue for penalized regression, and so we prefer creating one column for each level to maintain symmetry.
In some applications, there are category levels which we know exist but are so rare that they do not show up in our training dataset. In our little example above, the level "A" for X3 did not occur, even though it is possible. We should still include these levels in the training input matrix so that we don’t run into errors when testing or deploying the model.
makeX can handle this situation if a test set is available. You can pass the training and test set to makeX at the same time so that all the levels of factors present in the training and test sets are accommodated in the one-hot encoding. Let’s see this in action in the example below:
# training data
kable(df)| X1 | X2 | X3 |
|---|---|---|
| -0.3260365 | c | C |
| 0.5524619 | b | C |
| -0.6749438 | c | B |
| 0.2143595 | c | C |
| 0.3107692 | a | C |
# test data
set.seed(104)
n_test <- 3
X1t <- matrix(rnorm(n_test), nrow = n_test)
X2t <- sample(letters[1:3], n_test, replace = TRUE)
X3t <- sample(LETTERS[1:3], n_test, replace = TRUE)
dft <- data.frame(X1t, X2t, X3t)
colnames(dft) <- c("X1", "X2", "X3")
kable(dft)| X1 | X2 | X3 |
|---|---|---|
| -0.3465840 | a | B |
| 0.6276359 | b | B |
| 0.6437831 | c | A |
For the X3 variable, the training data only has levels "B" and "C" while the test data only has levels "A" and "B". By passing both training and test datasets to makeX, we get training and test data matrices with all 3 levels for X3.
makeX(df, dft)
#> $x
#> X1 X2a X2b X2c X3A X3B X3C
#> 1 -0.3260365 0 0 1 0 0 1
#> 2 0.5524619 0 1 0 0 0 1
#> 3 -0.6749438 0 0 1 0 1 0
#> 4 0.2143595 0 0 1 0 0 1
#> 5 0.3107692 1 0 0 0 0 1
#>
#> $xtest
#> X1 X2a X2b X2c X3A X3B X3C
#> 6 -0.3465840 1 0 0 0 1 0
#> 7 0.6276359 0 1 0 0 1 0
#> 8 0.6437831 0 0 1 1 0 0The makeX function can also help to impute missing data in the input matrix. Let’s create an example with some missing values denoted by NAs:
X1n <- X1; X1n[3,1] <- NA
X2n <- X2; X2n[1] <- NA
X3n <- X3; X3n[5] <- NA
dfn <- data.frame(X1n, X2n, X3n)
dfn
#> X1n X2n X3n
#> 1 -0.3260365 <NA> C
#> 2 0.5524619 b C
#> 3 NA c B
#> 4 0.2143595 c C
#> 5 0.3107692 a <NA>By default, makeX leaves NAs as is:
makeX(dfn)
#> X1n X2na X2nb X2nc X3nB X3nC
#> 1 -0.3260365 NA NA NA 0 1
#> 2 0.5524619 0 1 0 0 1
#> 3 NA 0 0 1 1 0
#> 4 0.2143595 0 0 1 0 1
#> 5 0.3107692 1 0 0 NA NATo impute missing values with column means, pass na.impute = TRUE to makeX:
makeX(dfn, na.impute = TRUE)
#> X1n X2na X2nb X2nc X3nB X3nC
#> 1 -0.3260365 0.25 0.25 0.5 0.00 1.00
#> 2 0.5524619 0.00 1.00 0.0 0.00 1.00
#> 3 0.1878885 0.00 0.00 1.0 1.00 0.00
#> 4 0.2143595 0.00 0.00 1.0 0.00 1.00
#> 5 0.3107692 1.00 0.00 0.0 0.25 0.75
#> attr(,"means")
#> X1n X2na X2nb X2nc X3nB X3nC
#> 0.1878885 0.2500000 0.2500000 0.5000000 0.2500000 0.7500000For example, the value of X3nC for row 5 is \(0.75 = (1 + 1 + 0 + 1) / 4\).
glmnet and cv.glmnet can accept the input matrix x as a sparse matrix (more specifically, as long as it inherits from the sparseMatrix class from the Matrix package). Sparse matrices allow us to store and compute efficiently with large matrices having only a few nonzero entries.
data(SparseExample)
x <- SparseExample$x
y <- SparseExample$y
class(x)
#> [1] "dgCMatrix"
#> attr(,"package")
#> [1] "Matrix"We can fit the elastic net model the same way as before:
fit <- glmnet(x, y)cv.glmnet works in the same way as well:
set.seed(1)
cvfit <- cv.glmnet(x, y)
plot(cvfit)
Sparse matrices can also be used for newx, the new input matrix in the predict method.
The makeX function introduced earlier in this chapter can also return a sparse matrix, simply pass sparse = TRUE as an argument:
makeX(df, sparse = TRUE)
#> 5 x 6 sparse Matrix of class "dgCMatrix"
#> X1 X2a X2b X2c X3B X3C
#> 1 -0.3260365 . . 1 . 1
#> 2 0.5524619 . 1 . . 1
#> 3 -0.6749438 . . 1 1 .
#> 4 0.2143595 . . 1 . 1
#> 5 0.3107692 1 . . . 1glmnet is designed for regularized linear models, so it is most often used for multivariate models. However, it is possible to fit a univariate model model with glmnet. As its first argument, glmnet expects an input matrix with more than one column (feature) – in fact, calling glmnet with a single-column matrix throws an error.
univariate.x <- rnorm(100)
y <- univariate.x + rnorm(100, sd = .1)
glmnet(univariate.x, y)
#> Error in `if (is.null(np) | (np[2] <= 1)) ...`:
#> ! argument is of length zeroTo get around this, we recommend modifying your matrix to include a column that is entirely \(0\).
multivariate.x <- cbind(0, univariate.x)
fit <- cv.glmnet(multivariate.x, y)This all-zero column is harmless: (1) in a preprocessing step, glmnet ignores any columns with zero variance and (2) any coefficient multiplied by \(0\) will still be \(0\). This is reflected in the fitted coefficients:
coef(fit)
#> 3 x 1 sparse Matrix of class "dgCMatrix"
#> lambda.1se
#> (Intercept) -0.0002191641
#> .
#> univariate.x 0.9588627183Just remember to do the same at prediction time!
new.univariate.x <- rnorm(100)
new.multivariate.x <- cbind(0, new.univariate.x)
head(predict(fit, new.multivariate.x))
#> lambda.1se
#> [1,] -1.0199288
#> [2,] 0.2358632
#> [3,] -0.2778093
#> [4,] -2.1719371
#> [5,] -1.3511133
#> [6,] 0.8781176