Regression Analysis for Finance: OLS, Inference, Robust Methods, and Regularization

Regression connects an outcome to explanatory variables. In finance, it can estimate market exposure, explain credit spreads, or predict returns. The central question is which relationship are we estimating, and what assumptions make the estimate and its uncertainty trustworthy? A fitted relationship alone establishes neither causality nor trading profitability.

Learning goals and notation

After this lesson, you should be able to derive OLS, interpret its geometry, distinguish estimation from inference, choose an appropriate response to problematic errors or predictors, and evaluate a financial regression without using future information.

Symbol	Dimensions / meaning
y	n × 1 vector of observed responses
X	n × p design matrix; each row is one observation
β	p × 1 vector of unknown coefficients
ε	n × 1 vector of unobserved model errors
β̂, ŷ, e	Estimated coefficients, fitted responses, and observed residuals
RSS	Residual sum of squares eᵀe

Throughout, p includes the intercept if one is fitted. A model with an intercept and three predictors has p=4. Classical residual variance estimation requires n>p and full column rank. Statistical statements below are conditional on X, or treat X as fixed, unless stated otherwise.

1. What “linear regression” actually means

The model is:

y=Xβ+ε, or yᵢ=xᵢᵀβ+εᵢ for each row i.

It is linear in the unknown coefficients, not necessarily in the original predictors. The model yᵢ=β₀+β₁xᵢ+β₂xᵢ²+εᵢ is a linear regression because x and x² are known columns of X. Interactions, Fourier components, indicator variables, and lagged observations can also be columns.

By contrast, β₁exp(β₂x) is nonlinear in its parameters. Transforming the response can change both the model and its interpretation; a log-response fit does not automatically give an unbiased level forecast after exponentiation.

Three questions should be kept distinct:

Description: how do variables co-move in this sample?
Prediction: does the relationship forecast unseen outcomes?
Causality: would changing a predictor change the outcome, holding relevant conditions fixed?

OLS provides a fitting rule. Causal interpretation requires identification assumptions or research design beyond obtaining a small residual sum of squares.

2. Deriving ordinary least squares

OLS chooses β to minimize the squared length of the residual:

Q(β)=(y−Xβ)ᵀ(y−Xβ).

Expand and differentiate:

Q(β)=yᵀy−2βᵀXᵀy+βᵀXᵀXβ,

∇Q(β)=−2Xᵀy+2XᵀXβ.

Setting the gradient to zero gives the normal equations:

Xᵀ(y−Xβ̂)=0.

If rank(X)=p, then XᵀX is positive definite and the unique minimizer is:

β̂=(XᵀX)⁻¹Xᵀy.

An exact duplicate predictor breaks uniqueness of the coefficients, but the projection ŷ onto the column space of X remains unique. A pseudoinverse selects a particular coefficient solution. Numerically, QR or SVD is usually preferable to explicitly forming (XᵀX)⁻¹; the displayed inverse is primarily an algebraic expression.

3. Projection geometry and the hat matrix

Define H=X(XᵀX)⁻¹Xᵀ. Then ŷ=Hy and e=(I−H)y. H is symmetric and idempotent: Hᵀ=H and H²=H. It projects y onto the column space of X.

The normal equations say Xᵀe=0. Consequently eᵀŷ=0: residuals are orthogonal to the fitted vector and every fitted predictor column. If the model has an intercept, Σᵢeᵢ=0 as well.

These are algebraic properties of the fit. They do not establish that the population errors have conditional mean zero, that an omitted variable is harmless, or that the model forecasts well.

For an intercept model, define TSS=Σᵢ(yᵢ−ȳ)². Orthogonality gives TSS=Σᵢ(ŷᵢ−ȳ)²+RSS and R²=1−RSS/TSS, provided TSS>0. In-sample R² cannot decrease when nested OLS models add columns; this is not evidence of better out-of-sample prediction. The usual centered decomposition needs care when no intercept is included.

4. A complete worked fit

Use five illustrative observations, not actual market data:

Observation	Predictor x	Response y
1	−2	−2
2	−1	−2
3	0	1
4	1	2
5	2	6

Fit an intercept and slope, so n=5, p=2. Here ȳ=1, x̄=0, Σxᵢ²=10, and Σxᵢyᵢ=20. Thus:

β̂₀=1; β̂₁=2; ŷ=1+2x.

The fitted vector is [−3,−1,1,3,5]ᵀ, and residuals are [1,−1,0,−1,1]ᵀ. Their sum is zero and their dot product with x is zero. RSS=4 and residual degrees of freedom are 5−2=3.

TSS=44, so R²=1−4/44≈90.91%. That high fit on five constructed observations does not validate a real-world model.

5. Which assumptions support which conclusions?

Assumptions	What they support
Full column rank	Unique OLS coefficients
Zero conditional error mean given X	Conditional unbiasedness of OLS
Conditional error covariance σ²I, plus zero conditional mean	Classical covariance formula and Gauss–Markov efficiency
Errors conditional on X are jointly Gaussian with covariance σ²I	Exact finite-sample normal, t, chi-square, and F distributions

Normality is not needed to calculate OLS or establish its unbiasedness. Conversely, Gaussian-looking residuals do not establish exogeneity. “Uncorrelated errors” also differs from “errors uncorrelated with predictors.”

From the model:

β̂−β=(XᵀX)⁻¹Xᵀε.

Therefore, under the classical moment assumptions:

E[β̂ | X]=β, and Cov(β̂ | X)=σ²(XᵀX)⁻¹.

If ε | X∼N(0,σ²I), a linear transformation of that Gaussian vector gives:

β̂ | X∼Nₚ(β,σ²(XᵀX)⁻¹).

For any p-vector t, the conditional MGF of β̂−β is exp{½σ²tᵀ(XᵀX)⁻¹t}, identifying that normal law. This uses the normal model’s MGF, rather than claiming arbitrary regression errors become exactly Gaussian.

6. Estimating error variance and testing coefficients

Since e=(I−H)ε, tr(I−H)=n−p, and I−H is idempotent:

E[RSS | X]=(n−p)σ².

Hence s²=RSS/(n−p) is unbiased under the classical covariance assumptions. Under Gaussian errors, additionally RSS/σ²∼χ²ₙ₋ₚ, independent of β̂. This independence produces exact t- and F-tests.

An individual coefficient

To test H₀:βⱼ=b₀, use:

SE(β̂ⱼ)=s√[(XᵀX)⁻¹]ⱼⱼ,

t=(β̂ⱼ−b₀)/SE(β̂ⱼ)∼tₙ₋ₚ under H₀.

A two-sided 95% confidence interval is β̂ⱼ ± t₀.₉₇₅,ₙ₋ₚ SE(β̂ⱼ). Its frequentist interpretation concerns the coverage of the procedure over repeated samples, not a 95% probability assigned to a fixed parameter after observing the interval.

For the worked fit, XᵀX=diag(5,10), s²=4/3, and SE(β̂₁)=√[(4/3)/10]≈0.3651. Testing slope zero gives t≈5.477 with 3 degrees of freedom. Using t₀.₉₇₅,₃≈3.182, the slope interval is approximately [0.838,3.162], conditional on the normal model.

A small p-value measures incompatibility with the null under the assumed model. It is not the probability that the null is true, and statistical significance does not establish economic importance.

Joint restrictions and nested models

Suppose a full model has p parameters and a nested restricted model has k<p. Fit both to the same response observations, obtaining RSS₁ and RSS₀ respectively. Under the Gaussian model and p−k independent valid restrictions:

F=[(RSS₀−RSS₁)/(p−k)]/[RSS₁/(n−p)]∼Fₚ₋ₖ,ₙ₋ₚ under H₀.

In the example, restricting the slope to zero gives the intercept-only fit, RSS₀=44, k=1. Therefore F=(44−4)/(4/3)=30. With one restriction, F=t²=30. This RSS test cannot directly compare arbitrary non-nested models.

Confidence intervals are narrower than prediction intervals

For a new design vector x₀, the fitted conditional mean is ŷ₀=x₀ᵀβ̂. Its classical standard error is s√[x₀ᵀ(XᵀX)⁻¹x₀]. An independent new response includes a new error, giving prediction standard error s√[1+x₀ᵀ(XᵀX)⁻¹x₀].

At x=0 in the example, the fitted mean is 1. Its standard error is √(4/15)≈0.5164, while a new response’s prediction standard error is √1.6≈1.2649. These formulas require the same error variance and model to apply at the new observation.

7. Gauss–Markov: what BLUE means

Under full column rank, E[ε | X]=0, and Cov(ε | X)=σ²I, OLS is the best linear unbiased estimator (BLUE). “Linear” means linear in y, and “best” means smallest covariance in the positive-semidefinite ordering among unbiased estimators of that form. It does not mean best among every possible nonlinear or biased estimator, or best for forecasting under all losses.

A short proof compares OLS L₀y, with L₀=(XᵀX)⁻¹Xᵀ, to an alternative Ly. Unbiasedness requires LX=I. Write L=L₀+D, so DX=0. Cross terms vanish and:

Cov(Ly | X)=σ²(XᵀX)⁻¹+σ²DDᵀ.

The added matrix is PSD. This establishes BLUE without assuming Gaussian errors.

8. GLS and robust standard errors solve different problems

If Cov(ε | X)=σ²Ω, with known positive-definite Ω, OLS remains unbiased when E[ε | X]=0, but its covariance becomes:

σ²(XᵀX)⁻¹XᵀΩX(XᵀX)⁻¹.

The usual σ²(XᵀX)⁻¹ formula is then generally wrong. Choose a whitening transformation L with LΩLᵀ=I, for example the symmetric Ω⁻¹ᐟ², and fit Ly=LXβ+Lε. This yields:

β̂_GLS=(XᵀΩ⁻¹X)⁻¹XᵀΩ⁻¹y.

For independent errors with variances σ²vᵢ, GLS becomes weighted least squares with weights proportional to 1/vᵢ, not 1/√vᵢ. The transformed observations are scaled by 1/√vᵢ. If Ω must be estimated, feasible GLS adds covariance-model uncertainty and need not retain exact known-Ω finite-sample properties.

An alternative is to keep the OLS coefficients and adjust their estimated uncertainty. Heteroscedasticity-consistent standard errors address unequal error variance; HAC estimates, such as Newey–West, also address suitable serial dependence. These rely on their own sampling and asymptotic conditions. They do not repair endogeneity, omitted-variable bias, structural breaks, or a misspecified conditional mean.

9. Gaussian maximum likelihood and M-estimation

Under the classical Gaussian error model, conditional log-likelihood is:

ℓ(β,σ²)=−(n/2)log(2πσ²)−RSS(β)/(2σ²).

Maximizing over β gives the OLS coefficients. At the optimum, σ̂²_ML=RSS/n. Its conditional expectation is [(n−p)/n]σ², so it is downward biased when p>0. The unbiased estimator divides by n−p. In the worked fit, σ̂²_ML=0.8, whereas s²=4/3.

An M-estimator minimizes a sum of losses, Σᵢρ(yᵢ−xᵢᵀβ), or a more general observation-specific criterion. The choice of loss changes what the regression targets.

Loss	Target / motivation	Limitation
Squared error e²	Conditional mean under suitable assumptions	Large residuals receive quadratic weight
Absolute error, abs(e)	Conditional median; less sensitive to large response residuals	Does not automatically protect against high-leverage predictor contamination
Huber loss	Quadratic near zero, linear in large residuals	Requires a scale and tuning choice
Quantile pinball loss	Conditional quantile τ	Tail quantiles need substantial information and stable specification

Least absolute deviations is often abbreviated LAD. Calling it “MAD regression” can be confusing because MAD commonly means median absolute deviation, a scale statistic.

For 0<τ<1, pinball loss is ρτ(e)=e[τ−1{e<0}]. At τ=0.5 it is half the absolute loss. Bounded or redescending losses can limit the effect of extreme residuals but may introduce nonconvex optimization; robust residual loss alone is not a universal solution to influential predictor values. Robust-loss inference also differs from classical Gaussian OLS inference.

10. Multicollinearity and regularization

If two predictors almost duplicate one another, individual coefficients can change sharply while fitted values remain similar. With the thin SVD X=UDVᵀ, positive singular values dⱼ determine the OLS coefficient directions:

β̂_OLS=Σⱼvⱼ(uⱼᵀy)/dⱼ.

Small dⱼ amplify noise; under classical errors, variance along vⱼ is σ²/dⱼ². Perfect collinearity causes rank deficiency. Adding more observations helps only if they provide useful independent design information.

Ridge: stabilize weakly identified directions

For the objective ‖y−Xβ‖²+λ‖β‖², λ>0:

β̂_ridge=(XᵀX+λI)⁻¹Xᵀy.

In singular directions, ridge uses dⱼ/(dⱼ²+λ) instead of 1/dⱼ, shrinking the OLS component by dⱼ²/(dⱼ²+λ). If d₁=10, d₂=1, and λ=4, the factors are 100/104≈0.9615 and 1/5=0.2. The weak direction shrinks much more.

The formula penalizes every coordinate. In common practice, center the response and predictors, penalize the slopes, then restore an unpenalized intercept. Equivalently use a penalty matrix with zero in the intercept position. Standardizing predictors makes the penalty comparable across measurement units; fit that scaling on training data only.

Ridge is biased but can lower mean squared estimation or prediction error. With known Gaussian noise variance σ² and independent slope priors N(0,τ²), its MAP coefficients have λ=σ²/τ² for the unscaled objective above. Other conventions, such as RSS/n or RSS/(2n), change the numeric mapping of λ. This is a penalized coefficient estimate, not automatic Bayesian uncertainty quantification.

Lasso: shrink and select

For ‖y−Xβ‖²+λΣⱼ|βⱼ|, the L₁ penalty can set coefficients exactly to zero. It encourages sparsity, but does not guarantee that the selected predictors are the truly informative ones; correlated columns can lead to unstable selections. Lasso generally has no ordinary inverse formula, and coefficient solutions need not be unique for every design.

Standardization and intercept handling matter here too. Ordinary OLS p-values applied after selecting variables on the same sample ignore selection uncertainty. Elastic net combines L₁ and L₂ penalties and can be useful when predictors form correlated groups.

Principal component regression: truncate directions

For a centered design, keep the first m singular components and define Zₘ=XVₘ=UₘDₘ. Regress the centered y on these orthogonal scores and map coefficients back:

β̂_PCR=VₘDₘ⁻¹Uₘᵀy, retaining only nonzero singular values.

PCR removes selected low-variance predictor directions; ridge shrinks directions continuously. The most variable directions in X need not be the directions most predictive of y. An omitted low-variance component can carry the signal. Both m and λ should be chosen using appropriate validation, rather than optimizing the final test sample.

For derivations and comparisons, see Hastie, Tibshirani, and Friedman’s Elements of Statistical Learning, especially its chapter on linear regression methods.

11. Diagnostics: investigate the observations and the model

Leverage is hᵢᵢ=xᵢᵀ(XᵀX)⁻¹xᵢ, the ith diagonal of H. It measures design-space extremeness, not the size of yᵢ. The leverage values sum to p, so their average is p/n.

For hᵢᵢ<1, an internally studentized residual is:

rᵢ=eᵢ/[s√(1−hᵢᵢ)].

An externally studentized residual replaces s with the residual scale estimated after deleting observation i. The two quantities have different reference distributions; neither should be treated as an ordinary raw residual.

Cook’s distance combines residual size and leverage:

Dᵢ=[eᵢ²/(p s²)]×hᵢᵢ/(1−hᵢᵢ)².

In the worked example, endpoint leverage is 1/5+2²/10=0.6, while center leverage is 0.2. For the first observation, e₁=1, so D₁=1.40625. This flags substantial influence in this small design, rather than proving the observation is wrong.

Diagnostic	What to investigate
Residuals versus fitted values	Curvature, changing variance, missing interactions
Residuals over time and autocorrelation	Serial dependence, volatility clustering, breaks
Q–Q plot	Deviations from the assumed error distribution
Leverage and Cook’s distance	Unusual predictor combinations and sensitivity to observations
Singular values or condition number	Weakly identified coefficient combinations
Rolling coefficients and held-out errors	Instability and generalization

Outliers, high leverage, and influence are related but distinct. Investigate data quality and compare sensitivity fits; automatic deletion can hide economically important crisis observations. Penn State’s influential-points lesson explains these distinctions; NIST’s diagnostic reference lists the associated statistics.

12. Financial applications and a validation workflow

CAPM regression

Using matched-period simple total returns and risk-free returns:

rᵢ,ₜ−r_f,ₜ=αᵢ+βᵢ(rₘ,ₜ−r_f,ₜ)+εₜ.

The slope estimates market exposure. The intercept is average excess return unexplained by this particular specification, not proof of skill. CAPM’s equilibrium restriction α=0 is a hypothesis to test under suitable assumptions. Monthly alpha is not automatically an annual compounded realized return; arithmetic annualization and compounding are different operations.

Sector ETF regressions

An illustrative contemporaneous model relates XLP total returns to SPY, MDY, QQQ, and DIA total returns. These broad-market predictors can overlap strongly, making individual slopes unstable. A large R² may mostly reflect common market exposure rather than several distinct economic effects.

A regression using today’s market returns explains today’s sector movement. A tradable forecast of tomorrow’s return must use predictors actually available at the decision time, such as appropriately lagged observations. Align currencies, trading calendars, corporate-action adjustments, and publication times. No empirical fit or performance claim is being made here.

Corporate credit spreads

Spreads might be modeled using rates, leverage, volatility, and liquidity proxies. Define whether the response is a level, change, or horizon-specific forecast. Persistent nonstationary levels can produce spurious regressions unless the time-series structure, such as cointegration, is justified. Differences can change the question and do not automatically cure every specification problem. Endogenous predictors may require an identification strategy rather than different standard errors.

A practical sequence

Define the outcome, units, horizon, and information available at the prediction time.
Align observations and establish a chronological training, validation, and final test split.
Fit imputations, transformations, scaling, and PCA only within each training window.
Fit a simple baseline, inspect rank and diagnostics, and choose an uncertainty estimator appropriate to the data.
Tune regularization using rolling or expanding validation windows. Overlapping future-return labels may require gaps or purging to prevent information overlap.
Evaluate the untouched test period against a relevant benchmark using forecast errors. If converting predictions into trades, also account for turnover, costs, financing, and drawdowns.
Report coefficient instability, specification choices, and the number of tried models. Repeatedly selecting the best backtest can overfit the research process itself.

Randomly shuffling a time series can leak future regimes or overlapping information into training. A high in-sample R², a low p-value, and a profitable future strategy are different claims that need different evidence.

13. Practice questions

A regression includes an intercept and four predictors, with n=60 observations. What are p and the residual degrees of freedom?
Is β₀+β₁x+β₂sin(x) linear regression? What about β₁exp(β₂x)?
A full model has RSS₁=80 and p=5; a nested restricted model has RSS₀=100 and k=3, using n=45 observations. Compute the classical F statistic.
Which exact identities in an OLS fit hold even when the error model is misspecified?
Does Gauss–Markov require Gaussian errors? Do exact classical finite-sample t-tests require more assumptions?
Does ridge shrink a singular direction with d=0.5 more or less than one with d=5, for the same positive λ?
Can PCR discard a useful predictive direction? Can LAD still be sensitive to high-leverage observations?
Does fitting XLP on same-day SPY returns establish a strategy for predicting tomorrow’s XLP return?

Answers

p=5 and n−p=55.
The first is linear in its coefficients. The second is nonlinear in β₂.
F=[(100−80)/2]/[80/40]=5, with 2 and 40 degrees of freedom under the null and classical normal model.
Xᵀe=0, eᵀŷ=0, and, with an intercept, Σeᵢ=0. They are consequences of fitting by OLS, not proof of population exogeneity.
No Gaussianity is required for BLUE. Exact classical t distributions additionally follow from the joint Gaussian spherical-error model stated above.
More: d²/(d²+λ) is smaller for the weaker singular direction.
Yes to both. Predictor variance is not predictive importance, and reducing response-outlier sensitivity does not eliminate leverage problems.
No. It estimates contemporaneous co-movement; a forecast needs a stated future target and predictors available before that outcome.