I've been looking into linear regression, and on the wikipedia page it says:

"In contrast, non-linear least squares problems generally must be solved by an iterative procedure"

This got me thinking more about OLS, and the differences between it and non-linear regression methods. Mores specifically, why equations that are non-linear in their parameters can't also be solved using the OLS assumption that $y=\beta x$ where $\beta =(X^TX)^{-1}X^Ty $.

So i guess my question is:

What is it about the process of solving OLS that requires the parameters to be linear? What would happen if they were non-linear and we tried to solve using OLS?

  • 2
    $\begingroup$ You can minimize $L(\beta)$ by solving $\nabla L(\beta) = 0$. Usually this is a nonlinear system of equations, which makes it difficult to solve. But if $L(\beta) = (1/2) \| X \beta - y \|^2$, then $\nabla L(\beta) = X^T (X \beta - y)$, so $\nabla L(\beta) = 0$ is equivalent to $X^T X \beta = X^T y$. This is a linear system of equations that can be solved with Gaussian elimination! That is the special thing that makes linear regression easy. $\endgroup$ – littleO Jun 1 at 5:48

What is it about the process of solving OLS that requires the parameters to be linear?

Because equations which are nonlinear in their parameters can't be written as $y=X\beta$. OLS estimates $\beta$ in the equation

$$ y = X\beta +\epsilon. $$

This is a linear relationship, so when we say that $\hat{\beta} = (X^\top X)^{-1}X^\top y$ is the optimal estimator of $\beta$, what we mean is that it's optimal in the sense that it minimizes $\|y - X\beta\|_2^2$. Minimizing $\|y - X\beta\|_2^2$ is only important if this objective is meaningful for your task; particularly, if the task isn't linear in these parameters, then the fit may be poor.

However, one reason that OLS is so flexible is that if you can find a way to represent your data in a linear way, then it is linear in the parameters, otherwise known as basis expansion.

A textbook example of a change of basis is using a polynomial basis, so you have $X_\text{polynomial} = [1, x, x^2, x^3, \dots, x^p]$. The model $X_\text{polynomial}\beta$ is linear in its parameters, but viewed as a function of $x$, it's a nonlinear polynomial.

What would happen if they were non-linear and we tried to solve using OLS?

It won't work very well!

This data's deterministic component is given by $$ y = \beta_0 + \beta_1 \sin (\beta_2 x + \beta_3) $$ which is not linear in $\beta$, the parameter vector to be estimated, because you can't write this in the form $y=X\beta$. I also add small, independent 0-mean Gaussian noise to each observation.

If we do the naive thing and assume that our output $y$ is a linear function of $x$, then we find a poor fit, in the sense that there is a large discrepancy between the estimated line (red) and the true function (blue). The model finds that the best linear approximation is a decreasing line, completely ignoring the sinusoidal behavior.

poor linear fit

One way to try to improve the fit is to re-express $x$. Since this looks like something sinusoidal, we might try a sine function. This gives the design matrix $X_\text{sine}=[1, \sin(x)]$. This give a flatter line, but it's still not a satisfying model. Even though the model and the desired function are both sine waves, we're implicitly using $\beta_0 + \beta_1 \sin(1 \times x + 0)$ to approximate $$ y = \beta_0 + \beta_1 \sin (\beta_2 x + \beta_3).$$ This is not a good approximation, because we've fixed $\beta_2=1$ and $\beta_3=0$, so the further the true values are from these assumed values, the poorer this approximation will be.

poor sine fit

What we really need is a way to recover all of the parameters in the function $$ y = \beta_0 + \beta_1 \sin (\beta_2 x + \beta_3), $$ but this is a nonlinear estimation task, so we need to use the appropriate tools to accommodate the nonlinearity of the $\beta$s. Nonlinear least squares is one method to achieve this, among many others.


N <- 1000

x <- runif(N, -pi, pi)
f <- function(x) pi + 2 * sin(4 * x) 
y <- f(x) + rnorm(N,sd=0.5)

model <- lm(y ~ x)

abline(model, col="red", lwd=2, lty="dashed")
lines(sort(x), f(sort(x)), lwd=2, col="blue")

model2 <- lm(y ~ sin(x) )

abline(model2, col="red", lwd=2, lty="dashed")
lines(sort(x), f(sort(x)), lwd=2, col="blue")
| cite | improve this answer | |
  • 1
    $\begingroup$ In this particular case, I’d change it to b1 sin (b3x) + b2 cos (b3 x), which gives exactly the same numbers but has one Non-linear parameter only. $\endgroup$ – gnasher729 Jun 1 at 14:29
  • $\begingroup$ That's a useful observation for a person carrying out NLLS on this toy problem, but it doesn't do much to explain why OLS is a poor substitute for NLLS. $\endgroup$ – Sycorax Jun 1 at 17:25

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.