Fitting Ridge Regression with glmnet

Once the response and predictor matrix are ready, fitting a ridge model in R is a short call to glmnet. The important work is making the modeling choices explicit: identify a numeric response, pass a numeric predictor matrix with matching rows, choose the Gaussian family, and set the penalty parameter. The returned object contains a fitted model at each penalty value supplied or generated.

This tutorial focuses on fitting and reading the main parts of that object. It does not treat a printed coefficient as if it were an ordinary least-squares estimate: with ridge, the fitted coefficients depend on the penalty, and with glmnet’s default standardization the penalty is applied on an internally standardized predictor scale.

Definition: For a numeric response $y$ and predictor matrix $X$, a Gaussian ridge fit in glmnet uses family = "gaussian" and alpha = 0. The parameter $\lambda$ controls the size of the ridge penalty. By default, the intercept is included and is not penalized; predictors are internally standardized when standardize = TRUE.

Inputs and a Basic glmnet Call

The response supplied as y should be numeric. The predictor input x should be a numeric matrix with one row per observation and one column per predictor. Row $i$ of x must describe the same observation as element $i$ of y. A mismatch can produce an error, but a more serious mismatch can silently fit the wrong response-predictor pairings.

If a design matrix was created with model.matrix, do not pass its intercept column as a predictor when using glmnet’s default intercept = TRUE. The model has its own intercept. A typical call is:

\texttt{glmnet(x, y, family = "gaussian", alpha = 0, lambda = chosen\_values)}

In R, the call might be written as follows. The particular values of lambda are chosen for the analysis; the example uses one positive value so the fitted object contains one penalized model.

Worked Example: Fit a Ridge Model with an Explicit Penalty

A hypothetical analyst models a numeric claim-cost measure using two numeric predictors. Use four observations, turn off standardization to make the arithmetic transparent, and fit at $\lambda=1$.

Step 1: Construct inputs. Let the predictor columns be $x_1=(-1,-1,1,1)$ and $x_2=(-1,1,-1,1)$, with response $y=(6,10,8,16)$. Each column has four entries, so the matrix has four rows and two columns, matching the four responses.

\texttt{x <- cbind(x1 = c(-1,-1,1,1), x2 = c(-1,1,-1,1))}\\ \texttt{y <- c(6,10,8,16)}\\ \texttt{fit <- glmnet::glmnet(x, y, family = "gaussian", alpha = 0,}\\ \texttt{\ \ lambda = 1, standardize = FALSE, intercept = TRUE)}

Step 2: Check the fit by hand. The mean response is $\bar y=(6+10+8+16)/4=10$. Both predictor columns sum to zero and are orthogonal: $x_1^\mathsf{T}x_1=x_2^\mathsf{T}x_2=4$ and $x_1^\mathsf{T}x_2=0$. The centered cross-products are $x_1^\mathsf{T}y=8$ and $x_2^\mathsf{T}y=12$, giving unpenalized slopes 2 and 3.

With $n=4$, $\lambda=1$, and standardize = FALSE, the ridge slopes are $8/(4+4)=1$ and $12/(4+4)=1.5$. Because the predictors are centered, the intercept is 10. Thus the expected intercept is $a_0=10$, and the coefficient column in fit$beta contains 1 and 1.5, up to numerical precision.

Step 3: Verify fit summaries. The fitted values are $7.5, 10.5, 9.5, 12.5$. Using residual = observed minus fitted, the residuals are $-1.5,-0.5,-1.5,3.5$, and the residual sum of squares is $2.25+0.25+2.25+12.25=17$. The null-model deviance is the sum of squared deviations from 10: $16+0+4+36=56$. For a Gaussian model the deviance ratio is $1-\text{RSS}/\text{null deviance}=1-17/56=39/56$, approximately 0.6964. The check is consistent with the expected fit and its summary.

The method is an explicit-$\lambda$ Gaussian ridge fit. A common mistake this example avoids is assuming the default predictor standardization while calculating coefficients for the supplied, unstandardized matrix.

What Objective Is glmnet Fitting?

For equal observation weights and a Gaussian response, glmnet minimizes a residual sum of squares term plus a penalty. With standardize = FALSE, the following objective is expressed directly in terms of the supplied predictors $x_{ij}$ and their slopes $\beta_j$:

\frac{1}{2n}\sum_{i=1}^{n}(y_i-\beta_0-\sum_{j=1}^{p}x_{ij}\beta_j)^2 +\frac{\lambda}{2}\sum_{j=1}^{p}\beta_j^2.

The intercept $\beta_0$ is not included in the penalty. With the default standardize = TRUE, the same form applies to internally standardized predictors $z_{ij}=(x_{ij}-\bar{x}_j)/s_j$ and their coefficients $\gamma_j$, not directly to the original-scale slopes. Here $s_j$ is the scale used for predictor $j$. Since $\gamma_j=s_j\beta_j$, the penalty written in terms of original-scale slopes is $\frac{\lambda}{2}\sum_j(s_j\beta_j)^2$. glmnet reports coefficients transformed back to the supplied predictor scale.

This distinction matters because changing predictor units changes a ridge penalty if predictors are not standardized. The penalty on a slope of one unit per dollar is not equivalent to the penalty on a slope of one unit per thousand dollars. Standardization makes the penalty act on comparable predictor scales, while standardize = FALSE leaves the supplied scales in effect.

Why the Ridge Slope Formula Has This Form

A short derivation shows where the ridge adjustment comes from. Consider centered predictors in a Gaussian model with standardize = FALSE, equal weights, and an unpenalized intercept. For simplicity, suppose the predictor columns are mutually orthogonal. Write $S_{jj}=x_j^\mathsf{T}x_j$ and $c_j=x_j^\mathsf{T}y$, where $y$ is centered. The terms of the objective that depend on slope $j$ are proportional to

\frac{1}{2n}\left(y^\mathsf{T}y-2\beta_j c_j+\beta_j^2 S_{jj}\right) +\frac{\lambda}{2}\beta_j^2.

Differentiate with respect to $\beta_j$ and set the derivative to zero. The derivative is $(-c_j+\beta_j S_{jj})/n+\lambda\beta_j$. Rearranging gives $\beta_j(S_{jj}+n\lambda)=c_j$, hence:

\widehat{\beta}_j=\frac{x_j^\mathsf{T}y}{x_j^\mathsf{T}x_j+n\lambda}.

At $\lambda=0$, this is the ordinary least-squares slope for that orthogonal predictor. For positive $\lambda$, the denominator is larger while the numerator is unchanged, so the absolute slope is reduced. The formula assumes centered, mutually orthogonal predictors and the stated objective; for correlated predictors, the slopes are determined jointly by the ridge normal equations rather than by separate one-predictor formulas.

Default Standardization and Original-Scale Coefficients

The next calculation makes the scale distinction concrete. Assume one predictor, an unpenalized intercept, equal weights, Gaussian errors, and glmnet’s default predictor standardization. The predictor scale is $s=\sqrt{\sum_i(x_i-\bar{x})^2/n}$. For one predictor, translating the standardized-scale ridge estimate back to the original scale gives a slope that can also be computed as

\widehat{\beta} = \frac{\sum_i(x_i-\bar{x})(y_i-\bar{y})} {\sum_i(x_i-\bar{x})^2+n\lambda s^2}, \qquad \widehat{\beta}_0=\bar{y}-\widehat{\beta}\bar{x}.

This expression is for the stated default-standardization setup. It is not the unstandardized formula with the penalty denominator simply added as $n\lambda$; the scale factor $s^2$ is essential.

Worked Example: Account for glmnet’s Default Standardization

A hypothetical pricing analyst fits $y=(1,2,3,4)$ using $x=(0,1,2,3)$. Fit a Gaussian ridge model at $\lambda=0.5$, using the default standardize = TRUE and an intercept.

Step 1: Find the center and scale. Both means are 2.5 and 1.5, respectively. The centered predictor values are $(-1.5,-0.5,0.5,1.5)$, so their sum of squares is $2.25+0.25+0.25+2.25=5$. Thus $s^2=5/4=1.25$.

Step 2: Compute the original-scale slope. The centered cross-product is also 5, since $y=x+1$. The slope is $5/(5+4(0.5)(1.25))=5/7.5=2/3$. The intercept is $2.5-(2/3)(1.5)=1.5$. These are the coefficients glmnet reports on the original predictor scale, subject to rounding.

Step 3: Check the fit. The fitted values are $1.5, 13/6, 17/6, 3.5$. The residuals are $-0.5,-1/6,1/6,0.5$, so the residual sum of squares is $1/4+1/36+1/36+1/4=5/9$. The null deviance is $5$, the sum of squared deviations of $y$ from 2.5. The deviance ratio is $1-(5/9)/5=8/9$, approximately 0.8889. The calculation confirms that the scale-weighted penalty produces the stated fit.

The method uses the default internal predictor scale and then expresses the fitted slope in original units. The common mistake avoided here is treating the penalty as $\lambda\beta^2/2$ on the original slope without its predictor scale factor.

Reading the Main Parts of a Fitted Object

A glmnet fit is an object containing information for the penalty values in the fit. When a single $\lambda$ is supplied, the coefficient matrix has one model column; when a sequence is supplied, it has one column per fitted penalty value. The main components to recognize are:

Component	What it represents	How to interpret it
`fit$a0`	Intercept for each fitted penalty value	With one value of lambda, the single entry is the fitted intercept.
`fit$beta`	Predictor coefficient matrix	Rows correspond to predictors and columns to fitted lambda values; coefficients are returned on the original predictor scale.
`fit$lambda`	Penalty values in the fitted path	Match a coefficient column or summary entry to its corresponding penalty value.
`fit$dev.ratio`	Fraction of null deviance explained at each fitted value	For Gaussian fits, it is based on the reduction in residual sum of squares relative to the null model.
`fit$df`	Number of nonzero predictor coefficients at each fitted value	Ridge usually shrinks coefficients without setting them exactly to zero, so this is not generally a count of variables selected by ridge.

The coefficient matrix is commonly stored as a sparse matrix, even when many or all entries are nonzero. Its dimensions can be checked with dim(fit$beta); the observation count can be inspected with fit$nobs. Object details can vary with model options and package version, so inspect names(fit) rather than assuming every possible component is present in every fit.

Worked Example: Read a Short Ridge Penalty Path

Use a new hypothetical dataset with $x=(0,1,2,3)$ and $y=(1,2,3,4)$. Fit three penalty values, with standardization turned off so the coefficient calculations use the supplied predictor scale.

\texttt{x <- matrix(c(0,1,2,3), ncol = 1, dimnames = list(NULL, "exposure"))}\\ \texttt{y <- c(1,2,3,4)}\\ \texttt{fit <- glmnet::glmnet(x, y, family = "gaussian", alpha = 0,}\\ \texttt{\ \ lambda = c(2,1,0.5), standardize = FALSE, intercept = TRUE)}

Step 1: Calculate the coefficients. Here $\bar{x}=1.5$, $\bar{y}=2.5$, $x^\mathsf{T}x$ after centering is 5, and the centered cross-product $x^\mathsf{T}y$ is 5. Therefore the slope is $5/(5+4\lambda)$, and the intercept is $2.5-1.5\widehat{\beta}$.

At $\lambda=2$, the slope is $5/13$ and the intercept is $25/13$. At $\lambda=1$, the slope is $5/9$ and the intercept is $5/3$. At $\lambda=0.5$, the slope is $5/7$ and the intercept is $10/7$. Thus fit$beta has one predictor row and three model columns, corresponding in order to the three supplied lambda values; fit$a0 has the three matching intercepts.

Step 2: Check the deviance summaries. The null deviance is 5. For a fitted slope $\widehat{\beta}$, the residual sum of squares is $5(1-\widehat{\beta})^2$, because the fitted line differs from the unpenalized line by a centered slope difference. At $\lambda=2$, RSS is $5(8/13)^2=320/169$ and the deviance ratio is $1-64/169=105/169$, approximately 0.6213. At $\lambda=1$, RSS is $5(4/9)^2=80/81$ and the ratio is $65/81$, approximately 0.8025. At $\lambda=0.5$, RSS is $5(2/7)^2=20/49$ and the ratio is $45/49$, approximately 0.9184. These checks agree with the corresponding model columns.

The method is to pair each column of the fitted object with its $\lambda$, then verify the Gaussian deviance ratio against the residual and null deviance. A common mistake avoided here is reading the first coefficient column as the smallest penalty simply because it appears first; always check fit$lambda and match columns to those values.

Common Input and Interpretation Traps

Misaligned rows. Sorting or filtering the predictor data without applying the same operation to the response changes which $y_i$ is paired with each row of $X$. Confirm the row order before fitting.
Non-numeric predictor input. glmnet expects a numeric matrix (or an appropriate sparse matrix), not a data frame containing character columns. Encode categorical predictors in a design matrix before fitting.
Duplicate intercept. Passing a column of ones from a design matrix while also using intercept = TRUE supplies a redundant intercept-like predictor. Usually omit that matrix column and retain glmnet’s intercept.
Wrong family or penalty setting. For a continuous Gaussian response and ridge, state both family = "gaussian" and alpha = 0; relying on defaults makes the intended model less clear.
Scale confusion. With default standardization, the penalty is applied to internally standardized coefficients, while the returned slopes are on the supplied predictor scale. Do not use an unstandardized objective to explain a default-standardized fit.
Reading deviance ratio as a validation score. It summarizes fit relative to the null model on the fitting data. By itself it does not establish how well the model predicts new observations.

Fitting correctly is a prerequisite to interpreting or using the model. A high deviance ratio on the training data is not, by itself, evidence of strong out-of-sample performance, and a nonzero ridge coefficient is not proof that a predictor is independently important. Those conclusions require further model assessment and careful interpretation.

Key takeaway: For a Gaussian ridge fit, align the numeric response and predictor matrix, specify family = "gaussian" and alpha = 0, and make the standardization and lambda choices visible in the call. Read a0, beta, lambda, and deviance summaries as corresponding parts of the fitted path, while keeping the internal penalty scale distinct from the reported coefficient scale.

Check Your Understanding

Work through each problem, then compare your result with the answer shown.

A Gaussian ridge model uses $x=(0,1,2)$, $y=(2,4,6)$, an intercept, $\lambda=1$, and standardize = FALSE. Calculate the fitted slope and intercept using the centered ridge formula. (Answer: slope $0.8$, intercept $3.2$.)
For $x=(0,2,4)$, $y=(1,3,5)$, and $n=3$, suppose glmnet uses its default standardization and $\lambda=1$. The centered predictor sum of squares is 8 and the scale satisfies $s^2=8/3$. Calculate the original-scale slope and intercept. (Answer: slope $0.5$, intercept $2$.)
For $x=(0,1,2,3)$, $y=(1,2,3,4)$, fit with an intercept, $\lambda=1$, and standardize = FALSE. Calculate the residual sum of squares and Gaussian deviance ratio. (Answer: RSS $80/81$, deviance ratio $65/81$, approximately $0.8025$.)
A fit uses two predictors and the lambda sequence $c(3,1,0.2)$. What are the dimensions of fit$beta, and what does its second column correspond to? (Answer: $2$ rows by $3$ columns; the second column contains the two predictor slopes for $\lambda=1$.)
A default-standardized fit has predictor scales $s_1=2$, $s_2=0.5$, original-scale slopes $\beta_1=1$, $\beta_2=4$, and $\lambda=0.25$. Calculate the penalty contribution in the original-scale expression, excluding the residual term and intercept. (Answer: $\frac{0.25}{2}[(2\cdot1)^2+(0.5\cdot4)^2]=1$.)

Fitting Ridge Regression with glmnet

What You'll Learn