# How to derive the ridge regression solution?

I am having some issues with the derivation of the solution for ridge regression.

I know the regression solution without the regularization term is given by:

$$\beta = (X^\top X)^{-1}X^\top y.$$

But after adding the L2 term $$\lambda\|\beta\|_2^2$$ to the cost function, how does the solution become:

$$\beta = (X^\top X + \lambda I)^{-1}X^\top y?$$

## 5 Answers

It suffices to modify the loss function by adding the penalty. In matrix terms, the initial quadratic loss function becomes $$(Y - X\beta)^{T}(Y-X\beta) + \lambda \beta^T\beta.$$ Deriving with respect to $\beta$ leads to the normal equation $$X^{T}Y = \left(X^{T}X + \lambda I\right)\beta$$ which leads to the Ridge estimator.

• How come the derivative of $\lambda \beta^T \beta$ is equal to $\lambda I \beta$ Commented Sep 4, 2013 at 16:21
• @user34790 It's not. It's equal to $2\lambda\beta$. But the 2 cancels with similar 2s on the other terms. Of course, the $I$ factor is like a factor of 1 in "regular" algebra, you can multiply it anywhere you like without changing anything.
– Bill
Commented Sep 4, 2013 at 18:42
• @bill: here you need the $I$ to get a matrix of the correct dimension so the addition works with $X^TX$: $\lambda$ is just a scalar Commented Aug 4, 2015 at 17:23

Let's build on what we know, which is that whenever the $$n\times p$$ model matrix is $$X$$, the response $$n$$-vector is $$y$$, and the parameter $$p$$-vector is $$\beta$$, the objective function

$$f(\beta) = (y - X\beta)^\prime(y - X\beta)$$

(which is the sum of squares of residuals) is minimized when $$\beta$$ solves the Normal equations

$$(X^\prime X)\beta = X^\prime y.$$

Ridge regression adds another term to the objective function (usually after standardizing all variables in order to put them on a common footing), asking to minimize

$$(y - X\beta)^\prime(y - X\beta) + \lambda \beta^\prime \beta$$

for some non-negative constant $$\lambda$$. It is the sum of squares of the residuals plus a multiple of the sum of squares of the coefficients themselves (making it obvious that it has a global minimum). Because $$\lambda\ge 0$$, it has a positive square root $$\nu^2 = \lambda$$.

Consider the matrix $$X$$ augmented with rows corresponding to $$\nu$$ times the $$p\times p$$ identity matrix $$I$$:

$$X_{*} = \pmatrix{X \\ \nu I}$$

When the vector $$y$$ is similarly extended with $$p$$ zeros at the end to $$y_{*}$$, the matrix product in the objective function adds $$p$$ additional terms of the form $$(0 - \nu \beta_i)^2 = \lambda \beta_i^2$$ to the original objective. Therefore

$$(y_{*} - X_{*}\beta)^\prime(y_{*} - X_{*}\beta) = (y - X\beta)^\prime(y - X\beta) + \lambda \beta^\prime \beta.$$

From the form of the left hand expression it is immediate that the Normal equations are

$$(X_{*}^\prime X_{*})\beta = X_{*}^\prime y_{*}.$$

Because we adjoined zeros to the end of $$y$$, the right hand side is the same as $$X^\prime y$$. On the left hand side $$\nu^2 I=\lambda I$$ is added to the original $$X^\prime X$$. Therefore the new Normal equations simplify to

$$(X^\prime X + \lambda I)\beta = X^\prime y.$$

Besides being conceptually economical--no new manipulations are needed to derive this result--it also is computationally economical: your software for doing ordinary least squares will also do ridge regression without any change whatsoever. (It nevertheless can be helpful in large problems to use software designed for this purpose, because it will exploit the special structure of $$X_{*}$$ to obtain results efficiently for a densely spaced interval of $$\lambda$$, enabling you to explore how the answers vary with $$\lambda$$.)

Another beauty of this way of looking at things is how it can help us understand ridge regression. When we want to really understand regression, it almost always helps to think of it geometrically: the columns of $$X$$ constitute $$p$$ vectors in a real vector space of dimension $$n$$. By adjoining $$\nu I$$ to $$X$$, thereby prolonging them from $$n$$-vectors to $$n+p$$-vectors, we are embedding $$\mathbb{R}^n$$ into a larger space $$\mathbb{R}^{n+p}$$ by including $$p$$ "imaginary", mutually orthogonal directions. The first column of $$X$$ is given a small imaginary component of size $$\nu$$, thereby lengthening it and moving it out of the space generated by the original $$p$$ columns. The second, third, ..., $$p^\text{th}$$ columns are similarly lengthened and moved out of the original space by the same amount $$\nu$$--but all in different new directions. Consequently, any collinearity present in the original columns will immediately be resolved. Moreover, the larger $$\nu$$ becomes, the more these new vectors approach the individual $$p$$ imaginary directions: they become more and more orthogonal. Consequently, the solution of the Normal equations will immediately become possible and it will rapidly become numerically stable as $$\nu$$ increases from $$0$$.

This description of the process suggests some novel and creative approaches to addressing the problems Ridge Regression was designed to handle. For instance, using any means whatsoever (such as the variance decomposition described by Belsley, Kuh, and Welsch in their 1980 book on Regression Diagnostics, Chapter 3), you might be able to identify subgroups of nearly collinear columns of $$X$$, where each subgroup is nearly orthogonal to any other. You only need adjoin as many rows to $$X$$ (and zeros to $$y$$) as there are elements in the largest group, dedicating one new "imaginary" dimension for displacing each element of a group away from its siblings: you don't need $$p$$ imaginary dimensions to do this.

A method of performing this was proposed by Rolf Sundberg some 30 years ago. See Continuum Regression and Ridge Regression. Journal of the Royal Statistical Society. Series B (Methodological) Vol. 55, No. 3 (1993), pp. 653-659.

• Whoa, this just blew my mind. Is there any discussion about what happens when this is generalized outside linear models, i.e. to glm's? The penalty shouldn't line up to be the same as ridge regression...but this interpretation implies that it would still be a potential useful estimator! Commented Mar 25, 2016 at 19:14
• @Cliff That's a very interesting suggestion. Since, however, GLM estimates depend in a more complicated way on $X$ and their estimators cannot usually be factored in the form $$\hat\beta = g(X)\cdot h(y)$$ as they are for OLS (where $g(X)=(X^\prime X)^{-1}X^\prime$ and $h(y)=y$), it may be difficult to establish a useful relationship between imposing a penalty function and modifying the columns of $X$. In particular, it's unclear how the values in $y$ would need to be augmented in order to make this work.
– whuber
Commented Mar 25, 2016 at 20:45
• Yes, it would take some thought to try establish what the penalty is, but I'm not so concerned about that. The idea of what $y_*$ to use is not generally easy either...except perhaps in the case of logistic regression, where we could add two $y_*$'s; one of 0's and one of 1's. This augmentation would then be a more general version of the "+2 binomial estimator" (there's a more proper name for this estimator that I'm blanking on, which is basically when you are estimating $p$ from a binomial distribution using the posterior mean as the estimate with a uniform prior on $p$). Commented Mar 25, 2016 at 21:15
• @whuber It might also be worth pointing out that this recipe for fitting ridge penalized regression models also extends to nonnegative least squares ridge penalized regression and I believe also to nonnegative (identity link) Poisson ridge penalized regression if one plugs this augmented covariate matrix & obs vector into nnls and nnpois (in the addreg package) or nnlm (using Kullback-Leibler divergence, in the NNLM package). For nonnegative fused ridge regression (stats.stackexchange.com/questions/327412/…) this didn't seem to work unfortunately Commented Apr 23, 2019 at 18:32
• The idea of looking at $X^TX + \lambda I$ as $X_*^T X_*$ is quite clever. Is there a particular way to think about matrix multiplication that makes this insight obvious? Commented Jul 29, 2021 at 17:52

The derivation includes matrix calculus, which can be quite tedious. We would like solve the following problem: $$\min_\beta (Y-\beta^T X)^T(Y-\beta^T X)+\lambda \beta^T \beta$$

Now note that $$\frac{\partial (Y-\beta^T X)^T (Y-\beta^T X)}{\partial \beta}=-2X^T(Y-\beta^T X)$$ and $$\frac{\partial \lambda \beta^T \beta}{\partial \beta}=2\lambda\beta.$$ Together we get to the first order condition $$X^TY = X^TX\beta + \lambda\beta.$$ Isolating $\beta$ yields the solution: $$\beta = (X^TX+ \lambda I )^{-1}X^T Y.$$

I have recently stumbled upon the same question in the context of P-Splines and as the concept is the same I want to give a more detailed answer on the derivation of the ridge estimator.

We start with a penalized criterion function that differs from the classic OLS-criterion function by its penalization term in the last summand:

$Criterion_{Ridge} = \sum_{i=1}^{n}(y_i-x_i^T\beta)^2 + \lambda \sum_{j=1}^p\beta_j^2$

where

• $p=$ the amount of covariables used in the model
• $x_i^T\beta =$ your standard linear predictor
• the first summand respresents the MSE (squared divergence of the prediction from the actual value) that we want to minimize as usual
• the second summand represents the penalization we apply on the coefficients. Here we are in the Ridge-context which implies a Euclidian Distance Measure and therefore the degree of 2 in the penalization term. In the case of a Lasso-Penalization we would apply a degree of 1 and yield a totally different estimator.

We can rewrite this criterion in matrix-notation and further break it down:

$Criterion_{Ridge} = (y-X\beta)^T(y-X\beta) + \lambda\beta^T\beta$

$= y^Ty - \beta^TX^Ty - y^TX\beta+ \beta^Tx^TX\beta + \lambda\beta^T\beta$

$= y^Ty - \beta^TX^Ty - \beta^TX^Ty + \beta^TX^TX\beta + \beta^T\lambda I\beta$ with $I$ being the identity matrix

$= y^Ty - 2\beta^TX^Ty + \beta^T(X^TX + \lambda I)\beta$

Now we search for the $\beta$ that minimizes our criterion. Amongst others we make use of the matrix differentiation rule $\frac{\partial x^TAx}{\partial x} = (A+A^T)x \overset{\text{A symmetric}}{=} 2Ax$ which we can apply here as $(X^TX + \lambda I) \in \mathbb{R}^{n \times n}$:

$\frac{\partial Criterion_{Ridge} }{\partial\beta} = -2X^Ty + 2(X^TX + \lambda I)\beta \overset{!}{=}0$

$(X^TX + \lambda I)\beta = X^Ty$

$\overset{\text{et voilà}}{\Rightarrow} \hat\beta = (X^TX + \lambda I)^{-1} X^Ty$

• @Jahn, can you please explain how $$y^TX\beta$$ became $$\beta ^TX^Ty$$? I think you just applied transpose on it, right. But, you cannot just apply transpose on one term without applying it on all equation. What am I missing here? Commented Jun 17, 2018 at 23:31
• @theateist A transposed scalar is the same scalar. Commented Jun 20, 2018 at 14:40

There are a few important things that are missing in the answers given.

1. The solution for $\beta$ is derived from the first-order necessary condition: $\frac{\partial f_{ridge}(\beta, \lambda)}{\partial \beta} = 0$ which yields $\beta = (X^TX+ \lambda I )^{-1}X^T Y$. But is this sufficient? That is, the solution is a global minimum only if $f_{ridge}(\beta, \lambda)$ is strictly convex. This can be shown to be true.

2. Another way to look at the problem is to see the equivalence between $f_{ridge}(\beta, \lambda)$ and $f_{OLS}(\beta) = (Y-\beta^T X)^T(Y-\beta^T X)$ constrained to $||\beta||^2_2 \leq t$. OLS stands for Ordinary Least Squares. From this perspective $f_{ridge}(\beta, \lambda)$ is just the Lagrangian function used to find the global minima of the convex objective function $f_{OLS}(\beta)$ constrained with convex function $||\beta||^2_2$.

A good explanation of these points and the derivation of $\beta$ can be found in these fine lecture notes: http://math.bu.edu/people/cgineste/classes/ma575/p/w14_1.pdf