--- title: "Exponential Population Growth" subtitle: "From the lake's water lilies to the differential equation" author: "Population Ecology Laboratory — Department of Biology" date: "`r format(Sys.Date(), '%d %B %Y')`" output: html_document: toc: true toc_float: true toc_depth: 3 number_sections: true theme: flatly highlight: tango code_folding: show df_print: paged pdf_document: toc: true number_sections: true latex_engine: xelatex --- ```{r setup, include=FALSE} knitr::opts_chunk$set( echo = TRUE, warning = FALSE, message = FALSE, fig.width = 7, fig.height = 4.5, fig.align = "center" ) # Packages we will need # install.packages(c("ggplot2", "deSolve")) # run once if you don't have them library(ggplot2) library(deSolve) theme_set(theme_minimal(base_size = 12)) ``` # The problem: a lake with water lilies Imagine that we arrive at a lake and count the water lilies. - **Day 0:** $N_0 = 100$ water lilies - **Day 1:** $N_1 = 110$ water lilies Returning the next day we observe a **10% increase**. Our goal is: > How, starting from this simple observation, do we arrive at the *exponential growth formula* by way of a differential equation? The journey has three stops: 1. **Observation → rate of change** (discrete thinking) 2. **Discrete → continuous** (we take the limit and write the differential equation $\frac{dN}{dt} = rN$) 3. **Solving the differential equation** (separation of variables → $N(t) = N_0 e^{rt}$) --- # From the observation to the rate of change ## The change in the population Over a time $\Delta t = 1$ day, the change in the population is: $$ \Delta N \;=\; N_1 - N_0 \;=\; 110 - 100 \;=\; 10 $$ ## The **per capita** rate $r$ We are not interested in the absolute change (10 water lilies), because it depends on how many water lilies we had at the start. We are interested in **how much the population changed per individual, per unit of time**: $$ r \;=\; \frac{\Delta N}{N \cdot \Delta t} \;=\; \frac{10}{100 \cdot 1} \;=\; 0.1 \;=\; 10\%/\text{day} $$ $r$ is called the **intrinsic rate of increase** of the population. It is the fundamental parameter of all classical population dynamics. ```{r per-capita-rate} # Computing r from the data N0 <- 100 N1 <- 110 dt <- 1 # days r <- (N1 - N0) / (N0 * dt) r ``` --- # From the discrete to the continuous ## Discrete equation If we assume that this rate stays constant, then at every time step $\Delta t$: $$ \Delta N \;=\; r \, N \, \Delta t \quad\Longleftrightarrow\quad \frac{\Delta N}{\Delta t} \;=\; r N $$ ### example values for r=0.1 (10%) | N | ΔΝ/Δt | | :-------- | :-------- | | 100 | 10 | | 200 | 20 | | 1000 | 100 | | 5000 | 500 | This is a **discrete** approximation: we measure the population at discrete steps (every day). ## The limit $\Delta t \to 0$ But the water lilies do not wait for dawn to grow — they grow *continuously*. Taking the limit as the time step becomes infinitesimal: $$ \lim_{\Delta t \to 0} \frac{\Delta N}{\Delta t} \;=\; \frac{dN}{dt} $$ We arrive at the **differential equation of exponential growth**: $$ \boxed{\;\frac{dN}{dt} \;=\; r\,N\;} $$ > **Intuitive interpretation:** The rate of change of the population at every moment is **proportional** to the current size of the population. The more water lilies there are, the faster new ones appear. --- # Solving the differential equation We will solve $\frac{dN}{dt} = rN$ using the method of **separation of variables**. ## Step 1 — Separate the variables $$ \frac{dN}{N} \;=\; r\,dt $$ (All the $N$ on the left-hand side, all the $t$ on the right.) ## Step 2 — Integrate both sides $$ \int \frac{dN}{N} \;=\; \int r\,dt \quad\Longrightarrow\quad \ln(N) \;=\; r\,t + C $$ where $C$ is a constant of integration. ## Step 3 — Solve for $N$ $$ N(t) \;=\; e^{rt + C} \;=\; e^{C} \cdot e^{rt} $$ ## Step 4 — Find the constant from the initial condition At $t = 0$ we have $N(0) = N_0$, so: $$ N_0 \;=\; e^{C} \cdot e^{0} \;=\; e^{C} \quad\Longrightarrow\quad e^{C} = N_0 $$ ## The final formula $$ \boxed{\;N(t) \;=\; N_0\, e^{r t}\;} $$ This is the **exponential growth formula**. --- # Implementation in R ## Analytical solution Let's plot $N(t) = 100 \cdot e^{0.1\,t}$ for 60 days: ```{r analytical-solution} N0 <- 100 r <- 0.1 t <- seq(0, 60, by = 0.1) N_analytical <- N0 * exp(r * t) df_analytical <- data.frame(t = t, N = N_analytical) ggplot(df_analytical, aes(x = t, y = N)) + geom_line(linewidth = 1, colour = "steelblue") + labs( title = "Exponential growth of the water lily population", subtitle = expression(paste(N(t) == N[0] * e^{r*t}, ", ", N[0] == 100, ", ", r == 0.1)), x = "Time (days)", y = "Population N(t)" ) ``` In 60 days the population reaches `r round(N0 * exp(r * 60))` water lilies — unrealistic for a finite lake! This is the **problem** of exponential growth that we will address later with the logistic model. ## Numerical solution with `deSolve` Very often in biology differential equations **cannot be solved analytically**. In those cases we use numerical integrators. For this particular model, which we already know analytically, this gives us an opportunity to check whether the numerical method agrees. ```{r numerical-solution} # 1. Define the differential equation as a function exponential_model <- function(t, state, parameters) { with(as.list(c(state, parameters)), { dN <- r * N # dN/dt = rN list(dN) }) } # 2. Parameters, initial condition, time parameters <- c(r = 0.1) state <- c(N = 100) times <- seq(0, 60, by = 0.1) # 3. Solve out <- ode(y = state, times = times, func = exponential_model, parms = parameters) out_df <- as.data.frame(out) head(out_df) ``` ## Comparison of the analytical and numerical solutions ```{r comparison-plot} df_compare <- data.frame( t = times, Analytical = N0 * exp(r * times), Numerical = out_df$N ) ggplot(df_compare, aes(x = t)) + geom_line(aes(y = Analytical), linewidth = 1.2, colour = "steelblue") + geom_point(aes(y = Numerical), data = df_compare[seq(1, nrow(df_compare), by = 50), ], colour = "tomato", size = 2) + labs( title = "Analytical (line) versus numerical (dots) solution", x = "Time (days)", y = "Population N(t)" ) # Maximum error max(abs(df_compare$Analytical - df_compare$Numerical)) ``` The error is negligible — the two solutions coincide. --- # A subtle but **very important** observation Notice something: we started by saying *"10% increase per day"*. That is a **discrete** statement. But then we wrote $r = 0.1$ inside a **continuous** differential equation. Are the two really equal? ## The discrete formula If every day the population is multiplied by $1.10$: $$ N(t) \;=\; N_0 \cdot (1+r)^{t} \;=\; 100 \cdot 1.1^{t} $$ ## Comparison of the two models ```{r discrete-vs-continuous} days <- 0:30 N_discrete <- 100 * (1 + 0.1)^days N_continuous <- 100 * exp(0.1 * days) df_dc <- data.frame( day = rep(days, 2), N = c(N_discrete, N_continuous), model = factor(rep(c("Discrete N0·(1+r)^t", "Continuous N0·exp(rt)"), each = length(days))) ) ggplot(df_dc, aes(x = day, y = N, colour = model)) + geom_line(linewidth = 1) + geom_point(size = 1.5) + labs( title = "Discrete versus continuous exponential growth", subtitle = "With the same value r = 0.1 the two models do **not** coincide", x = "Days", y = "Population", colour = "" ) + theme(legend.position = "bottom") ``` ```{r dc-table} data.frame( Day = c(1, 5, 10, 20, 30), Discrete = round(100 * 1.1^c(1, 5, 10, 20, 30), 1), Continuous = round(100 * exp(0.1 * c(1, 5, 10, 20, 30)), 1) ) ``` ## Why do they differ? In the **discrete** model, "10% increase/day" means: at the end of the day the population is 110% of the initial one. In the **continuous** model with $r = 0.1$, the *instantaneous* rate is 10%, but because the newly born water lilies **immediately produce new water lilies themselves** (continuous compounding!), at the end of the day we have: $$ N(1) \;=\; 100 \cdot e^{0.1} \;\approx\; 110.52 $$ that is, a little more than a 10% increase. ## How are they connected? If we want the continuous model to give **exactly** a 10% increase per day, we need: $$ e^{r_{\text{cont}}} = 1.1 \quad\Longrightarrow\quad r_{\text{cont}} \;=\; \ln(1.1) \;\approx\; 0.0953 $$ ```{r r-conversion} log(1.1) # in R, log() is the natural logarithm (ln) ``` > **In ecological practice**, most of the time we work with the continuous model and take $r$ to be the instantaneous rate. The difference is small for small $r$, but it is good to know that it **exists**. --- # The logarithmic axis: the biologist's diagnostic tool Taking the logarithm of the formula $N(t) = N_0 e^{rt}$: $$ \ln N(t) \;=\; \ln N_0 \;+\; r\,t $$ That is, $\ln N$ is a **linear** function of $t$ with slope equal to $r$. This is extremely useful in practice: > **If your data on a semi-log plot ($\ln N$ against $t$) form a straight line, then the population is growing exponentially**, and the slope of the line gives you $r$. ```{r log-plot} ggplot(df_analytical, aes(x = t, y = N)) + geom_line(linewidth = 1, colour = "steelblue") + scale_y_log10() + labs( title = "Exponential growth on a semi-log axis", subtitle = "The straight line indicates exponential growth", x = "Time (days)", y = "log10(N(t))" ) ``` ## Estimating $r$ from data with linear regression ```{r estimate-r} # Simulated data with a little noise set.seed(42) sim_days <- 0:30 sim_N <- 100 * exp(0.1 * sim_days) * exp(rnorm(length(sim_days), 0, 0.05)) # Linear regression log(N) ~ t fit <- lm(log(sim_N) ~ sim_days) summary(fit)$coefficients cat("Estimated r =", round(coef(fit)[2], 4), " (true r = 0.1)\n") cat("Estimated N0 =", round(exp(coef(fit)[1]), 2), " (true N0 = 100)\n") ``` --- # Doubling time A useful quantity: how long does it take for the population to double? $$ 2 N_0 = N_0 e^{r T_2} \;\Longleftrightarrow\; T_2 = \frac{\ln 2}{r} $$ ```{r doubling-time} T2 <- log(2) / r cat("Doubling time of the water lilies:", round(T2, 2), "days\n") ``` So the water lilies double roughly every **7 days**. --- # Exercises 1. **Play with $r$.** Plot on the same graph the exponential growth for $r = 0.05$, $r = 0.1$ and $r = 0.2$. How does the doubling time change? 2. **Negative $r$.** What happens if $r < 0$? Plot the case $r = -0.05$ and interpret it biologically. 3. **Estimation from real data.** You are given the following measurements: ```r day <- c(0, 3, 7, 14, 21, 28) N <- c(50, 73, 109, 218, 437, 870) ``` Estimate $r$ with linear regression on $\ln N$ and predict the population at 60 days. 4. **The mathematical danger.** With $r = 0.1$ and $N_0 = 100$ water lilies, in how many days will the water lilies exceed the number of humans on Earth ($\approx 8 \times 10^{9}$)? What does this tell us about the *limitations* of the exponential model? --- # Next step: the law returns to reality The exponential model **cannot** hold forever. In a real lake there are constraints: space, light, nutrients. At some point the population will *saturate*. In the next lesson we will modify the differential equation so that it includes the **carrying capacity** $K$ of the environment: $$ \frac{dN}{dt} \;=\; r N \left(1 - \frac{N}{K}\right) $$ This is the **logistic equation** — and it is the subject of the next laboratory. --- # R session information ```{r session-info} sessionInfo() ```