structure(function (param = parameter_template)
{
if (is.data.frame(param)) {
param_lower <- structure(param$lower, names = param$switch)
param_upper <- structure(param$upper, names = param$switch)
param <- structure(param$value, names = param$switch)
}
else {
param_lower <- lapply(param, function(x) NA)
param_upper <- lapply(param, function(x) NA)
}
stopifnot("retro_years" %in% names(param))
stopifnot("fish.Linf" %in% names(param))
stopifnot("fish.K" %in% names(param))
stopifnot("fish.t0" %in% names(param))
stopifnot("fish.lencv" %in% names(param))
stopifnot("fish.init.scalar" %in% names(param))
stopifnot("fish.init.1" %in% names(param))
stopifnot("fish.init.2" %in% names(param))
stopifnot("fish.init.3" %in% names(param))
stopifnot("fish.init.4" %in% names(param))
stopifnot("fish.init.5" %in% names(param))
stopifnot("fish.init.6" %in% names(param))
stopifnot("fish.init.7" %in% names(param))
stopifnot("fish.init.8" %in% names(param))
stopifnot("fish.init.9" %in% names(param))
stopifnot("fish.init.10" %in% names(param))
pt.fish.init <- list(`1` = param[["fish.init.1"]], `2` = param[["fish.init.2"]], `3` = param[["fish.init.3"]], `4` = param[["fish.init.4"]], `5` = param[["fish.init.5"]], `6` = param[["fish.init.6"]], `7` = param[["fish.init.7"]], `8` = param[["fish.init.8"]], `9` = param[["fish.init.9"]], `10` = param[["fish.init.10"]])
stopifnot("fish.M.1" %in% names(param))
stopifnot("fish.M.2" %in% names(param))
stopifnot("fish.M.3" %in% names(param))
stopifnot("fish.M.4" %in% names(param))
stopifnot("fish.M.5" %in% names(param))
stopifnot("fish.M.6" %in% names(param))
stopifnot("fish.M.7" %in% names(param))
stopifnot("fish.M.8" %in% names(param))
stopifnot("fish.M.9" %in% names(param))
stopifnot("fish.M.10" %in% names(param))
pt.fish.M <- list(`1` = param[["fish.M.1"]], `2` = param[["fish.M.2"]], `3` = param[["fish.M.3"]], `4` = param[["fish.M.4"]], `5` = param[["fish.M.5"]], `6` = param[["fish.M.6"]], `7` = param[["fish.M.7"]], `8` = param[["fish.M.8"]], `9` = param[["fish.M.9"]], `10` = param[["fish.M.10"]])
stopifnot("init.F" %in% names(param))
stopifnot("recage" %in% names(param))
stopifnot("fish.walpha" %in% names(param))
stopifnot("fish.wbeta" %in% names(param))
stopifnot("report_detail" %in% names(param))
stopifnot("fish.comm.alpha" %in% names(param))
stopifnot("fish.comm.l50" %in% names(param))
stopifnot("fish.bbin" %in% names(param))
stopifnot("fish.rec.proj" %in% names(param))
stopifnot("fish.rec.1990" %in% names(param))
stopifnot("fish.rec.1991" %in% names(param))
stopifnot("fish.rec.1992" %in% names(param))
stopifnot("fish.rec.1993" %in% names(param))
stopifnot("fish.rec.1994" %in% names(param))
stopifnot("fish.rec.1995" %in% names(param))
stopifnot("fish.rec.1996" %in% names(param))
stopifnot("fish.rec.1997" %in% names(param))
stopifnot("fish.rec.1998" %in% names(param))
stopifnot("fish.rec.1999" %in% names(param))
stopifnot("fish.rec.2000" %in% names(param))
pt.fish.rec <- list(`1990` = param[["fish.rec.1990"]], `1991` = param[["fish.rec.1991"]], `1992` = param[["fish.rec.1992"]], `1993` = param[["fish.rec.1993"]], `1994` = param[["fish.rec.1994"]], `1995` = param[["fish.rec.1995"]], `1996` = param[["fish.rec.1996"]], `1997` = param[["fish.rec.1997"]], `1998` = param[["fish.rec.1998"]], `1999` = param[["fish.rec.1999"]], `2000` = param[["fish.rec.2000"]])
stopifnot("fish.rec.scalar" %in% names(param))
stopifnot("adist_surveyindices_log_acoustic_dist_weight" %in% names(param))
stopifnot("cdist_sumofsquares_comm_ldist_weight" %in% names(param))
assert_msg <- function(expr, message) {
if (isFALSE(expr)) {
warning(message)
return(TRUE)
}
return(FALSE)
}
dif_pmax <- function(a, b, scale) {
logspace_add <- function(a, b) pmax(a, b) + log1p(exp(pmin(a, b) - pmax(a, b)))
b <- as.vector(b)
logspace_add(a * scale, b * scale)/scale
}
avoid_zero <- function(a) {
dif_pmax(a, 0, 1000)
}
nvl <- function(...) {
for (i in seq_len(...length())) if (!is.null(...elt(i)))
return(...elt(i))
return(NULL)
}
normalize_vec <- function(a) {
a/avoid_zero(sum(a))
}
as_numeric_arr <- function(x) x
REPORT <- function(var) {
var_name <- as.character(sys.call()[[2]])
attr(nll, var_name) <<- if (var_name == "nll")
as.vector(var)
else var
}
g3_cast_vector <- function(x) x
intlookup_getdefault <- function(lookup, key, def) {
if (is.environment(lookup)) {
out <- lookup[[as.character(key)]]
}
else {
out <- if (key < 1)
NULL
else lookup[key][[1]]
}
return(if (is.null(out)) def else out)
}
nonconform_add <- function(base_ar, extra_ar) {
base_ar + as.vector(extra_ar)
}
dif_pmin <- function(a, b, scale) {
dif_pmax(a, b, -scale)
}
nonconform_mult <- function(base_ar, extra_ar) {
base_ar * as.vector(extra_ar)
}
growth_bbinom <- function(delt_l, binn, beta) {
alpha <- (beta * delt_l)/(binn - delt_l)
x <- 0:binn
na <- length(alpha)
n <- length(x) - 1
alpha <- rep(alpha, n + 1)
x <- rep(x, each = na)
val <- exp(lgamma(n + 1) + lgamma(alpha + beta) + lgamma(n - x + beta) + lgamma(x + alpha) - lgamma(n - x + 1) - lgamma(x + 1) - lgamma(n + alpha + beta) - lgamma(beta) - lgamma(alpha))
dim(val) <- c(na, n + 1)
return(val)
}
g3a_grow_vec_rotate <- function(vec, a) {
out <- vapply(seq_len(a), function(i) vec[i:(i + length(vec) - 1)], numeric(length(vec)))
out[is.na(out)] <- vec[length(vec)]
out
}
g3a_grow_vec_extrude <- function(vec, a) {
array(vec, dim = c(length(vec), a))
}
g3a_grow_matrix_wgt <- function(delta_w) {
na <- dim(delta_w)[[1]]
n <- dim(delta_w)[[2]] - 1
wgt.matrix <- array(0, c(na, na))
for (lg in 1:na) {
if (lg == na) {
wgt.matrix[lg, lg:na] <- delta_w[lg, 1:(na - lg + 1)]
}
else if (lg + n > na) {
wgt.matrix[lg, lg:na] <- delta_w[lg, 1:(na - lg + 1)]
}
else {
wgt.matrix[lg, lg:(n + lg)] <- delta_w[lg, ]
}
}
return(wgt.matrix)
}
g3a_grow_matrix_len <- function(delta_l) {
na <- dim(delta_l)[[1]]
n <- dim(delta_l)[[2]] - 1
growth.matrix <- array(0, c(na, na))
for (lg in 1:na) {
if (lg == na) {
growth.matrix[na, na] <- sum(delta_l[lg, ])
}
else if (lg + n > na) {
growth.matrix[lg, lg:(na - 1)] <- delta_l[lg, 1:(na - lg)]
growth.matrix[lg, na] <- sum(delta_l[lg, (na - lg + 1):(n + 1)])
}
else {
growth.matrix[lg, lg:(n + lg)] <- delta_l[lg, ]
}
}
return(growth.matrix)
}
g3a_grow_apply <- function(growth.matrix, wgt.matrix, input_num, input_wgt) {
na <- dim(growth.matrix)[[1]]
growth.matrix <- growth.matrix * as.vector(input_num)
wgt.matrix <- growth.matrix * (wgt.matrix + as.vector(input_wgt))
growth.matrix.sum <- colSums(growth.matrix)
return(array(c(growth.matrix.sum, colSums(wgt.matrix)/avoid_zero(growth.matrix.sum)), dim = c(na, 2)))
}
ratio_add_vec <- function(orig_vec, orig_amount, new_vec, new_amount) {
(orig_vec * orig_amount + new_vec * new_amount)/avoid_zero(orig_amount + new_amount)
}
surveyindices_linreg <- function(N, I, fixed_alpha, fixed_beta) {
meanI <- mean(I)
meanN <- mean(N)
beta <- if (is.nan(fixed_beta))
sum((I - meanI) * (N - meanN))/avoid_zero(sum((N - meanN)^2))
else fixed_beta
alpha <- if (is.nan(fixed_alpha))
meanI - beta * meanN
else fixed_alpha
return(c(alpha, beta))
}
cur_time <- -1L
stopifnot("project_years" %in% names(param))
project_years <- param[["project_years"]]
cur_year <- 0L
start_year <- 1990L
step_lengths <- 12L
step_count <- length(step_lengths)
cur_year_projection <- FALSE
end_year <- 2000L
cur_step <- 0L
cur_step_size <- step_lengths[[1]]/12
cur_step_final <- FALSE
fish__minage <- 1L
fish__maxage <- 10L
fish__area <- 1L
fish__num <- array(0, dim = c(length = 6L, area = 1L, age = 10L), dimnames = list(length = c("50:60", "60:70", "70:80", "80:90", "90:100", "100:Inf"), area = "all", age = c("age1", "age2", "age3", "age4", "age5", "age6", "age7", "age8", "age9", "age10")))
fish__wgt <- array(1, dim = c(length = 6L, area = 1L, age = 10L), dimnames = list(length = c("50:60", "60:70", "70:80", "80:90", "90:100", "100:Inf"), area = "all", age = c("age1", "age2", "age3", "age4", "age5", "age6", "age7", "age8", "age9", "age10")))
as_integer <- as.integer
retro_years <- param[["retro_years"]]
total_steps <- length(step_lengths) * (end_year - as_integer(retro_years) - start_year + as_integer(project_years)) + length(step_lengths) - 1L
dstart_fish__num <- array(0, dim = c(length = 6L, area = 1L, age = 10L, time = as_integer(total_steps + 1)), dimnames = list(length = c("50:60", "60:70", "70:80", "80:90", "90:100", "100:Inf"), area = "all", age = c("age1", "age2", "age3", "age4", "age5", "age6", "age7", "age8", "age9", "age10"), time = attributes(gen_dimnames(param))[["time"]]))
dstart_fish__wgt <- array(1, dim = c(length = 6L, area = 1L, age = 10L, time = as_integer(total_steps + 1)), dimnames = list(length = c("50:60", "60:70", "70:80", "80:90", "90:100", "100:Inf"), area = "all", age = c("age1", "age2", "age3", "age4", "age5", "age6", "age7", "age8", "age9", "age10"), time = attributes(gen_dimnames(param))[["time"]]))
suit_fish_comm__report <- array(NA, dim = c(length = 6L), dimnames = list(length = c("50:60", "60:70", "70:80", "80:90", "90:100", "100:Inf")))
adist_surveyindices_log_acoustic_dist_model__wgt <- array(0, dim = c(length = 1L, time = 11L, area = 1L), dimnames = list(length = "0:Inf", time = c("1990", "1991", "1992", "1993", "1994", "1995", "1996", "1997", "1998", "1999", "2000"), area = "all"))
cdist_sumofsquares_comm_ldist_model__wgt <- array(0, dim = c(length = 5L, time = 11L, area = 1L), dimnames = list(length = c("50:60", "60:70", "70:80", "80:90", "90:Inf"), time = c("1990", "1991", "1992", "1993", "1994", "1995", "1996", "1997", "1998", "1999", "2000"), area = "all"))
detail_fish__predby_comm <- array(NA, dim = c(length = 6L, area = 1L, age = 10L, time = as_integer(total_steps + 1)), dimnames = list(length = c("50:60", "60:70", "70:80", "80:90", "90:100", "100:Inf"), area = "all", age = c("age1", "age2", "age3", "age4", "age5", "age6", "age7", "age8", "age9", "age10"), time = attributes(gen_dimnames(param))[["time"]]))
detail_fish__renewalnum <- array(0, dim = c(length = 6L, area = 1L, age = 10L, time = as_integer(total_steps + 1)), dimnames = list(length = c("50:60", "60:70", "70:80", "80:90", "90:100", "100:Inf"), area = "all", age = c("age1", "age2", "age3", "age4", "age5", "age6", "age7", "age8", "age9", "age10"), time = attributes(gen_dimnames(param))[["time"]]))
detail_fish_comm__cons <- array(NA, dim = c(length = 6L, area = 1L, age = 10L, time = as_integer(total_steps + 1)), dimnames = list(length = c("50:60", "60:70", "70:80", "80:90", "90:100", "100:Inf"), area = "all", age = c("age1", "age2", "age3", "age4", "age5", "age6", "age7", "age8", "age9", "age10"), time = attributes(gen_dimnames(param))[["time"]]))
detail_fish_comm__suit <- array(NA, dim = c(length = 6L, area = 1L, age = 10L, time = as_integer(total_steps + 1)), dimnames = list(length = c("50:60", "60:70", "70:80", "80:90", "90:100", "100:Inf"), area = "all", age = c("age1", "age2", "age3", "age4", "age5", "age6", "age7", "age8", "age9", "age10"), time = attributes(gen_dimnames(param))[["time"]]))
nll <- 0
nll_adist_surveyindices_log_acoustic_dist__weight <- array(0, dim = c(time = as_integer(total_steps + 1L)), dimnames = list(time = attributes(gen_dimnames(param))[["time"]]))
nll_adist_surveyindices_log_acoustic_dist__wgt <- array(0, dim = c(time = as_integer(total_steps + 1L)), dimnames = list(time = attributes(gen_dimnames(param))[["time"]]))
nll_cdist_sumofsquares_comm_ldist__weight <- array(0, dim = c(time = as_integer(total_steps + 1L)), dimnames = list(time = attributes(gen_dimnames(param))[["time"]]))
nll_cdist_sumofsquares_comm_ldist__wgt <- array(0, dim = c(time = as_integer(total_steps + 1L)), dimnames = list(time = attributes(gen_dimnames(param))[["time"]]))
nll_understocking__wgt <- array(0, dim = c(time = as_integer(total_steps + 1L)), dimnames = list(time = attributes(gen_dimnames(param))[["time"]]))
fish__totalpredate <- array(NA, dim = c(length = 6L, area = 1L, age = 10L), dimnames = list(length = c("50:60", "60:70", "70:80", "80:90", "90:100", "100:Inf"), area = "all", age = c("age1", "age2", "age3", "age4", "age5", "age6", "age7", "age8", "age9", "age10")))
comm__totalsuit <- array(NA, dim = c(area = 1L), dimnames = list(area = "all"))
fish_comm__suit <- array(NA, dim = c(length = 6L, area = 1L, age = 10L), dimnames = list(length = c("50:60", "60:70", "70:80", "80:90", "90:100", "100:Inf"), area = "all", age = c("age1", "age2", "age3", "age4", "age5", "age6", "age7", "age8", "age9", "age10")))
comm__area <- 1L
fish_comm__cons <- array(NA, dim = c(length = 6L, area = 1L, age = 10L), dimnames = list(length = c("50:60", "60:70", "70:80", "80:90", "90:100", "100:Inf"), area = "all", age = c("age1", "age2", "age3", "age4", "age5", "age6", "age7", "age8", "age9", "age10")))
intlookup_zip <- function(keys, values) {
if (min(keys) > 0 && max(keys) < 1e+05) {
out <- list()
out[as.integer(keys)] <- as.list(values)
}
else {
out <- as.list(values)
names(out) <- keys
out <- as.environment(out)
}
attr(out, "key_var") <- deparse(sys.call()[[2]])
attr(out, "value_var") <- deparse(sys.call()[[3]])
return(out)
}
comm_landings <- intlookup_zip(comm_landings_keys, comm_landings_values)
fish__consratio <- array(NA, dim = c(length = 6L, area = 1L, age = 10L), dimnames = list(length = c("50:60", "60:70", "70:80", "80:90", "90:100", "100:Inf"), area = "all", age = c("age1", "age2", "age3", "age4", "age5", "age6", "age7", "age8", "age9", "age10")))
fish__overconsumption <- structure(0, desc = "Total overconsumption of fish")
fish__consconv <- array(NA, dim = c(length = 6L, area = 1L, age = 10L), dimnames = list(length = c("50:60", "60:70", "70:80", "80:90", "90:100", "100:Inf"), area = "all", age = c("age1", "age2", "age3", "age4", "age5", "age6", "age7", "age8", "age9", "age10")))
fish__predby_comm <- array(NA, dim = c(length = 6L, area = 1L, age = 10L), dimnames = list(length = c("50:60", "60:70", "70:80", "80:90", "90:100", "100:Inf"), area = "all", age = c("age1", "age2", "age3", "age4", "age5", "age6", "age7", "age8", "age9", "age10")))
fish__growth_lastcalc <- -1L
fish__growth_l <- array(NA, dim = c(length = 6L, delta = 6L), dimnames = list(length = c("50:60", "60:70", "70:80", "80:90", "90:100", "100:Inf"), delta = c("0", "1", "2", "3", "4", "5")))
fish__plusdl <- 10
fish__growth_w <- array(NA, dim = c(length = 6L, delta = 6L), dimnames = list(length = c("50:60", "60:70", "70:80", "80:90", "90:100", "100:Inf"), delta = c("0", "1", "2", "3", "4", "5")))
fish__prevtotal <- 0
fish__renewalnum <- array(0, dim = c(length = 6L, area = 1L, age = 10L), dimnames = list(length = c("50:60", "60:70", "70:80", "80:90", "90:100", "100:Inf"), area = "all", age = c("age1", "age2", "age3", "age4", "age5", "age6", "age7", "age8", "age9", "age10")))
fish__renewalwgt <- array(0, dim = c(length = 6L, area = 1L, age = 10L), dimnames = list(length = c("50:60", "60:70", "70:80", "80:90", "90:100", "100:Inf"), area = "all", age = c("age1", "age2", "age3", "age4", "age5", "age6", "age7", "age8", "age9", "age10")))
adist_surveyindices_log_acoustic_dist_model__area <- 1L
adist_surveyindices_log_acoustic_dist_model__times <- intlookup_zip(adist_surveyindices_log_acoustic_dist_model__times_keys, adist_surveyindices_log_acoustic_dist_model__times_values)
fish_adist_surveyindices_log_acoustic_dist_model_lgmatrix <- array(1, dim = c(1L, 6L), dimnames = NULL)
adist_surveyindices_log_acoustic_dist_obs__area <- 1L
cdist_sumofsquares_comm_ldist_model__area <- 1L
cdist_sumofsquares_comm_ldist_model__times <- intlookup_zip(cdist_sumofsquares_comm_ldist_model__times_keys, cdist_sumofsquares_comm_ldist_model__times_values)
cdist_sumofsquares_comm_ldist_obs__area <- 1L
cdist_sumofsquares_comm_ldist_obs__times <- intlookup_zip(cdist_sumofsquares_comm_ldist_obs__times_keys, cdist_sumofsquares_comm_ldist_obs__times_values)
g3l_understocking_total <- 0
nll_understocking__weight <- array(0, dim = c(time = as_integer(total_steps + 1L)), dimnames = list(time = attributes(gen_dimnames(param))[["time"]]))
while (TRUE) {
{
comment("g3a_time: Start of time period")
cur_time <- cur_time + 1L
if (cur_time == 0 && assert_msg(param[["retro_years"]] >= 0, "retro_years must be >= 0"))
return(NaN)
if (cur_time == 0 && assert_msg(project_years >= 0, "project_years must be >= 0"))
return(NaN)
cur_year <- start_year + (cur_time%/%step_count)
cur_year_projection <- cur_year > end_year - param[["retro_years"]]
cur_step <- (cur_time%%step_count) + 1L
cur_step_size <- step_lengths[[cur_step]]/12
cur_step_final <- cur_step == step_count
}
{
comment("g3a_initialconditions for fish")
for (age in seq(fish__minage, fish__maxage, by = 1)) if (cur_time == 0L) {
fish__age_idx <- age - fish__minage + 1L
area <- fish__area
fish__area_idx <- (1L)
ren_dnorm <- dnorm(fish__midlen, (param[["fish.Linf"]] * (1 - exp(-1 * param[["fish.K"]] * ((age - cur_step_size) - param[["fish.t0"]])))), avoid_zero(((param[["fish.Linf"]] * (1 - exp(-1 * param[["fish.K"]] * ((age - cur_step_size) - param[["fish.t0"]])))) * param[["fish.lencv"]])))
factor <- (param[["fish.init.scalar"]] * nvl(pt.fish.init[[paste(age, sep = ".")]], {
warning("No value found in g3_param_table fish.init, ifmissing not specified")
NaN
}) * exp(-1 * (nvl(pt.fish.M[[paste(age, sep = ".")]], {
warning("No value found in g3_param_table fish.M, ifmissing not specified")
NaN
}) + param[["init.F"]]) * (age - param[["recage"]])))
{
fish__num[, fish__area_idx, fish__age_idx] <- normalize_vec(ren_dnorm) * 10000 * factor
fish__wgt[, fish__area_idx, fish__age_idx] <- param[["fish.walpha"]] * fish__midlen^param[["fish.wbeta"]]
}
}
}
if ((cur_time <= total_steps && param[["report_detail"]] == 1))
dstart_fish__num[, , , cur_time + 1] <- as_numeric_arr(fish__num)
if ((cur_time <= total_steps && param[["report_detail"]] == 1))
dstart_fish__wgt[, , , cur_time + 1] <- as_numeric_arr(fish__wgt)
if (cur_time == 0L) {
suit_fish_comm__report[] <- 1/(1 + exp(-param[["fish.comm.alpha"]] * (fish__midlen - param[["fish.comm.l50"]])))
REPORT(suit_fish_comm__report)
}
if (reporting_enabled > 0L && cur_time > total_steps)
REPORT(adist_surveyindices_log_acoustic_dist_model__params)
if (reporting_enabled > 0L && cur_time > total_steps)
REPORT(adist_surveyindices_log_acoustic_dist_model__wgt)
if (reporting_enabled > 0L && cur_time > total_steps)
REPORT(adist_surveyindices_log_acoustic_dist_obs__wgt)
if (reporting_enabled > 0L && cur_time > total_steps)
REPORT(cdist_sumofsquares_comm_ldist_model__wgt)
if (reporting_enabled > 0L && cur_time > total_steps)
REPORT(cdist_sumofsquares_comm_ldist_obs__wgt)
if (reporting_enabled > 0L && cur_time > total_steps)
REPORT(detail_fish__predby_comm)
if (reporting_enabled > 0L && cur_time > total_steps)
REPORT(detail_fish__renewalnum)
if (reporting_enabled > 0L && cur_time > total_steps)
REPORT(detail_fish_comm__cons)
if (reporting_enabled > 0L && cur_time > total_steps)
REPORT(detail_fish_comm__suit)
if (reporting_enabled > 0L && cur_time > total_steps)
REPORT(dstart_fish__num)
if (reporting_enabled > 0L && cur_time > total_steps)
REPORT(dstart_fish__wgt)
if (reporting_enabled > 0L && cur_time > total_steps)
REPORT(nll)
if (reporting_enabled > 0L && cur_time > total_steps)
REPORT(nll_adist_surveyindices_log_acoustic_dist__weight)
if (reporting_enabled > 0L && cur_time > total_steps)
REPORT(nll_adist_surveyindices_log_acoustic_dist__wgt)
if (reporting_enabled > 0L && cur_time > total_steps)
REPORT(nll_cdist_sumofsquares_comm_ldist__weight)
if (reporting_enabled > 0L && cur_time > total_steps)
REPORT(nll_cdist_sumofsquares_comm_ldist__wgt)
if (reporting_enabled > 0L && cur_time > total_steps)
REPORT(nll_understocking__wgt)
if (reporting_enabled > 0L && cur_time > total_steps)
REPORT(step_lengths)
{
if (cur_time > total_steps)
return(nll)
}
fish__totalpredate[] <- 0
comm__totalsuit[] <- 0
{
suitability <- g3_cast_vector(suit_fish_comm__report[])
{
comment("g3a_predate for comm predating fish")
fish_comm__suit[] <- 0
for (age in seq(fish__minage, fish__maxage, by = 1)) {
fish__age_idx <- age - fish__minage + 1L
area <- fish__area
fish__area_idx <- (1L)
if (area == comm__area) {
catchability <- (suitability * fish__num[, fish__area_idx, fish__age_idx])
comm__area_idx <- (1L)
predator_area <- area
{
comment("Collect all suitable fish biomass for comm")
fish_comm__suit[, fish__area_idx, fish__age_idx] <- catchability
comm__totalsuit[comm__area_idx] <- comm__totalsuit[comm__area_idx] + sum(fish_comm__suit[, fish__area_idx, fish__age_idx])
}
}
}
}
}
{
comment("Scale comm catch of fish by total expected catch")
fish_comm__cons[] <- 0
for (age in seq(fish__minage, fish__maxage, by = 1)) {
fish__age_idx <- age - fish__minage + 1L
area <- fish__area
fish__area_idx <- (1L)
if (area == comm__area) {
comm__area_idx <- (1L)
predator_area <- area
total_predsuit <- comm__totalsuit[comm__area_idx]
fish_comm__cons[, fish__area_idx, fish__age_idx] <- fish_comm__suit[, fish__area_idx, fish__age_idx] * ((if (area != 1L)
0
else intlookup_getdefault(comm_landings, cur_year, 0))/total_predsuit) * fish__wgt[, fish__area_idx, fish__age_idx]
}
}
fish__totalpredate[] <- nonconform_add(fish__totalpredate[], fish_comm__cons[, , ])
}
{
comment("Calculate fish overconsumption coefficient")
comment("Apply overconsumption to fish")
fish__consratio <- fish__totalpredate/avoid_zero(fish__num * fish__wgt)
fish__consratio <- dif_pmin(fish__consratio, 0.95, 1000)
fish__overconsumption <- sum(fish__totalpredate)
fish__consconv <- 1/avoid_zero(fish__totalpredate)
fish__totalpredate <- (fish__num * fish__wgt) * fish__consratio
fish__overconsumption <- fish__overconsumption - sum(fish__totalpredate)
fish__consconv <- fish__consconv * fish__totalpredate
fish__num <- fish__num * (1 - fish__consratio)
}
{
comment("Apply overconsumption to fish_comm__cons")
fish_comm__cons <- nonconform_mult(fish_comm__cons, fish__consconv)
}
{
fish__predby_comm[] <- 0
fish__predby_comm[] <- nonconform_add(fish__predby_comm, fish_comm__cons[, , ])
}
{
comment("Natural mortality for fish")
for (age in seq(fish__minage, fish__maxage, by = 1)) {
fish__age_idx <- age - fish__minage + 1L
area <- fish__area
fish__area_idx <- (1L)
fish__num[, fish__area_idx, fish__age_idx] <- fish__num[, fish__area_idx, fish__age_idx] * exp(-(nvl(pt.fish.M[[paste(age, sep = ".")]], {
warning("No value found in g3_param_table fish.M, ifmissing not specified")
NaN
})) * cur_step_size)
}
}
{
growth_delta_l <- if (fish__growth_lastcalc == floor(cur_step_size * 12L))
fish__growth_l
else (fish__growth_l[] <- growth_bbinom(avoid_zero(avoid_zero((param[["fish.Linf"]] - fish__midlen) * (1 - exp(-(param[["fish.K"]]) * cur_step_size)))/fish__plusdl), 5L, avoid_zero(param[["fish.bbin"]])))
growth_delta_w <- if (fish__growth_lastcalc == floor(cur_step_size * 12L))
fish__growth_w
else (fish__growth_w[] <- (g3a_grow_vec_rotate(fish__midlen^param[["fish.wbeta"]], 5L + 1) - g3a_grow_vec_extrude(fish__midlen^param[["fish.wbeta"]], 5L + 1)) * param[["fish.walpha"]])
growthmat_w <- g3a_grow_matrix_wgt(growth_delta_w)
growthmat_l <- g3a_grow_matrix_len(growth_delta_l)
{
comment("g3a_grow for fish")
for (age in seq(fish__minage, fish__maxage, by = 1)) {
fish__age_idx <- age - fish__minage + 1L
area <- fish__area
fish__area_idx <- (1L)
growthresult <- g3a_grow_apply(growthmat_l, growthmat_w, fish__num[, fish__area_idx, fish__age_idx], fish__wgt[, fish__area_idx, fish__age_idx])
{
if (FALSE)
fish__prevtotal <- sum(fish__num[, fish__area_idx, fish__age_idx])
comment("Update fish using delta matrices")
fish__num[, fish__area_idx, fish__age_idx] <- growthresult[, (1)]
fish__wgt[, fish__area_idx, fish__age_idx] <- growthresult[, (2)]
if (FALSE)
assert_msg(~abs(fish__prevtotal - sum(fish__num[, fish__area_idx, fish__age_idx])) < 1e-04, "g3a_growmature: fish__num totals are not the same before and after growth")
}
}
fish__growth_lastcalc <- floor(cur_step_size * 12L)
}
}
{
factor <- (nvl(pt.fish.rec[[paste(cur_year, sep = ".")]], param[["fish.rec.proj"]]) * param[["fish.rec.scalar"]])
{
comment("g3a_renewal for fish")
for (age in seq(fish__minage, fish__maxage, by = 1)) if (age == fish__minage && cur_step == 1) {
fish__age_idx <- age - fish__minage + 1L
area <- fish__area
fish__area_idx <- (1L)
ren_dnorm <- dnorm(fish__midlen, (param[["fish.Linf"]] * (1 - exp(-1 * param[["fish.K"]] * (age - param[["fish.t0"]])))), avoid_zero(((param[["fish.Linf"]] * (1 - exp(-1 * param[["fish.K"]] * (age - param[["fish.t0"]])))) * param[["fish.lencv"]])))
{
fish__renewalnum[, fish__area_idx, fish__age_idx] <- normalize_vec(ren_dnorm) * 10000 * factor
fish__renewalwgt[, fish__area_idx, fish__age_idx] <- param[["fish.walpha"]] * fish__midlen^param[["fish.wbeta"]]
comment("Add result to fish")
fish__wgt[, fish__area_idx, fish__age_idx] <- ratio_add_vec(fish__wgt[, fish__area_idx, fish__age_idx], fish__num[, fish__area_idx, fish__age_idx], fish__renewalwgt[, fish__area_idx, fish__age_idx], fish__renewalnum[, fish__area_idx, fish__age_idx])
fish__num[, fish__area_idx, fish__age_idx] <- fish__num[, fish__area_idx, fish__age_idx] + fish__renewalnum[, fish__area_idx, fish__age_idx]
}
}
}
}
if (param[["adist_surveyindices_log_acoustic_dist_weight"]] > 0) {
comment("g3l_abundancedistribution_surveyindices_log: Collect abundance from fish for adist_surveyindices_log_acoustic_dist")
for (age in seq(fish__minage, fish__maxage, by = 1)) {
fish__age_idx <- age - fish__minage + 1L
area <- fish__area
fish__area_idx <- (1L)
if (area == adist_surveyindices_log_acoustic_dist_model__area) {
adist_surveyindices_log_acoustic_dist_model__area_idx <- (1L)
adist_surveyindices_log_acoustic_dist_model__time_idx <- intlookup_getdefault(adist_surveyindices_log_acoustic_dist_model__times, (cur_year * 100L + cur_step * 0L), -1L)
if (adist_surveyindices_log_acoustic_dist_model__time_idx >= (1L)) {
comment("Convert fish to wgt")
adist_surveyindices_log_acoustic_dist_model__wgt[, adist_surveyindices_log_acoustic_dist_model__time_idx, adist_surveyindices_log_acoustic_dist_model__area_idx] <- adist_surveyindices_log_acoustic_dist_model__wgt[, adist_surveyindices_log_acoustic_dist_model__time_idx, adist_surveyindices_log_acoustic_dist_model__area_idx] + as.vector(fish_adist_surveyindices_log_acoustic_dist_model_lgmatrix %*% (fish__num[, fish__area_idx, fish__age_idx] * fish__wgt[, fish__area_idx, fish__age_idx]))
}
}
}
}
{
adist_surveyindices_log_acoustic_dist_model__max_time_idx <- (11L)
{
comment("g3l_abundancedistribution_surveyindices_log: Compare adist_surveyindices_log_acoustic_dist_model to adist_surveyindices_log_acoustic_dist_obs")
if (param[["adist_surveyindices_log_acoustic_dist_weight"]] > 0 && cur_step_final) {
area <- adist_surveyindices_log_acoustic_dist_model__area
adist_surveyindices_log_acoustic_dist_model__area_idx <- (1L)
adist_surveyindices_log_acoustic_dist_model__time_idx <- intlookup_getdefault(adist_surveyindices_log_acoustic_dist_model__times, (cur_year * 100L + cur_step * 0L), -1L)
if (adist_surveyindices_log_acoustic_dist_model__time_idx >= (1L))
if (area == adist_surveyindices_log_acoustic_dist_obs__area) {
adist_surveyindices_log_acoustic_dist_obs__area_idx <- (1L)
{
adist_surveyindices_log_acoustic_dist_model__params <- if (adist_surveyindices_log_acoustic_dist_model__time_idx != adist_surveyindices_log_acoustic_dist_model__max_time_idx)
adist_surveyindices_log_acoustic_dist_model__params
else surveyindices_linreg(log(avoid_zero(adist_surveyindices_log_acoustic_dist_model__wgt[, , adist_surveyindices_log_acoustic_dist_model__area_idx])), log(avoid_zero(adist_surveyindices_log_acoustic_dist_obs__wgt[, , adist_surveyindices_log_acoustic_dist_obs__area_idx])), NaN, 1)
{
cur_cdist_nll <- if (adist_surveyindices_log_acoustic_dist_model__time_idx != adist_surveyindices_log_acoustic_dist_model__max_time_idx)
0
else sum((adist_surveyindices_log_acoustic_dist_model__params[[1]] + adist_surveyindices_log_acoustic_dist_model__params[[2]] * log(avoid_zero(adist_surveyindices_log_acoustic_dist_model__wgt[, , adist_surveyindices_log_acoustic_dist_model__area_idx])) - log(avoid_zero(adist_surveyindices_log_acoustic_dist_obs__wgt[, , adist_surveyindices_log_acoustic_dist_obs__area_idx])))^2)
{
nll <- nll + param[["adist_surveyindices_log_acoustic_dist_weight"]] * cur_cdist_nll
nll_adist_surveyindices_log_acoustic_dist__wgt[cur_time + 1L] <- nll_adist_surveyindices_log_acoustic_dist__wgt[cur_time + 1L] + cur_cdist_nll
nll_adist_surveyindices_log_acoustic_dist__weight[cur_time + 1L] <- param[["adist_surveyindices_log_acoustic_dist_weight"]]
}
}
}
}
}
}
}
if (param[["cdist_sumofsquares_comm_ldist_weight"]] > 0) {
comment("g3l_catchdistribution_sumofsquares: Collect catch from comm/fish for cdist_sumofsquares_comm_ldist")
for (age in seq(fish__minage, fish__maxage, by = 1)) {
fish__age_idx <- age - fish__minage + 1L
area <- fish__area
fish__area_idx <- (1L)
if (area == cdist_sumofsquares_comm_ldist_model__area) {
cdist_sumofsquares_comm_ldist_model__area_idx <- (1L)
cdist_sumofsquares_comm_ldist_model__time_idx <- intlookup_getdefault(cdist_sumofsquares_comm_ldist_model__times, (cur_year * 100L + cur_step * 0L), -1L)
if (cdist_sumofsquares_comm_ldist_model__time_idx >= (1L)) {
comment("Convert fish_comm to wgt")
cdist_sumofsquares_comm_ldist_model__wgt[, cdist_sumofsquares_comm_ldist_model__time_idx, cdist_sumofsquares_comm_ldist_model__area_idx] <- cdist_sumofsquares_comm_ldist_model__wgt[, cdist_sumofsquares_comm_ldist_model__time_idx, cdist_sumofsquares_comm_ldist_model__area_idx] + as.vector(fish_comm_cdist_sumofsquares_comm_ldist_model_lgmatrix %*% fish_comm__cons[, fish__area_idx, fish__age_idx])
}
}
}
}
{
comment("g3l_catchdistribution_sumofsquares: Compare cdist_sumofsquares_comm_ldist_model to cdist_sumofsquares_comm_ldist_obs")
if (param[["cdist_sumofsquares_comm_ldist_weight"]] > 0 && cur_step_final) {
area <- cdist_sumofsquares_comm_ldist_model__area
cdist_sumofsquares_comm_ldist_model__area_idx <- (1L)
cdist_sumofsquares_comm_ldist_model__time_idx <- intlookup_getdefault(cdist_sumofsquares_comm_ldist_model__times, (cur_year * 100L + cur_step * 0L), -1L)
if (cdist_sumofsquares_comm_ldist_model__time_idx >= (1L))
if (area == cdist_sumofsquares_comm_ldist_obs__area) {
cdist_sumofsquares_comm_ldist_model__sstotal <- avoid_zero(sum(cdist_sumofsquares_comm_ldist_model__wgt[, cdist_sumofsquares_comm_ldist_model__time_idx, cdist_sumofsquares_comm_ldist_model__area_idx]))
{
cdist_sumofsquares_comm_ldist_obs__area_idx <- (1L)
cdist_sumofsquares_comm_ldist_obs__time_idx <- intlookup_getdefault(cdist_sumofsquares_comm_ldist_obs__times, (cur_year * 100L + cur_step * 0L), -1L)
if (cdist_sumofsquares_comm_ldist_obs__time_idx >= (1L)) {
cdist_sumofsquares_comm_ldist_obs__sstotal <- avoid_zero(sum(cdist_sumofsquares_comm_ldist_obs__wgt[, cdist_sumofsquares_comm_ldist_obs__time_idx, cdist_sumofsquares_comm_ldist_obs__area_idx]))
cur_cdist_nll <- sum((((cdist_sumofsquares_comm_ldist_model__wgt[, cdist_sumofsquares_comm_ldist_model__time_idx, cdist_sumofsquares_comm_ldist_model__area_idx]/cdist_sumofsquares_comm_ldist_model__sstotal) - (cdist_sumofsquares_comm_ldist_obs__wgt[, cdist_sumofsquares_comm_ldist_obs__time_idx, cdist_sumofsquares_comm_ldist_obs__area_idx]/cdist_sumofsquares_comm_ldist_obs__sstotal))^2))
{
nll <- nll + param[["cdist_sumofsquares_comm_ldist_weight"]] * cur_cdist_nll
nll_cdist_sumofsquares_comm_ldist__wgt[cur_time + 1L] <- nll_cdist_sumofsquares_comm_ldist__wgt[cur_time + 1L] + cur_cdist_nll
nll_cdist_sumofsquares_comm_ldist__weight[cur_time + 1L] <- param[["cdist_sumofsquares_comm_ldist_weight"]]
}
}
}
}
}
}
{
comment("Reset understocking total")
g3l_understocking_total <- 0
}
{
comment("g3l_understocking for fish")
comment("Add understocking from fish as biomass to nll")
g3l_understocking_total <- g3l_understocking_total + fish__overconsumption
}
{
comment("g3l_understocking: Combine and add to nll")
g3l_understocking_total <- g3l_understocking_total^2
nll <- nll + 1e+08 * g3l_understocking_total
nll_understocking__wgt[cur_time + 1L] <- nll_understocking__wgt[cur_time + 1L] + g3l_understocking_total
nll_understocking__weight[cur_time + 1L] <- 1e+08
}
if (param[["report_detail"]] == 1)
detail_fish__predby_comm[, , , cur_time + 1] <- as_numeric_arr(fish__predby_comm)
if (param[["report_detail"]] == 1)
detail_fish__renewalnum[, , , cur_time + 1] <- as_numeric_arr(fish__renewalnum)
if (param[["report_detail"]] == 1)
detail_fish_comm__cons[, , , cur_time + 1] <- as_numeric_arr(fish_comm__cons)
if (param[["report_detail"]] == 1)
detail_fish_comm__suit[, , , cur_time + 1] <- as_numeric_arr(fish_comm__suit)
if (cur_step_final) {
comment("g3a_age for fish")
for (age in seq(fish__maxage, fish__minage, by = -1)) {
fish__age_idx <- age - fish__minage + 1L
{
comment("Check stock has remained finite for this step")
if (age == fish__maxage) {
comment("Oldest fish is a plus-group, combine with younger individuals")
fish__wgt[, , fish__age_idx] <- ratio_add_vec(fish__wgt[, , fish__age_idx], fish__num[, , fish__age_idx], fish__wgt[, , fish__age_idx - 1], fish__num[, , fish__age_idx - 1])
fish__num[, , fish__age_idx] <- fish__num[, , fish__age_idx] + fish__num[, , fish__age_idx - 1]
}
else if (age == fish__minage) {
comment("Empty youngest fish age-group")
fish__num[, , fish__age_idx] <- 0
}
else {
comment("Move fish age-group to next one up")
fish__num[, , fish__age_idx] <- fish__num[, , fish__age_idx - 1]
fish__wgt[, , fish__age_idx] <- fish__wgt[, , fish__age_idx - 1]
}
}
}
}
}
}, class = c("g3_r", "function"), parameter_template = list(retro_years = 0, fish.Linf = 1, fish.K = 1, fish.t0 = 0, fish.lencv = 0.1, fish.init.scalar = 1, fish.init.1 = 1, fish.init.2 = 1, fish.init.3 = 1, fish.init.4 = 1, fish.init.5 = 1, fish.init.6 = 1, fish.init.7 = 1, fish.init.8 = 1, fish.init.9 = 1, fish.init.10 = 1, fish.M.1 = 0, fish.M.2 = 0, fish.M.3 = 0, fish.M.4 = 0, fish.M.5 = 0, fish.M.6 = 0, fish.M.7 = 0, fish.M.8 = 0, fish.M.9 = 0, fish.M.10 = 0, init.F = 0, recage = 0, fish.walpha = 0,
fish.wbeta = 0, report_detail = 1L, fish.comm.alpha = 0, fish.comm.l50 = 0, fish.bbin = 0, fish.rec.proj = 0, fish.rec.1990 = 0, fish.rec.1991 = 0, fish.rec.1992 = 0, fish.rec.1993 = 0, fish.rec.1994 = 0, fish.rec.1995 = 0, fish.rec.1996 = 0, fish.rec.1997 = 0, fish.rec.1998 = 0, fish.rec.1999 = 0, fish.rec.2000 = 0, fish.rec.scalar = 0, adist_surveyindices_log_acoustic_dist_weight = 1, cdist_sumofsquares_comm_ldist_weight = 1, project_years = 0))
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.