I'm looking to fit a model with about 1k-40k variables and up to a few million observations.

Can anyone with a bit more experience speak to its performance for larger datasets?

It looks like I can use R's sparse matrix support to represent the data in a much smaller form (>90% of variables will have zero values for any given observation). However, the few examples I've found that use glmnet with a sparse matrix appear to also have very few observations.

I'm left a bit unsure whether to prep my data for glmnet if it won't be able to converge to a solution in a reasonable amount of time, so any advice would be appreciated.


1 Answer 1


glmnet is very strong in this respect. Of course it depends on the situation but to give you an impression of its performance, I made an example with 3 mio lines and a sparse model matrix of dimension 40'001 (one factor with 40k levels and one normally distributed variable). It takes about half a minute to run on my normal laptop to find the regularization path of length 44.


n <- 3e6 # 3 mio

y <- sample(0:1, n, replace = TRUE) 
x1 <- rnorm(n)
x2 <- sample(1:40000, n, replace = TRUE)

X <- sparse.model.matrix(~ x1 + factor(x2))
dim(X)  # 3000000   40001

# Runs 25.32 seconds on normal laptop
system.time(fit <- glmnet(y = y, x = X, family = "binomial", alpha = 0.5))

Unfortunately, glmnet does not return standard errors etc. for its coefficients, even in the unregularized mode. (Using regularization invalidates quantities like this.)

An alternative is the h2o R package that connects to h2o.ai backend, a small JAVA program designed for high performance machine learning. Their elastic net GLM implementation is able to provide p-values, standard errors etc. for your coefficients as long as you use the slow "IRWLS" optimizer and no regularization. If you are not interested in p values etc., it is as fast as glmnet. You can fit data sets of any size as long as it fits into memory.

h2o.init(max_mem_size = "4G", nthreads = -1)

# 20 seconds
train <- as.h2o(data.frame(y = y, x1 = x1, x2 = factor(x2)))

# still at 2% done after 1h
system.time(fit_h2o <- h2o.glm(x = c("x1", "x2"), y = "y", 
                               training_frame = train, 
                               family = "binomial",
                               compute_p_values = TRUE,
                               lambda = 0,
                               alpha = 0,
                               solver = "IRLSM"))

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.