Least Squares Curve Fitting

In engineering and science, you often collect data and need to find a mathematical model that best describes the underlying relationship. Least squares curve fitting is the foundational technique for this task, finding the model parameters that minimize the sum of the squared differences between the observed data and the model's predictions. This method provides an objective, optimal fit that is computationally tractable and forms the backbone of regression analysis, calibration, and predictive modeling across every engineering discipline, from signal processing to structural analysis.

The Foundational Idea: Minimizing Squared Residuals

At its core, least squares fitting is an optimization problem. You start with a set of $N$ data points $(x_i,y_i)$ and a proposed model—a function $f(x,\beta )$ that predicts $y$ based on $x$ and a set of unknown parameters $\beta$. The difference between the observed value $y_i$ and the predicted value $f(x_i,\beta )$ is called the residual, $r_i=y_i-f(x_i,\beta )$.

The goal is to find the parameter values $\beta$ that make the model predictions as close as possible to the data. The least squares criterion achieves this by minimizing the sum of squared residuals (SSR):

$$SSR(\beta )=\sum _{i=1}^N{r}_i^2=\sum _{i=1}^N[y_i-f(x_i,\beta ){]}^2$$

Squaring the residuals has two key effects: it ensures all errors are positive (so a negative residual doesn't cancel a positive one), and it penalizes larger errors more severely. This yields a unique solution for many common model forms and has desirable statistical properties if the errors in the data are random and normally distributed.

Linear Least Squares and the Normal Equations

The most common application is linear least squares, where the model is a linear combination of the parameters. A quintessential example is fitting a straight line, $y={\beta }_0+{\beta }_{1}x$. Here, $f(x,\beta )={\beta }_0+{\beta }_{1}x$, and the parameters ${\beta }_0$ (intercept) and ${\beta }_1$ (slope) are unknown.

To find the optimal ${\beta }_0$ and ${\beta }_1$, you set the partial derivatives of the SSR with respect to each parameter to zero. This process leads to a system of normal equations. For the line example, the normal equations are:

$$\begin{array}{rl}N{\beta }_0+\left(\sum x_i\right){\beta }_1& =\sum y_i\\ \left(\sum x_i\right){\beta }_0+\left(\sum x_i^2\right){\beta }_1& =\sum x_i{y}_i\end{array}$$

Solving this 2x2 system gives you the familiar formulas for the slope and intercept. This framework extends elegantly to any model linear in its parameters, such as a polynomial. Polynomial regression fits a model like $y={\beta }_0+{\beta }_{1}x+{\beta }_2{x}^2+...+{\beta }_m{x}^m$. Although the model is nonlinear in $x$, it is linear in the parameters ${\beta }_0,{\beta }_1,...,{\beta }_m$, so the same least squares principle applies.

You can express any linear least squares problem in matrix form. Define the design matrix $\mathbf{X}$, where each row corresponds to a data point and each column to a basis function (e.g., for a quadratic, $1$, $x$, $x^2$). Let $\mathbf{y}$ be the vector of observed $y$-values. The parameter vector $\beta$ that minimizes $SSR=||\mathbf{y}-\mathbf{X}\beta |{|}^2$ is given by the solution to the matrix normal equation:

$$({\mathbf{X}}^T\mathbf{X})\beta ={\mathbf{X}}^T\mathbf{y}$$

Solving for $\beta$ gives $\beta =({\mathbf{X}}^T\mathbf{X}{)}^{-1}{\mathbf{X}}^T\mathbf{y}$. This is a powerful result, as it allows you to fit complex linear models (including multivariable ones) using straightforward matrix algebra.

Analyzing Residuals for Model Adequacy

Finding the best-fit parameters is only half the job. You must then assess whether the model is adequate. Residual analysis is the primary tool for this. After fitting, calculate the residuals $r_i=y_i-{\hat{y}}_i$, where ${\hat{y}}_i$ is the fitted value.

Plotting residuals versus the independent variable $x$ (or versus the fitted values $\hat{y}$) is crucial. A good fit will show residuals randomly scattered around zero with constant variance (homoscedasticity). Patterns in this plot—like a curve, a funnel shape, or systematic streaks—indicate model inadequacy. A curved pattern suggests a higher-order term (e.g., a quadratic) is missing. A funnel shape indicates non-constant variance, violating a key least squares assumption.

You should also check a histogram or a normal probability plot of the residuals. The least squares method is most statistically efficient if the errors (and thus the residuals) are approximately normally distributed. Significant skew or outliers can suggest problems with the data or the model. Outliers, points with exceptionally large residuals, can disproportionately influence the fit; it's essential to investigate whether they are data errors or indicate a need for a different model.

Nonlinear Least Squares and the Gauss-Newton Method

Many engineering models are inherently nonlinear in their parameters. Examples include exponential decay ($y={\beta }_1{e}^{{\beta }_{2}x}$) or a sinusoidal model ($y={\beta }_1\sin ({\beta }_{2}x+{\beta }_3)$). Here, the sum of squared residuals landscape is not a simple parabola, and no closed-form solution like the normal equations exists.

Nonlinear least squares relies on iterative, numerical methods, with the Gauss-Newton method being a standard approach. The strategy is to approximate the nonlinear model with a linear one at each iteration. You start with an initial guess for the parameters ${\beta }^{(0)}$. The nonlinear model $f(x_i,\beta )$ is then linearized using a first-order Taylor series expansion around the current guess.

This linearization creates a design matrix $\mathbf{J}$, called the Jacobian, which contains the partial derivatives of the model with respect to each parameter, evaluated at the current guess and each data point. This turns the problem into a linear least squares problem for a small parameter update $\delta$:

$$\mathbf{J}\delta \approx \mathbf{y}-\mathbf{f}({\beta }^{(k)})$$

You solve this linear system for $\delta$ using methods for linear least squares, then update the parameters: ${\beta }^{(k+1)}={\beta }^{(k)}+\delta$. This process repeats until the parameter changes become negligibly small, indicating convergence to a local minimum of the SSR. The success of Gauss-Newton depends heavily on a good initial guess; poor guesses can lead to convergence to a suboptimal local minimum or failure to converge at all.

Common Pitfalls

1. Overfitting with High-Order Polynomials: It's tempting to keep adding polynomial terms to make the curve pass through more points. A model like a 10th-order polynomial can fit a 12-point dataset nearly perfectly (SSR near zero), but it will oscillate wildly between points and fail to predict new data. This is overfitting. Always prefer the simplest model that adequately describes the trend. Use residual analysis and consider the physical plausibility of the model.

1. Ignoring Residual Patterns: Declaring the job done after obtaining the parameter values is a major error. A low SSR does not automatically mean a good model. You must examine residual plots. Failing to detect a curved pattern means you might miss a systematic bias in your predictions, rendering your model unreliable for its intended use.

1. Misapplying Linear Least Squares to Nonlinear Problems: You cannot directly use the linear normal equations $({\mathbf{X}}^T\mathbf{X})\beta ={\mathbf{X}}^T\mathbf{y}$ for models nonlinear in their parameters (e.g., $y={\beta }_1/(x+{\beta }_2)$). Attempting to do so will yield incorrect, meaningless parameters. Recognize the structure of your model and apply nonlinear least squares methods like Gauss-Newton when required.

1. Forgetting to Assess the Conditioning of ${\mathbf{X}}^T\mathbf{X}$: In linear least squares, solving the normal equations requires inverting ${\mathbf{X}}^T\mathbf{X}$. If your basis functions (e.g., $1,x,x^2,...$) are highly correlated for your data, this matrix can be ill-conditioned (nearly singular). This makes the solution for $\beta$ extremely sensitive to tiny errors or noise in the data, leading to numerically unstable and unreliable parameters. Using polynomial bases with centered and scaled $x$-values can help mitigate this.

Summary

• Least squares fitting is the process of optimizing a model's parameters by minimizing the sum of squared residuals between the model predictions and the observed data.
• For models linear in their parameters (like straight lines or polynomials), the optimal parameters are found analytically by forming and solving the normal equations, often expressed in the compact matrix form $\beta =({\mathbf{X}}^T\mathbf{X}{)}^{-1}{\mathbf{X}}^T\mathbf{y}$.
• Residual analysis—plotting and examining the differences between data and the fit—is non-negotiable for diagnosing model adequacy, checking assumptions, and detecting outliers or missing model terms.
• Models that are nonlinear in their parameters require iterative numerical solutions, such as the Gauss-Newton method, which repeatedly linearizes the problem to converge on an optimal parameter set.
• Key practical dangers include overfitting complex models to limited data, ignoring diagnostic residual patterns, and failing to account for numerical ill-conditioning in the solution process.