0


pop.ss <- nls(MRDRSLT ~ SSlogis(TIME, phi1, phi2, phi3), data = testing, na.action = na.omit)


theta <- coef(pop.ss)  #extracting coefficients
plot(MRDRSLT ~ TIME, data = testing, main = "Logistic Growth Model", 
     xlab = "Time", ylab = "MRD")  # Census data
curve(theta[1]/(1 + exp(-(x - theta[2])/theta[3])), add = T, col = "green")  # Fitted model


summary(pop.ss)

Formula: MRDRSLT ~ SSlogis(TIME, phi1, phi2, phi3)

Parameters:
      Estimate Std. Error t value Pr(>|t|)   
phi1 9.618e-05  6.935e-06  13.869  0.00516 **
phi2 2.480e+02  1.512e+01  16.403  0.00370 **
phi3 2.896e+01  2.392e+01   1.211  0.34960   
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.197e-05 on 2 degrees of freedom

Number of iterations to convergence: 1 
Achieved convergence tolerance: 1.687e-06

I want to able to calculate the doubling time of my variable of interest. How do I pull this from the coefficients. https://rdrr.io/r/stats/SSlogis.html

Asym/(1+exp((xmid-input)/scal))

In expontial growth models it is easy to calculate as it is 1/rate.

t_doubling <- (1 / mu) * log10(2)
Bobby
  • 13
  • 3

2 Answers2

2

This might be better for CrossValidated. Doubling time is not an intrinsic characteristic of logistic curves, because it depends on the starting value (e.g. any logistic curve that starts above K/2, where K is the asymptotic value, will never double, because it can never increase beyond K ...). If the starting value is very small then the formula for the doubling time is the same as for the exponential curve (your formula is wrong: the doubling time is the natural log of 2/rate (log(2)/mu, not log10(2)/mu):

x0*exp(mu*T_dbl) = 2*x0
exp(mu*T_dbl) = 2
mu*T_dbl = log(2)
T_dbl = log(2)/mu

In the parameterization given above (which uses a scale parameter rather than a rate parameter) this would be log(2)*scal (but again note that it is only applicable (approximately) to logistic curves starting at a low density and not projected beyond a level of about Asym*0.3 [beyond which the logistic curve looks more and more linear rather than exponential])

Ben Bolker
  • 211,554
  • 25
  • 370
  • 453
2

In general, the doubling time is only constant for an exponential curve and for the logistic and other curves we can only calculate an instantaneous doubling time which varies from point to point as opposed to being a single value.

The situation is similar to that of the slope which is only constant for a line. For other (differentiable) curves we can only calculate the slope of the tangent at each point and that varies from point to point as opposed to being a single value.

In general the instantaneous doubling time at point t equals

1/(log2 f(t))'

i.e. the reciprocal of the derivative of log base 2 of f(t) with respect to t evaluated at t. Here f(t) is the logistic or other curve. Note that because log2(x) = log(x) / log(2) and because of the rule for evaluating the derivative of the log of a function the instantaneous doubling time equals the following where log is log base e and f is short for f(t) and f' is short for the derivative of f(t) with respect to t evaluated at t.

log(2) * f / f'

Logistic curve

In the case of the logistic it satisfies a well known differential equation which can help us calculate the doubling time.

f = SSlogis(Time, Asym, xmid, scal) # Asym / (1 + exp((xmid - x)/scal))
f' = (1/scal) * f * (1 - f / Asym)
doubling time = log(2) * f / f' = log(2) * scal / (1 - f/Asym)

In R consider the Chickwt.1 data example used in example(SSlogis). With that we can carry out the R calculations:

# example data
Chick.1 <- ChickWeight[ChickWeight$Chick == 1, ] # ChickWeight comes with R
Asym <- 368; xmid <- 14; scal <- 6
Time <- Chick.1$Time

f <- SSlogis(Time, Asym, xmid, scal)

fderiv <- (1/scal) * f * (1-f/Asym)
log(2) * f / fderiv # doubling time

# or equivalently
log(2) * scal / (1 - f/Asym) # doubling time

plot(log(2) * scal / (1 - f/Asym) ~ Time, type = "o", pch = 20, 
  ylab = "Instantaneous Doubling Time")

This is the instantaneous doubling time of f at each time point in the vector Time.

In particular, consider three situations:

  • when f is small relative to Asym the denominator of log(2) * scal / (1 - f/Asym) is approximately 1 so the doubling time is about log(2) * scal. This expression is the same as for the exponential (see next section).
  • at the inflection point of the logistic f = Asym/2 so the doubling time is 2 * log(2) * scal.
  • when f is near Asym the denominator of log(2) * scal / (1 - f/Asym) is close to 0 so the doubling time approaches infinity.

(continued after plot)

screenshot

Exponential curve

As another example, the definition and differential equation of the exponential curve are:

f = Asym * exp((xmid - Time)/scal)
f' = (1/scal) * f

so without even using R, given the cancellations in the numerator and denominator, we can calculate the doubling time as:

log(2) * f / f'
= log(2) * scal

Fitting

We can estimate the parameters of the logistic using nls and SSlogis for the Chick.1 data:

fm.nls <- nls(weight ~ SSlogis(Time, Asym, xmid, scal), Chick.1)
fm.nls  # displays coefficients
log(2) * coef(fm.nls)[["scal"]]  # doubling time in initial portion of curve

In the case of this data it actually appears exponential. Plotting log2(weight) vs. Time gives approximately a straight line if it is exponential and that is the case here. In that case the doubling time is the reciprocal of the slope of the line of best fit.

fm.lm <- lm(log2(weight) ~ Time, Chick.1)
1/coef(fm.lm)[2]  # doubling time

# seems like a straight line
plot(log2(weight) ~ Time, Chick.1)   
abline(fm.lm)

screenshot

G. Grothendieck
  • 254,981
  • 17
  • 203
  • 341