mth-335-project-1/project.R

141 lines
4.8 KiB
R

####################
# Global Constants #
####################
# age at which birds in captivity begin breeding (become adults)
captive_begin_breeding <- 7
# age at which birds in the wild begin breeding (become adults)
wild_begin_breeding <- 8
# if the captive population exceeds this amount, then the entire year one
# population will be released
max_in_captivity <- 150
# carrying capacity in the wild
wild_carrying_capacity <- 400
# percentage of children which will make it to year 1 in captivity,
# per mating pair
captive_productivity_rate <- 0.9
# percentage of children which will make it to year 1 in captivity,
# per mating pair (with double clutching)
captive_productivity_rate_dc <- 1.8
# percentage of children which will make it to year 1, per mating pair
wild_productivity_rate <- 0.3729
################################
# Global Simulation Parameters #
################################
# number of years to simulate
simulation_years <- 50
######################
# Derived Parameters #
######################
# rate at which captive immature birds will become adults
# an approximation which assumes an equal distribution of immature bird ages
captive_immature_to_adult_rate <- 1 / (captive_begin_breeding - 2)
# rate at which wild immature birds will become adults
# an approximation which assumes an equal distribution of immature bird ages
wild_immature_to_adult_rate <- 1 / (wild_begin_breeding - 2)
# Simulates the Condor population over time using the given parameters.
# * `rd1` - death rate of year one birds
# * `rdI` - death rate of immature birds
# * `rdA` - death rate of adult birds
simulate_condor_population <- function(rd1, rdI, rdA) {
# captive population groups
c_Y0 <- c(30)
c_Y1 <- c(30)
c_I <- c(30)
c_A <- c(30)
# wild population groups
w_Y0 <- c(100)
w_Y1 <- c(100)
w_I <- c(100)
w_A <- c(100)
years <- seq(1, simulation_years)
for (t in 2:simulation_years) {
# model captive population changes (year 1, immature, and adult)
c_Y1[t] <- c_Y0[t - 1] - # last year's population
rd1 * c_Y0[t - 1] # outflow to deaths
c_I[t] <- c_I[t - 1] - # last year's population
rdI * c_I[t - 1] - # outflow to deaths
captive_immature_to_adult_rate * c_I[t - 1] + # outflow to adult
c_Y1[t - 1] # inflow from year 1
c_A[t] <- c_A[t - 1] - # last year's population
rdA * c_A[t - 1] + # outflow to deaths
captive_immature_to_adult_rate * c_I[t - 1] # inflow from immature
# model new births in captive population
c_Y0[t] <- captive_productivity_rate * (c_A[t] / 2)
# sum of all population groups
last_captive_population <- c_Y0[t - 1] + c_Y1[t - 1] + c_I[t - 1] + c_A[t - 1]
released_to_wild <- if (last_captive_population > max_in_captivity) {
c_Y1[t]
} else {
0.8 * c_Y1[t]
}
# subtract the released birds from the captive year one population
c_Y1[t] <- c_Y1[t] - released_to_wild
# carrying capacity death rate
last_wild_population <- w_Y0[t - 1] + w_Y1[t - 1] + c_I[t - 1] + c_A[t - 1]
approx_wild_growth_rate <- (wild_productivity_rate * (w_A[t - 1] / 2) + released_to_wild) / last_wild_population
rd_carrying <- approx_wild_growth_rate * (last_wild_population / wild_carrying_capacity) # r * P / M
w_Y1[t] <- (w_Y0[t - 1] + # inflow from year zero
released_to_wild) * # inflow from captive
(1 - rd1) * # outflow to deaths (death rate)
(1 - rd_carrying) # outflow to deaths (carrying capacity)
w_I[t] <- w_I[t - 1] + # last year's population
w_Y1[t - 1] - # inflow from year one
w_I[t - 1] * wild_immature_to_adult_rate - # outflow to adult
w_I[t - 1] * rdI - # outflow to deaths (death rate)
w_I[t - 1] * rd_carrying # outflow to deaths (carrying capacity)
w_A[t] <- w_A[t - 1] + # last year's population
w_I[t - 1] * wild_immature_to_adult_rate - # inflow from immature
w_A[t - 1] * rdA - # outflow to deaths (death rate)
w_A[t - 1] * rd_carrying # outflow to deaths (carrying capacity)
w_Y0[t] <- wild_productivity_rate * (w_A[t] / 2) *
(1 - rd_carrying) # deaths due to carrying capacity
}
captive_population <- c_Y0 + c_Y1 + c_I + c_A
plot(
years,
c_Y0,
ylim = c(0, 150),
col = "red",
main = "Condor Population in Captivity",
xlab = "Time (years)",
ylab = "Population"
)
points(years, c_Y1, col = "orange")
points(years, c_I, col = "darkgreen")
points(years, c_A, col = "blue")
points(years, captive_population)
wild_population <- w_Y0 + w_Y1 + w_I + w_A
plot(
years,
w_Y0,
ylim = c(0, wild_carrying_capacity + 100),
col = "red",
main = "Condor Population in the Wild",
xlab = "Time (years)",
ylab = "Population"
)
points(years, w_Y1, col = "orange")
points(years, w_I, col = "darkgreen")
points(years, w_A, col = "blue")
points(years, wild_population)
}