Newer
Older
## NOTE : Demographic models are POST-BREEDING
##~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
## Model 1: No Density-Dependence (noDD), No Demographic Stochasticity ~~
##~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
#' Demographic model without Density-Dependence and without Demographic Stochasticity
#' In addition to natural survivals and fecundities, this demographic model includes
#' a specific harvest / fatality parameter.
#' NOTE : This is a POST-BREEDING demographic model.
#'
#' @param N1 a vector of population sizes for each age class at time t1
#' @param s a vector of survival probabilities for each age class
#' @param f a vector of fecundity values for each age class
#' @param h a number. The harvest or fatality rate
#' @param DD_params density-dependence parameters. Not used in this model.
#'
#' @return a vector of population sizes for each age class at time t2
#' @export
#'
#' @examples
#' s <- c(0.5, 0.7, 0.8, 0.95)
#' f <- c(0, 0, 0.05, 0.55)
#' N1 <- c(50, 60, 75, 100)
#' h <- 0.05
#' M1_noDD_noDemoStoch(N1, s, f, h)
#'
M1_noDD_noDemoStoch <- function(N1, s, f, h, DD_params = NULL){

thierrychambert
committed
## M1_noDD_noDemoStoch
# Build the LESLIE matrix
A <- build_Leslie(s = s, f = f)
# Apply the LESLIE matrix calculation at t+1
N2 <- A%*%N1*(1-h)
return(N2)
} # END FUNCTION
################################################################################
##~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
## Model 2: No Density-Dependence (noDD), WITH Demographic Stochasticity ~~
##~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
#' Demographic model without Density-Dependence, but INCLUDING Demographic Stochasticity
#' In addition to natural survivals and fecundities, this demographic model includes
#' a specific harvest / fatality parameter.
#' NOTE : This is a POST-BREEDING demographic model.
#'
#' @param N1 a vector of population sizes for each age class at time t1
#' @param s a vector of survival probabilities for each age class
#' @param f a vector of fecundity values for each age class
#' @param h a number. The harvest or fatality rate
#' @param DD_params density-dependence parameters. Not used in this model.
#'
#' @return a vector of population sizes for each age class at time t2
#' @export
#'
#' @import popbio
#' @import magrittr
#'
#' @examples
#' s <- c(0.5, 0.7, 0.8, 0.95)
#' f <- c(0, 0, 0.05, 0.55)
#' N1 <- c(50, 60, 75, 100)
#' h <- 0.05
#' M2_noDD_WithDemoStoch(N1, s, f, h)
#'
M2_noDD_WithDemoStoch <- function(N1, s, f, h, DD_params = NULL){

thierrychambert
committed
## M2_noDD_WithDemoStoch
nac = length(s)
# Survivors (to "natural mortality" (s) and Wind Turbine Fatalities (1-h))
S <- rbinom(nac, N1, (1-h)*s)
N2 <- c(rep(0, nac-1), tail(S,1)) + c(0, head(S,-1))
# Births
N2[1] <- sum(rpois(nac, f*N2))
return(N2)
} # END FUNCTION
################################################################################
##==============================================================================
## ADD DENSITY DEPENDENCE ==
##==============================================================================
##~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
## Model 3: WITH Density-Dependence (DD), No Demographic Stochasticity ~~
##~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
#' Demographic model with Density-Dependence, but without Demographic Stochasticity
#' In addition to natural survivals and fecundities, this demographic model includes
#' a specific harvest / fatality parameter.
#' NOTE : This is a POST-BREEDING demographic model.
#'
#' @param N1 a vector of population sizes for each age class at time t1
#' @param s a vector of survival probabilities for each age class
#' @param f a vector of fecundity values for each age class
#' @param h a number. The harvest or fatality rate
#' @param DD_params a list containing the 3 parameters required to model density-dependence :
#' rMAX (maximum population intinsic rate of increase: lambda_max - 1),
#' K (carrying capacity), and
#' theta (shape of DD relationshp)
#'
#' @return a vector of population sizes for each age class at time t2
#' @export
#'
#' @import popbio
#' @import magrittr
#'
#' @examples
#' s <- c(0.5, 0.7, 0.8, 0.95)
#' f <- c(0, 0, 0.05, 0.55)
#' N1 <- c(50, 60, 75, 100)
#' h <- 0.05
#' DD_params <- list(rMAX = 0.15, K = 1200, theta = 1)
#' M3_WithDD_noDemoStoch(N1, s, f, h,DD_params = DD_params)

thierrychambert
committed
M3_WithDD_noDemoStoch <- function(N1, s, f, h, DD_params){

thierrychambert
committed
## M3_WithDD_noDemoStoch
rMAX <- DD_params$rMAX
K <- DD_params$K
theta <- DD_params$theta
# Apply density dependence effect
lam_Nt <- 1 + rMAX*(1-(sum(N1)/K)^theta)
# Calibrate vital rates to match lam_Nt

thierrychambert
committed
A <- build_Leslie(s = s, f = f)
diff_rel_lam <- (lam_Nt - lambda(A))/lambda(A)
d <- match_lam_delta(diff_rel_lam = diff_rel_lam, s=s, f=f)
vr_Nt <- c(s,f) + d
s_Nt <- head(vr_Nt, length(s)) %>% sapply(min, 0.999)

thierrychambert
committed
## Check if approximation is close enough to desired lambda
if( abs((lambda(build_Leslie(s = s_Nt, f = f_Nt)) - lam_Nt) / lam_Nt) > 0.005 ){
# If difference is too large : Use optimisation function for better calibration
inits <- c(tail(vr_Nt, length(f)), head(vr_Nt, length(s)) %>% sapply(min, 0.999))
inits <- inits[inits != 0]
vr_calib <- calibrate_params(inits = inits, f = f_Nt, s = s_Nt, lam0 = lam_Nt)
s_Nt <- head(vr_calib, length(s_Nt))
f_Nt <- tail(vr_calib, length(f_Nt))
} # if
# Build the LESLIE matrix
A_Nt <- build_Leslie(s = s_Nt, f = f_Nt)
# Apply the LESLIE matrix calculation at t+1
N2 <- A_Nt%*%N1*(1-h)
return(N2)
} # END FUNCTION
################################################################################
##~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
## Model 4: With Density-Dependence (DD), WITH Demographic Stochasticity ~~
##~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
#' Demographic model with Density-Dependence and INCLUDING Demographic Stochasticity
#' In addition to natural survivals and fecundities, this demographic model includes
#' a specific harvest / fatality parameter.
#' NOTE : This is a POST-BREEDING demographic model.
#'
#' @param N1 a vector of population sizes for each age class at time t1
#' @param s a vector of survival probabilities for each age class
#' @param f a vector of fecundity values for each age class
#' @param h a number. The harvest or fatality rate
#' @param DD_params density-dependence parameters. Not used in this model.
#'
#' @import popbio
#' @import magrittr
#'
#' @return a vector of population sizes for each age class at time t2
#' @export
#'
#' @examples
#' s <- c(0.5, 0.7, 0.8, 0.95)
#' f <- c(0, 0, 0.05, 0.55)
#' N1 <- c(50, 60, 75, 100)
#' h <- 0.05
#' DD_params <- list(rMAX = 0.15, K = 1200, theta = 1)
#' M4_WithDD_WithDemoStoch(N1, s, f, h, DD_params = DD_params)
#'
M4_WithDD_WithDemoStoch <- function(N1, s, f, h, DD_params){

thierrychambert
committed
## M4_WithDD_WithDemoStoch
rMAX <- DD_params$rMAX
K <- DD_params$K
theta <- DD_params$theta
# Apply density dependence effect
lam_Nt <- 1 + rMAX*(1-(sum(N1)/K)^theta)
# Calibrate vital rates to match lam_Nt
A <- build_Leslie(s = s, f = f)
diff_rel_lam <- (lam_Nt - lambda(A))/lambda(A)
d <- match_lam_delta(diff_rel_lam = diff_rel_lam, s=s, f=f)
vr_Nt <- c(s,f) + d
s_Nt <- head(vr_Nt, length(s)) %>% sapply(min, 0.999)
f_Nt <- tail(vr_Nt, length(f))

thierrychambert
committed
## Check if approximation is close enough to desired lambda
if( abs((lambda(build_Leslie(s = s_Nt, f = f_Nt)) - lam_Nt) / lam_Nt) > 0.005 ){
# If difference is too large : Use optimisation function for better calibration
inits <- c(tail(vr_Nt, length(f)), head(vr_Nt, length(s)) %>% sapply(min, 0.999))
inits <- inits[inits != 0]
vr_calib <- calibrate_params(inits = inits, f = f_Nt, s = s_Nt, lam0 = lam_Nt)
s_Nt <- head(vr_calib, length(s_Nt))
f_Nt <- tail(vr_calib, length(f_Nt))
} # if
# Number of age classes
nac = length(s)
# Survivors (to "natural mortality" (s) and Wind Turbine Fatalities (1-h))
S <- rbinom(nac, N1, (1-h)*s_Nt)
N2 <- c(rep(0, nac-1), tail(S,1)) + c(0, head(S,-1))
# Births
N2[1] <- sum(rpois(nac, f_Nt*N2))
return(N2)
} # END FUNCTION
################################################################################