4
$\begingroup$

We -- my friend -- had asked a question here, which was not completely related to Mathematica, however, in a hope that someone familiar with the topic could help. There are still some aspects which we would like to know.

We have a dataset as:

data = {{0, 0.046}, {40, 0.111}, {80, 0.291}, {120, 0.808}, {160, 1.742}, {200, 3.319}, {240, 5.017}, {280, 5.503}, {320, 5.897}};

where the first components are time; and, with the corresponding standard deviations as:

sd = {0.003, 0.012, 0.023, 0.056, 0.083, 0.216, 0.526, 0.366, 0.313};

We would like to obtain the growth rate, which is the slope of line in the exponential phase of the growth.

To this end, our first approach is to employ nonlinear regression, that is, to fit the log of logistic function to the log of data, which in this case, the logistic rate approximates the slope of line in the exponential phase of growth. Thus, we have:

data2 = data;
data2[[All, 2]] = Log[data2[[All, 2]]];
nlm2 = NonlinearModelFit[data2, Log[a/(1 + E^(-k (t - t0)))], {{a, 1}, {t0, 50}, {k, 0.1}}, t];
nlm2["ParameterConfidenceIntervalTable"]

which $k$, i.e., the logistic rate, reads as: $k = 0.025$, with the confidence interval of $0.024 < k < 0.026$.

Now, we would like to see how to approximate the growth rate using linear regression. In this paper, it suggests that, first to fit a polynomial to the growth curve data, in order to obtain the maximum growth rate $t_{\text{max}}$, by differentiating the fitted polynomial ($t_{\text{max}}$ is where the maximum of the derivative occurs). Then, to fit $a t + b$ to the log of the growth curve data in the neighborhood of $t_{\text{max}}$.

Now, our question is: it seems Mathematica does not return sensible result if we want to fit a polynomial, let's say, degree $6$ to the data:

nlm3 = NonlinearModelFit[data, 
a t^6 + b t^5 + c t^4 + d t^3 + e t^2 + f t + g, {{a, 1}, {b, 1}, {c, 1}, {d, 1}, {e, 1}, {f, 1}, {g, 1}}, t];
nlm3["ParameterConfidenceIntervalTable"]

which $a$, $b$, $c$, $d$, and $e$ are effectively $0$.

Any help is appreciated how to move forward from here!

$\endgroup$
3
  • $\begingroup$ The most natural way of fitting a growth rate is to do a linear fit on the logarithms of the data points. $\endgroup$
    – Roman
    Nov 29, 2021 at 20:35
  • $\begingroup$ @Dave What are the organisms? What's wrong with your $k=0.025$ estimate? $\endgroup$
    – Chris K
    Nov 29, 2021 at 20:39
  • $\begingroup$ The equation you use ($\frac{a}{e^{-k (t-\text{t0})}+1}$) isn't linear "in the exponential phase" or pretty much anywhere between 0 and 350 time units so I don't know why you'd want to fit a model that you know from the start won't fit. Consider the plot of the fitted slope: slope = D[a/(1 + E^(-k (t - t0))), t] /. nlm2["BestFitParameters"]; Plot[slope, {t, 0, Max[data[[All, 1]]]}]. $\endgroup$
    – JimB
    Nov 30, 2021 at 4:23

2 Answers 2

2
$\begingroup$

First, plot your data ($n$ vs $t$) on a log scale:

ListLogPlot[data]

enter image description here

As Roman suggests, a linear fit on log scale is the growth rate you're interested in. It doesn't look like there's any lag phase at the beginning (the slope doesn't increase as the cells get ready to grow - maybe they were already preconditioned?). However we see the declining growth rate with density (density-dependence) that we want to avoid.

I've seen people just do a linear fit for the first data points that look linear (judged by eye). I'd say use the first four points.

ldata = Transpose[{data[[All, 1]], Log[data[[All, 2]]]}]
LinearModelFit[ldata[[1 ;; 4]], {1, t}, t]

enter image description here

$k=0.0239$ seems close enough to your other estimate. Of course you can vary the number of points included to test the robustness of this estimate. There's probably a better way to do it (@JimB might know) -- maybe that's your nonlinear fit :)

$\endgroup$
4
  • $\begingroup$ @Dave While the referee you mentioned might be a bit rigid, "judging by eye" or the "I'll know it when I see it" method doesn't seem appropriate for a scientific journal (unless it's in the speculation section - which I think should be a legitimate section). $\endgroup$
    – JimB
    Nov 29, 2021 at 23:21
  • $\begingroup$ @JimB That's how you can tell I'm a theoretical ecologist, not a statistical ecologist! Anyhow, the nonlinear approach is probably better (although it assumes that growth is logistic, which it looks like here but isn't necessarily the case). What strikes me is that all estimates are so similar, so that it doesn't really matter. $\endgroup$
    – Chris K
    Nov 29, 2021 at 23:45
  • $\begingroup$ This isn't the site to be pushing R, but check out this package and also this article if you need a reference. $\endgroup$
    – Chris K
    Nov 29, 2021 at 23:51
  • $\begingroup$ @ChrisK I already knew you couldn't be a statistical ecologist because I've never seen you use the religious terms AIC and AICc. $\endgroup$
    – JimB
    Nov 30, 2021 at 3:36
3
$\begingroup$

Note, a polynomial fit is a linear fit. Therefore, you need a different function "LinearModelFit":

data = {{0, 0.046}, {40, 0.111}, {80, 0.291}, {120, 0.808}, {160, 
    1.742}, {200, 3.319}, {240, 5.017}, {280, 5.503}, {320, 5.897}};
fit = LinearModelFit[data, Table[t^i, {i, 0, 6}], t]
Plot[fit[t], {t, 0, 320}, Epilog -> Point[data]]

enter image description here

$\endgroup$
1
  • 1
    $\begingroup$ <<1>> means that in the display there is a term not displayed to save space. You may display the full poly by e.g. writing fit[x] $\endgroup$ Nov 30, 2021 at 9:00

Your Answer

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