Ridge regression
View on Wikipedia| Part of a series on |
| Regression analysis |
|---|
| Models |
| Estimation |
| Background |
Ridge regression (also known as Tikhonov regularization, named for Andrey Tikhonov) is a method of estimating the coefficients of multiple-regression models in scenarios where the independent variables are highly correlated.[1] It has been used in many fields including econometrics, chemistry, and engineering.[2] It is a method of regularization of ill-posed problems.[a] It is particularly useful to mitigate the problem of multicollinearity in linear regression, which commonly occurs in models with large numbers of parameters.[3] In general, the method provides improved efficiency in parameter estimation problems in exchange for a tolerable amount of bias (see bias–variance tradeoff).[4]
The theory was first introduced by Hoerl and Kennard in 1970 in their Technometrics papers "Ridge regressions: biased estimation of nonorthogonal problems" and "Ridge regressions: applications in nonorthogonal problems".[5][6][1]
Ridge regression was developed as a possible solution to the imprecision of least square estimators when linear regression models have some multicollinear (highly correlated) independent variables—by creating a ridge regression estimator (RR). This provides a more precise ridge parameters estimate, as its variance and mean square estimator are often smaller than the least square estimators previously derived.[7][2]
Overview
[edit]In the ordinary least squares solution of
the problem of a near-singular moment matrix is alleviated by adding positive elements to the diagonals, thereby decreasing its condition number. Compared to the ordinary least squares estimator, the simple ridge estimator has an extra term in the denominator: where is the regressand, is the design matrix, is the identity matrix, and the ridge parameter serves as the constant shifting the diagonals of the moment matrix.[8] It can be shown that this estimator is the solution to the least squares problem subject to the constraint , which can be expressed as a Lagrangian minimization: which shows that is nothing but the Lagrange multiplier of the constraint.[9] In fact, there is a one-to-one relationship between and and since, in practice, we do not know , we define heuristically or find it via additional data-fitting strategies, see Determination of the Tikhonov factor.
Note that, when , in which case the constraint is non-binding, the ridge estimator reduces to ordinary least squares. A more general approach to Tikhonov regularization is discussed below.
History
[edit]Tikhonov regularization was invented independently in many different contexts. It became widely known through its application to integral equations in the works of Andrey Tikhonov[10][11][12][13][14] and David L. Phillips.[15] Some authors use the term Tikhonov–Phillips regularization. The finite-dimensional case was expounded by Arthur E. Hoerl, who took a statistical approach,[16] and by Manus Foster, who interpreted this method as a Wiener–Kolmogorov (Kriging) filter.[17] Following Hoerl, it is known in the statistical literature as ridge regression,[18] named after ridge analysis ("ridge" refers to the path from the constrained maximum).[19]
Tikhonov regularization
[edit]Suppose that for a known real matrix and vector , we wish to find a vector such that where and may be of different sizes and may be non-square.
The standard approach is ordinary least squares linear regression.[clarification needed] However, if no satisfies the equation or more than one does—that is, the solution is not unique—the problem is said to be ill posed. In such cases, ordinary least squares estimation leads to an overdetermined, or more often an underdetermined system of equations. Most real-world phenomena have the effect of low-pass filters[clarification needed] in the forward direction where maps to . Therefore, in solving the inverse-problem, the inverse mapping operates as a high-pass filter that has the undesirable tendency of amplifying noise (eigenvalues / singular values are largest in the reverse mapping where they were smallest in the forward mapping). In addition, ordinary least squares implicitly nullifies every element of the reconstructed version of that is in the null-space of , rather than allowing for a model to be used as a prior for . Ordinary least squares seeks to minimize the sum of squared residuals, which can be compactly written as where is the Euclidean norm.
In order to give preference to a particular solution with desirable properties, a regularization term can be included in this minimization: where and , for some suitably chosen Tikhonov matrix . In many cases, this matrix is chosen as a scalar multiple of the identity matrix (), giving preference to solutions with smaller norms; this is known as L2 regularization.[20] In other cases, high-pass operators (e.g., a difference operator or a weighted Fourier operator) may be used to enforce smoothness if the underlying vector is believed to be mostly continuous. This regularization improves the conditioning of the problem, thus enabling a direct numerical solution. Treating it as an ordinary least squares problem with augmented matrices and , the solution is The effect of regularization may be varied by the scale of the matrix . For this reduces to the unregularized least-squares solution, provided that (ATA)−1 exists. Note that in case of a complex matrix , as usual the transpose has to be replaced by the Hermitian transpose .
L2 regularization is used in many contexts aside from linear regression, such as classification with logistic regression or support vector machines,[21] and matrix factorization.[22]
Application to existing fit results
[edit]Since Tikhonov Regularization simply adds a quadratic term to the objective function in optimization problems, it is possible to do so after the unregularised optimisation has taken place. E.g., if the above problem with yields the solution , the solution in the presence of can be expressed as: with the "regularisation matrix" .
If the parameter fit comes with a covariance matrix of the estimated parameter uncertainties , then the regularisation matrix will be and the regularised result will have a new covariance
In the context of arbitrary likelihood fits, this is valid, as long as the quadratic approximation of the likelihood function is valid. This means that, as long as the perturbation from the unregularised result is small, one can regularise any result that is presented as a best fit point with a covariance matrix. No detailed knowledge of the underlying likelihood function is needed.[23]
Generalized Tikhonov regularization
[edit]For general multivariate normal distributions for and the data error, one can apply a transformation of the variables to reduce to the case above. Equivalently, one can seek an to minimize where we have used to stand for the weighted norm squared (compare with the Mahalanobis distance). In the Bayesian interpretation is the inverse covariance matrix of , is the expected value of , and is the inverse covariance matrix of .
The Tikhonov matrix is not explicitly included because the corresponding regularization term reduces to above with and . For normal regularization where , the Tikhonov matrix then appears in the Cholesky factorization and is considered a whitening filter.
This generalized problem has an optimal solution which can be written explicitly using the formula
Lavrentyev regularization
[edit]In some situations, one can avoid using the transpose , as proposed by Mikhail Lavrentyev.[24] For example, if is symmetric positive definite, i.e. , so is its inverse , which can thus be used to set up the weighted norm squared in the generalized Tikhonov regularization, leading to minimizing or, equivalently up to a constant term,
This minimization problem has an optimal solution which can be written explicitly using the formula which is nothing but the solution of the generalized Tikhonov problem where
The Lavrentyev regularization, if applicable, is advantageous to the original Tikhonov regularization, since the Lavrentyev matrix can be better conditioned, i.e., have a smaller condition number, compared to the Tikhonov matrix
Regularization in Hilbert space
[edit]Typically discrete linear ill-conditioned problems result from discretization of integral equations, and one can formulate a Tikhonov regularization in the original infinite-dimensional context. In the above we can interpret as a compact operator on Hilbert spaces, and and as elements in the domain and range of . The operator is then a self-adjoint bounded invertible operator.
Relation to singular-value decomposition and Wiener filter
[edit]With , this least-squares solution can be analyzed in a special way using the singular-value decomposition. Given the singular value decomposition with singular values , the Tikhonov regularized solution can be expressed as where has diagonal values and is zero elsewhere. This demonstrates the effect of the Tikhonov parameter on the condition number of the regularized problem. For the generalized case, a similar representation can be derived using a generalized singular-value decomposition.[25]
Finally, it is related to the Wiener filter: where the Wiener weights are and is the rank of .
Determination of the Tikhonov factor
[edit]The optimal regularization parameter is usually unknown and often in practical problems is determined by an ad hoc method. A possible approach relies on the Bayesian interpretation described below. Other approaches include the discrepancy principle, cross-validation, L-curve method,[26] restricted maximum likelihood and unbiased predictive risk estimator. Grace Wahba proved that the optimal parameter, in the sense of leave-one-out cross-validation minimizes[27][28] where is the residual sum of squares, and is the effective number of degrees of freedom.
Using the previous SVD decomposition, we can simplify the above expression: and
Relation to probabilistic formulation
[edit]The probabilistic formulation of an inverse problem introduces (when all uncertainties are Gaussian) a covariance matrix representing the a priori uncertainties on the model parameters, and a covariance matrix representing the uncertainties on the observed parameters.[29] In the special case when these two matrices are diagonal and isotropic, and , and, in this case, the equations of inverse theory reduce to the equations above, with .[30][31]
Bayesian interpretation
[edit]Although at first the choice of the solution to this regularized problem may look artificial, and indeed the matrix seems rather arbitrary, the process can be justified from a Bayesian point of view.[32] Note that for an ill-posed problem one must necessarily introduce some additional assumptions in order to get a unique solution. Statistically, the prior probability distribution of is sometimes taken to be a multivariate normal distribution.[33] For simplicity here, the following assumptions are made: the means are zero; their components are independent; the components have the same standard deviation . The data are also subject to errors, and the errors in are also assumed to be independent with zero mean and standard deviation . Under these assumptions the Tikhonov-regularized solution is the most probable solution given the data and the a priori distribution of , according to Bayes' theorem.[34]
If the assumption of normality is replaced by assumptions of homoscedasticity and uncorrelatedness of errors, and if one still assumes zero mean, then the Gauss–Markov theorem entails that the solution is the minimal unbiased linear estimator.[35]
See also
[edit]- LASSO estimator is another regularization method in statistics.
- Elastic net regularization
- Matrix regularization
- L-curve
Notes
[edit]- ^ In statistics, the method is known as ridge regression, in machine learning it and its modifications are known as weight decay, and with multiple independent discoveries, it is also variously known as the Tikhonov–Miller method, the Phillips–Twomey method, the constrained linear inversion method, L2 regularization, and the method of linear regularization. It is related to the Levenberg–Marquardt algorithm for non-linear least-squares problems.
References
[edit]- ^ a b Hilt, Donald E.; Seegrist, Donald W. (1977). Ridge, a computer program for calculating ridge regression estimates. doi:10.5962/bhl.title.68934.[page needed]
- ^ a b Gruber, Marvin (1998). Improving Efficiency by Shrinkage: The James--Stein and Ridge Regression Estimators. CRC Press. p. 2. ISBN 978-0-8247-0156-7.
- ^ Kennedy, Peter (2003). A Guide to Econometrics (Fifth ed.). Cambridge: The MIT Press. pp. 205–206. ISBN 0-262-61183-X.
- ^ Gruber, Marvin (1998). Improving Efficiency by Shrinkage: The James–Stein and Ridge Regression Estimators. Boca Raton: CRC Press. pp. 7–15. ISBN 0-8247-0156-9.
- ^ Hoerl, Arthur E.; Kennard, Robert W. (1970). "Ridge Regression: Biased Estimation for Nonorthogonal Problems". Technometrics. 12 (1): 55–67. doi:10.2307/1267351. JSTOR 1267351.
- ^ Hoerl, Arthur E.; Kennard, Robert W. (1970). "Ridge Regression: Applications to Nonorthogonal Problems". Technometrics. 12 (1): 69–82. doi:10.2307/1267352. JSTOR 1267352.
- ^ Jolliffe, I. T. (2006). Principal Component Analysis. Springer Science & Business Media. p. 178. ISBN 978-0-387-22440-4.
- ^ For the choice of in practice, see Khalaf, Ghadban; Shukur, Ghazi (2005). "Choosing Ridge Parameter for Regression Problems". Communications in Statistics – Theory and Methods. 34 (5): 1177–1182. doi:10.1081/STA-200056836. S2CID 122983724.
- ^ van Wieringen, Wessel (2021-05-31). "Lecture notes on ridge regression". arXiv:1509.09169 [stat.ME].
- ^ Tikhonov, Andrey Nikolayevich (1943). "Об устойчивости обратных задач" [On the stability of inverse problems]. Doklady Akademii Nauk SSSR. 39 (5): 195–198. Archived from the original on 2005-02-27.
- ^ Tikhonov, A. N. (1963). "О решении некорректно поставленных задач и методе регуляризации". Doklady Akademii Nauk SSSR. 151: 501–504.. Translated in "Solution of incorrectly formulated problems and the regularization method". Soviet Mathematics. 4: 1035–1038.
- ^ Tikhonov, A. N.; V. Y. Arsenin (1977). Solution of Ill-posed Problems. Washington: Winston & Sons. ISBN 0-470-99124-0.
- ^ Tikhonov, Andrey Nikolayevich; Goncharsky, A.; Stepanov, V. V.; Yagola, Anatolij Grigorevic (30 June 1995). Numerical Methods for the Solution of Ill-Posed Problems. Netherlands: Springer Netherlands. ISBN 0-7923-3583-X. Retrieved 9 August 2018.
- ^ Tikhonov, Andrey Nikolaevich; Leonov, Aleksandr S.; Yagola, Anatolij Grigorevic (1998). Nonlinear ill-posed problems. London: Chapman & Hall. ISBN 0-412-78660-5. Retrieved 9 August 2018.
- ^ Phillips, D. L. (1962). "A Technique for the Numerical Solution of Certain Integral Equations of the First Kind". Journal of the ACM. 9: 84–97. doi:10.1145/321105.321114. S2CID 35368397.
- ^ Hoerl, Arthur E. (1962). "Application of Ridge Analysis to Regression Problems". Chemical Engineering Progress. 58 (3): 54–59.
- ^ Foster, M. (1961). "An Application of the Wiener-Kolmogorov Smoothing Theory to Matrix Inversion". Journal of the Society for Industrial and Applied Mathematics. 9 (3): 387–392. doi:10.1137/0109031.
- ^ Hoerl, A. E.; R. W. Kennard (1970). "Ridge regression: Biased estimation for nonorthogonal problems". Technometrics. 12 (1): 55–67. doi:10.1080/00401706.1970.10488634.
- ^ Hoerl, Roger W. (2020-10-01). "Ridge Regression: A Historical Context". Technometrics. 62 (4): 420–425. doi:10.1080/00401706.2020.1742207. ISSN 0040-1706.
- ^ Ng, Andrew Y. (2004). Feature selection, L1 vs. L2 regularization, and rotational invariance (PDF). Proc. ICML.
- ^ R.-E. Fan; K.-W. Chang; C.-J. Hsieh; X.-R. Wang; C.-J. Lin (2008). "LIBLINEAR: A library for large linear classification". Journal of Machine Learning Research. 9: 1871–1874.
- ^ Guan, Naiyang; Tao, Dacheng; Luo, Zhigang; Yuan, Bo (2012). "Online nonnegative matrix factorization with robust stochastic approximation". IEEE Transactions on Neural Networks and Learning Systems. 23 (7): 1087–1099. Bibcode:2012ITNNL..23.1087G. doi:10.1109/TNNLS.2012.2197827. PMID 24807135. S2CID 8755408.
- ^ Koch, Lukas (2022). "Post-hoc regularisation of unfolded cross-section measurements". Journal of Instrumentation. 17 (10) 10021. arXiv:2207.02125. Bibcode:2022JInst..17P0021K. doi:10.1088/1748-0221/17/10/P10021.
- ^ Lavrentiev, M. M. (1967). Some Improperly Posed Problems of Mathematical Physics. New York: Springer.
- ^ Hansen, Per Christian (Jan 1, 1998). Rank-Deficient and Discrete Ill-Posed Problems: Numerical Aspects of Linear Inversion (1st ed.). Philadelphia, USA: SIAM. ISBN 978-0-89871-403-6.
- ^ P. C. Hansen, "The L-curve and its use in the numerical treatment of inverse problems", [1]
- ^ Wahba, G. (1990). "Spline Models for Observational Data". CBMS-NSF Regional Conference Series in Applied Mathematics. Society for Industrial and Applied Mathematics. Bibcode:1990smod.conf.....W.
- ^ Golub, G.; Heath, M.; Wahba, G. (1979). "Generalized cross-validation as a method for choosing a good ridge parameter" (PDF). Technometrics. 21 (2): 215–223. doi:10.1080/00401706.1979.10489751.
- ^ Tarantola, Albert (2005). Inverse Problem Theory and Methods for Model Parameter Estimation (1st ed.). Philadelphia: Society for Industrial and Applied Mathematics (SIAM). ISBN 0-89871-792-2. Retrieved 9 August 2018.
- ^ Huang, Yunfei.; et al. (2019). "Traction force microscopy with optimized regularization and automated Bayesian parameter selection for comparing cells". Scientific Reports. 9 (1) 539: 537. arXiv:1810.05848. Bibcode:2019NatSR...9..539H. doi:10.1038/s41598-018-36896-x. PMC 6345967. PMID 30679578.
- ^ Huang, Yunfei; Gompper, Gerhard; Sabass, Benedikt (2020). "A Bayesian traction force microscopy method with automated denoising in a user-friendly software package". Computer Physics Communications. 256 107313. arXiv:2005.01377. Bibcode:2020CoPhC.25607313H. doi:10.1016/j.cpc.2020.107313.
- ^ Greenberg, Edward; Webster, Charles E. Jr. (1983). Advanced Econometrics: A Bridge to the Literature. New York: John Wiley & Sons. pp. 207–213. ISBN 0-471-09077-8.
- ^ Huang, Yunfei.; et al. (2019). "Traction force microscopy with optimized regularization and automated Bayesian parameter selection for comparing cells". Scientific Reports. 9 (1) 539: 537. arXiv:1810.05848. Bibcode:2019NatSR...9..539H. doi:10.1038/s41598-018-36896-x. PMC 6345967. PMID 30679578.
- ^ Vogel, Curtis R. (2002). Computational methods for inverse problems. Philadelphia: Society for Industrial and Applied Mathematics. ISBN 0-89871-550-4.
- ^ Amemiya, Takeshi (1985). Advanced Econometrics. Harvard University Press. pp. 60–61. ISBN 0-674-00560-0.
Further reading
[edit]- Gruber, Marvin (1998). Improving Efficiency by Shrinkage: The James–Stein and Ridge Regression Estimators. Boca Raton: CRC Press. ISBN 0-8247-0156-9.
- Kress, Rainer (1998). "Tikhonov Regularization". Numerical Analysis. New York: Springer. pp. 86–90. ISBN 0-387-98408-9.
- Press, W. H.; Teukolsky, S. A.; Vetterling, W. T.; Flannery, B. P. (2007). "Section 19.5. Linear Regularization Methods". Numerical Recipes: The Art of Scientific Computing (3rd ed.). New York: Cambridge University Press. ISBN 978-0-521-88068-8.
- Saleh, A. K. Md. Ehsanes; Arashi, Mohammad; Kibria, B. M. Golam (2019). Theory of Ridge Regression Estimation with Applications. New York: John Wiley & Sons. ISBN 978-1-118-64461-4.
- Taddy, Matt (2019). "Regularization". Business Data Science: Combining Machine Learning and Economics to Optimize, Automate, and Accelerate Business Decisions. New York: McGraw-Hill. pp. 69–104. ISBN 978-1-260-45277-8.
Ridge regression
View on GrokipediaIntroduction
Definition and Purpose
Ridge regression is a biased estimation technique for linear regression models that addresses issues arising from multicollinearity among predictor variables by incorporating an L2 penalty term into the objective function. Specifically, it minimizes the sum of squared residuals between observed and predicted values, augmented by a regularization term , where represents the vector of regression coefficients and is a tuning parameter controlling the strength of the penalty.[2] This formulation shrinks the coefficients toward zero, preventing extreme values that can occur in the presence of highly correlated predictors.[1] The primary purpose of ridge regression is to stabilize coefficient estimates and improve prediction accuracy in ill-posed problems, such as those with multicollinearity, where ordinary least squares (OLS) estimates can exhibit high variance and instability due to near-singular design matrices. By introducing a small amount of bias, ridge regression reduces the variance of the estimates, leading to a better bias-variance trade-off and more reliable out-of-sample predictions, particularly when predictors are intercorrelated.[5] For instance, in a real estate pricing model where square footage and number of rooms are highly correlated predictors of house price, OLS might produce unstable coefficients sensitive to small data changes, whereas ridge regression shrinks these coefficients proportionally, yielding more consistent estimates across similar datasets.[5] Geometrically, ridge regression can be understood as finding the coefficients that minimize the residual sum of squares subject to a constraint on the Euclidean norm of , which traces out a circular boundary in the standardized parameter space. The solutions lie at the points where this circle intersects the elliptical contours of the residual sum of squares, which are elongated due to multicollinearity; this intersection shrinks coefficients toward the origin without setting any to exactly zero, unlike some other regularization methods.[6]Comparison to Ordinary Least Squares
Ridge regression introduces a deliberate bias into the coefficient estimates through the L2 penalty term, which shrinks the magnitudes of the coefficients toward zero, in contrast to ordinary least squares (OLS), which yields unbiased estimates that minimize the residual sum of squares without regularization. This bias reduces the variance of the estimates, particularly when predictors exhibit multicollinearity or in high-dimensional settings where the number of features approaches the sample size. The bias-variance decomposition reveals that the total expected prediction error, or mean squared error (MSE), is the sum of irreducible error, bias squared, and variance; ridge regression trades a modest increase in bias for a substantial decrease in variance, often yielding lower overall MSE than OLS in scenarios prone to overfitting or instability.[2] In the presence of multicollinearity, where independent variables are highly correlated, OLS estimates become unstable and exhibit large standard errors, as the design matrix becomes ill-conditioned with small eigenvalues, amplifying the impact of noise on coefficient values. Ridge regression addresses this by stabilizing the estimates through uniform shrinkage across all coefficients, preventing extreme values and improving the reliability of predictions without eliminating variables entirely. This makes ridge particularly advantageous for datasets with near-linear dependencies among predictors, where OLS might produce coefficients that are difficult to interpret or overly sensitive to minor data perturbations.[2] Under classical asymptotic assumptions with a fixed number of parameters p and sample size n approaching infinity, OLS is a consistent estimator that converges in probability to the true parameters and is asymptotically efficient, attaining the minimum variance among unbiased estimators. Ridge regression with a fixed regularization parameter λ > 0 introduces bias in finite samples but is asymptotically unbiased and consistent, similar to OLS, attaining asymptotic normality. In finite-sample multicollinear contexts, its variance reduction typically leads to superior MSE performance despite the bias. Simulations highlight these contrasts effectively; for instance, in generated datasets with severe multicollinearity (e.g., a condition number exceeding 10^3 for the correlation matrix), ridge regression with λ tuned via cross-validation substantially outperforms OLS by reducing MSE and enhancing out-of-sample prediction accuracy. These results underscore ridge's practical edge over OLS when collinearity inflates variance, though OLS remains preferable in well-conditioned, low-dimensional orthogonal designs.[7]Mathematical Formulation
Linear Model Setup
In the standard setup for linear regression, the observed response vector , an vector, is expressed as , where is the design matrix containing the predictor variables, is the vector of unknown regression coefficients, and represents the error term.[8] The errors are typically assumed to be independently and identically distributed as , implying zero mean, constant variance , and independence across observations. The foundational assumptions of this model include linearity of the response in the parameters, independence of the errors, homoscedasticity (equal variance of errors), and no perfect multicollinearity among the predictors (ensuring is invertible).[9] While normality of errors is often invoked for statistical inference under the full Gaussian linear model, ridge regression applications relax this requirement, focusing instead on bias-variance trade-offs without relying on distributional assumptions for consistency.[3] A frequent practical violation of these assumptions is multicollinearity, where predictors exhibit high linear correlations, inflating the variance of coefficient estimates and rendering the model sensitive to minor data perturbations.[10] To facilitate analysis, particularly in contexts addressing multicollinearity, the predictors in and the response are conventionally centered by subtracting their respective means, yielding and with zero column and vector means, respectively; this centering simplifies the ridge solutions by eliminating the need for an intercept term in centered data.[3] Scaling the centered variables—dividing each column of by its standard deviation and similarly for —is also common to ensure comparable magnitudes across coefficients, though not strictly required for the model setup.[11] This configuration becomes ill-posed when is nearly singular due to multicollinearity, leading to unstable inverses in estimation procedures and highly variable coefficient estimates that poorly generalize beyond the sample.[12]Ridge Estimator Derivation
The ridge regression estimator addresses multicollinearity in linear regression by introducing a penalty term to the least squares objective function. Specifically, it seeks to minimize the sum of the squared residuals and the squared Euclidean norm of the coefficient vector, scaled by a positive regularization parameter . This formulation is given byHistorical Context
Early Developments
In the 1950s, econometric analysis encountered persistent issues with multicollinearity, where highly correlated explanatory variables in economic datasets—such as those involving time series or cross-sectional observations—produced near-singular information matrices, rendering ordinary least squares (OLS) estimates highly unstable and sensitive to small perturbations in the data. This challenge was highlighted in foundational econometric works, including the 1953 Cowles Commission monograph Studies in Econometric Method, edited by William C. Hood and Tjalling C. Koopmans, which discussed the impact of multicollinearity on coefficient precision and inference in economic modeling, a term coined by Ragnar Frisch in 1934. Independently, in the field of inverse problems, Andrey Tikhonov developed similar regularization techniques in 1943, later known as Tikhonov regularization, which is mathematically equivalent to ridge regression. Preceding the formalization of ridge regression, early shrinkage concepts emerged in statistical theory through Charles Stein's 1956 demonstration of the inadmissibility of the standard maximum likelihood estimator for the mean of a multivariate normal distribution under squared error loss, suggesting that intentionally biased estimators could achieve lower overall risk in high dimensions—a principle that indirectly influenced later biased regression techniques, though not directly applied to linear models at the time. Practical motivations for ridge methods arose in fields like chemistry and engineering, where OLS often failed due to multicollinearity in experimental data; for instance, in chemical process optimization, correlated factors such as temperature, pressure, and concentration in mixture experiments led to ill-conditioned design matrices, causing extreme variance in parameter estimates and unreliable predictions. Arthur E. Hoerl addressed this in his 1962 work on ridge analysis, initially developed for response surface methodology in chemical engineering to trace paths of steepest ascent or descent along constrained ridges, stabilizing interpretations when the Hessian matrix was nearly singular. These applications underscored the need for estimators that traded unbiasedness for reduced variance in real-world datasets with inherent correlations. The explicit introduction of ridge regression occurred in the 1970 paper by Arthur E. Hoerl and Robert W. Kennard, who proposed it as a deliberate bias-introducing approach to enhance estimator stability in multicollinear settings, building on ridge analysis by augmenting the diagonal of the cross-product matrix to mitigate the effects of nonorthogonality while preserving predictive accuracy.[2]Key Milestones and Contributors
Ridge regression emerged as a response to challenges posed by multicollinearity in linear regression models, building on earlier recognition of instability in ordinary least squares estimates. The foundational contributions came from Arthur E. Hoerl and Robert W. Kennard, who published a series of influential papers in the early 1970s that established ridge regression as a practical tool for biased estimation. Their seminal 1970 article, "Ridge Regression: Biased Estimation for Nonorthogonal Problems," introduced the ridge estimator as a method to stabilize coefficients by adding a penalty term, supported by simulation studies showing reduced variance and improved mean squared error in multicollinear settings. A companion paper that same year, "Ridge Regression: Applications to Nonorthogonal Problems," applied the technique to real datasets, including chemical process data, and demonstrated its empirical benefits through ridge traces—graphical tools for selecting the shrinkage parameter. Subsequent works by Hoerl and Kennard, such as their 1975 paper on simulations and 1976 exploration of generalized ridge estimators, further validated the approach across diverse scenarios, solidifying its role as a standard method for ill-conditioned regression problems. Concurrent with Hoerl and Kennard's efforts, Donald W. Marquardt contributed to the theoretical and computational foundations in his 1970 paper, "Generalized Inverses, Ridge Regression, Biased Linear Estimation, and Nonlinear Estimation." Marquardt highlighted the connections between ridge regression and generalized inverses, emphasizing computational efficiency for solving the augmented normal equations and properties like bias-variance trade-offs in high-dimensional settings.[14] His work provided early insights into implementation challenges, such as numerical stability, which were crucial for practical adoption. The 1970s saw significant debate over the admissibility of ridge estimators, particularly their introduction of bias, which some statisticians argued violated classical principles of unbiasedness and invariance under data transformations. Critics like G. Smith and R. Campbell, in their 1980 critique, questioned the data-dependent selection of the ridge parameter and potential for inconsistent interpretations across reparameterizations of the model. Hoerl and colleagues responded in subsequent publications, including a 1986 paper, by focusing on mean squared error (MSE) as the primary criterion, where ridge regression consistently outperformed ordinary least squares in prediction accuracy and parameter stability under multicollinearity, as evidenced by extensive Monte Carlo simulations. This MSE-based resolution shifted emphasis from strict admissibility to pragmatic performance, helping to mainstream the technique. Later formalizations advanced the field's rigor. In 1981, Hrishikesh D. Vinod and Aman Ullah published the textbook Recent Advances in Regression Methods, which systematically derived ridge estimators, discussed optimality under various loss functions, and integrated them into broader shrinkage estimation frameworks, serving as a key reference for theoretical developments.[15] Key milestones in adoption occurred in the 1980s with ridge regression's integration into major statistical software packages, notably the RIDGE option in SAS PROC REG introduced around 1980, which automated ridge trace plots and parameter selection, enabling routine use in applied research across industries like engineering and economics. By the 1990s, extensions to generalized linear models emerged, adapting ridge penalties to non-normal responses such as logistic and Poisson regression; for instance, E. C. Malthouse's 1999 work applied ridge methods to scoring models in direct marketing, demonstrating improved stability in generalized settings with binary outcomes. These advancements broadened ridge regression's applicability beyond linear models, influencing modern regularization techniques.Regularization Frameworks
Tikhonov Regularization Equivalence
Tikhonov regularization addresses ill-posed inverse problems by solving the optimization problemLavrentyev Regularization Variant
Lavrentyev regularization, developed by Mikhail M. Lavrentyev in 1967, represents a foundational approach to addressing ill-posed inverse problems, particularly in the context of Fredholm integral equations of the first kind arising in physics.[16] The method transforms the original ill-posed equation $ Ax = b $, where $ A $ is a compact self-adjoint positive operator on a Hilbert space, into the well-posed regularized equation $ (A + \lambda I)x_{\lambda} = b $ for a positive parameter $ \lambda > 0 $. This solution $ x_{\lambda} $ offers stability by damping the amplification of data errors in the inverse process. Historically, Lavrentyev's formulation provided one of the earliest systematic solutions to ill-posed problems in mathematical physics, such as recovering potentials or sources from integral data, well before the statistical adaptation of similar ideas in ridge regression during the 1970s. It emphasized the need for regularization to ensure continuous dependence of solutions on input data, building on concepts of improper problems introduced in Soviet mathematical literature. In contrast to the more general Tikhonov regularization, which incorporates an arbitrary positive definite matrix $ L $ in the penalty term $ |Lx|^2 $, Lavrentyev's variant assumes $ L = I $ and directly perturbs the forward operator $ A $ itself, making it particularly suited to self-adjoint operators without requiring a separate design of the regularization matrix.[17] This simplicity restricts its direct application to finite-dimensional regression but excels in abstract operator settings over infinite-dimensional function spaces.[18] Today, Lavrentyev regularization remains relevant in numerical analysis for solving discretized versions of integral equations, where the finite-dimensional approximations yield systems analogous to ridge regression, facilitating stable computations in applications like geophysics and imaging.[19]Theoretical Connections
Hilbert Space Perspective
In the Hilbert space perspective, ridge regression generalizes to Tikhonov regularization for solving ill-posed operator equations in a Hilbert space , where the goal is to find minimizing for a bounded linear operator between Hilbert spaces and , observed data , and regularization parameter . The solution satisfies the normal equation , where denotes the adjoint operator and is the identity on . This framework extends the finite-dimensional ridge estimator to infinite-dimensional settings, enabling the treatment of continuous models such as partial differential equations (PDEs) or functional data analysis, where finite-dimensional ridge regression arises as a discretization of the infinite-dimensional problem. It provides stability against noise and ill-posedness inherent in such operators, which often have compact or unbounded inverses. Convergence theory establishes that, as the regularization parameter and the noise level (with ), the regularized solutions converge to the true minimizer of the original problem under suitable source conditions, such as belonging to the range of for some . These conditions ensure not only convergence but also optimal rates, balancing bias and variance in the infinite-dimensional context. A representative application is the inverse heat conduction problem, where Tikhonov regularization in Hilbert spaces recovers the initial temperature distribution from final measurements, stabilizing the inherently ill-posed backward heat equation by incorporating the -norm penalty on the solution.[20]Singular Value Decomposition Link
The singular value decomposition (SVD) of the design matrix is given by , where and are orthogonal matrices, and is a rectangular diagonal matrix containing the singular values along its main diagonal. Substituting this decomposition into the ridge estimator yields the explicit form , or equivalently in terms of the singular values, the -th principal component coefficient is shrunk by the factor relative to the ordinary least squares estimates in the principal component basis.[21] This shrinkage mechanism damps the contributions from principal components associated with small singular values , which correspond to directions in the feature space dominated by noise or multicollinearity, thereby stabilizing the estimator in those ill-conditioned subspaces while preserving signal in directions with large . For , the factors reduce to the ordinary least squares solution, but as increases, the damping becomes more pronounced for smaller , effectively filtering out high-variance noise components.[21] In signal processing, ridge regression corresponds to a Wiener filter, which is the optimal linear estimator minimizing mean squared error under additive noise assumptions, where controls the trade-off between signal fidelity and noise suppression based on the signal-to-noise ratio in each spectral component. The filter function mirrors the Wiener form, attenuating components where noise dominates (small ). Computationally, the SVD enables efficient evaluation of the ridge estimator, particularly when , by avoiding direct inversion of the potentially ill-conditioned matrix; instead, the diagonal structure of allows element-wise operations, reducing complexity from to dominated by the SVD itself.[22] This approach is especially beneficial for high-dimensional problems, as it leverages stable numerical libraries for SVD and permits rapid recomputation for multiple values.[23]Parameter Estimation
Criteria for Tikhonov Factor Selection
The Tikhonov factor λ in ridge regression governs the balance between fidelity to the observed data and the penalty on the magnitude of the coefficient vector, thereby controlling the bias-variance trade-off in the estimator. For small values of λ approaching zero, the ridge estimator converges to the ordinary least squares solution, retaining high variance in the presence of multicollinearity. Conversely, as λ increases to large values, the coefficients are increasingly shrunk toward zero, which diminishes variance but introduces greater bias. Analytical methods for selecting λ focus on minimizing estimates of prediction error or mean squared error without relying on resampling. One prominent approach is the generalized cross-validation (GCV) criterion, which provides an approximately unbiased estimate of the prediction mean squared error. The GCV function is defined asCross-Validation Approaches
Cross-validation approaches provide data-driven methods for selecting the regularization parameter λ in ridge regression, enabling the estimation of prediction error without relying on independent test data. These techniques involve partitioning the dataset into subsets, fitting the model on training portions while evaluating performance on held-out portions, and choosing the λ that minimizes the average prediction error across partitions. This process helps balance the bias-variance tradeoff inherent in ridge regression by penalizing excessive model complexity. k-fold cross-validation is a widely used resampling method for tuning λ, where the dataset is randomly divided into k equally sized folds. For each candidate λ in a predefined grid (typically spanning several orders of magnitude, such as from 10^{-5} to 10^{5}), the model is trained on k-1 folds and its prediction error—often measured by mean squared error—is computed on the remaining fold; this is repeated for all folds, and the average error determines the performance for that λ. The optimal λ is then the one yielding the lowest cross-validation error, after which the final model is refit on the full dataset. Leave-one-out cross-validation serves as a special case of k-fold where k equals the sample size n, providing an nearly unbiased estimate of prediction error but at higher computational expense, as it requires n separate model fits. To address the computational demands of repeated refitting, variants like generalized cross-validation (GCV) offer efficient approximations. GCV estimates the leave-one-out prediction error using a single fit of the ridge model, leveraging the trace of the hat matrix to compute an effective degrees-of-freedom adjustment without needing to refit for each observation; it is particularly advantageous for ridge regression due to the closed-form solution for the ridge estimator. Random subsampling, also known as Monte Carlo cross-validation, accelerates the process by repeatedly drawing random subsets (e.g., 80% of the data for training and 20% for validation) over fewer iterations than full k-fold, trading some precision for speed in large datasets. Practical considerations arise in implementing these approaches, particularly with small sample sizes where cross-validation can introduce bias in error estimates due to high variance in fold assignments, potentially leading to overly optimistic λ selections. For high-dimensional settings with large n or p (number of predictors), the computational cost of grid search combined with k-fold or GCV can be substantial, often necessitating parallelization or coarser grids to maintain feasibility.[26][26]Probabilistic Interpretations
Relation to Maximum Likelihood
In the standard linear regression model $ y = X\beta + \epsilon $, where the errors $ \epsilon $ are independent and identically distributed as Gaussian with mean zero and variance $ \sigma^2 $, the ordinary least squares (OLS) estimator coincides with the maximum likelihood estimator (MLE) of the parameter vector $ \beta $.[27] Ridge regression emerges as the solution to a constrained version of this MLE problem, where an additional restriction bounds the magnitude of the parameters: $ \min_\beta |y - X\beta|^2 $ subject to $ \beta^T \beta \leq t $ for some positive constant $ t $. This formulation addresses instability in the unconstrained OLS estimates arising from multicollinearity in $ X $, by implicitly assuming the true $ \beta $ lies within a spherical region of radius $ \sqrt{t} $.[28] To derive the ridge estimator, introduce a Lagrange multiplier $ \lambda \geq 0 $ for the constraint, yielding the LagrangianBayesian Prior Formulation
Ridge regression can be interpreted within a Bayesian framework by modeling the regression coefficients with a zero-mean Gaussian prior, , where is the noise variance, is the regularization parameter, and is the identity matrix. This prior assumes independence among coefficients and favors values close to zero, with the variance controlling the degree of shrinkage. The likelihood is specified as multivariate normal, , for observations and predictors.[30] Given the conjugate Gaussian prior and likelihood, the posterior distribution of is also Gaussian, with a closed-form expression for the mean that exactly matches the ridge estimator:This equivalence arises because maximizing the log-posterior is equivalent to minimizing the ridge regression objective, , up to scaling by . The derivation follows directly from the properties of the multivariate normal distribution under this setup.[30][31] The hyperparameter relates to the prior variance as , where is the variance of the Gaussian prior on each coefficient; smaller (larger ) induces stronger shrinkage toward zero.[30] This Bayesian perspective offers advantages such as the ability to derive credible intervals for from the posterior covariance matrix, providing uncertainty quantification beyond point estimates. It also enables extensions to hierarchical Bayesian models, where (or ) is treated as a hyperparameter and estimated via the marginal likelihood of the data, facilitating empirical Bayes approaches for regularization parameter selection.[30]