Every symbol, named
Single feature
- y
- Target. The number you want to predict. Continuous — price, temperature, score.
- x
- Feature. The input. Also called predictor, covariate, independent variable.
- β₀
- Intercept / bias. Prediction when every feature is zero. Often meaningless alone — it anchors the line vertically.
- β₁
- Coefficient / weight. Expected change in y per one-unit increase in x, everything else held fixed.
- ε
- Irreducible error. Everything the features cannot explain. Assumed random, mean zero.
- ŷ
- Prediction. What the fitted model outputs — the line's height at that x.
Many features — matrix form
Past one input, stacking into matrices is the only sane notation — and it's exactly how the code runs.
- X
- Design matrix, n × (p+1). One row per observation, one column per feature, plus a leading column of 1s so the intercept falls out of the same multiplication.
- β
- Column vector of p+1 weights — the only thing being learned.
- n, p
- Observations, features. When p ≥ n, ordinary least squares breaks down and you must regularize.
What "best fit" optimizes
Fitting is a minimization problem. Define a number that measures how wrong the model is, then search for the weights that make it smallest. Here that number is the sum of squared residuals.
Cost function
Same minimum, different scale. MSE is the average; RMSE is its square root and reads in the units of y.
Why squared
- Signs cancel otherwise. Raw errors sum to roughly zero for any sensible line.
- Smooth and convex. Differentiable everywhere, one global minimum, no local traps.
- Closed-form solution exists. Absolute error gives you none.
- Maximum likelihood under Gaussian noise. If errors are normal, least squares is the most likely explanation.
The loss surface
Plot cost against the two parameters and you get a bowl. Convex, one bottom. Every fitting method is a different way of reaching it.
Closed form vs. gradient descent
Linear regression is unusual: it has an exact algebraic answer. Most models don't. Learn both routes — the iterative one generalizes to everything else in machine learning.
Normal equation — solve directly
Set the derivative of SSE to zero and solve. One shot, no tuning, exact.
- No learning rate, no iterations, no convergence worry
- Roughly O(p³) — painful past a few thousand features
- Fails when XᵀX is singular: perfectly correlated features, or p > n
- Real libraries use QR or SVD rather than literally inverting
Gradient descent — walk downhill
Start anywhere, step against the gradient until the cost stops improving.
- Scales to millions of rows and features; works in mini-batches and streaming
- Needs a learning rate α — too big diverges, too small crawls
- Requires feature scaling, or the bowl becomes a ravine and progress zig-zags
- The same algorithm that trains every neural network
| Variant | Update uses | Character |
|---|---|---|
| Batch GD | All n rows per step | Smooth, stable, slow on big data |
| Stochastic GD | One row per step | Noisy path, very fast, needs decaying α |
| Mini-batch GD | 32–512 rows per step | The practical default everywhere |
How good is the fit?
| Metric | Definition | Read it as |
|---|---|---|
| R² | 1 − SSE / SST | Fraction of variance explained. 1 = perfect, 0 = no better than predicting the mean, negative = worse than the mean. |
| Adjusted R² | Penalizes extra features | Use when comparing models with different feature counts — plain R² never decreases when you add a feature, even a random one. |
| RMSE | √(SSE / n) | Typical error in the units of y. Punishes large misses. The usual headline number. |
| MAE | mean |y − ŷ| | Typical error, outlier-tolerant. Report next to RMSE — a big gap between them signals outliers. |
| MAPE | mean |y − ŷ| / |y| | Percentage error. Convenient for business, explodes when y is near zero. |
| p-values, CIs | Per-coefficient inference | Is this weight distinguishable from zero? Statistics cares deeply; pure prediction work often ignores them. |
Residual plots say more than R²
Plot residuals against predicted values. You want structureless noise. Any pattern is the model telling you what it failed to capture.
Constant spread, centred on zero. Nothing left to model.
Error grows with the prediction. Try log(y), or weighted least squares.
A straight line can't reach the bend. Add polynomial or interaction terms.
Five conditions, and what breaks without them
These matter most when you want to interpret coefficients or quote confidence intervals. For raw prediction accuracy, violations are less fatal — but they still usually mean a better model exists.
| Assumption | Meaning | Detect | Fix |
|---|---|---|---|
| Linearity | The relationship really is a weighted sum | Residuals vs fitted; partial plots | Transform features, add polynomial or interaction terms, switch model |
| Independence | Observations don't influence each other | Durbin–Watson; residuals in time order | Time-series models, mixed effects, clustered errors |
| Homoscedasticity | Error spread is constant across predictions | Funnel shape; Breusch–Pagan | Log-transform y, weighted least squares, robust standard errors |
| Normality of ε | Errors are Gaussian | Q–Q plot of residuals | Only affects inference; large n makes it largely moot |
| No multicollinearity | Features aren't near-duplicates | Correlation matrix; VIF above ~5–10 | Drop or combine features, use PCA, use Ridge |
Shrinking the weights on purpose
Add a penalty on coefficient size to the cost. The model trades a little training accuracy for weights that don't swing wildly — usually better performance on data it hasn't seen. The bias–variance tradeoff, made into a dial.
Ridge · L2
Shrinks all coefficients smoothly toward zero, never exactly to zero. Handles correlated features gracefully by splitting weight between them. Keeps the closed form and fixes singular XᵀX.
Lasso · L1
Drives some coefficients to exactly zero — automatic feature selection. Sparse, readable models. With correlated features it arbitrarily picks one and discards the rest.
Elastic Net
Both penalties. Sparsity from L1, stability from L2. The default when you have many correlated features and want selection anyway.
Why L1 zeroes things out
The penalty defines a region the coefficients must stay inside. Cost contours grow outward from the unpenalized solution until they touch that region — and where they touch is the answer.
L2's region is a circle, so contact happens almost anywhere. L1's is a diamond with corners on the axes — and expanding ellipses hit corners first. A corner means a coefficient of exactly zero.
Choosing λ: cross-validate over a logarithmic grid. λ = 0 is plain OLS; as λ → ∞ every coefficient collapses to zero and the model predicts the mean. The sweet spot is empirical.
Underfit, right, overfit
Same data, three models. The middle one is the goal — and you can only identify it with data the model never trained on.
Too rigid to follow the signal. Poor on training and test. Add features or complexity.
Follows the trend, ignores the noise. Training and test scores are close and both decent.
Memorized the noise. Near-perfect on training, bad on test. Regularize or simplify.
Diagnostic shortcut: a large gap between training and validation score is a variance problem — regularize, get more data, cut features. Both scores poor is a bias problem — more expressive model, better features.
Getting data ready
Scaling
Standardize (mean 0, sd 1) or min–max normalize. Mandatory for gradient descent and any regularized fit; irrelevant for the plain normal equation. Fit the scaler on training data only, then apply it to the rest.
Categorical features
One-hot encode, dropping one level to avoid the dummy-variable trap. For high-cardinality categories consider target encoding — carefully, with cross-fitting, or it leaks.
Missing values
Impute with median or a model, and add a binary "was missing" indicator — absence often carries signal. Dropping rows quietly biases the sample.
Transforming the target
Right-skewed targets (income, price, counts) often behave far better as log(y). Exponentiate predictions back, and note that coefficients then read as approximate percentage changes.
Outliers & leverage
Squared loss makes single points powerful. A high-leverage point — extreme in x — can pivot the whole line. Check Cook's distance; investigate before deleting.
Leakage
The most common way a great model turns out worthless. If a feature encodes information unavailable at prediction time, R² looks spectacular and the model fails in production.
Extensions and close relatives
| Model | What changes | Use when |
|---|---|---|
| Polynomial | Add x², x³ … as features | Curved relationship, still fit by least squares |
| Interactions | Add x₁·x₂ terms | One feature's effect depends on another's value |
| Splines / GAM | Piecewise smooth basis functions | Flexible curves with interpretability intact |
| Logistic | Sigmoid output, log-loss | Binary classification — same linear core |
| Poisson | Log link, count distribution | Counts and rates |
| Huber | Squared near zero, absolute in the tails | Outliers you can't remove |
| Quantile | Pinball loss | You want the median or the 90th percentile, not the mean |
| Bayesian linear | Priors on β, posterior distribution | Small data, or you need honest uncertainty |
| PLS / PCR | Regress on components | Many highly correlated features |
Its role in machine learning
Why it stays relevant
- The baseline. Fit it first, always. If a gradient-boosted ensemble can't beat it meaningfully, the complexity isn't earning its keep.
- Interpretable by construction. Each coefficient is a sentence: one more bedroom, £14,000 more. Regulated domains often require exactly this.
- Cheap. Trains in milliseconds, predicts with a dot product, deploys anywhere, needs no GPU.
- The atom of deep learning. A single neuron is a linear regression; an activation function is what makes stacking them worthwhile.
- A teaching device. Loss functions, gradients, regularization, bias–variance, cross-validation — every ML concept appears here in its simplest form.
When to reach for something else
- Genuinely non-linear relationship of unknown shape → gradient boosting, random forests
- Complicated interactions you can't enumerate → tree ensembles find these automatically
- Unstructured input: images, audio, text → neural networks
- Very high dimension with strong correlation → regularized variants, or reduce dimensions first
- Calibrated class probabilities → logistic regression or a probabilistic classifier
The whole workflow, minimally
# scikit-learn: pipeline keeps scaling inside cross-validation, so no leakage from sklearn.pipeline import make_pipeline from sklearn.preprocessing import StandardScaler from sklearn.linear_model import LinearRegression, RidgeCV from sklearn.model_selection import train_test_split, cross_val_score from sklearn.metrics import mean_squared_error, r2_score import numpy as np X_tr, X_te, y_tr, y_te = train_test_split(X, y, test_size=0.2, random_state=42) # baseline ols = LinearRegression().fit(X_tr, y_tr) # regularized, with lambda chosen by cross-validation model = make_pipeline( StandardScaler(), RidgeCV(alphas=np.logspace(-3, 3, 25)), ).fit(X_tr, y_tr) pred = model.predict(X_te) print("RMSE", mean_squared_error(y_te, pred, squared=False)) print("R2 ", r2_score(y_te, pred)) print("CV ", cross_val_score(model, X_tr, y_tr, cv=5, scoring="r2").mean()) # read the model: coefficients compare only because features were scaled coefs = model[-1].coef_ for name, w in sorted(zip(feature_names, coefs), key=lambda t: -abs(t[1])): print(f"{name:>22} {w:+.3f}")
Use statsmodels.OLS when you want the inference table — standard errors, t-statistics, p-values, confidence intervals, diagnostics. scikit-learn is built for prediction; statsmodels is built for explanation.
Mistakes that look like success
An order that works
| Step | Do this | You'll know it when |
|---|---|---|
| 1 | Fit a one-feature model by hand: compute slope and intercept from the formulas | The numbers match what sklearn returns |
| 2 | Implement least squares in NumPy — normal equation first, then gradient descent | Both converge to the same β |
| 3 | Plot cost against iterations for three learning rates | You recognize divergence, crawling, and healthy descent on sight |
| 4 | Take a real dataset. Full loop: EDA, split, scale, fit, residual plots, metrics | Your residual plot is structureless |
| 5 | Add Ridge and Lasso, sweep λ, plot the coefficient paths | You can watch coefficients shrink and, under L1, hit zero |
| 6 | Deliberately break an assumption and observe the damage | You can diagnose the problem from the residuals alone |
| 7 | Move to logistic regression | You can explain what changed: the link function and the loss, not the linear core |
The one-paragraph summary
Linear regression predicts a continuous number as a weighted sum of inputs. Fitting means choosing weights that minimize squared error, either by solving a matrix equation directly or by walking downhill on the loss surface. Evaluate on held-out data with RMSE and R², and always look at the residuals. Add an L1 or L2 penalty when you have many features or limited data. Its assumptions govern how far you can trust the coefficients as explanations. It is simultaneously the simplest useful model in machine learning, the baseline everything else must beat, and the single neuron that deep learning is built from.