Least squares
Least squares estimates model parameters by making predicted values close to observed values according to squared residual distance. It is one of the most widely used fitting methods in mathematical biology.
Observed and predicted values
Suppose observations are
\[y_1,y_2,\ldots,y_n\]at times or experimental conditions
\[t_1,t_2,\ldots,t_n.\]Let the corresponding model predictions be
\[\mu_i(\boldsymbol\theta)=f(t_i;\boldsymbol\theta).\]Residuals
The residual for observation \(i\) is
\[\boxed{r_i(\boldsymbol\theta)=y_i-\mu_i(\boldsymbol\theta)}.\]If \(r_i>0\), the observation lies above the model prediction. If \(r_i<0\), it lies below the prediction.
A good model should not merely have small residuals; the residual pattern should also be compatible with the assumed observation model.
Residual sum of squares
Ordinary least squares uses
\[\boxed{S(\boldsymbol\theta)=\sum_{i=1}^{n}r_i(\boldsymbol\theta)^2}.\]The estimate is
\[\boxed{\hat{\boldsymbol\theta}=\arg\min_{\boldsymbol\theta}S(\boldsymbol\theta)}.\]Squaring prevents positive and negative residuals from cancelling and penalises large discrepancies more strongly than small ones.
Why not minimise the sum of residuals?
If we used
\[\sum_i r_i,\]positive and negative errors could cancel. A model with large errors above and below the data could therefore appear artificially good.
Squaring produces non-negative contributions:
\[r_i^2\ge0.\]Geometric interpretation
Collect the residuals into a vector
\[\mathbf r=(r_1,\ldots,r_n)^T.\]Then
\[\boxed{S(\boldsymbol\theta)=\mathbf r^T\mathbf r=\|\mathbf r\|_2^2}.\]Least squares therefore chooses the model prediction vector with the smallest Euclidean squared distance from the observed data vector.
A simple straight-line example
Suppose
\[y_i=a+bt_i+\varepsilon_i.\]The residual is
\[r_i=y_i-a-bt_i,\]and least squares minimises
\[S(a,b)=\sum_i(y_i-a-bt_i)^2.\]Because the parameters enter linearly, this problem has an analytic solution.
Linear least squares in matrix form
For the linear model
\[\mathbf y=X\boldsymbol\beta+\boldsymbol\varepsilon,\]ordinary least squares minimises
\[\|\mathbf y-X\boldsymbol\beta\|_2^2.\]If \(X\) has full column rank, the estimator satisfies the normal equations
\[X^TX\hat{\boldsymbol\beta}=X^T\mathbf y.\]Formally,
\[\boxed{\hat{\boldsymbol\beta}=(X^TX)^{-1}X^T\mathbf y}.\]In numerical computation, directly forming the inverse is usually unnecessary; QR or related stable linear-algebra methods are preferred.
Nonlinear least squares
Most mechanistic biological models are nonlinear in their parameters. For example, an ODE model may produce
\[\mu_i(\boldsymbol\theta)=h(\mathbf x(t_i;\boldsymbol\theta)).\]There is generally no closed-form expression for \(\hat{\boldsymbol\theta}\), so numerical optimisation is required.
Least squares with an ODE model
For each trial parameter vector, one typically solves
\[\frac{d\mathbf x}{dt}=\mathbf f(\mathbf x,\boldsymbol\theta),\]evaluates the solution at observation times, constructs residuals and computes
\[S(\boldsymbol\theta)=\sum_i\left[y_i-h(\mathbf x(t_i;\boldsymbol\theta))\right]^2.\]The optimiser then proposes another parameter vector.
Statistical interpretation
Suppose the observation model is
\[\boxed{Y_i=\mu_i(\boldsymbol\theta)+\varepsilon_i},\]with independent errors
\[\varepsilon_i\sim N(0,\sigma^2).\]Then
\[Y_i\mid\boldsymbol\theta\sim N(\mu_i(\boldsymbol\theta),\sigma^2).\]Under these assumptions, minimising the residual sum of squares is equivalent to maximising the likelihood with respect to \(\boldsymbol\theta\).
Why Gaussian errors lead to squares
The Gaussian density contributes a factor proportional to
\[\exp\left[-\frac{(y_i-\mu_i)^2}{2\sigma^2}\right].\]For independent observations, the log-likelihood contains
\[-\frac{1}{2\sigma^2}\sum_i(y_i-\mu_i)^2.\]Thus maximising the likelihood is equivalent to minimising
\[\sum_i(y_i-\mu_i)^2.\]Assumptions behind ordinary least squares
For its standard statistical interpretation, ordinary least squares commonly assumes that residual errors are approximately independent, centred at zero and have constant variance.
Normality is additionally required for the exact Gaussian likelihood interpretation and many classical finite-sample inferential results.
Homoscedasticity
Constant error variance means
\[\operatorname{Var}(Y_i\mid\boldsymbol\theta)=\sigma^2\]for all observations.
This is called homoscedasticity.
Many biological datasets instead have variance that increases with the mean.
Heteroscedasticity
If
\[\operatorname{Var}(Y_i\mid\boldsymbol\theta)=\sigma_i^2\]varies between observations, the errors are heteroscedastic.
Ordinary least squares then gives the same numerical importance to high-precision and low-precision observations.
Weighted least squares
When observation variances are known or modelled, use
\[\boxed{S_w(\boldsymbol\theta)=\sum_{i=1}^{n}w_i r_i^2}.\]If errors are independent Gaussian with known variances \(\sigma_i^2\), the natural weights are
\[\boxed{w_i=\frac{1}{\sigma_i^2}}.\]More precise observations therefore receive greater weight.
Standardised residuals
With known standard deviations, define
\[z_i=\frac{y_i-\mu_i}{\sigma_i}.\]Then weighted least squares becomes
\[S_w=\sum_i z_i^2.\]This shows that weighting measures discrepancies relative to their expected uncertainty.
Correlated residuals
Time-series, spatial and repeated-measures data may contain correlated errors. Let the residual covariance matrix be \(\Sigma\).
A generalised least-squares criterion is
\[\boxed{S_G(\boldsymbol\theta)=\mathbf r^T\Sigma^{-1}\mathbf r}.\]Ordinary least squares corresponds to the special case \(\Sigma=\sigma^2I\).
Why independence matters
If neighbouring residuals are positively correlated, treating them as independent exaggerates the amount of independent information in the dataset.
This can make parameter uncertainty appear smaller than it really is.
Relative-error models
When measurement error grows approximately in proportion to the magnitude of the observation, an additive constant-variance model may be inappropriate.
One possibility is a multiplicative model such as
\[Y_i=\mu_i e^{\varepsilon_i}.\]Taking logarithms gives
\[\log Y_i=\log\mu_i+\varepsilon_i.\]Least squares on the log scale therefore corresponds to a different observation model from least squares on the original scale.
Fitting multiple variables
Suppose a model is fitted simultaneously to susceptible and infectious measurements. If one variable is measured in millions and another in hundreds, unweighted squared errors can be dominated by the larger numerical scale.
Weights or an explicit multivariate observation model may be needed so that the objective reflects measurement uncertainty rather than arbitrary units.
Counts are not automatically Gaussian
Disease cases, cell counts and ecological counts are discrete. Their variance may depend strongly on their mean.
For such data, Poisson, negative-binomial or other likelihood models can be more appropriate than ordinary least squares.
The next lesson develops likelihood-based fitting.
Residual plots
After fitting, residuals should be inspected against time, fitted values and relevant covariates.
A reasonable residual pattern should not show unexplained systematic trends, changing spread or long runs of the same sign.
Systematic residual patterns
If residuals are repeatedly positive during one phase and negative during another, the problem may not be parameter values alone.
The model may be missing a biological mechanism, intervention change, delay or observation feature.
Residual mean
A residual mean close to zero is useful but insufficient. Positive and negative structured errors can average to zero while the model remains systematically wrong.
Residual structure matters more than a single summary statistic.
Root mean squared error
A common descriptive measure is
\[\boxed{\operatorname{RMSE}=\sqrt{\frac1n\sum_i r_i^2}}.\]RMSE has the same units as the observed variable.
It describes fit magnitude but does not by itself determine whether the statistical assumptions are valid.
Residual sum of squares always decreases with flexibility
Adding additional free parameters generally cannot increase the minimum residual sum of squares, because the larger model can often reproduce the smaller model as a special case.
Therefore a lower training RSS does not automatically mean a scientifically better model.
Outliers and squared loss
Because residuals are squared, large residuals receive strong influence:
\[10^2=100,\qquad1^2=1.\]This makes least squares sensitive to extreme observations.
Outliers should first be investigated scientifically rather than automatically deleted.
Robust alternatives
When heavy-tailed observation errors are plausible, alternative loss functions or probability models can reduce the influence of extreme residuals.
The choice should represent the expected observation process rather than merely force the model to fit more closely.
Numerical solver error
For mechanistic ODE models, residuals should primarily reflect mismatch between model and observations, not poor numerical integration.
Solver tolerances should therefore be sufficiently accurate relative to the measurement precision.
Parameter uncertainty
The location of the least-squares minimum gives a point estimate, but the shape of the objective surface around the minimum contains information about parameter uncertainty and correlation.
A flat direction indicates that substantially different parameter values produce similar residual sums of squares.
Least squares and identifiability
A very small minimum RSS does not imply that every parameter is identifiable.
Many parameter combinations can sometimes generate almost identical fitted curves. This issue is examined explicitly in the identifiability lesson.
When least squares is appropriate
Ordinary least squares is most natural when the measured response is continuous and additive errors have approximately constant variance and weak dependence.
Weighted or generalised least squares extends this framework when variances differ or errors are correlated.
When the data-generating distribution is intrinsically discrete or strongly non-Gaussian, an explicit likelihood may be preferable.
Transition to maximum likelihood
Least squares can be derived from one particular observation model: independent Gaussian errors with constant variance. The next lesson generalises this idea by specifying a probability distribution appropriate to the observed data and estimating parameters through likelihood.