0
$\begingroup$

I'm attempting to fit a mixed effects model to some data with several predictors. However, one of the predictors has a distinctly non-linear relationship with the response variable. After some visualization and model comparison in the drc package, there's decent agreement using a Weibull function (specifically fct W2.3).

I now want to include this non-linear relationship in my lme model, along with the other predictor variables. I'm familiar with the I() for simple arithmetic transformations, but how do I implement a more complex, 3-parameter Weibull distribution for a single predictor variable using the nlme package?

Sample data (actual is ~1300 rows):

dput(df)
df <- structure(list(Ctot = c(198.7458618, 155.6201223, 59.87797928, 
105.5318751, 61.62373459, 38.13153444, 49.34262602, 73.83036104, 
68.04493868, 121.5958957, 102.4991762, 101.3973061, 40.00028914, 
76.83397851, 90.01915178, 79.44320778, 79.44320778, 109.9208446, 
78.50240298, 67.63455447), Dcm = c(22.987, 42.0878, 18.288, 30.734, 
23.1648, 15.875, 17.78, 18.288, 26.543, 25.019, 21.1328, 34.2138, 
21.59, 20.066, 25.781, 23.368, 23.368, 23.368, 26.797, 23.495
), slope_deg = c(24.22774532, 4.57392126, 26.56505118, 0, 5.710593137, 
9.648045316, 1.145762838, 0, 5.710593137, 1.718358002, 2.290610043, 
11.30993247, 1.718358002, 0, 24.22774532, 4.57392126, 4.57392126, 
1.145762838, 17.74467163, 11.30993247), elev_m = c(427, 631, 
1116, 636, 296, 233, 508, 738, 655, 645, 651, 807, 143, 102, 
811, 261, 261, 513, 838, 385), PLOT = c("21471", "54421", "51904", 
"25850", "1010", "44138", "11691", "18226", "20391", "29547", 
"20467", "32394", "12389", "3882", "16438", "52477", "52477", 
"19013", "15014", "21827")), row.names = c(NA, -20L), class = "data.frame")

Model

nlme1r <- lme(Ctot ~ I(0+(d-0)/exp(-exp(b(log(Dcm - e))))) + 
slope_deg + elev_m, 
random = ~1 | PLOT, data=df)
  • How do I get nlme to estimate parameters d, b and e? I've included start=c(b=1.1, d=250, e=75) but obviously this fails.

  • Is lme the wrong function and I should switch to nlme for this? Or something else entirely?

  • I've had success (model converges with reasonable parameter estimates) using poly(Dcm, 3). However, I'm not sure this method provides a good approximation of the relationship between Ctot and Dcm. Should I just be happy that poly() fits the data and run with it?

I wasn't sure if this belongs on CV or SO, apologies if it's misplaced.

$\endgroup$

Your Answer

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

Browse other questions tagged or ask your own question.