This is my first question on here, so please bear with me (and please let me know how to improve my questions in the future).

I have a dataset with two independent (continuous) variables: temperature (1, 3, 6 degrees) and nutrient concentration (~ 9-1400) and one dependent variable (size), with replicates. I'd like to analyze the interactive effects of nutrients and temperature on size.

I was able to get a non-linear fit for size vs nutrient at 1, 3 and 6 degree separately (using nlMS package), but I need to do a non-linear mixed model to get the interactive effects of temperature and nutrients on size. I have no idea where to begin, or what package to use. Any help would be appreciated!

temp <- c(1,1,1,3,3,3,6,6,6,6,6,1,1,1,3,3,3,6,6,6,6,6,6,6,1,1,1,3,3,3,6,6,6,6,6,6,6,1,1,1,3,3,3,6,6,6,6,6,6,6,1,1,1,3,3,3,6,6,6,6,6,6,6,6,1,3,3,3,6,6,6,6,6,6,1,1,1,3,3,3,6,6,6,6,6,6,1,1,1,3,3,3,6,6,6,6,6,6,6)

nutrient<- c(41.922282,41.922282,41.922282,37.23794,37.23794,37.23794,31.662541,31.662541,31.662541,31.662541,31.662541,279.720746,279.720746,279.720746,248.465109,248.465109,248.465109,211.264016,211.264016,211.264016,211.264016,211.264016,211.264016,211.264016,27.946784,27.946784,27.946784,24.824046,24.824046,24.824046,21.1073,21.1073,21.1073,21.1073,21.1073,21.1073,21.1073,55.899183,55.899183,55.899183,49.65308,49.65308,49.65308,42.218842,42.218842,42.218842,42.218842,42.218842,42.218842,42.218842,97.838316,97.838316,97.838316,86.905988,86.905988,86.905988,73.89411,73.89411,73.89411,73.89411,73.89411,73.89411,73.89411,73.89411,10.479381,9.308428,9.308428,9.308428,7.914737,7.914737,7.914737,7.914737,7.914737,7.914737,1404.251692,1404.251692,1404.251692,1247.342413,1247.342413,1247.342413,1060.585803,1060.585803,1060.585803,1060.585803,1060.585803,1060.585803,139.790093,139.790093,139.790093,124.170127,124.170127,124.170127,105.578928,105.578928,105.578928,105.578928,105.578928,105.578928,105.578928)

size <- c (0.089477044,0.07255548,0.08559649,0.148534783,0.152138103,0.147775201,0.150692377,0.145413964,0.150045565,0.139712503,0.157196943,0.156799728,0.138026014,0.148072372,0.243282066,0.253498973,0.236256636,0.310144261,0.281107498,0.318315753,0.34032368,0.28983498,0.346447123,0.290460835,0.023320566,0.009979184,0.033718296,0.090433546,0.091796231,0.096347009,0.157691429,0.112364272,0.155331758,0.163524305,0.154034127,0.12710671,0.169722436,0.065031661,0.086993854,0.058166144,0.164571185,0.168336541,0.169147685,0.22096527,0.230259712,0.216570836,0.246243675,0.220842629,0.209989795,0.223320508,0.118446307,0.143604328,0.109049472,0.222921061,0.225832377,0.221248684,0.197760991,0.366515371,0.367449366,0.21394256,0.282823072,0.364060778,0.285209082,0.198528104,0.025291159,0.073610059,0.060269992,0.0557738,0.049028256,0.088362402,0.047770389,0.057592633,0.094855429,0.060442358,0.154164802,0.165561525,0.164140252,0.243798695,0.242919888,0.256692927,0.309521766,0.323766393,0.294349209,0.31684002,0.316469317,0.317274182,0.121243793,0.142287948,0.11136973,0.219106574,0.227289561,0.222398572,0.32907817,0.287993773,0.297453928,0.315788911,0.318406081,0.328165403,0.29831968)
  • 1
    $\begingroup$ A linear or non-linear model has nothing to do with the underlying relationship you use to connect your variables. It is about the (error) functions behind the models. Or why do you have to use nonlinear models? The most present packages are 'lme4' and "'nlme'. $\endgroup$
    – Ben
    Commented Jan 23, 2020 at 6:11

2 Answers 2


A good first approach to this would be to use a linear mixed effects model, with random intercepts for replicate (provided that you have enough replicates). Using the usual notation adopted in lme4 it would look like:

size ~ temp * nutrient + (1| replicate)

You can allow for non-linear associations by adding non-linear terms (eg quadratic) or splines. If you really need to model non-linearity that cannot be accomplished like this, then transformations such as log or square root can be used, and if this is not suficient then nlme can fit non linear mixed effects models.

  • $\begingroup$ Thank you for the reply. I have 3 replicates at each treatment. I did the following ``` linear <- lmer(size~ nutrient*temp +(1|replicate)) log <- lmer (size~ log(nutrient+)*temp +(1|replicate)) sqrt <- lmer (size~ sqrt(nutrient)*temp +(1|replicate)) squared <- lmer (size^2 ~(nutrient)*temp +(1|replicate)) anova(linear, log, sqrt, squared) #this showed that the squared transformation had the lowest AIC. #I'm just not sure how to interpret the summary output from the squared model ``` $\endgroup$
    – Loay Jabre
    Commented Jan 23, 2020 at 11:42
  • $\begingroup$ @LoayJabre You are welcome. However, 3 replicates is not really sufficient for a random intercept. I would model it as a fixed effect instead. As for interpretation, that is a different question. Please post a new question about that, and make sure to include the output from the model. And if this answers this question please consider marking it as the accepted answer. $\endgroup$ Commented Jan 23, 2020 at 12:16

As a first cut at modeling the data I found that "size = a * pow(temp, b) + c * pow(nutrient, d)" is a simple equation that gave an OK fit with parameters a = 2.0951264467116359E-01, b = 3.0824180572670512E-01, c = -7.6822508590637906E-01, and d = -4.6064198833929815E-01 yielding R-squared = 0.865 and RMSE = 0.0345.




(here "absolute error" means "not relative error")


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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