中文 | Español | English | português | Turkish

knitr::opts_chunk$set(echo = TRUE)

library(ggplot2)
library(tidyr)
library(purrr)
library(dplyr)
library(deSolve)
library(ecoevoapps)
library(patchwork)
library(RColorBrewer)
library(kableExtra)

This app simulates the dynamics of two consumer species ($P_1$ and $P_2$) that both require the same resource species ($H$) to grow. This allows us to implicitly model resource competition between the two consumer species, since an increase in $P_1$ results in a decreased growth rate of $P_2$ because of $P_1$'s suppression of resource $H$.

Even this fairly straightforward scenario can be modeled in many ways. E.g. we can consider a resource species $H$ that, in the absence of either consumer, grows either exponentially (with no intrinsic limit to population growth), or logistically.

The models can also consider different relationships between the density of resources and rate at which they are attacked by consumers. This relationship is called the functional response, and we consider two types of relationships. First, we consider a "Type 1 functional response", in which the rate at which consumers consume resources is a linear function of the resource density. The second relationship is the "Type II functional response", in which the amount of resources consumed is a saturating function of resource density. The idea of the Type II functional response was developed by C.S. Holling in the 1959 paper "Some characteristics of simple types of predation and parasitism".

This system of equations describes the most complex version of this model:

[ \begin{align} \frac{dH}{dt} &= rH \bigl(1-qH\bigr) - \frac{a_1HP_1}{1+a_1T_{h1}H} - \frac{a_2HP_2}{1+a_2T_{h2}H} \ \ \frac{dP_1}{dt} &= e_{1} \frac{a_1HP_1}{1+a_1T_{h1}H} - d_1P_1 \ \ \frac{dP2}{dt} &= e_2 \frac{a_2HP_2}{1+a_2T_{h2}H} - d_2P_2 \end{align} ]

pars_vars <- c("$H$", 
               "$P_i$", 
               "$r$", 
               "$q$", 
               "$a_i$", 
               "$T_{h_i}$", 
               "$e_i$", 
               "$d_i$")
descriptions <- c("Population size of the prey",
                 "Population size of predator $i$",
                 "Per capita growth rate of prey",
                 "Per capita prey self-limitation",
                 "Attack rate of predator $i$",
                 "Handling time of predator $i$",
                 "Conversion efficiency of predator $i$",
                 "Death rate of predator $i$")
param_df <- data.frame(pars_vars, descriptions)
kable(x = param_df, format = "html", 
      col.names = c("Parameter/Variable", "Description")) %>%
  kable_styling(full_width = FALSE, 
                bootstrap_options = c("striped", "hover", "condensed"),
                position = "center")
wellPanel(
  fluidRow(
    column(4, 
           ### Ask users for parameter values ----
           ## start with Prey params
           numericInput("H", label = "Initial population size of Prey", 
                        min = 1, value = 30),
           sliderInput("r", label = "Per capita growth rate of Prey (r)", 
                       min = .0001, max = 1.0, value = .2),
           sliderInput("q", label = "Prey self-limitation rate (q)",
                       min = 0, max = 0.10, step = 0.0001, value = .0066),
           br(),
           hr(),
           numericInput("t", label = "Timesteps", min = 1, value = 1000),
           checkboxGroupInput("fxnrep", 
                              label = "Display plot of functional responses?",
                              choices = c("Yes" = "Yes"),
                              selected = NULL)

    ),

    column(4,
           ## now pred1 params + initial conditions
           numericInput("P1", label = "Initial population size of Predator 1", 
                        min = 1, value = 25),
           sliderInput("a1", label = "Predator 1 attack rate (a1)", 
                       min = .001, max = 1.0, value = .02),
           sliderInput("T_h1", label = "Predator 1 handling time (Th_1)", 
                       min = 0, max = 1.0, value = 0.1),
           sliderInput("e1", label = "Predator 1 conversion efficiency (e1)",
                       min = .001, max = 1.0, value = 0.4),
           sliderInput("d1", label = "Per capita death rate of Predator 1 (d1)",
                       min = .0001, max = 1.0, value = .1)
    ),

    column(4,
           ## now pred2 params + initial conditions
           numericInput("P2", label = "Initial population size of Predator 2",
                        min = 1, value = 25),
           sliderInput("a2", label = "Predator 2 attack rate (a2)",
                       min = .001, max = 1.0, value = .02),
           sliderInput("T_h2", label = "Predator 2 handling time (Th_2)", 
                       min = 0, max = 1.0, value = 0.1),
           sliderInput("e2", label = "Predator 2 conversion efficiency (e2)",
                       min = .001, max = 1.0, value = 0.39),
           sliderInput("d2", label = "Per capita death rate of Predator 2 (d2)", 
                       min = .0001, max = 1.0, value = .1)
    )
  ))


fluidRow(
  column(12,
         renderPlot(plots_to_print())
  ))


# Set the initial population sizes

init <- reactive({c(H = input$H , P1 = input$P1, P2 = input$P2)})

# Set the parameter values

# description of parameters:
# r = per capita growth rate (prey)
# a = attack rate 
# T_h = handling time
# e = conversion efficiency
# d = predator death rate 

params <- reactive({c(r = input$r, q = input$q,
                      a1 = input$a1, T_h1 = input$T_h1, 
                      e1 = input$e1, d1 = input$d1, 
                      a2 = input$a2, T_h2 = input$T_h2, 
                      e2 = input$e2, d2 = input$d2)})

# Time over which to simulate model dynamics
time <- reactive({seq(0,input$t, by = .1)})

# Simulate the dynamics
out <- reactive({ run_biotic_comp_model(time = time(), 
                                        init = init(), 
                                        params = params())
})

# Make a plot of abundance over time
abund_plot <- reactive({
  plot_biotic_comp_time(out())
})

## functional responses of the two consumers
plot_fxnreps <- reactive({
    wrap_plots(plot_functional_responses(params()), ncol = 1)
})

# Make a list of plots and print out plots based on which ones were requested ----
plots_to_print <- reactive({
  if("Yes" %in% input$fxnrep){
    wrap_plots(abund_plot(), plot_fxnreps(), design = "AAB")
  } else {
    wrap_plots(abund_plot(), ncol = 1)
  }

})

References


suppressWarnings(ecoevoapps::print_app_footer())


gauravsk/ecoevoapps documentation built on July 9, 2024, 9:37 p.m.