Linear Regression: Modeling and Gradient Descent

This page develops the prediction model and iterative fitting. The companion ordinary least-squares note explains the closed-form solution, rank conditions, statistical uncertainty, and coefficient penalties.

Start with two dots and a straight line. The slope tilts the line; the intercept moves it up or down. Fitting means adjusting those two numbers to reduce prediction errors. Try it before reading the matrix notation.

Try it yourself

Teach a line to fit two dots

The green dots are observations; the blue line is your prediction. Move the sliders to shorten the dashed errors, then let a gradient step do it for you. MSE is the average of the squared errors: smaller means a closer fit to these points.

ObservationsPrediction

Intercept b · move up or down

0.000

Slope w · change the tilt

0.000

Learning rate · step size0.010.11

MSE 5.0000

ŷ = 0.000x + (0.000)∇L [b, w] = [-4.000, -3.000]

Try this: reset, take one step at 0.1, then reset and try 1. A bigger step can make the error worse. Can you fit both dots yourself?

The controls stay between −3 and 4. A step outside that range is disabled; lower the learning rate or reset. The chart clips lines outside its vertical window. Playback stops after 60 steps, near zero error, or when this panel is hidden. Reduced-motion mode keeps manual steps available.

Open the code · NumPy / JAX

The interaction above runs in JavaScript. These Python examples use the same data and equations, starting from the defaults—not your current controls. JAX can compute the loss gradient with jax.grad. JAX autodiff guide

In a new local project, save the example as demo.py. CPU setup: uv init ai-lab → cd ai-lab → uv add numpy jax → uv run demo.py

NumPy

import numpy as np

x = np.array([0., 1.])
y = np.array([1., 3.])
p = np.array([0., 0.])  # [b, w]

def loss(p):
    return np.mean((p[1] * x + p[0] - y) ** 2)

def gradient(p):
    error = p[1] * x + p[0] - y
    return np.array([2 * error.mean(), 2 * (error * x).mean()])

p = p - 0.1 * gradient(p)
print(p, loss(p))  # [0.4, 0.3], 2.825

JAX

import jax
import jax.numpy as jnp

x = jnp.array([0., 1.])
y = jnp.array([1., 3.])
p = jnp.array([0., 0.])  # [b, w]

def loss(p):
    return jnp.mean((p[1] * x + p[0] - y) ** 2)

gradient = jax.grad(loss)
p = p - 0.1 * gradient(p)
print(p, loss(p))  # approximately [0.4, 0.3], 2.825

Linear regression predicts a numerical target as an affine function of the features:

$$\hat{y} = \mathbf{x}^{\top}\mathbf{w} + b.$$

For $n$ observations, collect the rows into a design matrix $\mathbf{X}$. If the intercept is represented by an added column of ones, the model becomes

$$\hat{\mathbf{y}} = \mathbf{X}\boldsymbol{\theta}.$$

The word linear refers to the parameters. Features may include transformations such as interactions or polynomial terms while the model remains linear in $\boldsymbol{\theta}$.

Squared-Error Objective

Ordinary least squares minimizes the residual sum of squares, or equivalently mean squared error:

$$\hat{\boldsymbol{\theta}} = \arg\min_{\boldsymbol{\theta}} \frac{1}{n}\left\|\mathbf{X}\boldsymbol{\theta}-\mathbf{y}\right\|_2^2.$$

Gradient descent reaches the same least-squares objective through successive updates. It is especially useful for large datasets and models embedded in a wider differentiable system.

From Residuals to an Update

With $d$ features and an intercept, $X$ has shape $n\times(d+1)$ and $\theta$ has $d+1$ entries. Differentiating the MSE above gives

$$\nabla_\theta L=\frac{2}{n}X^\top(X\theta-y), \qquad \theta\leftarrow\theta-\eta\nabla_\theta L.$$

The learning rate $\eta>0$ sets the update size. For points $(x,y)=(0,1),(1,3)$, start at $(b,w)=(0,0)$. Prediction-minus-target errors $X\theta-y$ are $[-1,-3]$, so the gradient in $(b,w)$ order is $[-4,-3]$. With $\eta=0.1$, the new parameters are $(0.4,0.3)$ and MSE falls from $5$ to $2.825$. The exact fit is $(b,w)=(1,2)$. Inference only evaluates $wx+b$; it does not use the target or update weights.

Scale and step size

Changing a feature from meters to millimeters multiplies its numerical values by 1,000. Unless weights and step size are adjusted, its contribution to the gradient can become much larger. Fit centering and scaling on training data, then reuse those statistics at prediction time. Scaling changes the optimization geometry; it does not create evidence that a linear relationship is appropriate.

The training-process note derives a learning-rate stability bound for a small quadratic example and connects this update to minibatches, momentum, and validation. Reducing training loss is an optimization result; compare held-out error against a training-mean baseline before claiming useful prediction.

Interpretation Boundaries

• A coefficient is conditional on the included features and their scaling.
• Good prediction does not establish a causal relationship.
• Statistical inference needs assumptions beyond minimizing squared error; prediction and inference are different goals.
• Collinearity can make coefficients unstable even when predictions remain adequate.
• Extrapolation outside the observed feature range relies heavily on the chosen functional form.
• Residual plots and out-of-sample error reveal failures that training loss cannot.

Linear regression is valuable both as a model and as a baseline. A more complex method should justify itself against this simpler alternative.

For a full worked fit, uncertainty assumptions and regularization, see Ordinary Least Squares. That note uses half-MSE, so its loss gradient is half the MSE gradient shown here; learning-rate and penalty conventions must be compared accordingly.

Continue with Softmax Regression for mutually exclusive classification. The main external reference is Dive into Deep Learning: Linear Neural Networks for Regression.

My reading · Opened notes are remembered in this browser.