Logistic Growth with Nutrient Limitation

Set Models Parameter Below

Population Over Time

P[i+1] = P[i] + growth * (1 - P[i]/capacity - 1/lifetime) * P[i]


Logistic growth is widely used in biology and ecology to describe populations growth under environmental constraints. Essentially this is a modification of exponential growth with a growth rate modified (reduced) as the population approaches limits set by nutrients or other environmental conditions.

Exponential growth is written: \[ \frac{dP}{dt} = g P \]

where P is population and g is the growth rate (fractional increase in poulation per unit time).

Next we modify the growth rate to decrease as the population approaches a limit imposed by the environment, for example by a limiting nutrient:

\[ \frac{dP}{dt} = g(1-\frac{P}{K}) P \].

It’s helpful to think of \(g(1-\frac{P}{K})\) as the modified or constrained growth rate, which is the probability of reproduction for each individual in the population per tinmestep.

Finally, we add death. Here mortality is implemented as an average lifetime \(L\), so that in each time step the probability of an individual dying is \(\frac{g}{L}\).

The final model as implemented here is then

\[ \frac{dP}{dt} = g(1-\frac{P}{K}-\frac{1}{L}) P \].

The program that does the calculation is very simple, and the code that controls this website is surprisingly simple too! It is all written in the programming language R, using a web programming package called shiny. You can read all about it on the “Website Code” tab to the right.

How this website works (including all the code!)

This website is controlled using the R package “shiny.” There are three important components:

  • A script (logistic.R) that integrates the logistic growth model
  • The user interface (ui.R) that controls the sliders and passes inputs to the model
  • The “server” (server.R) that responds to user interface events in the web browser

Scroll down or click links in the list above to read all about it!


logistic <- function(growth=0.1, lifetime=10, capacity=100, nTime=100){

  # Set up an array to hold the values of population
  P <- array(0,nTime)

  # Set the initial population
  P[1] <- 1.

  # Loop through nTime steps, updating the population each time 
  # and storing the values in the array P for later plotting

  for (i in 1:(nTime-1)) {
    P[i+1] <- P[i] + growth * (1 - P[i]/capacity - 1/lifetime) * P[i]

  # Return to the main program, passing back the timeseries of P   

This is the user interface that actually controls this website



# Define user interface for logistic growth model

  #  Application title
  headerPanel("Logistic Growth with Nutrient Limitation"),

  # Sidebar on the left with controls for the user to specify atmospheric properties


    h4('Set Models Parameter Below'),            

    # Set drainage times for each pool
    sliderInput("growth", "Maximum Growth Rate (percent)", 
                min = 1, max = 100, value = 10, step= 1),
    sliderInput("lifetime", "Average Lifetime (time steps)", 
                min = 1, max = 50, value = 10, step= 1),
    sliderInput("capacity", "Nutrient-Limited Maximum Population", 
                min = 1, max = 1000, value = 100, step=1),
    sliderInput("nTime", "Number of time steps", 
                min = 10, max = 500, value = 100, step=10)

  # Show a tabset in the main panel of the browser that displays model output

      # First tab shows plots of model output
               h4('Population Over Time'), 
               p('  P[i+1] = P[i] + growth * (1 - P[i]/capacity - 1/lifetime) * P[i]'),

      # Second tab displays a brief model description
      tabPanel("Model Description", 

      # Third tab displays website code
      tabPanel("Website code", 



# Source required R scripts

shinyServer(function(input, output) {

  #  plot of model output
  output$population.plot <- renderPlot(
      plot(logistic(growth=input$growth/100, lifetime=input$lifetime, 
                    capacity=input$capacity, nTime=input$nTime),
           type='l', lwd=10, col='darkgreen', ylab='', xlab='time', 
           cex.lab=2, cex.axis=2)