Patched aspatial_siri function

aspatial_siri function threw errors when abundance threshold was zero; this has now been fixed.

R/aspatial_siri_seasons.R

@@ -21,7 +21,7 @@
 #'     \item{\code{stages}}{Number of life cycle stages.}
 #'     \item{\code{compartments}}{Number of disease compartments (e.g., 3 for a
 #'     SIR model).}
-#'     \item{\code{abundance_threshold}}{A quasi-extinction threshold below
+#'     \item{\code{abundance_threshold}}{A quasi-extinction threshold at
 #'     which a population becomes extinct.}
 #'     \item{\code{mortality}}{A vector of mortality rates, one for each
 #'     combination of stages and compartments.}
@@ -251,7 +251,7 @@ siri_model_winter <- function(inputs) {
 #' susceptible/recovered stage in the season in question.
 #' @param fecundity Only necessary when `season = "breeding"` (see below).
 #' Default NULL. A single numeric with the daily fecundity of adults.
-#' @param abundance_threshold A quasi-extinction threshold below which a
+#' @param abundance_threshold A quasi-extinction threshold at which a
 #' population becomes extinct.
 #' @param carrying_capacity A single numeric that indicates the carrying
 #' capacity of the population in this season.
@@ -282,87 +282,118 @@ aspatial_siri <- function(initial_pop,
 
   N <- pmin(total_pop, carrying_capacity)
 
-  if (N < abundance_threshold) {
+  if (N <= abundance_threshold) {
     # Fill everything with 0
     final_state <- rep(0, 8)
   } else {
-
-  dd_multiplier <- 1.0 + N / carrying_capacity
-
-  dd_mortality <- mortality * pmin(dd_multiplier, 1.0)
-
-  new_juv <- 0.0
-
-  if (season == "breeding") {
-    adults <- Sa + I1a + Ra + I2a
-    births <- rpois(adults, fecundity * (1.0 - N / carrying_capacity))
-    new_juv <- sum(births)
-  }
-
-  infection1_juv_raw <- rbinom(1, Sj * (I1j + I2j + I1a + I2a), transmission[1])
-  infection1_juv <- pmin(infection1_juv_raw, Sj)
-
-  infection1_adult_raw <- rbinom(1, Sa * (I1j + I2j + I1a + I2a), transmission[2])
-  infection1_adult <- pmin(infection1_adult_raw, Sa)
-
-  susceptible_adult_death <- rbinom(1, Sa - infection1_adult, dd_mortality[2])
-
-  if (season == "breeding") {
-    susceptible_juvenile_death <- rbinom(1, Sj + new_juv - infection1_juv, dd_mortality[1])
-  } else {
-    susceptible_juvenile_death <- rbinom(1, Sj - infection1_juv, dd_mortality[1])
-  }
-
-  infected1_juvenile_death <- rbinom(1, I1j + infection1_juv, dd_mortality[3])
-
-  infected1_adult_death <- rbinom(1, I1a + infection1_adult, dd_mortality[4])
-
-  recovery1_juv_raw <- rbinom(1, I1j + infection1_juv - infected1_juvenile_death, recovery[3])
-  recovery1_juv <- pmin(recovery1_juv_raw, I1j + infection1_juv - infected1_juvenile_death)
-
-  recovery1_adult_raw <- rbinom(1, I1a + infection1_adult - infected1_adult_death, recovery[4])
-  recovery1_adult <- pmin(recovery1_adult_raw, I1a + infection1_adult - infected1_adult_death)
-
-  infection2_juv_raw <- rbinom(1, (Rj + recovery1_juv) * (I1j + I2j + I1a + I2a), transmission[3])
-  infection2_juv <- pmin(infection2_juv_raw, Rj + infection1_juv)
-
-  infection2_adult_raw <- rbinom(1, (Ra + recovery1_adult) * (I1j + I2j + I1a + I2a), transmission[4])
-  infection2_adult <- pmin(infection2_adult_raw, Ra + recovery1_adult)
-
-  infected2_juvenile_death <- rbinom(1, I2j + infection2_juv, dd_mortality[7])
-
-  infected2_adult_death <- rbinom(1, I2a + infection2_adult, dd_mortality[8])
-
-  recovery2_juv_raw <- rbinom(1, I2j + infection2_juv - infected2_juvenile_death, recovery[3])
-  recovery2_juv <- pmin(recovery2_juv_raw, I2j + infection2_juv - infected2_juvenile_death)
-
-  recovery2_adult_raw <- rbinom(1, I2a + infection2_adult - infected2_adult_death, recovery[4])
-  recovery2_adult <- pmin(recovery2_adult_raw, I2a + infection2_adult - infected2_adult_death)
-
-  recovered_juvenile_death <- rbinom(1, Rj + recovery1_juv + recovery2_juv - infection2_juv, dd_mortality[5])
-
-  recovered_adult_death <- rbinom(1, Ra + recovery1_adult + recovery2_adult - infection2_adult, dd_mortality[6])
-
-  # Update state for the next time step
-  if (season == "breeding") {
-    final_state <- c(Sj + new_juv - infection1_juv - susceptible_juvenile_death,
-                     Sa - infection1_adult - susceptible_adult_death,
-                     I1j + infection1_juv - recovery1_juv - infected1_juvenile_death,
-                     I1a + infection1_adult - recovery1_adult - infected1_adult_death,
-                     Rj + recovery1_juv + recovery2_juv - infection2_juv - recovered_juvenile_death,
-                     Ra + recovery1_adult + recovery2_adult - infection2_adult - recovered_adult_death,
-                     I2j + infection2_juv - recovery2_juv - infected2_juvenile_death,
-                     I2a + infection2_adult - recovery2_adult - infected2_adult_death)
-  } else {
-    final_state <- c(Sj - infection1_juv - susceptible_juvenile_death,
-                     Sa - infection1_adult - susceptible_adult_death,
-                     I1j + infection1_juv - recovery1_juv - infected1_juvenile_death,
-                     I1a + infection1_adult - recovery1_adult - infected1_adult_death,
-                     Rj + recovery1_juv + recovery2_juv - infection2_juv - recovered_juvenile_death,
-                     Ra + recovery1_adult + recovery2_adult - infection2_adult - recovered_adult_death,
-                     I2j + infection2_juv - recovery2_juv - infected2_juvenile_death,
-                     I2a + infection2_adult - recovery2_adult - infected2_adult_death)
-  }
+    dd_multiplier <- 1.0 + N / carrying_capacity
+
+    dd_mortality <- mortality * pmin(dd_multiplier, 1.0)
+
+    new_juv <- 0.0
+
+    if (season == "breeding") {
+      adults <- Sa + I1a + Ra + I2a
+      births <-
+        rpois(adults, fecundity * (1.0 - N / carrying_capacity))
+      new_juv <- sum(births)
+    }
+
+    infection1_juv_raw <-
+      rbinom(1, Sj * (I1j + I2j + I1a + I2a), transmission[1])
+    infection1_juv <- pmin(infection1_juv_raw, Sj)
+
+    infection1_adult_raw <-
+      rbinom(1, Sa * (I1j + I2j + I1a + I2a), transmission[2])
+    infection1_adult <- pmin(infection1_adult_raw, Sa)
+
+    susceptible_adult_death <-
+      rbinom(1, Sa - infection1_adult, dd_mortality[2])
+
+    if (season == "breeding") {
+      susceptible_juvenile_death <-
+        rbinom(1, Sj + new_juv - infection1_juv, dd_mortality[1])
+    } else {
+      susceptible_juvenile_death <-
+        rbinom(1, Sj - infection1_juv, dd_mortality[1])
+    }
+
+    infected1_juvenile_death <-
+      rbinom(1, I1j + infection1_juv, dd_mortality[3])
+
+    infected1_adult_death <-
+      rbinom(1, I1a + infection1_adult, dd_mortality[4])
+
+    recovery1_juv_raw <-
+      rbinom(1, I1j + infection1_juv - infected1_juvenile_death, recovery[3])
+    recovery1_juv <-
+      pmin(recovery1_juv_raw,
+           I1j + infection1_juv - infected1_juvenile_death)
+
+    recovery1_adult_raw <-
+      rbinom(1, I1a + infection1_adult - infected1_adult_death, recovery[4])
+    recovery1_adult <-
+      pmin(recovery1_adult_raw,
+           I1a + infection1_adult - infected1_adult_death)
+
+    infection2_juv_raw <-
+      rbinom(1, (Rj + recovery1_juv) * (I1j + I2j + I1a + I2a), transmission[3])
+    infection2_juv <- pmin(infection2_juv_raw, Rj + infection1_juv)
+
+    infection2_adult_raw <-
+      rbinom(1, (Ra + recovery1_adult) * (I1j + I2j + I1a + I2a), transmission[4])
+    infection2_adult <- pmin(infection2_adult_raw, Ra + recovery1_adult)
+
+    infected2_juvenile_death <-
+      rbinom(1, I2j + infection2_juv, dd_mortality[7])
+
+    infected2_adult_death <-
+      rbinom(1, I2a + infection2_adult, dd_mortality[8])
+
+    recovery2_juv_raw <-
+      rbinom(1, I2j + infection2_juv - infected2_juvenile_death, recovery[3])
+    recovery2_juv <- pmin(recovery2_juv_raw,
+           I2j + infection2_juv - infected2_juvenile_death)
+
+    recovery2_adult_raw <-
+      rbinom(1, I2a + infection2_adult - infected2_adult_death, recovery[4])
+    recovery2_adult <-
+      pmin(recovery2_adult_raw,
+           I2a + infection2_adult - infected2_adult_death)
+
+    recovered_juvenile_death <- rbinom(1,
+             Rj + recovery1_juv + recovery2_juv - infection2_juv,
+             dd_mortality[5])
+
+    recovered_adult_death <- rbinom(1,
+             Ra + recovery1_adult + recovery2_adult - infection2_adult,
+             dd_mortality[6])
+
+    # Update state for the next time step
+    if (season == "breeding") {
+      final_state <-
+        c(
+          Sj + new_juv - infection1_juv - susceptible_juvenile_death,
+          Sa - infection1_adult - susceptible_adult_death,
+          I1j + infection1_juv - recovery1_juv - infected1_juvenile_death,
+          I1a + infection1_adult - recovery1_adult - infected1_adult_death,
+          Rj + recovery1_juv + recovery2_juv - infection2_juv - recovered_juvenile_death,
+          Ra + recovery1_adult + recovery2_adult - infection2_adult - recovered_adult_death,
+          I2j + infection2_juv - recovery2_juv - infected2_juvenile_death,
+          I2a + infection2_adult - recovery2_adult - infected2_adult_death
+        )
+    } else {
+      final_state <- c(
+        Sj - infection1_juv - susceptible_juvenile_death,
+        Sa - infection1_adult - susceptible_adult_death,
+        I1j + infection1_juv - recovery1_juv - infected1_juvenile_death,
+        I1a + infection1_adult - recovery1_adult - infected1_adult_death,
+        Rj + recovery1_juv + recovery2_juv - infection2_juv - recovered_juvenile_death,
+        Ra + recovery1_adult + recovery2_adult - infection2_adult - recovered_adult_death,
+        I2j + infection2_juv - recovery2_juv - infected2_juvenile_death,
+        I2a + infection2_adult - recovery2_adult - infected2_adult_death
+      )
+    }
   }
 
   return(final_state)

R/check_aspatial_siri_inputs.R

@@ -15,7 +15,7 @@
 #'     \item{\code{stages}}{Number of life cycle stages.}
 #'     \item{\code{compartments}}{Number of disease compartments (e.g., 3 for a
 #'     SIR model).}
-#'     \item{\code{abundance_threshold}}{A quasi-extinction threshold below
+#'     \item{\code{abundance_threshold}}{A quasi-extinction threshold at
 #'     which a population becomes extinct.}
 #'     \item{\code{mortality}}{A vector of mortality rates, one for each
 #'     combination of stages and compartments.}

R/check_simulator_inputs.R

@@ -71,7 +71,7 @@
 #'  \code{stages * compartments}. A list of vectors may be provided if this
 #'  varies by season. If no fecundity mask is provided, then it is assumed that
 #'  all stages and compartments reproduce.}
-#'  \item{\code{abundance_threshold}}{A quasi-extinction threshold below which a
+#'  \item{\code{abundance_threshold}}{A quasi-extinction threshold at which a
 #'  population becomes extinct. Default: 0.}
 #'  \item{\code{demographic_stochasticity}}{Boolean for choosing demographic
 #'  stochasticity for transition, dispersal, and/or other processes (default is

R/disease_simulator.R

@@ -77,7 +77,7 @@
 #'  \code{stages * compartments}. A list of vectors may be provided if this
 #'  varies by season. If no fecundity mask is provided, then it is assumed that
 #'  all stages and compartments reproduce.}
-#'  \item{\code{abundance_threshold}}{A quasi-extinction threshold below which a
+#'  \item{\code{abundance_threshold}}{A quasi-extinction threshold at which a
 #'  population becomes extinct. Default: 0.}
 #'  \item{\code{demographic_stochasticity}}{Boolean for choosing demographic
 #'  stochasticity for transition, dispersal, and/or other processes (default is

R/disease_transformation.R

@@ -20,7 +20,7 @@
 #'     demographic stochasticity for the transformation.}
 #'     \item{\code{density_stages}}{Array of booleans or numeric (0,1) for each
 #'     stage to indicate which stages are affected by density.}
-#'     \item{\code{abundance_threshold}}{A quasi-extinction threshold below
+#'     \item{\code{abundance_threshold}}{A quasi-extinction threshold at
 #'     which a population becomes extinct.}
 #'     \item{\code{mortality}}{A vector of mortality rates, one for each
 #'     combination of stages and compartments.}