2
$\begingroup$

We have time series of plant growth (length vs time) and we can reasonably model these with two approaches:

  1. With a nonlinear model length(time)=length0*exp(r*time) where length0 is lenght at initial time and r0 is a good initial guess for r (relative growth rate), like in Reference: https://doi.org/10.1111/j.2041-210X.2011.00155.x
    This approach uses length as a proxy for biomass and models the growth of the plant as a whole through time.
length0 = 0.8532
r0 = 0.01
model1 <- nls(length~length0*exp(r*time),start = list(r=r0))
  1. With a piecewise linear model with the package 'segmented', where psi is a guess for the breakpoint. This approach assumes constant absolute growth rate per stage and two growth stages.
library(segmented)
model2 <- segmented(lm(length~time), seg = ~time, psi=200)

Plot: Red - nonlinear, Blue - piecewise linearenter image description here

Is it rigorous to select one of this two models with AIC (lower is better)?

AIC(model1)
AIC(model2)

If not with AIC, how?

The data of this example is:

time <-c(1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34,35,36,37,38,39,40,41,42,43,44,45,46,47,48,49,50,51,52,53,54,55,56,57,58,59,60,61,62,63,64,65,66,67,68,69,70,71,72,73,74,75,76,77,78,79,80,81,82,83,84,85,86,87,88,89,90,91,92,93,94,95,96,97,98,99,100,101,102,103,104,105,106,107,108,109,110,111,112,113,114,115,116,117,118,119,120,121,122,123,124,125,126,127,128,129,130,131,132,133,134,135,136,137,138,139,140,141,142,143,144,145,146,147,148,149,150,151,152,153,154,155,156,157,158,159,160,161,162,163,164,165,166,167,168,169,170,171,172,173,174,175,176,177,178,179,180,181,182,183,184,185,186,187,188,189,190,191,192,193,194,195,196,197,198,199,200,201,202,203,204,205,206,207,208,209,210,211,212,213,214,215,216,217,218,219,220,221,222,223,224,225,226,227,228,229,230,231,232,233,234,235,236,237,238,239,240,241,242,243,244,245,246,247,248,249,250,251,252,253,254,255,256,257,258,259,260,261,262,263,264,265,266,267,268,269,270,271,272,273,274,275,276,277,278,279,280,281,282,283,284,285,286,287,288,289,290,291,292,293,294,295,296,297,298,299,300,301,302,303,304,305,306,307,308,309,310,311,312)
length <-c(0.8532,0.8532,0.8532,0.8532,0.8532,0.8532,0.8532,0.8532,0.8916,0.9,0.9,0.9,0.9,0.9,0.9,0.9,0.9372,0.9,0.9,0.9492,0.9492,0.9492,0.9492,0.9492,0.9984,0.9984,0.984,0.9984,0.984,0.9984,0.9984,0.9984,1.05,1.05,1.05,1.05,1.0812,1.0812,1.0812,1.0812,1.05,1.05,1.1016,1.1016,1.1016,1.1016,1.1016,1.1316,1.1316,1.1832,1.1832,1.1832,1.1832,1.1832,1.1544,1.1544,1.1544,1.2072,1.2072,1.2072,1.2072,1.236,1.2888,1.2888,1.2888,1.2888,1.236,1.2612,1.2612,1.2612,1.2612,1.2612,1.3164,1.3164,1.3416,1.3416,1.3416,1.3416,1.3956,1.3416,1.3704,1.3704,1.3704,1.3704,1.3704,1.4256,1.4256,1.4256,1.4256,1.4496,1.4496,1.4496,1.4496,1.482,1.482,1.482,1.482,1.482,1.482,1.482,1.5384,1.5384,1.5384,1.5384,1.5384,1.5384,1.5384,1.5384,1.5948,1.5948,1.5948,1.5948,1.5948,1.5948,1.5948,1.5948,1.6512,1.6512,1.6512,1.6512,1.6512,1.6512,1.6512,1.6512,1.6512,1.7076,1.7076,1.7076,1.7076,1.7076,1.7076,1.7076,1.7076,1.7076,1.7652,1.728,1.728,1.8408,1.7844,1.7844,1.7844,1.7652,1.8216,1.8216,1.8216,1.8792,1.8792,1.8792,1.8792,1.8792,1.9212,1.9368,1.9368,1.9788,1.9944,1.9788,1.9944,1.9944,2.052,2.052,2.1108,2.1108,2.1108,2.1684,2.1684,2.1684,2.1684,2.226,2.226,2.226,2.2848,2.2848,2.2848,2.2848,2.3436,2.3436,2.3436,2.4012,2.4012,2.4012,2.46,2.46,2.46,2.46,2.5188,2.5188,2.5188,2.5776,2.5776,2.5776,2.6244,2.6364,2.6832,2.6832,2.6832,2.6952,2.742,2.742,2.802,2.802,2.802,2.8608,2.8608,2.9196,2.9196,2.9196,2.9784,2.9784,2.9784,3.0384,3.0384,3.0972,3.0972,3.0972,3.0972,3.1572,3.1572,3.1572,3.216,3.216,3.2748,3.2676,3.2748,3.2748,3.3348,3.3264,3.3864,3.3864,3.4452,3.4452,3.5052,3.5052,3.5052,3.5652,3.624,3.624,3.624,3.684,3.6912,3.6912,3.744,3.8028,3.8028,3.87,3.87,3.8628,3.9228,3.9228,3.9828,3.9888,4.0416,4.056,4.1016,4.1076,4.1616,4.1676,4.2204,4.2276,4.2864,4.2936,4.3464,4.3464,4.4388,4.4064,4.4664,4.4664,4.4724,4.5324,4.5672,4.5996,4.626,4.7352,4.7352,4.7184,4.7784,4.794,4.8456,4.8372,4.8972,4.9728,5.016,5.0232,5.0832,5.142,5.0916,5.1348,5.1948,5.202,5.2692,5.2548,5.3292,5.3136,5.3736,5.3736,5.4336,5.4336,5.4996,5.4336,5.5068,5.4924,5.5752,5.6856,5.694,5.6184,5.7384,5.7972,5.7972,5.7972,5.8572,5.8716,5.9316,5.9172)
$\endgroup$
5
  • 2
    $\begingroup$ IMHO, you use the nonlinear model because you have some scientific theory to support that particular model. Ordinary linear regression is "just" a Taylor series type approximation to some "scientific" model. What makes your question so interesting is that the breakpoint is almost a competing model, although you could have used an approach to help you find the breakpoint (rather than guess at 200). This is just a personal view, but I think I'd choose based on the science rather than the model fit. I would however do a lot of residual checking....! $\endgroup$ Commented Nov 25, 2019 at 14:51
  • $\begingroup$ @Paul The one reason we are interested in the piecewise model is that 'it is possible' that two growth phases exist and the exponential model does not give us any information about them, while the piecewise model gives us two slopes (growth rates) and one breakpoint (growth phase shift). Regarding the initial guess, we can consider an unbiased approach instead of a first guess. We will check the residuals. Thanks very much for your comment. $\endgroup$
    – GeorgeP
    Commented Nov 25, 2019 at 15:20
  • $\begingroup$ @PaulHewson I have frequently estimate precisely that parameter ($\theta$ in the below) using NLLS regression with the hinge function as follows: $y_{i} = \beta_{0} + \beta_{x}x_{i} + \beta_{\text{c}}\max (x_{i}-\theta,0)_{i} + \varepsilon_{i}$. This models a linear regression using the first two terms, and then $\beta_{\text{c}}$ change in the effect of $x$ at the value $\theta$. The $\max()$ variable takes the value 0 below $\theta$, and increases by 1 $x$-unit thereafter. $\endgroup$
    – Alexis
    Commented Nov 26, 2019 at 3:36
  • $\begingroup$ @GeorgeP Do you have a test set, which you can use to make predictions? $\endgroup$ Commented Nov 26, 2019 at 9:16
  • $\begingroup$ @user2974951 sorry no, we only have experimental data from one experiment and all are very similar to the example above. $\endgroup$
    – GeorgeP
    Commented Nov 27, 2019 at 18:13

1 Answer 1

2
$\begingroup$

The issue with the approach you described is that the breakpoint is a parameter, which appears to have been chosen by hand to fit the data. But, it's treated as a known value when calculating the AIC. This fails to account for the uncertainty in estimating the breakpoint from the data. As a result, the estimated AIC will be overoptimistically biased, and unfairly favor this model more than it should.

If the breakpoint truly is a known parameter (i.e. constrained by theory), then disregard the above. Otherwise, the proper approach is to estimate the breakpoint as part of the fitting algorithm (e.g. see the comment from @alexis above). And, then account for this during model selection (whether AIC or otherwise).

AIC and its cousin BIC are certainly established model selection criteria, but their validity rests on asymptotic assumptions. It's not particularly clear when these assumptions are or aren't justified for real world datasets (especially small ones). Their main advantage lies in their computational efficiency, but that's not really a limiting factor in your case. In general, I have greater trust in cross validation.

The fact that you have time series data raises another important point. Standard model selection procedures are often designed for i.i.d. data, and can lose their validity when dependencies are present (as is often the case for time series, spatial data, etc.). You must be very careful to perform model selection properly in this setting. For example, there are specialized forms of cross validation designed for time series and other dependent data. This can also be true for statistical inference (e.g. p-values and confidence intervals on model parameters).

$\endgroup$
3
  • $\begingroup$ Please see my comment to @PaulHewson above for a nonlinear least squares approach to precisely the estimation you mention, $\endgroup$
    – Alexis
    Commented Nov 26, 2019 at 3:38
  • $\begingroup$ Thanks for the reply @user20160. Regarding the issue of the break-point being a parameter: it's my understanding that the number supplied to the function 'segmented' in the form of 'psi = 200' it's just a starting value/best guess and it's also subjected to fitting. It's a problem for us that this break-point it's not supported nor undermined by theory. $\endgroup$
    – GeorgeP
    Commented Nov 27, 2019 at 18:07
  • $\begingroup$ @user20160 Regarding the issue of using AIC as selection criteria for time series: I can only say that it is common practice as a selection criteria between nonlinear models of plant growth, although I can't judge if it is rigorous. $\endgroup$
    – GeorgeP
    Commented Nov 27, 2019 at 18:18

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.