Modelos no lineales: crecimiento poblacional
Vamos a estudiar dos modelos no-lineales que se usan para describir
el crecimiento poblacional de organismos, y vamos a usar el ejemplo de
crecimiento bacteriano.
Función Gompertz
La forma general de la función Gompertz, es la siguiente: \[y = ae^{-bc^x}\]
Crecimiento Logístico
Esta es la forma general de la función de crecimiento logístico:
\[N_t = \frac{K}{1 + (\frac{K -
N_0}{N_0})e^{-rt}}\]
Ambas son no-lineales, y antes debemos encontrar valores iniciales de
los parámetros.
Datos
Utilizaremos los datos de crecimiento de un cultivo de Salmonella
sp. (ver proyecto de modelo no-lineal).
library(readxl)
gbac <- read.csv("bacteria19.csv")
Gráfica
library(ggplot2)
ggplot(gbac, aes(x = Time, y = Logc)) + geom_point(alpha=0.7) +
ylim(0,12)
Método 1: Parámetros del modelo Gompertz modificado
Una forma de la ecuación re-parametrizada (para usar parámetros más
relacionados al proceso de crecimiento) del modelo Gompertz es la
siguiente:
\[y = y_0 + (y_{max} -
y_0)*exp(-exp(k*(lag - x)/(y_{max} - y_0) + 1))\] Para utilizar
el comando nls más fácilmente, primero creamos una
función R para la ecuación Gompertz.
Gompertz <- function(x, y0, ymax, k, lag){
result <- y0 + (ymax -y0)*exp(-exp(k*(lag-x)/(ymax-y0) + 1) )
return(result)
}
Luego obtenemos unos valores iniciales para los parámetros, a partir
de la gráfico y/o los datos.
- y0 es el valor inicial de \(y\), cuando \(x =
0\), es decir el número de bacterias (o absorbancia)
inicial.
- ymax es el valor máximo de \(y\) (cuando \(x
\rightarrow \infty\)).
- lag es el tiempo que transcurre antes que
se inicie el crecimiento bacteriano (“fase lag”).
- k está relacionado con la velocidad de
crecimiento máxima; puede obtenerse como un estimado de la pendiente
máxima en la gráfica de puntos (el valor real usualmente es más
alto).
y0 = 4.5
ymax = 10.2
lag = 5
k = 0.42
Finalmente se calculan los parámetros del modelo ajustado usando
nls.
# modelo
Gomp1 <- nls(Logc ~ Gompertz(Time, y0, ymax, k, lag),
data = gbac,
start = list(y0=y0, ymax=ymax, k=k, lag=lag))
# parámetros y estadisticas
summary(Gomp1)
##
## Formula: Logc ~ Gompertz(Time, y0, ymax, k, lag)
##
## Parameters:
## Estimate Std. Error t value Pr(>|t|)
## y0 4.43112 0.07977 55.55 < 2e-16 ***
## ymax 10.35959 0.06175 167.78 < 2e-16 ***
## k 0.87926 0.03772 23.31 < 2e-16 ***
## lag 7.66629 0.50661 15.13 8.96e-14 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1736 on 24 degrees of freedom
##
## Number of iterations to convergence: 7
## Achieved convergence tolerance: 1.979e-06
Gráfica de puntos y del
modelo
Se grafican los puntos originales y valores de la predicción con la
función y los parámetros calculados.
time.levels <- seq(0, 100, by=0.1)
pred <- predict(Gomp1, list(x=time.levels) )
plot(Logc ~ Time, data=gbac, xlim=c(0,100), ylim=c(0, 12),
xlab="time", ylab="logCFU")
lines(gbac$Time, pred,
lty = 1, col = "red")
Referencia
Tjørve,
K.M.C., Tjørve, E., 2017. The use of Gompertz models in growth analyses,
and new Gompertz-model approach: An addition to the Unified-Richards
family. PLoS One 12, e0178691.
https://doi.org/10.1371/journal.pone.0178691
Método 2: Paquete especializado para curvas de crecimiento
bacteriano con modelo logístico
Cálculo
de los parámetros del modelo y del análisis de crecimiento
bacteriano.
Con este método no es necesario estimar valores iniciales, y se puede
aplicar a muchas curvas de crecimiento a la vez.
model.bac <- SummarizeGrowth(gbac$Time, gbac$Logc, bg_correct = "none")
model.bac$vals
## k k_se k_p n0 n0_se n0_p
## 10.534 0.171 7e-29 3.491 0.179 1e-16
##
## r r_se r_p sigma df t_mid
## 0.094 0.007 9e-13 0.437 25 7.497
##
## t_gen auc_l auc_e
## 7.406 922.81 915.104
## Nonlinear regression model
## model: n ~ k/(1 + ((k - n0)/n0) * exp(-r * t))
## data: d
## k n0 r
## 10.53397 3.49142 0.09359
## residual sum-of-squares: 4.769
##
## Number of iterations to convergence: 8
## Achieved convergence tolerance: 1.49e-08
Significado e
interpretación de los parámetros
A continuación la lista de parámetros poblacionales que se extraen
con el paquete growthcurver y que tienen un significado en el
análisis del crecimiento de microorganismos:
- k es la capacidad de acarreo (\(K\)) de la población en el medio y
condiciones de crecimiento. k_se es el error estándar y
k_p el nivel de significancia del parámetro (la probabilidad de
error al rechazar que el parámetro es igual a 0).
- n0 es el tamaño de la población inicial de
bacterias (\(N_0\)). n0_se es
el error estándar y n0_p el nivel de significancia del
parámetro (la probabilidad de error al rechazar que el parámetro es
igual a 0).
- r es la tasa de crecimiento poblacional
(\(r\)). r_se es el error
estándar y r_p el nivel de significancia del parámetro (la
probabilidad de error al rechazar que el parámetro es igual a 0). No
necesariamente igual a la k del modelo Gompertz.
- sigma es el error estándar
residual, una medida de la bondad de ajuste de los parámetros
del modelo a los datos (mejor cuanto más pequeño sea).
- df son los grados de libertad, calculados
como el número de datos (\(n\)) menos
el número de parámetros del modelo (\(n -
3\)).
- t_mid es el tiempo cuando la población
alcanza la mitad de la capacidad de acarreo (\(\frac{1}{2}K\)).
- t_gen es el tiempo de generación más
rápida o tiempo de duplicación de la población. Se puede calcular como
\(\frac{ln2}{r}\).
- auc_l es el área bajo la curva del modelo
logístico y sirve de medida integradora de \(K, r\) y \(N_0\).
- auc_e es el área bajo la curva empírica,
es decir usando los datos medidos del crecimiento bacteriano.
Gráfica de puntos y del
modelo
plot(model.bac, xlab="time", ylab="logCFU")
Referencias
Sprouffske,
K., Wagner, A. Growthcurver: an R package for obtaining interpretable
metrics from microbial growth curves. BMC Bioinformatics 17, 172 (2016).
https://doi.org/10.1186/s12859-016-1016-7.
Sprouffske,
K., 2020. Using Growthcurver. CRAN Repository