# How does one perform multiple non-linear regression?

I performed an experiment where I took the heights of plants and measured a number of environmental conditions (air temp, soil temp, lux, air humidity, soil pH, wind) for each of those plants. I want to then determine the effect that these conditions had on plant height. I have used multiple linear regression many times before in R using the form of:

LinearModelTemp<-lm(Height~ AirTemp+SoilTemp+lux+..., data=data)


but I would not expect the response of plants to all these conditions to be linear. How do you perform multiple non-linear regression?

So far the options I have found are non-linear least squares and segmented linear regression. For non-linear least squares I would have to set the parameters of the curve and I have no prior ideas for what these are. Furthermore, I am not aware of being able to perform multiple regression using this format. The other option is segmented linear regression, but again I would have to choose where my breakpoints are. I have some idea for where I would expect these to be but surely choosing these breakpoints will have a significant impact on form of the final regression model.

Is there a standard protocol I can follow when performing this or are there other options that I have not considered?

• If you know the curve then non-linear regression finds the optimal parameters. You only need to set the initial values for optimisation method. – mpiktas Aug 12 '15 at 10:42
• Thank you for your help. The problem is I am doing multiple regression as the parameters were non-constant for each plant. Therefore I don't know how I would know what the curve would be as I can't graph it to get an idea. – PaulBarr Aug 12 '15 at 10:46
• See Distinction between linear and nonlinear model & Why is polynomial regression considered a special case of multiple linear regression?. It's not clear that allowing for curvilinear relationships between predictors & the response in a linear model isn't what you want. – Scortchi - Reinstate Monica Aug 12 '15 at 11:25
• Thankyou. It might well be that allowing curvilinear relationships is ok for what I am trying to achieve. I have little experience in this and therefore want to ensure that what I do is theoretically correct. – PaulBarr Aug 12 '15 at 11:50

There are a few methods that do what you want, which is to allow functional forms to be flexible. Probably the best one for your case here however is the additive model (or generalized additive model if your response isn't continuous).

The AM has the form $$y = \alpha + X'\beta + \displaystyle\sum_m f_m(Z_m) + \epsilon$$

where $X$ are variables that you represent linearly as in a parametric model, and the $Z$ are terms that you allow to be flexible.

AMs don't remove all functional form assumptions. You still have that your terms have an additive effect on the response, which might not be true (interactions can help with some violations of this). And when they are fit the usual way, they assume that the functions $f_m$ are smooth (so, like, no sharp thresholds).

They are usually fit with a method similar to what you call "segmented linear regression", but modified. One uses way, way more breakpoints than one thinks one needs, and then lets the data tell them how to "penalize" the coefficients on those breakpoints, so that the function doesn't get too wiggly (or too straight).

Specifically, if the response is continuous, the parameters ($\theta$ $\equiv$ $\beta$ along with the splines that make up each of the $f_z$) can be estimated by $$\theta = \left(W'W+\lambda D\right)^{-1}W'y$$ where $W$ is $X$ along with the data gotten from the basis expansion (the segmenting) of $Z$. $D$ is a square matrix of dimension $p$ (the number of variables including the spline terms). In a really simple AM, it is zero on al the off-diagonals, and zero for the $X$'s and $Z$'s, but one for all the spline terms. $\lambda$ is a vector with different values in the entries associated with the sets of splined terms for each $Z$. When the entries of $\lambda$ get bigger, the splined terms shrink, and the function $f(Z)$ reduces to the main effect $Z$ (because in this simple example there is no penalty on $Z$ by the way $D$ is constructed).

Choosing $\lambda$ from data is a little more involved, as is much else about AMs and GAMs. Good books include those by Ruppert, Wand and Carrol ("Semiparametric Regression" -- dated, but easy to read), Wood (introduction to GAMs with R), and the chapter in Hastie & Tibshirani's Intro to Statistical Learning.

While I think there are a few R packages that do GAMs, I've had good luck with gam in mgcv.

Try the following as an example of the kind of output you'd get:

library(mgcv)
m = gam(Sepal.Width~s(Sepal.Length)+s(Petal.Length)+s(Petal.Width)+Species,data=iris)
summary(m)
plot(m,pages=1)


All the above said: GAMs can be kind of a rabbit hole. If you don't have the time and inclination to learn all this stuff, see if your needs aren't served by including a few polynomial terms.

• Thank you for your comprehensive answer. With this form of smoothing would it allow for a function to increase and then decrease, rather than increase at a decreasing rate? I.e. with temperature I would expect there to be an optimum temperature for plant height and then the height to decrease either side of this optimum. – PaulBarr Aug 12 '15 at 11:48
• In general, yes. In practice, such things can be hard to find sometimes with small datasets and weaker signals. But a polynomial will also find something that increases and then decreases. – generic_user Aug 12 '15 at 11:49
• Part of this research is to indicate that the dataset does need to be much larger and comprehensive to tease out these responses so its not too much of a problem if the signals are missed out. This may be a whole new question but when you say include polynomial terms how would I do that? Thank you for your help, I have found there to be sparse information on this topic at a level that I can understand - I am a biologist not a statistician. – PaulBarr Aug 12 '15 at 11:52
• en.wikipedia.org/wiki/Polynomial_regression m = lm(y~poly(x,2,raw=T)) is the same thing as x2 = x^2 and then lm(y~x+x2) – generic_user Aug 12 '15 at 11:54
• @PaulBarr: See Frank's Harrell's Regression Modelling Strategies site for a useful approach: estimating the total degrees of freedom you can use in the regression without overfitting, & divvying them up between predictors based on subject-matter knowledge or strength of relationship with the response. – Scortchi - Reinstate Monica Aug 12 '15 at 12:13

The answer from @generic_user should be very helpful. I just want to suggest ways to minimize the number of variables for which you include truly non-linear terms, and to discuss a bit when to use a nonlinear model function versus polynomial regression.

The relations between the dependent variable and each of the independent variables don't have to be exactly linear for linear regression to work. Changes in the dependent variable with changes in each of the independent variables just have to be well enough represented by linear relations over the range of interest. This is like doing a Taylor expansion of a function and only keeping terms up to first order: sometimes a linear approximation is good enough.

You can do linear regression with non-linearly transformed variables. For example, you already list $pH$ as an independent variable, which can be used directly in linear regression even though it is related logarithmically to the underlying hydrogen-ion activity. Such a non-linear transformation doesn't rule out a linear regression, as the "linear" part of the name only refers to the relation (in this case) between changes in Height and changes in $pH$, regardless of whether you consider $pH$ to have been obtained by a non-linear transformation. That's why regression against polynomials, as suggested by some here, is also linear regression. So for variables which you expect to have a monotonic relation to Height (over the range of your data), an appropriate transformation of an independent variable will typically allow linear regression.

Variables with a peak or a trough in the relation to Height (as you suggest may be the case for temperature) require additional attention, if such a peak or trough occurs within the range of your data. Including polynomial terms in your linear regression can capture such relations (analogous to higher-order terms in a Taylor expansion), but they will be hard to interpret in terms of biological processes. I suspect that in your field there are some general formulas for relations of plant height to temperature, which include as parameters the temperature for peak height and how quickly height drops off with temperature on either side of the peak. If you find such formulas, your analysis can be set up in a way that identifies the peak and drop-off parameter values that best match your data. That may be easier to interpret biologically than a generic polynomial regression.