Appendix D — Comparing glmnet and scikit-learn

scikit-learn (also known as sklearn) an open-source Python library for a vast array of machine learning tools. scikit-learn has its own implementation of the elastic net. This note demonstrates how to match the outputs of fitting the elastic net in scikit-learn and in adelie.

Both implementations use the same parametrization of the elastic net, but unfortunately give the hyperparameters different names. alpha in glmnet and adelie corresponds to l1_ratio in scikit-learn, and somewhat confusingly lambda in glmnet and adelie corresponds to alpha in scikit-learn.

Let’s load the relevant packages and create a simulated dataset:

import adelie as ad
import numpy as np
from sklearn.linear_model import ElasticNet

n = 100      # number of observations
p = 10    # number of features

np.random.seed(0)
# adelie works fastest with column-major (Fortran-ordered) matrices, and warns
# if given a row-major (C-ordered) one, which is numpy's default.
X = np.asfortranarray(np.random.normal(0, 1, (n, p)))
y = X[:, -1] * np.random.normal(0, 1) + np.random.normal(0, 1, n)

D.1 Model-fitting with adelie

First, we use ad.grpnet to fit the elastic net model with (lambda, alpha) = (0.3, 0.9):

state = ad.grpnet(
    X=X,
    glm=ad.glm.gaussian(y=y),
    lmda_path = [0.1],
    alpha = 0.9,
)

  0%|          | 0/1 [00:00:00<?, ?it/s]

100%|██████████| 1/1 [00:00:00<00:00:00, 828.29it/s] [dev:37.0%]

100%|██████████| 1/1 [00:00:00<00:00:00, 689.80it/s] [dev:37.0%]

The coefficients and intercepts are stored in the returned state object. In particular, the coefficients are stored in a sparse matrix format.

print(state.betas)
print(state.intercepts)

<Compressed Sparse Row sparse matrix of dtype 'float64'
    with 5 stored elements and shape (1, 10)>
  Coords    Values
  (0, 0)    0.09329161369646884
  (0, 1)    -0.07735457264972584
  (0, 2)    0.036164819426146876
  (0, 6)    0.11936608257545449
  (0, 9)    0.558516180791759
[-0.00370799]

D.2 Model-fitting with sklearn

The code below fits the elastic model in sklearn. To match the hyperparameters for the adelie example, we set alpha = 0.3 and l1_ratio = 0.9.

regr = ElasticNet(alpha = 0.1, l1_ratio = 0.9, random_state=0)
regr.fit(X, y)
print(regr.coef_)
print(regr.intercept_)

[ 0.09329161 -0.07735457  0.03616482 -0.         -0.         -0.
  0.11936608  0.          0.          0.55851618]
-0.0037079903203665054

D.3 Comparing the results

The code below compares the intercept and coefficients from the two methods up to 6 decimal places. You will see that the coefficients agree.

# Convert sparse matrix to dense array
adelie_coefs = state.betas.toarray().flatten()
sklearn_coefs = np.array(regr.coef_)

# Stack them side by side and compare
comparison_matrix = np.vstack([adelie_coefs, sklearn_coefs]).T
intercepts_row = np.array([[state.intercepts[0], regr.intercept_]])
comparison_matrix_with_intercepts = np.vstack([intercepts_row, comparison_matrix])

print("Coef index |    adelie |   sklearn")
for i, (b, r) in enumerate(comparison_matrix_with_intercepts):
    print(f"{i:10} | {b:9.6f} | {r:9.6f}")

Coef index |    adelie |   sklearn
         0 | -0.003708 | -0.003708
         1 |  0.093292 |  0.093292
         2 | -0.077355 | -0.077355
         3 |  0.036165 |  0.036165
         4 |  0.000000 | -0.000000
         5 |  0.000000 | -0.000000
         6 |  0.000000 | -0.000000
         7 |  0.119366 |  0.119366
         8 |  0.000000 |  0.000000
         9 |  0.000000 |  0.000000
        10 |  0.558516 |  0.558516