From 10d03d4c6242170136f7f4fa0bba2d45335f1753 Mon Sep 17 00:00:00 2001 From: ruthmccabe Date: Fri, 10 Jan 2025 10:58:35 +0000 Subject: [PATCH 01/26] children_ind_raw -> is_child; adult_ind_raw -> 1 - is_child --- inst/odin/model-targeted-vax.R | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/inst/odin/model-targeted-vax.R b/inst/odin/model-targeted-vax.R index 36c543f..90d61fb 100644 --- a/inst/odin/model-targeted-vax.R +++ b/inst/odin/model-targeted-vax.R @@ -12,10 +12,10 @@ ## below are indicators for children vs adult groups ## use of _raw indicates that we are in prop form for the boundary category ## of 15 - 19 -children_ind_raw <- parameter() -dim(children_ind_raw) <- c(n_group) -adults_ind_raw <- parameter() -dim(adults_ind_raw) <- c(n_group) +is_child <- parameter() +dim(is_child) <- c(n_group) +# adults_ind_raw <- parameter() +# dim(adults_ind_raw) <- c(n_group) ### setup the number of daily doses ## daily_doses_value has dimensions n_vax x n_time @@ -80,7 +80,7 @@ dim(target_met_children_t) <- c(n_group, n_vax) ## children but leaving it in in case we do expand this to 2 doses in future target_met_children_t[, ] <- 0 target_met_children_t[, 3] <- - ((sum(N[i, 3:4]) * children_ind_raw[i]) > + ((sum(N[i, 3:4]) * is_child[i]) > prioritisation_strategy_children[ i, prioritisation_step_1st_dose_children] * sum(N[i, ])) @@ -93,14 +93,14 @@ dim(target_met_adults_t) <- c(n_group, n_vax) ## for this in the 1st dose target target_met_adults_t[, ] <- 0 target_met_adults_t[, 3] <- - ((sum(N[i, 3:4]) * adults_ind_raw[i]) > + ((sum(N[i, 3:4]) * (1 - is_child[i])) > prioritisation_strategy_adults[i, prioritisation_step_1st_dose_adults] * sum(N[i, ])) ## 2nd doses target_met_adults_t[, 4] <- - ((sum(N[i, 4]) * adults_ind_raw[i]) > + ((sum(N[i, 4]) * (1 - is_child[i])) > prioritisation_strategy_adults[i, prioritisation_step_2nd_dose_adults] * sum(N[i, ])) From 2ade9b553b1a27c9df8b7c66b4564744d5b54a88 Mon Sep 17 00:00:00 2001 From: ruthmccabe Date: Fri, 10 Jan 2025 11:00:12 +0000 Subject: [PATCH 02/26] children_ind_raw -> is_child --- R/parameters.R | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/R/parameters.R b/R/parameters.R index 32e787c..3eada1f 100644 --- a/R/parameters.R +++ b/R/parameters.R @@ -475,8 +475,7 @@ parameters_fixed <- function(region, initial_infections, use_ve_D = FALSE, overr daily_doses_children_time = 1, daily_doses_adults_value = matrix(0, nrow = n_vax, ncol = 1), daily_doses_adults_time = 1, - children_ind_raw = group_bins$children, - adults_ind_raw = group_bins$adults) + is_child = group_bins$children) # Ensure overridden parameters are passed as a list if (!is.list(overrides)) { From ca5a6ff432e77a37a94db3bcfacf26789c2272e3 Mon Sep 17 00:00:00 2001 From: ruthmccabe Date: Fri, 10 Jan 2025 11:15:10 +0000 Subject: [PATCH 03/26] update tests children_ind_raw -> is_child --- tests/testthat/test-model-targeted-vax.R | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/tests/testthat/test-model-targeted-vax.R b/tests/testthat/test-model-targeted-vax.R index 9a3d67e..aa1ad90 100644 --- a/tests/testthat/test-model-targeted-vax.R +++ b/tests/testthat/test-model-targeted-vax.R @@ -365,7 +365,7 @@ test_that("vaccines are only given in the prioritised groups", { ## identify which child groups aren't prioritised for vaccination in the first step idx_novax_children <- (pars$prioritisation_strategy_children[, 1] == 0) & - (pars$children_ind_raw > 0) + (pars$is_child > 0) idx_vax <- c(idx$vax$one_dose, idx$vax$two_dose) if(all(res["prioritisation_step_1st_dose_children", , ] == 1)){ @@ -376,7 +376,7 @@ test_that("vaccines are only given in the prioritised groups", { ## repeat above for adults, including first and second doses idx_novax_adults <- (pars$prioritisation_strategy_adults[, 1] == 0) & - (pars$adults_ind_raw > 0) + ((1 - pars$is_child) > 0) if(all(res["prioritisation_step_1st_dose_adults",,] == 1)){ expect_true(all(y$N[idx_novax_adults, idx_vax, , ] == 0)) @@ -442,7 +442,7 @@ test_that("2nd doses are not given to children", { expect_true(all(res["total_vax_2nddose",,max(t)]>0)) # children_idx - idx_children <- (rep(ceiling(pars$children_ind_raw),pars$n_vax)) *seq(1:(pars$n_group*pars$n_vax)) + idx_children <- (rep(ceiling(pars$is_child),pars$n_vax)) *seq(1:(pars$n_group*pars$n_vax)) # 2nd dose compartments idx_children <- idx_children[which(idx_children!=0&idx_children>=3*pars$n_group)] @@ -530,7 +530,7 @@ test_that("no 2nd doses are given if no 1st doses are given", { res <- dust2::dust_system_simulate(sys, t) rownames(res) <- names(unlist(dust2::dust_unpack_index(sys))) - idx_adults <- ceiling(pars$adults_ind_raw)*seq(1:(pars$n_group*pars$n_vax)) + idx_adults <- ceiling(1-pars$is_child)*seq(1:(pars$n_group*pars$n_vax)) idx_adults <- idx_adults[which(idx_adults!=0&idx_adults>2*pars$n_group)] idx_adults_vax <- paste0(rep("N", each = length(idx_adults)), @@ -567,7 +567,7 @@ test_that("vaccination of children still occurs if adult vaccination is turned o expect_true(all(res["total_vax_2nddose",,max(t)]==0)) # check that the doses aren't in the adult groups - idx_adults <- ceiling(pars$adults_ind_raw)*seq(1:(pars$n_group*pars$n_vax)) + idx_adults <- ceiling(1-pars$is_child)*seq(1:(pars$n_group*pars$n_vax)) idx_adults <- idx_adults[which(idx_adults!=0&idx_adults>2*pars$n_group)] idx_adults_vax <- paste0(rep("N", each = length(idx_adults)), @@ -602,7 +602,7 @@ test_that("vaccination of adults still occurs if child vaccination is turned off expect_true(all(res["total_vax_2nddose",,max(t)]>0)) # check that the doses aren't in the child groups - idx_children <- ceiling(pars$children_ind_raw)*seq(1:(pars$n_group*pars$n_vax)) + idx_children <- ceiling(pars$is_child)*seq(1:(pars$n_group*pars$n_vax)) idx_children <- idx_children[which(idx_children!=0&idx_children>2*pars$n_group)] idx_children_vax <- paste0(rep("N", each = length(idx_children)), From 764866b59385aa4e52741ca43dc9e80b48d1a2f7 Mon Sep 17 00:00:00 2001 From: ruthmccabe Date: Fri, 10 Jan 2025 11:15:22 +0000 Subject: [PATCH 04/26] regenerate model --- R/dust.R | 12 ++++++------ inst/dust/model-targeted-vax.cpp | 29 +++++++++++------------------ src/model-targeted-vax.cpp | 31 ++++++++++++------------------- 3 files changed, 29 insertions(+), 43 deletions(-) diff --git a/R/dust.R b/R/dust.R index ad816bb..b9f5c13 100644 --- a/R/dust.R +++ b/R/dust.R @@ -1,4 +1,4 @@ -## Generated by dust2 (version 0.3.15) - do not edit +## Generated by dust2 (version 0.3.11) - do not edit model_targeted_vax <- structure( function() get("model_targeted_vax"), class = "dust_system_generator", @@ -6,11 +6,11 @@ model_targeted_vax <- structure( package = "mpoxseir", path = NULL, parameters = data.frame( - name = c("children_ind_raw", "adults_ind_raw", "daily_doses_children_value", "daily_doses_children_time", "daily_doses_adults_value", "daily_doses_adults_time", "N_prioritisation_steps_children", "N_prioritisation_steps_adults", "prioritisation_strategy_children", "prioritisation_strategy_adults", "m_gen_pop", "m_sex", "S0", "Ea0", "Eb0", "Ir0", "Id0", "R0", "D0", "beta_h", "beta_s", "beta_z", "beta_hcw", "gamma_E", "gamma_Ir", "gamma_Id", "CFR", "ve_T", "ve_I", "n_vax", "n_group", "exp_noise", "alpha_cases", "alpha_cases_00_04", "alpha_cases_05_14", "alpha_cases_15_plus", "alpha_deaths", "alpha_deaths_00_04", "alpha_deaths_05_14", "alpha_deaths_15_plus"), - type = c("real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "int", "int", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "int", "int", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type"), - constant = c(FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE), - required = c(TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, FALSE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE), - rank = c(1L, 1L, 2L, 1L, 2L, 1L, 0L, 0L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 0L, 0L, 1L, 0L, 0L, 0L, 0L, 2L, 1L, 2L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L)), + name = c("is_child", "daily_doses_children_value", "daily_doses_children_time", "daily_doses_adults_value", "daily_doses_adults_time", "N_prioritisation_steps_children", "N_prioritisation_steps_adults", "prioritisation_strategy_children", "prioritisation_strategy_adults", "m_gen_pop", "m_sex", "S0", "Ea0", "Eb0", "Ir0", "Id0", "R0", "D0", "beta_h", "beta_s", "beta_z", "beta_hcw", "gamma_E", "gamma_Ir", "gamma_Id", "CFR", "ve_T", "ve_I", "n_vax", "n_group", "exp_noise", "alpha_cases", "alpha_cases_00_04", "alpha_cases_05_14", "alpha_cases_15_plus", "alpha_deaths", "alpha_deaths_00_04", "alpha_deaths_05_14", "alpha_deaths_15_plus"), + type = c("real_type", "real_type", "real_type", "real_type", "real_type", "int", "int", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "int", "int", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type", "real_type"), + constant = c(FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE), + required = c(TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, FALSE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE), + rank = c(1L, 2L, 1L, 2L, 1L, 0L, 0L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 0L, 0L, 1L, 0L, 0L, 0L, 0L, 2L, 1L, 2L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L)), properties = list( time_type = "discrete", has_compare = TRUE, diff --git a/inst/dust/model-targeted-vax.cpp b/inst/dust/model-targeted-vax.cpp index cc8ac58..e740c43 100644 --- a/inst/dust/model-targeted-vax.cpp +++ b/inst/dust/model-targeted-vax.cpp @@ -3,8 +3,7 @@ // [[dust2::class(model_targeted_vax)]] // [[dust2::time_type(discrete)]] // [[dust2::has_compare()]] -// [[dust2::parameter(children_ind_raw, type = "real_type", rank = 1, required = TRUE, constant = FALSE)]] -// [[dust2::parameter(adults_ind_raw, type = "real_type", rank = 1, required = TRUE, constant = FALSE)]] +// [[dust2::parameter(is_child, type = "real_type", rank = 1, required = TRUE, constant = FALSE)]] // [[dust2::parameter(daily_doses_children_value, type = "real_type", rank = 2, required = TRUE, constant = FALSE)]] // [[dust2::parameter(daily_doses_children_time, type = "real_type", rank = 1, required = TRUE, constant = FALSE)]] // [[dust2::parameter(daily_doses_adults_value, type = "real_type", rank = 2, required = TRUE, constant = FALSE)]] @@ -58,8 +57,7 @@ class model_targeted_vax { } offset; } odin; struct dim_type { - dust2::array::dimensions<1> children_ind_raw; - dust2::array::dimensions<1> adults_ind_raw; + dust2::array::dimensions<1> is_child; dust2::array::dimensions<2> daily_doses_children_value; dust2::array::dimensions<1> daily_doses_children_time; dust2::array::dimensions<2> daily_doses_adults_value; @@ -171,8 +169,7 @@ class model_targeted_vax { std::vector daily_doses_children_time; std::vector daily_doses_adults_value; std::vector daily_doses_adults_time; - std::vector children_ind_raw; - std::vector adults_ind_raw; + std::vector is_child; dust2::interpolate::InterpolateConstantArray interpolate_daily_doses_children_t; dust2::interpolate::InterpolateConstantArray interpolate_daily_doses_adults_t; std::vector prioritisation_strategy_children; @@ -304,8 +301,7 @@ class model_targeted_vax { const real_type alpha_deaths_00_04 = dust2::r::read_real(parameters, "alpha_deaths_00_04"); const real_type alpha_deaths_05_14 = dust2::r::read_real(parameters, "alpha_deaths_05_14"); const real_type alpha_deaths_15_plus = dust2::r::read_real(parameters, "alpha_deaths_15_plus"); - dim.children_ind_raw.set({static_cast(n_group)}); - dim.adults_ind_raw.set({static_cast(n_group)}); + dim.is_child.set({static_cast(n_group)}); std::vector daily_doses_children_value(dim.daily_doses_children_value.size); dust2::r::read_real_array(parameters, dim.daily_doses_children_value, daily_doses_children_value.data(), "daily_doses_children_value", true); std::vector daily_doses_children_time(dim.daily_doses_children_time.size); @@ -397,10 +393,8 @@ class model_targeted_vax { dim.delta_Ea_n_vaccination.set({static_cast(n_group), static_cast(n_vax)}); dim.delta_Eb_n_vaccination.set({static_cast(n_group), static_cast(n_vax)}); dim.delta_R_n_vaccination.set({static_cast(n_group), static_cast(n_vax)}); - std::vector children_ind_raw(dim.children_ind_raw.size); - dust2::r::read_real_array(parameters, dim.children_ind_raw, children_ind_raw.data(), "children_ind_raw", true); - std::vector adults_ind_raw(dim.adults_ind_raw.size); - dust2::r::read_real_array(parameters, dim.adults_ind_raw, adults_ind_raw.data(), "adults_ind_raw", true); + std::vector is_child(dim.is_child.size); + dust2::r::read_real_array(parameters, dim.is_child, is_child.data(), "is_child", true); const auto interpolate_daily_doses_children_t = dust2::interpolate::InterpolateConstantArray(daily_doses_children_time, daily_doses_children_value, dim.daily_doses_children_t, "daily_doses_children_time", "daily_doses_children_value"); const auto interpolate_daily_doses_adults_t = dust2::interpolate::InterpolateConstantArray(daily_doses_adults_time, daily_doses_adults_value, dim.daily_doses_adults_t, "daily_doses_adults_time", "daily_doses_adults_value"); std::vector prioritisation_strategy_children(dim.prioritisation_strategy_children.size); @@ -554,7 +548,7 @@ class model_targeted_vax { {"cases_cumulative_by_age", std::vector(dim.cases_cumulative_by_age.dim.begin(), dim.cases_cumulative_by_age.dim.end())} }; odin.packing.state.copy_offset(odin.offset.state.begin()); - return shared_state{odin, dim, N_prioritisation_steps_children, N_prioritisation_steps_adults, beta_h, beta_s, beta_hcw, gamma_E, gamma_Ir, gamma_Id, n_vax, n_group, exp_noise, alpha_cases, alpha_cases_00_04, alpha_cases_05_14, alpha_cases_15_plus, alpha_deaths, alpha_deaths_00_04, alpha_deaths_05_14, alpha_deaths_15_plus, daily_doses_children_value, daily_doses_children_time, daily_doses_adults_value, daily_doses_adults_time, children_ind_raw, adults_ind_raw, interpolate_daily_doses_children_t, interpolate_daily_doses_adults_t, prioritisation_strategy_children, prioritisation_strategy_adults, m_gen_pop, m_sex, S0, Ea0, Eb0, Ir0, Id0, R0, D0, beta_z, CFR, ve_T, ve_I, lambda_z}; + return shared_state{odin, dim, N_prioritisation_steps_children, N_prioritisation_steps_adults, beta_h, beta_s, beta_hcw, gamma_E, gamma_Ir, gamma_Id, n_vax, n_group, exp_noise, alpha_cases, alpha_cases_00_04, alpha_cases_05_14, alpha_cases_15_plus, alpha_deaths, alpha_deaths_00_04, alpha_deaths_05_14, alpha_deaths_15_plus, daily_doses_children_value, daily_doses_children_time, daily_doses_adults_value, daily_doses_adults_time, is_child, interpolate_daily_doses_children_t, interpolate_daily_doses_adults_t, prioritisation_strategy_children, prioritisation_strategy_adults, m_gen_pop, m_sex, S0, Ea0, Eb0, Ir0, Id0, R0, D0, beta_z, CFR, ve_T, ve_I, lambda_z}; } static internal_state build_internal(const shared_state& shared) { std::vector n_IrR(shared.dim.n_IrR.size); @@ -662,8 +656,7 @@ class model_targeted_vax { dust2::r::read_real_array(parameters, shared.dim.daily_doses_children_time, shared.daily_doses_children_time.data(), "daily_doses_children_time", false); dust2::r::read_real_array(parameters, shared.dim.daily_doses_adults_value, shared.daily_doses_adults_value.data(), "daily_doses_adults_value", false); dust2::r::read_real_array(parameters, shared.dim.daily_doses_adults_time, shared.daily_doses_adults_time.data(), "daily_doses_adults_time", false); - dust2::r::read_real_array(parameters, shared.dim.children_ind_raw, shared.children_ind_raw.data(), "children_ind_raw", false); - dust2::r::read_real_array(parameters, shared.dim.adults_ind_raw, shared.adults_ind_raw.data(), "adults_ind_raw", false); + dust2::r::read_real_array(parameters, shared.dim.is_child, shared.is_child.data(), "is_child", false); const auto interpolate_daily_doses_children_t = dust2::interpolate::InterpolateConstantArray(shared.daily_doses_children_time, shared.daily_doses_children_value, shared.dim.daily_doses_children_t, "daily_doses_children_time", "daily_doses_children_value"); const auto interpolate_daily_doses_adults_t = dust2::interpolate::InterpolateConstantArray(shared.daily_doses_adults_time, shared.daily_doses_adults_value, shared.dim.daily_doses_adults_t, "daily_doses_adults_time", "daily_doses_adults_value"); dust2::r::read_real_array(parameters, shared.dim.prioritisation_strategy_children, shared.prioritisation_strategy_children.data(), "prioritisation_strategy_children", false); @@ -957,7 +950,7 @@ class model_targeted_vax { } } for (size_t i = 1; i <= shared.dim.target_met_children_t.dim[0]; ++i) { - internal.target_met_children_t[i - 1 + 2 * shared.dim.target_met_children_t.mult[1]] = ((dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {2, 3}) * shared.children_ind_raw[i - 1]) > shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]] * dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1})); + internal.target_met_children_t[i - 1 + 2 * shared.dim.target_met_children_t.mult[1]] = ((dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {2, 3}) * shared.is_child[i - 1]) > shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]] * dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1})); } for (size_t i = 1; i <= shared.dim.target_met_adults_t.dim[0]; ++i) { for (size_t j = 1; j <= shared.dim.target_met_adults_t.dim[1]; ++j) { @@ -965,10 +958,10 @@ class model_targeted_vax { } } for (size_t i = 1; i <= shared.dim.target_met_adults_t.dim[0]; ++i) { - internal.target_met_adults_t[i - 1 + 2 * shared.dim.target_met_adults_t.mult[1]] = ((dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {2, 3}) * shared.adults_ind_raw[i - 1]) > shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_1st_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]] * dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1})); + internal.target_met_adults_t[i - 1 + 2 * shared.dim.target_met_adults_t.mult[1]] = ((dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {2, 3}) * (1 - shared.is_child[i - 1])) > shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_1st_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]] * dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1})); } for (size_t i = 1; i <= shared.dim.target_met_adults_t.dim[0]; ++i) { - internal.target_met_adults_t[i - 1 + 3 * shared.dim.target_met_adults_t.mult[1]] = ((dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {3, 3}) * shared.adults_ind_raw[i - 1]) > shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]] * dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1})); + internal.target_met_adults_t[i - 1 + 3 * shared.dim.target_met_adults_t.mult[1]] = ((dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {3, 3}) * (1 - shared.is_child[i - 1])) > shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]] * dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1})); } for (size_t i = 1; i <= shared.dim.coverage_achieved_1st_dose_children.size; ++i) { internal.coverage_achieved_1st_dose_children[i - 1] = monty::math::ceil(shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]]); diff --git a/src/model-targeted-vax.cpp b/src/model-targeted-vax.cpp index e7306fe..ca2d0b3 100644 --- a/src/model-targeted-vax.cpp +++ b/src/model-targeted-vax.cpp @@ -1,12 +1,11 @@ -// Generated by dust2 (version 0.3.15) - do not edit +// Generated by dust2 (version 0.3.11) - do not edit // Generated by odin2 (version 0.3.17) - do not edit #include // [[dust2::class(model_targeted_vax)]] // [[dust2::time_type(discrete)]] // [[dust2::has_compare()]] -// [[dust2::parameter(children_ind_raw, type = "real_type", rank = 1, required = TRUE, constant = FALSE)]] -// [[dust2::parameter(adults_ind_raw, type = "real_type", rank = 1, required = TRUE, constant = FALSE)]] +// [[dust2::parameter(is_child, type = "real_type", rank = 1, required = TRUE, constant = FALSE)]] // [[dust2::parameter(daily_doses_children_value, type = "real_type", rank = 2, required = TRUE, constant = FALSE)]] // [[dust2::parameter(daily_doses_children_time, type = "real_type", rank = 1, required = TRUE, constant = FALSE)]] // [[dust2::parameter(daily_doses_adults_value, type = "real_type", rank = 2, required = TRUE, constant = FALSE)]] @@ -60,8 +59,7 @@ class model_targeted_vax { } offset; } odin; struct dim_type { - dust2::array::dimensions<1> children_ind_raw; - dust2::array::dimensions<1> adults_ind_raw; + dust2::array::dimensions<1> is_child; dust2::array::dimensions<2> daily_doses_children_value; dust2::array::dimensions<1> daily_doses_children_time; dust2::array::dimensions<2> daily_doses_adults_value; @@ -173,8 +171,7 @@ class model_targeted_vax { std::vector daily_doses_children_time; std::vector daily_doses_adults_value; std::vector daily_doses_adults_time; - std::vector children_ind_raw; - std::vector adults_ind_raw; + std::vector is_child; dust2::interpolate::InterpolateConstantArray interpolate_daily_doses_children_t; dust2::interpolate::InterpolateConstantArray interpolate_daily_doses_adults_t; std::vector prioritisation_strategy_children; @@ -306,8 +303,7 @@ class model_targeted_vax { const real_type alpha_deaths_00_04 = dust2::r::read_real(parameters, "alpha_deaths_00_04"); const real_type alpha_deaths_05_14 = dust2::r::read_real(parameters, "alpha_deaths_05_14"); const real_type alpha_deaths_15_plus = dust2::r::read_real(parameters, "alpha_deaths_15_plus"); - dim.children_ind_raw.set({static_cast(n_group)}); - dim.adults_ind_raw.set({static_cast(n_group)}); + dim.is_child.set({static_cast(n_group)}); std::vector daily_doses_children_value(dim.daily_doses_children_value.size); dust2::r::read_real_array(parameters, dim.daily_doses_children_value, daily_doses_children_value.data(), "daily_doses_children_value", true); std::vector daily_doses_children_time(dim.daily_doses_children_time.size); @@ -399,10 +395,8 @@ class model_targeted_vax { dim.delta_Ea_n_vaccination.set({static_cast(n_group), static_cast(n_vax)}); dim.delta_Eb_n_vaccination.set({static_cast(n_group), static_cast(n_vax)}); dim.delta_R_n_vaccination.set({static_cast(n_group), static_cast(n_vax)}); - std::vector children_ind_raw(dim.children_ind_raw.size); - dust2::r::read_real_array(parameters, dim.children_ind_raw, children_ind_raw.data(), "children_ind_raw", true); - std::vector adults_ind_raw(dim.adults_ind_raw.size); - dust2::r::read_real_array(parameters, dim.adults_ind_raw, adults_ind_raw.data(), "adults_ind_raw", true); + std::vector is_child(dim.is_child.size); + dust2::r::read_real_array(parameters, dim.is_child, is_child.data(), "is_child", true); const auto interpolate_daily_doses_children_t = dust2::interpolate::InterpolateConstantArray(daily_doses_children_time, daily_doses_children_value, dim.daily_doses_children_t, "daily_doses_children_time", "daily_doses_children_value"); const auto interpolate_daily_doses_adults_t = dust2::interpolate::InterpolateConstantArray(daily_doses_adults_time, daily_doses_adults_value, dim.daily_doses_adults_t, "daily_doses_adults_time", "daily_doses_adults_value"); std::vector prioritisation_strategy_children(dim.prioritisation_strategy_children.size); @@ -556,7 +550,7 @@ class model_targeted_vax { {"cases_cumulative_by_age", std::vector(dim.cases_cumulative_by_age.dim.begin(), dim.cases_cumulative_by_age.dim.end())} }; odin.packing.state.copy_offset(odin.offset.state.begin()); - return shared_state{odin, dim, N_prioritisation_steps_children, N_prioritisation_steps_adults, beta_h, beta_s, beta_hcw, gamma_E, gamma_Ir, gamma_Id, n_vax, n_group, exp_noise, alpha_cases, alpha_cases_00_04, alpha_cases_05_14, alpha_cases_15_plus, alpha_deaths, alpha_deaths_00_04, alpha_deaths_05_14, alpha_deaths_15_plus, daily_doses_children_value, daily_doses_children_time, daily_doses_adults_value, daily_doses_adults_time, children_ind_raw, adults_ind_raw, interpolate_daily_doses_children_t, interpolate_daily_doses_adults_t, prioritisation_strategy_children, prioritisation_strategy_adults, m_gen_pop, m_sex, S0, Ea0, Eb0, Ir0, Id0, R0, D0, beta_z, CFR, ve_T, ve_I, lambda_z}; + return shared_state{odin, dim, N_prioritisation_steps_children, N_prioritisation_steps_adults, beta_h, beta_s, beta_hcw, gamma_E, gamma_Ir, gamma_Id, n_vax, n_group, exp_noise, alpha_cases, alpha_cases_00_04, alpha_cases_05_14, alpha_cases_15_plus, alpha_deaths, alpha_deaths_00_04, alpha_deaths_05_14, alpha_deaths_15_plus, daily_doses_children_value, daily_doses_children_time, daily_doses_adults_value, daily_doses_adults_time, is_child, interpolate_daily_doses_children_t, interpolate_daily_doses_adults_t, prioritisation_strategy_children, prioritisation_strategy_adults, m_gen_pop, m_sex, S0, Ea0, Eb0, Ir0, Id0, R0, D0, beta_z, CFR, ve_T, ve_I, lambda_z}; } static internal_state build_internal(const shared_state& shared) { std::vector n_IrR(shared.dim.n_IrR.size); @@ -664,8 +658,7 @@ class model_targeted_vax { dust2::r::read_real_array(parameters, shared.dim.daily_doses_children_time, shared.daily_doses_children_time.data(), "daily_doses_children_time", false); dust2::r::read_real_array(parameters, shared.dim.daily_doses_adults_value, shared.daily_doses_adults_value.data(), "daily_doses_adults_value", false); dust2::r::read_real_array(parameters, shared.dim.daily_doses_adults_time, shared.daily_doses_adults_time.data(), "daily_doses_adults_time", false); - dust2::r::read_real_array(parameters, shared.dim.children_ind_raw, shared.children_ind_raw.data(), "children_ind_raw", false); - dust2::r::read_real_array(parameters, shared.dim.adults_ind_raw, shared.adults_ind_raw.data(), "adults_ind_raw", false); + dust2::r::read_real_array(parameters, shared.dim.is_child, shared.is_child.data(), "is_child", false); const auto interpolate_daily_doses_children_t = dust2::interpolate::InterpolateConstantArray(shared.daily_doses_children_time, shared.daily_doses_children_value, shared.dim.daily_doses_children_t, "daily_doses_children_time", "daily_doses_children_value"); const auto interpolate_daily_doses_adults_t = dust2::interpolate::InterpolateConstantArray(shared.daily_doses_adults_time, shared.daily_doses_adults_value, shared.dim.daily_doses_adults_t, "daily_doses_adults_time", "daily_doses_adults_value"); dust2::r::read_real_array(parameters, shared.dim.prioritisation_strategy_children, shared.prioritisation_strategy_children.data(), "prioritisation_strategy_children", false); @@ -959,7 +952,7 @@ class model_targeted_vax { } } for (size_t i = 1; i <= shared.dim.target_met_children_t.dim[0]; ++i) { - internal.target_met_children_t[i - 1 + 2 * shared.dim.target_met_children_t.mult[1]] = ((dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {2, 3}) * shared.children_ind_raw[i - 1]) > shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]] * dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1})); + internal.target_met_children_t[i - 1 + 2 * shared.dim.target_met_children_t.mult[1]] = ((dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {2, 3}) * shared.is_child[i - 1]) > shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]] * dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1})); } for (size_t i = 1; i <= shared.dim.target_met_adults_t.dim[0]; ++i) { for (size_t j = 1; j <= shared.dim.target_met_adults_t.dim[1]; ++j) { @@ -967,10 +960,10 @@ class model_targeted_vax { } } for (size_t i = 1; i <= shared.dim.target_met_adults_t.dim[0]; ++i) { - internal.target_met_adults_t[i - 1 + 2 * shared.dim.target_met_adults_t.mult[1]] = ((dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {2, 3}) * shared.adults_ind_raw[i - 1]) > shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_1st_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]] * dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1})); + internal.target_met_adults_t[i - 1 + 2 * shared.dim.target_met_adults_t.mult[1]] = ((dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {2, 3}) * (1 - shared.is_child[i - 1])) > shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_1st_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]] * dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1})); } for (size_t i = 1; i <= shared.dim.target_met_adults_t.dim[0]; ++i) { - internal.target_met_adults_t[i - 1 + 3 * shared.dim.target_met_adults_t.mult[1]] = ((dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {3, 3}) * shared.adults_ind_raw[i - 1]) > shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]] * dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1})); + internal.target_met_adults_t[i - 1 + 3 * shared.dim.target_met_adults_t.mult[1]] = ((dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {3, 3}) * (1 - shared.is_child[i - 1])) > shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]] * dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1})); } for (size_t i = 1; i <= shared.dim.coverage_achieved_1st_dose_children.size; ++i) { internal.coverage_achieved_1st_dose_children[i - 1] = monty::math::ceil(shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]]); From ba3f9bd20ea43b4db1e1070591e4abd1521e0ac4 Mon Sep 17 00:00:00 2001 From: ruthmccabe Date: Fri, 10 Jan 2025 11:39:44 +0000 Subject: [PATCH 05/26] FOI tidy - instead of having ve_I included in each of the 4 components separately this is now moved to the line where they are combined --- inst/odin/model-targeted-vax.R | 14 +++++--------- 1 file changed, 5 insertions(+), 9 deletions(-) diff --git a/inst/odin/model-targeted-vax.R b/inst/odin/model-targeted-vax.R index 90d61fb..96930ff 100644 --- a/inst/odin/model-targeted-vax.R +++ b/inst/odin/model-targeted-vax.R @@ -10,12 +10,8 @@ ## Emergency use of the MVA-BN has not been granted for children (u18s) in DRC ## and so we split vaccination into two components depending on age ## below are indicators for children vs adult groups -## use of _raw indicates that we are in prop form for the boundary category -## of 15 - 19 is_child <- parameter() dim(is_child) <- c(n_group) -# adults_ind_raw <- parameter() -# dim(adults_ind_raw) <- c(n_group) ### setup the number of daily doses ## daily_doses_value has dimensions n_vax x n_time @@ -642,15 +638,15 @@ s_ij_gen_pop[, ] <- m_gen_pop[i, j] * prop_infectious[j] # as above but for the sexual contacts only s_ij_sex[, ] <- m_sex[i, j] * prop_infectious[j] -lambda_hh[, ] <- beta_h * sum(s_ij_gen_pop[i, ]) * (1 - ve_I[i, j]) -lambda_s[, ] <- beta_s * sum(s_ij_sex[i, ]) * (1 - ve_I[i, j]) +lambda_hh[, ] <- beta_h * sum(s_ij_gen_pop[i, ]) +lambda_s[, ] <- beta_s * sum(s_ij_sex[i, ]) # additional foi in HCW only (i = 20) homogeneous from infected as assumed equally # likely to attend hospital lambda_hc[, ] <- - if (i == 20) beta_hcw * sum(I_infectious[, ]) * (1 - ve_I[i, j]) else 0 -lambda_z[, ] <- beta_z[i] * (1 - ve_I[i, j]) + if (i == 20) beta_hcw * sum(I_infectious[, ]) else 0 +lambda_z[, ] <- beta_z[i] -lambda[, ] <- lambda_hh[i, j] + lambda_s[i, j] + lambda_hc[i, j] + lambda_z[i, j] +lambda[, ] <- (lambda_hh[i, j] + lambda_s[i, j] + lambda_hc[i, j] + lambda_z[i, j]) * (1 - ve_I[i, j]) ## Draws from binomial distributions for numbers changing between compartments # accounting for vaccination: From eb210acc6d5212e33fc911aee37ed1421babe7ce Mon Sep 17 00:00:00 2001 From: ruthmccabe Date: Fri, 10 Jan 2025 11:55:57 +0000 Subject: [PATCH 06/26] coverage_achieved -> coverage_target as I think more intuitive --- inst/odin/model-targeted-vax.R | 19 +++++++++---------- 1 file changed, 9 insertions(+), 10 deletions(-) diff --git a/inst/odin/model-targeted-vax.R b/inst/odin/model-targeted-vax.R index 96930ff..a802835 100644 --- a/inst/odin/model-targeted-vax.R +++ b/inst/odin/model-targeted-vax.R @@ -100,7 +100,6 @@ target_met_adults_t[, 4] <- prioritisation_strategy_adults[i, prioritisation_step_2nd_dose_adults] * sum(N[i, ])) - ## prioritisation step proposal to account for the fact that this would update ## every single time step if we vaccinate quickly enough (unlikely but would ## like it to be done properly) @@ -108,22 +107,22 @@ target_met_adults_t[, 4] <- ## the target goal for the current prioritisation strategy ## with ceiling basically saying if there is any target in that group to account ## for it -coverage_achieved_1st_dose_children[] <- ceiling( +coverage_target_1st_dose_children[] <- ceiling( prioritisation_strategy_children[i, prioritisation_step_1st_dose_children]) -dim(coverage_achieved_1st_dose_children) <- c(n_group) +dim(coverage_target_1st_dose_children) <- c(n_group) -coverage_achieved_1st_dose_adults[] <- ceiling( +coverage_target_1st_dose_adults[] <- ceiling( prioritisation_strategy_adults[i, prioritisation_step_1st_dose_adults]) -dim(coverage_achieved_1st_dose_adults) <- c(n_group) +dim(coverage_target_1st_dose_adults) <- c(n_group) -coverage_achieved_2nd_dose_adults[] <- ceiling( +coverage_target_2nd_dose_adults[] <- ceiling( prioritisation_strategy_adults[i, prioritisation_step_2nd_dose_adults]) -dim(coverage_achieved_2nd_dose_adults) <- c(n_group) +dim(coverage_target_2nd_dose_adults) <- c(n_group) ## children prioritisation_step_1st_dose_children_proposal <- if (sum(target_met_children_t[, 3]) == - sum(coverage_achieved_1st_dose_children[])) + sum(coverage_target_1st_dose_children[])) prioritisation_step_1st_dose_children + 1 else prioritisation_step_1st_dose_children @@ -136,7 +135,7 @@ update(prioritisation_step_1st_dose_children) <- ## adults prioritisation_step_1st_dose_adults_proposal <- - if (sum(target_met_adults_t[, 3]) == sum(coverage_achieved_1st_dose_adults[])) + if (sum(target_met_adults_t[, 3]) == sum(coverage_target_1st_dose_adults[])) prioritisation_step_1st_dose_adults + 1 else prioritisation_step_1st_dose_adults @@ -147,7 +146,7 @@ update(prioritisation_step_1st_dose_adults) <- prioritisation_step_1st_dose_adults_proposal prioritisation_step_2nd_dose_adults_proposal <- - if (sum(target_met_adults_t[, 4]) == sum(coverage_achieved_2nd_dose_adults[])) + if (sum(target_met_adults_t[, 4]) == sum(coverage_target_2nd_dose_adults[])) prioritisation_step_2nd_dose_adults + 1 else prioritisation_step_2nd_dose_adults From 87764c9be67875c04d8e496875e2a2e7ad9dce46 Mon Sep 17 00:00:00 2001 From: ruthmccabe Date: Fri, 10 Jan 2025 13:32:08 +0000 Subject: [PATCH 07/26] simplify target_met_children_t (can be vector not matrix) --- inst/odin/model-targeted-vax.R | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/inst/odin/model-targeted-vax.R b/inst/odin/model-targeted-vax.R index a802835..c42a436 100644 --- a/inst/odin/model-targeted-vax.R +++ b/inst/odin/model-targeted-vax.R @@ -68,14 +68,14 @@ dim(prioritisation_strategy_adults) <- c(n_group, N_prioritisation_steps_adults) ## neither children or adults, for children col 4 (2nd dose) also doesn't matter ## children -dim(target_met_children_t) <- c(n_group, n_vax) +dim(target_met_children_t) <- c(n_group) ## 1st doses ## if you have a 2nd dose this implies you also have had a 1st dose so account ## for this in the 1st dose target - this now isn't strictly relevant for ## children but leaving it in in case we do expand this to 2 doses in future -target_met_children_t[, ] <- 0 -target_met_children_t[, 3] <- +target_met_children_t[] <- 0 +target_met_children_t[i] <- ((sum(N[i, 3:4]) * is_child[i]) > prioritisation_strategy_children[ i, prioritisation_step_1st_dose_children] * sum(N[i, ])) From 2d8ee57100a0af6436c43c453dda0c9e39252c57 Mon Sep 17 00:00:00 2001 From: ruthmccabe Date: Fri, 10 Jan 2025 13:32:50 +0000 Subject: [PATCH 08/26] simplify target_met_children_t --- inst/odin/model-targeted-vax.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/inst/odin/model-targeted-vax.R b/inst/odin/model-targeted-vax.R index c42a436..b96b084 100644 --- a/inst/odin/model-targeted-vax.R +++ b/inst/odin/model-targeted-vax.R @@ -121,7 +121,7 @@ dim(coverage_target_2nd_dose_adults) <- c(n_group) ## children prioritisation_step_1st_dose_children_proposal <- - if (sum(target_met_children_t[, 3]) == + if (sum(target_met_children_t[]) == sum(coverage_target_1st_dose_children[])) prioritisation_step_1st_dose_children + 1 else prioritisation_step_1st_dose_children From 4e240f89c92faf9c0f4b8692ab47c62ca7f4cffa Mon Sep 17 00:00:00 2001 From: ruthmccabe Date: Fri, 10 Jan 2025 14:03:08 +0000 Subject: [PATCH 09/26] undo FOI tidy which it did NOT like; regenerate model --- inst/dust/model-targeted-vax.cpp | 60 +++++++++++++++----------------- inst/odin/model-targeted-vax.R | 29 ++++++++++----- src/model-targeted-vax.cpp | 60 +++++++++++++++----------------- 3 files changed, 79 insertions(+), 70 deletions(-) diff --git a/inst/dust/model-targeted-vax.cpp b/inst/dust/model-targeted-vax.cpp index e740c43..4284ced 100644 --- a/inst/dust/model-targeted-vax.cpp +++ b/inst/dust/model-targeted-vax.cpp @@ -66,11 +66,11 @@ class model_targeted_vax { dust2::array::dimensions<1> daily_doses_adults_t; dust2::array::dimensions<2> prioritisation_strategy_children; dust2::array::dimensions<2> prioritisation_strategy_adults; - dust2::array::dimensions<2> target_met_children_t; + dust2::array::dimensions<1> target_met_children_t; dust2::array::dimensions<2> target_met_adults_t; - dust2::array::dimensions<1> coverage_achieved_1st_dose_children; - dust2::array::dimensions<1> coverage_achieved_1st_dose_adults; - dust2::array::dimensions<1> coverage_achieved_2nd_dose_adults; + dust2::array::dimensions<1> coverage_target_1st_dose_children; + dust2::array::dimensions<1> coverage_target_1st_dose_adults; + dust2::array::dimensions<1> coverage_target_2nd_dose_adults; dust2::array::dimensions<1> n_eligible_for_dose1_children; dust2::array::dimensions<1> n_eligible_for_dose1_adults; dust2::array::dimensions<1> n_eligible_for_dose2_adults; @@ -194,9 +194,9 @@ class model_targeted_vax { std::vector n_IdD; std::vector target_met_children_t; std::vector target_met_adults_t; - std::vector coverage_achieved_1st_dose_children; - std::vector coverage_achieved_1st_dose_adults; - std::vector coverage_achieved_2nd_dose_adults; + std::vector coverage_target_1st_dose_children; + std::vector coverage_target_1st_dose_adults; + std::vector coverage_target_2nd_dose_adults; std::vector n_eligible_for_dose1_children; std::vector n_eligible_for_dose1_adults; std::vector n_eligible_for_dose2_adults; @@ -314,11 +314,11 @@ class model_targeted_vax { dim.daily_doses_adults_t.set({static_cast(n_vax)}); dim.prioritisation_strategy_children.set({static_cast(n_group), static_cast(N_prioritisation_steps_children)}); dim.prioritisation_strategy_adults.set({static_cast(n_group), static_cast(N_prioritisation_steps_adults)}); - dim.target_met_children_t.set({static_cast(n_group), static_cast(n_vax)}); + dim.target_met_children_t.set({static_cast(n_group)}); dim.target_met_adults_t.set({static_cast(n_group), static_cast(n_vax)}); - dim.coverage_achieved_1st_dose_children.set({static_cast(n_group)}); - dim.coverage_achieved_1st_dose_adults.set({static_cast(n_group)}); - dim.coverage_achieved_2nd_dose_adults.set({static_cast(n_group)}); + dim.coverage_target_1st_dose_children.set({static_cast(n_group)}); + dim.coverage_target_1st_dose_adults.set({static_cast(n_group)}); + dim.coverage_target_2nd_dose_adults.set({static_cast(n_group)}); dim.n_eligible_for_dose1_children.set({static_cast(n_group)}); dim.n_eligible_for_dose1_adults.set({static_cast(n_group)}); dim.n_eligible_for_dose2_adults.set({static_cast(n_group)}); @@ -555,9 +555,9 @@ class model_targeted_vax { std::vector n_IdD(shared.dim.n_IdD.size); std::vector target_met_children_t(shared.dim.target_met_children_t.size); std::vector target_met_adults_t(shared.dim.target_met_adults_t.size); - std::vector coverage_achieved_1st_dose_children(shared.dim.coverage_achieved_1st_dose_children.size); - std::vector coverage_achieved_1st_dose_adults(shared.dim.coverage_achieved_1st_dose_adults.size); - std::vector coverage_achieved_2nd_dose_adults(shared.dim.coverage_achieved_2nd_dose_adults.size); + std::vector coverage_target_1st_dose_children(shared.dim.coverage_target_1st_dose_children.size); + std::vector coverage_target_1st_dose_adults(shared.dim.coverage_target_1st_dose_adults.size); + std::vector coverage_target_2nd_dose_adults(shared.dim.coverage_target_2nd_dose_adults.size); std::vector n_eligible_for_dose1_children(shared.dim.n_eligible_for_dose1_children.size); std::vector n_eligible_for_dose1_adults(shared.dim.n_eligible_for_dose1_adults.size); std::vector n_eligible_for_dose2_adults(shared.dim.n_eligible_for_dose2_adults.size); @@ -617,7 +617,7 @@ class model_targeted_vax { std::vector new_N(shared.dim.N.size); std::vector n_SEa_hc(shared.dim.n_SEa_hc.size); std::vector n_SEa_z(shared.dim.n_SEa_z.size); - return internal_state{n_IrR, n_IdD, target_met_children_t, target_met_adults_t, coverage_achieved_1st_dose_children, coverage_achieved_1st_dose_adults, coverage_achieved_2nd_dose_adults, n_eligible_for_dose1_children, n_eligible_for_dose1_adults, n_eligible_for_dose2_adults, I_infectious, delta_R, delta_D, daily_doses_children_t, daily_doses_adults_t, n_vaccination_t_S_children, n_vaccination_t_S_adults, n_vaccination_t_Ea_children, n_vaccination_t_Ea_adults, n_vaccination_t_Eb_children, n_vaccination_t_Eb_adults, n_vaccination_t_R_children, n_vaccination_t_R_adults, new_D, prop_infectious, lambda_hc, n_vaccination_t_S, n_vaccination_t_Ea, n_vaccination_t_Eb, n_vaccination_t_R, s_ij_gen_pop, s_ij_sex, delta_S_n_vaccination, delta_Ea_n_vaccination, delta_Eb_n_vaccination, delta_R_n_vaccination, n_vaccination_t, lambda_hh, lambda_s, new_R, lambda, n_EaEb, n_EbI, p_SE, p_hh, p_s, p_hc, n_EbId, delta_Eb, new_Eb, n_SEa, n_EbIr, delta_Id, new_S, new_Id, n_SEa_hh, delta_Ea, delta_Ir, new_Ea, new_Ir, n_SEa_s, new_E, new_I, new_N, n_SEa_hc, n_SEa_z}; + return internal_state{n_IrR, n_IdD, target_met_children_t, target_met_adults_t, coverage_target_1st_dose_children, coverage_target_1st_dose_adults, coverage_target_2nd_dose_adults, n_eligible_for_dose1_children, n_eligible_for_dose1_adults, n_eligible_for_dose2_adults, I_infectious, delta_R, delta_D, daily_doses_children_t, daily_doses_adults_t, n_vaccination_t_S_children, n_vaccination_t_S_adults, n_vaccination_t_Ea_children, n_vaccination_t_Ea_adults, n_vaccination_t_Eb_children, n_vaccination_t_Eb_adults, n_vaccination_t_R_children, n_vaccination_t_R_adults, new_D, prop_infectious, lambda_hc, n_vaccination_t_S, n_vaccination_t_Ea, n_vaccination_t_Eb, n_vaccination_t_R, s_ij_gen_pop, s_ij_sex, delta_S_n_vaccination, delta_Ea_n_vaccination, delta_Eb_n_vaccination, delta_R_n_vaccination, n_vaccination_t, lambda_hh, lambda_s, new_R, lambda, n_EaEb, n_EbI, p_SE, p_hh, p_s, p_hc, n_EbId, delta_Eb, new_Eb, n_SEa, n_EbIr, delta_Id, new_S, new_Id, n_SEa_hh, delta_Ea, delta_Ir, new_Ea, new_Ir, n_SEa_s, new_E, new_I, new_N, n_SEa_hc, n_SEa_z}; } static data_type build_data(cpp11::list data, const shared_state& shared) { auto cases = dust2::r::read_real(data, "cases", NA_REAL); @@ -944,13 +944,11 @@ class model_targeted_vax { internal.n_IdD[i - 1 + (j - 1) * shared.dim.n_IdD.mult[1]] = monty::random::binomial(rng_state, Id[i - 1 + (j - 1) * shared.dim.Id.mult[1]], p_IdD); } } - for (size_t i = 1; i <= shared.dim.target_met_children_t.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.target_met_children_t.dim[1]; ++j) { - internal.target_met_children_t[i - 1 + (j - 1) * shared.dim.target_met_children_t.mult[1]] = 0; - } + for (size_t i = 1; i <= shared.dim.target_met_children_t.size; ++i) { + internal.target_met_children_t[i - 1] = 0; } - for (size_t i = 1; i <= shared.dim.target_met_children_t.dim[0]; ++i) { - internal.target_met_children_t[i - 1 + 2 * shared.dim.target_met_children_t.mult[1]] = ((dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {2, 3}) * shared.is_child[i - 1]) > shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]] * dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1})); + for (size_t i = 1; i <= shared.dim.target_met_children_t.size; ++i) { + internal.target_met_children_t[i - 1] = ((dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {2, 3}) * shared.is_child[i - 1]) > shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]] * dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1})); } for (size_t i = 1; i <= shared.dim.target_met_adults_t.dim[0]; ++i) { for (size_t j = 1; j <= shared.dim.target_met_adults_t.dim[1]; ++j) { @@ -963,14 +961,14 @@ class model_targeted_vax { for (size_t i = 1; i <= shared.dim.target_met_adults_t.dim[0]; ++i) { internal.target_met_adults_t[i - 1 + 3 * shared.dim.target_met_adults_t.mult[1]] = ((dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {3, 3}) * (1 - shared.is_child[i - 1])) > shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]] * dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1})); } - for (size_t i = 1; i <= shared.dim.coverage_achieved_1st_dose_children.size; ++i) { - internal.coverage_achieved_1st_dose_children[i - 1] = monty::math::ceil(shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]]); + for (size_t i = 1; i <= shared.dim.coverage_target_1st_dose_children.size; ++i) { + internal.coverage_target_1st_dose_children[i - 1] = monty::math::ceil(shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]]); } - for (size_t i = 1; i <= shared.dim.coverage_achieved_1st_dose_adults.size; ++i) { - internal.coverage_achieved_1st_dose_adults[i - 1] = monty::math::ceil(shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_1st_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]); + for (size_t i = 1; i <= shared.dim.coverage_target_1st_dose_adults.size; ++i) { + internal.coverage_target_1st_dose_adults[i - 1] = monty::math::ceil(shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_1st_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]); } - for (size_t i = 1; i <= shared.dim.coverage_achieved_2nd_dose_adults.size; ++i) { - internal.coverage_achieved_2nd_dose_adults[i - 1] = monty::math::ceil(shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]); + for (size_t i = 1; i <= shared.dim.coverage_target_2nd_dose_adults.size; ++i) { + internal.coverage_target_2nd_dose_adults[i - 1] = monty::math::ceil(shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]); } for (size_t i = 1; i <= shared.dim.n_eligible_for_dose1_children.size; ++i) { internal.n_eligible_for_dose1_children[i - 1] = (S[i - 1 + shared.dim.S.mult[1]] + Ea[i - 1 + shared.dim.Ea.mult[1]] + Eb[i - 1 + shared.dim.Eb.mult[1]] + R[i - 1 + shared.dim.R.mult[1]]) * shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]]; @@ -1004,9 +1002,9 @@ class model_targeted_vax { } shared.interpolate_daily_doses_children_t.eval(time, internal.daily_doses_children_t); shared.interpolate_daily_doses_adults_t.eval(time, internal.daily_doses_adults_t); - const real_type prioritisation_step_1st_dose_children_proposal = (dust2::array::sum(internal.target_met_children_t.data(), shared.dim.target_met_children_t, {0, shared.dim.target_met_children_t.dim[0] - 1}, {2, 2}) == dust2::array::sum(internal.coverage_achieved_1st_dose_children.data(), shared.dim.coverage_achieved_1st_dose_children) ? prioritisation_step_1st_dose_children + 1 : prioritisation_step_1st_dose_children); - const real_type prioritisation_step_1st_dose_adults_proposal = (dust2::array::sum(internal.target_met_adults_t.data(), shared.dim.target_met_adults_t, {0, shared.dim.target_met_adults_t.dim[0] - 1}, {2, 2}) == dust2::array::sum(internal.coverage_achieved_1st_dose_adults.data(), shared.dim.coverage_achieved_1st_dose_adults) ? prioritisation_step_1st_dose_adults + 1 : prioritisation_step_1st_dose_adults); - const real_type prioritisation_step_2nd_dose_adults_proposal = (dust2::array::sum(internal.target_met_adults_t.data(), shared.dim.target_met_adults_t, {0, shared.dim.target_met_adults_t.dim[0] - 1}, {3, 3}) == dust2::array::sum(internal.coverage_achieved_2nd_dose_adults.data(), shared.dim.coverage_achieved_2nd_dose_adults) ? prioritisation_step_2nd_dose_adults + 1 : prioritisation_step_2nd_dose_adults); + const real_type prioritisation_step_1st_dose_children_proposal = (dust2::array::sum(internal.target_met_children_t.data(), shared.dim.target_met_children_t) == dust2::array::sum(internal.coverage_target_1st_dose_children.data(), shared.dim.coverage_target_1st_dose_children) ? prioritisation_step_1st_dose_children + 1 : prioritisation_step_1st_dose_children); + const real_type prioritisation_step_1st_dose_adults_proposal = (dust2::array::sum(internal.target_met_adults_t.data(), shared.dim.target_met_adults_t, {0, shared.dim.target_met_adults_t.dim[0] - 1}, {2, 2}) == dust2::array::sum(internal.coverage_target_1st_dose_adults.data(), shared.dim.coverage_target_1st_dose_adults) ? prioritisation_step_1st_dose_adults + 1 : prioritisation_step_1st_dose_adults); + const real_type prioritisation_step_2nd_dose_adults_proposal = (dust2::array::sum(internal.target_met_adults_t.data(), shared.dim.target_met_adults_t, {0, shared.dim.target_met_adults_t.dim[0] - 1}, {3, 3}) == dust2::array::sum(internal.coverage_target_2nd_dose_adults.data(), shared.dim.coverage_target_2nd_dose_adults) ? prioritisation_step_2nd_dose_adults + 1 : prioritisation_step_2nd_dose_adults); for (size_t i = 1; i <= shared.dim.n_vaccination_t_S_children.size; ++i) { internal.n_vaccination_t_S_children[i - 1] = 0; } @@ -1068,7 +1066,7 @@ class model_targeted_vax { } for (size_t i = 1; i <= shared.dim.lambda_hc.dim[0]; ++i) { for (size_t j = 1; j <= shared.dim.lambda_hc.dim[1]; ++j) { - internal.lambda_hc[i - 1 + (j - 1) * shared.dim.lambda_hc.mult[1]] = (i == 20 ? shared.beta_hcw * dust2::array::sum(internal.I_infectious.data(), shared.dim.I_infectious) * (1 - shared.ve_I[i - 1 + (j - 1) * shared.dim.ve_I.mult[1]]) : 0); + internal.lambda_hc[i - 1 + (j - 1) * shared.dim.lambda_hc.mult[1]] = (i == 20 ? shared.beta_hcw * dust2::array::sum(internal.I_infectious.data(), shared.dim.I_infectious) / dust2::array::sum(N, shared.dim.N) * (1 - shared.ve_I[i - 1 + (j - 1) * shared.dim.ve_I.mult[1]]) : 0); } } for (size_t i = 1; i <= shared.dim.n_vaccination_t_S.dim[0]; ++i) { diff --git a/inst/odin/model-targeted-vax.R b/inst/odin/model-targeted-vax.R index b96b084..8af0c16 100644 --- a/inst/odin/model-targeted-vax.R +++ b/inst/odin/model-targeted-vax.R @@ -75,13 +75,14 @@ dim(target_met_children_t) <- c(n_group) ## for this in the 1st dose target - this now isn't strictly relevant for ## children but leaving it in in case we do expand this to 2 doses in future target_met_children_t[] <- 0 -target_met_children_t[i] <- +target_met_children_t[] <- ((sum(N[i, 3:4]) * is_child[i]) > prioritisation_strategy_children[ i, prioritisation_step_1st_dose_children] * sum(N[i, ])) ## adults +## keep this as a vector because it means the columns match (E.g. 3 and 3, 4 and 4) dim(target_met_adults_t) <- c(n_group, n_vax) ## 1st doses @@ -637,15 +638,27 @@ s_ij_gen_pop[, ] <- m_gen_pop[i, j] * prop_infectious[j] # as above but for the sexual contacts only s_ij_sex[, ] <- m_sex[i, j] * prop_infectious[j] -lambda_hh[, ] <- beta_h * sum(s_ij_gen_pop[i, ]) -lambda_s[, ] <- beta_s * sum(s_ij_sex[i, ]) + +lambda_hh[, ] <- beta_h * sum(s_ij_gen_pop[i, ]) * (1 - ve_I[i, j]) +lambda_s[, ] <- beta_s * sum(s_ij_sex[i, ]) * (1 - ve_I[i, j]) # additional foi in HCW only (i = 20) homogeneous from infected as assumed equally # likely to attend hospital -lambda_hc[, ] <- - if (i == 20) beta_hcw * sum(I_infectious[, ]) else 0 -lambda_z[, ] <- beta_z[i] - -lambda[, ] <- (lambda_hh[i, j] + lambda_s[i, j] + lambda_hc[i, j] + lambda_z[i, j]) * (1 - ve_I[i, j]) +lambda_hc[, ] <- + if (i == 20) beta_hcw * sum(I_infectious) / sum(N) * (1 - ve_I[i, j]) else 0 +lambda_z[, ] <- beta_z[i] * (1 - ve_I[i, j]) + +lambda[, ] <- lambda_hh[i, j] + lambda_s[i, j] + lambda_hc[i, j] + lambda_z[i, j] + +#### did not like this +# lambda_hh[, ] <- beta_h * sum(s_ij_gen_pop[i, ]) +# lambda_s[, ] <- beta_s * sum(s_ij_sex[i, ]) +# # additional foi in HCW only (i = 20) homogeneous from infected as assumed equally +# # likely to attend hospital +# lambda_hc[, ] <- +# if (i == 20) beta_hcw * sum(I_infectious[, ]) else 0 +# lambda_z[, ] <- beta_z[i] +# +# lambda[, ] <- (lambda_hh[i, j] + lambda_s[i, j] + lambda_hc[i, j] + lambda_z[i, j]) * (1 - ve_I[i, j]) ## Draws from binomial distributions for numbers changing between compartments # accounting for vaccination: diff --git a/src/model-targeted-vax.cpp b/src/model-targeted-vax.cpp index ca2d0b3..5668e49 100644 --- a/src/model-targeted-vax.cpp +++ b/src/model-targeted-vax.cpp @@ -68,11 +68,11 @@ class model_targeted_vax { dust2::array::dimensions<1> daily_doses_adults_t; dust2::array::dimensions<2> prioritisation_strategy_children; dust2::array::dimensions<2> prioritisation_strategy_adults; - dust2::array::dimensions<2> target_met_children_t; + dust2::array::dimensions<1> target_met_children_t; dust2::array::dimensions<2> target_met_adults_t; - dust2::array::dimensions<1> coverage_achieved_1st_dose_children; - dust2::array::dimensions<1> coverage_achieved_1st_dose_adults; - dust2::array::dimensions<1> coverage_achieved_2nd_dose_adults; + dust2::array::dimensions<1> coverage_target_1st_dose_children; + dust2::array::dimensions<1> coverage_target_1st_dose_adults; + dust2::array::dimensions<1> coverage_target_2nd_dose_adults; dust2::array::dimensions<1> n_eligible_for_dose1_children; dust2::array::dimensions<1> n_eligible_for_dose1_adults; dust2::array::dimensions<1> n_eligible_for_dose2_adults; @@ -196,9 +196,9 @@ class model_targeted_vax { std::vector n_IdD; std::vector target_met_children_t; std::vector target_met_adults_t; - std::vector coverage_achieved_1st_dose_children; - std::vector coverage_achieved_1st_dose_adults; - std::vector coverage_achieved_2nd_dose_adults; + std::vector coverage_target_1st_dose_children; + std::vector coverage_target_1st_dose_adults; + std::vector coverage_target_2nd_dose_adults; std::vector n_eligible_for_dose1_children; std::vector n_eligible_for_dose1_adults; std::vector n_eligible_for_dose2_adults; @@ -316,11 +316,11 @@ class model_targeted_vax { dim.daily_doses_adults_t.set({static_cast(n_vax)}); dim.prioritisation_strategy_children.set({static_cast(n_group), static_cast(N_prioritisation_steps_children)}); dim.prioritisation_strategy_adults.set({static_cast(n_group), static_cast(N_prioritisation_steps_adults)}); - dim.target_met_children_t.set({static_cast(n_group), static_cast(n_vax)}); + dim.target_met_children_t.set({static_cast(n_group)}); dim.target_met_adults_t.set({static_cast(n_group), static_cast(n_vax)}); - dim.coverage_achieved_1st_dose_children.set({static_cast(n_group)}); - dim.coverage_achieved_1st_dose_adults.set({static_cast(n_group)}); - dim.coverage_achieved_2nd_dose_adults.set({static_cast(n_group)}); + dim.coverage_target_1st_dose_children.set({static_cast(n_group)}); + dim.coverage_target_1st_dose_adults.set({static_cast(n_group)}); + dim.coverage_target_2nd_dose_adults.set({static_cast(n_group)}); dim.n_eligible_for_dose1_children.set({static_cast(n_group)}); dim.n_eligible_for_dose1_adults.set({static_cast(n_group)}); dim.n_eligible_for_dose2_adults.set({static_cast(n_group)}); @@ -557,9 +557,9 @@ class model_targeted_vax { std::vector n_IdD(shared.dim.n_IdD.size); std::vector target_met_children_t(shared.dim.target_met_children_t.size); std::vector target_met_adults_t(shared.dim.target_met_adults_t.size); - std::vector coverage_achieved_1st_dose_children(shared.dim.coverage_achieved_1st_dose_children.size); - std::vector coverage_achieved_1st_dose_adults(shared.dim.coverage_achieved_1st_dose_adults.size); - std::vector coverage_achieved_2nd_dose_adults(shared.dim.coverage_achieved_2nd_dose_adults.size); + std::vector coverage_target_1st_dose_children(shared.dim.coverage_target_1st_dose_children.size); + std::vector coverage_target_1st_dose_adults(shared.dim.coverage_target_1st_dose_adults.size); + std::vector coverage_target_2nd_dose_adults(shared.dim.coverage_target_2nd_dose_adults.size); std::vector n_eligible_for_dose1_children(shared.dim.n_eligible_for_dose1_children.size); std::vector n_eligible_for_dose1_adults(shared.dim.n_eligible_for_dose1_adults.size); std::vector n_eligible_for_dose2_adults(shared.dim.n_eligible_for_dose2_adults.size); @@ -619,7 +619,7 @@ class model_targeted_vax { std::vector new_N(shared.dim.N.size); std::vector n_SEa_hc(shared.dim.n_SEa_hc.size); std::vector n_SEa_z(shared.dim.n_SEa_z.size); - return internal_state{n_IrR, n_IdD, target_met_children_t, target_met_adults_t, coverage_achieved_1st_dose_children, coverage_achieved_1st_dose_adults, coverage_achieved_2nd_dose_adults, n_eligible_for_dose1_children, n_eligible_for_dose1_adults, n_eligible_for_dose2_adults, I_infectious, delta_R, delta_D, daily_doses_children_t, daily_doses_adults_t, n_vaccination_t_S_children, n_vaccination_t_S_adults, n_vaccination_t_Ea_children, n_vaccination_t_Ea_adults, n_vaccination_t_Eb_children, n_vaccination_t_Eb_adults, n_vaccination_t_R_children, n_vaccination_t_R_adults, new_D, prop_infectious, lambda_hc, n_vaccination_t_S, n_vaccination_t_Ea, n_vaccination_t_Eb, n_vaccination_t_R, s_ij_gen_pop, s_ij_sex, delta_S_n_vaccination, delta_Ea_n_vaccination, delta_Eb_n_vaccination, delta_R_n_vaccination, n_vaccination_t, lambda_hh, lambda_s, new_R, lambda, n_EaEb, n_EbI, p_SE, p_hh, p_s, p_hc, n_EbId, delta_Eb, new_Eb, n_SEa, n_EbIr, delta_Id, new_S, new_Id, n_SEa_hh, delta_Ea, delta_Ir, new_Ea, new_Ir, n_SEa_s, new_E, new_I, new_N, n_SEa_hc, n_SEa_z}; + return internal_state{n_IrR, n_IdD, target_met_children_t, target_met_adults_t, coverage_target_1st_dose_children, coverage_target_1st_dose_adults, coverage_target_2nd_dose_adults, n_eligible_for_dose1_children, n_eligible_for_dose1_adults, n_eligible_for_dose2_adults, I_infectious, delta_R, delta_D, daily_doses_children_t, daily_doses_adults_t, n_vaccination_t_S_children, n_vaccination_t_S_adults, n_vaccination_t_Ea_children, n_vaccination_t_Ea_adults, n_vaccination_t_Eb_children, n_vaccination_t_Eb_adults, n_vaccination_t_R_children, n_vaccination_t_R_adults, new_D, prop_infectious, lambda_hc, n_vaccination_t_S, n_vaccination_t_Ea, n_vaccination_t_Eb, n_vaccination_t_R, s_ij_gen_pop, s_ij_sex, delta_S_n_vaccination, delta_Ea_n_vaccination, delta_Eb_n_vaccination, delta_R_n_vaccination, n_vaccination_t, lambda_hh, lambda_s, new_R, lambda, n_EaEb, n_EbI, p_SE, p_hh, p_s, p_hc, n_EbId, delta_Eb, new_Eb, n_SEa, n_EbIr, delta_Id, new_S, new_Id, n_SEa_hh, delta_Ea, delta_Ir, new_Ea, new_Ir, n_SEa_s, new_E, new_I, new_N, n_SEa_hc, n_SEa_z}; } static data_type build_data(cpp11::list data, const shared_state& shared) { auto cases = dust2::r::read_real(data, "cases", NA_REAL); @@ -946,13 +946,11 @@ class model_targeted_vax { internal.n_IdD[i - 1 + (j - 1) * shared.dim.n_IdD.mult[1]] = monty::random::binomial(rng_state, Id[i - 1 + (j - 1) * shared.dim.Id.mult[1]], p_IdD); } } - for (size_t i = 1; i <= shared.dim.target_met_children_t.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.target_met_children_t.dim[1]; ++j) { - internal.target_met_children_t[i - 1 + (j - 1) * shared.dim.target_met_children_t.mult[1]] = 0; - } + for (size_t i = 1; i <= shared.dim.target_met_children_t.size; ++i) { + internal.target_met_children_t[i - 1] = 0; } - for (size_t i = 1; i <= shared.dim.target_met_children_t.dim[0]; ++i) { - internal.target_met_children_t[i - 1 + 2 * shared.dim.target_met_children_t.mult[1]] = ((dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {2, 3}) * shared.is_child[i - 1]) > shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]] * dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1})); + for (size_t i = 1; i <= shared.dim.target_met_children_t.size; ++i) { + internal.target_met_children_t[i - 1] = ((dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {2, 3}) * shared.is_child[i - 1]) > shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]] * dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1})); } for (size_t i = 1; i <= shared.dim.target_met_adults_t.dim[0]; ++i) { for (size_t j = 1; j <= shared.dim.target_met_adults_t.dim[1]; ++j) { @@ -965,14 +963,14 @@ class model_targeted_vax { for (size_t i = 1; i <= shared.dim.target_met_adults_t.dim[0]; ++i) { internal.target_met_adults_t[i - 1 + 3 * shared.dim.target_met_adults_t.mult[1]] = ((dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {3, 3}) * (1 - shared.is_child[i - 1])) > shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]] * dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1})); } - for (size_t i = 1; i <= shared.dim.coverage_achieved_1st_dose_children.size; ++i) { - internal.coverage_achieved_1st_dose_children[i - 1] = monty::math::ceil(shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]]); + for (size_t i = 1; i <= shared.dim.coverage_target_1st_dose_children.size; ++i) { + internal.coverage_target_1st_dose_children[i - 1] = monty::math::ceil(shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]]); } - for (size_t i = 1; i <= shared.dim.coverage_achieved_1st_dose_adults.size; ++i) { - internal.coverage_achieved_1st_dose_adults[i - 1] = monty::math::ceil(shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_1st_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]); + for (size_t i = 1; i <= shared.dim.coverage_target_1st_dose_adults.size; ++i) { + internal.coverage_target_1st_dose_adults[i - 1] = monty::math::ceil(shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_1st_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]); } - for (size_t i = 1; i <= shared.dim.coverage_achieved_2nd_dose_adults.size; ++i) { - internal.coverage_achieved_2nd_dose_adults[i - 1] = monty::math::ceil(shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]); + for (size_t i = 1; i <= shared.dim.coverage_target_2nd_dose_adults.size; ++i) { + internal.coverage_target_2nd_dose_adults[i - 1] = monty::math::ceil(shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]); } for (size_t i = 1; i <= shared.dim.n_eligible_for_dose1_children.size; ++i) { internal.n_eligible_for_dose1_children[i - 1] = (S[i - 1 + shared.dim.S.mult[1]] + Ea[i - 1 + shared.dim.Ea.mult[1]] + Eb[i - 1 + shared.dim.Eb.mult[1]] + R[i - 1 + shared.dim.R.mult[1]]) * shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]]; @@ -1006,9 +1004,9 @@ class model_targeted_vax { } shared.interpolate_daily_doses_children_t.eval(time, internal.daily_doses_children_t); shared.interpolate_daily_doses_adults_t.eval(time, internal.daily_doses_adults_t); - const real_type prioritisation_step_1st_dose_children_proposal = (dust2::array::sum(internal.target_met_children_t.data(), shared.dim.target_met_children_t, {0, shared.dim.target_met_children_t.dim[0] - 1}, {2, 2}) == dust2::array::sum(internal.coverage_achieved_1st_dose_children.data(), shared.dim.coverage_achieved_1st_dose_children) ? prioritisation_step_1st_dose_children + 1 : prioritisation_step_1st_dose_children); - const real_type prioritisation_step_1st_dose_adults_proposal = (dust2::array::sum(internal.target_met_adults_t.data(), shared.dim.target_met_adults_t, {0, shared.dim.target_met_adults_t.dim[0] - 1}, {2, 2}) == dust2::array::sum(internal.coverage_achieved_1st_dose_adults.data(), shared.dim.coverage_achieved_1st_dose_adults) ? prioritisation_step_1st_dose_adults + 1 : prioritisation_step_1st_dose_adults); - const real_type prioritisation_step_2nd_dose_adults_proposal = (dust2::array::sum(internal.target_met_adults_t.data(), shared.dim.target_met_adults_t, {0, shared.dim.target_met_adults_t.dim[0] - 1}, {3, 3}) == dust2::array::sum(internal.coverage_achieved_2nd_dose_adults.data(), shared.dim.coverage_achieved_2nd_dose_adults) ? prioritisation_step_2nd_dose_adults + 1 : prioritisation_step_2nd_dose_adults); + const real_type prioritisation_step_1st_dose_children_proposal = (dust2::array::sum(internal.target_met_children_t.data(), shared.dim.target_met_children_t) == dust2::array::sum(internal.coverage_target_1st_dose_children.data(), shared.dim.coverage_target_1st_dose_children) ? prioritisation_step_1st_dose_children + 1 : prioritisation_step_1st_dose_children); + const real_type prioritisation_step_1st_dose_adults_proposal = (dust2::array::sum(internal.target_met_adults_t.data(), shared.dim.target_met_adults_t, {0, shared.dim.target_met_adults_t.dim[0] - 1}, {2, 2}) == dust2::array::sum(internal.coverage_target_1st_dose_adults.data(), shared.dim.coverage_target_1st_dose_adults) ? prioritisation_step_1st_dose_adults + 1 : prioritisation_step_1st_dose_adults); + const real_type prioritisation_step_2nd_dose_adults_proposal = (dust2::array::sum(internal.target_met_adults_t.data(), shared.dim.target_met_adults_t, {0, shared.dim.target_met_adults_t.dim[0] - 1}, {3, 3}) == dust2::array::sum(internal.coverage_target_2nd_dose_adults.data(), shared.dim.coverage_target_2nd_dose_adults) ? prioritisation_step_2nd_dose_adults + 1 : prioritisation_step_2nd_dose_adults); for (size_t i = 1; i <= shared.dim.n_vaccination_t_S_children.size; ++i) { internal.n_vaccination_t_S_children[i - 1] = 0; } @@ -1070,7 +1068,7 @@ class model_targeted_vax { } for (size_t i = 1; i <= shared.dim.lambda_hc.dim[0]; ++i) { for (size_t j = 1; j <= shared.dim.lambda_hc.dim[1]; ++j) { - internal.lambda_hc[i - 1 + (j - 1) * shared.dim.lambda_hc.mult[1]] = (i == 20 ? shared.beta_hcw * dust2::array::sum(internal.I_infectious.data(), shared.dim.I_infectious) * (1 - shared.ve_I[i - 1 + (j - 1) * shared.dim.ve_I.mult[1]]) : 0); + internal.lambda_hc[i - 1 + (j - 1) * shared.dim.lambda_hc.mult[1]] = (i == 20 ? shared.beta_hcw * dust2::array::sum(internal.I_infectious.data(), shared.dim.I_infectious) / dust2::array::sum(N, shared.dim.N) * (1 - shared.ve_I[i - 1 + (j - 1) * shared.dim.ve_I.mult[1]]) : 0); } } for (size_t i = 1; i <= shared.dim.n_vaccination_t_S.dim[0]; ++i) { From b119a1cdeb12ae816c92a474b60263ac00e0ca16 Mon Sep 17 00:00:00 2001 From: ruthmccabe Date: Fri, 10 Jan 2025 18:25:10 +0000 Subject: [PATCH 10/26] rewriting of vaccination components based on 1) eligibility 2) prioritisation and 3) target reached or not --- inst/odin/model-targeted-vax.R | 208 ++++++++++++++++++--------------- 1 file changed, 114 insertions(+), 94 deletions(-) diff --git a/inst/odin/model-targeted-vax.R b/inst/odin/model-targeted-vax.R index 8af0c16..dc71e14 100644 --- a/inst/odin/model-targeted-vax.R +++ b/inst/odin/model-targeted-vax.R @@ -164,51 +164,87 @@ initial(prioritisation_step_1st_dose_children) <- 1 initial(prioritisation_step_1st_dose_adults) <- 1 initial(prioritisation_step_2nd_dose_adults) <- 1 -## now that we know the step we are on, we know who is eligible per age group to +## now that we know the step we are on, we know who is eligible per age group to ## be vaccinated - -## who can get a children's dose -n_eligible_for_dose1_children[] <- (S[i, 2] + Ea[i, 2] + Eb[i, 2] + R[i, 2]) * - prioritisation_strategy_children[i, prioritisation_step_1st_dose_children] -dim(n_eligible_for_dose1_children) <- c(n_group) -## who can get a 1st dose -## now we have to look in the preceding j class (e.g. who in unvaccinated j = 1 -## can get a 1st dose and move to j = 2) as these are the people who will be -## eligible -n_eligible_for_dose1_adults[] <- (S[i, 2] + Ea[i, 2] + Eb[i, 2] + R[i, 2]) * - prioritisation_strategy_adults[i, prioritisation_step_1st_dose_adults] -dim(n_eligible_for_dose1_adults) <- c(n_group) -## who can get a 2nd dose -n_eligible_for_dose2_adults[] <- (S[i, 3] + Ea[i, 3] + Eb[i, 3] + R[i, 3]) * - prioritisation_strategy_adults[i, prioritisation_step_2nd_dose_adults] -dim(n_eligible_for_dose2_adults) <- c(n_group) +#### Notes from meeting 8 Jan +### if target_met == 1 then vax is 0; +## if n_eligb == 0 then 0; +## else (daily_doses * S / n_elig) * ceiling(prioiritisation) (remove prioritisation) +## give_dose_X X = {dose_1_children,dose_1_adult,dose_2_adult} length n_group +## is target met or elgiible or prioiritsation step - single indicator + +## old eligibility +# ## isolate the relevant classes +# n_eligible_for_dose1_children[] <- (S[i, 2] + Ea[i, 2] + Eb[i, 2] + R[i, 2]) +# dim(n_eligible_for_dose1_children) <- c(n_group) +# ## who can get a 1st dose +# ## now we have to look in the preceding j class (e.g. who in unvaccinated j = 1 +# ## can get a 1st dose and move to j = 2) as these are the people who will be +# ## eligible +# n_eligible_for_dose1_adults[] <- (S[i, 2] + Ea[i, 2] + Eb[i, 2] + R[i, 2]) +# dim(n_eligible_for_dose1_adults) <- c(n_group) +# ## who can get a 2nd dose +# n_eligible_for_dose2_adults[] <- (S[i, 3] + Ea[i, 3] + Eb[i, 3] + R[i, 3]) +# dim(n_eligible_for_dose2_adults) <- c(n_group) + +# ## isolate the relevant compartments who can be vaccinated (relevant groups will be highlighted further down) +# comps_eligible_dose1[] <- (S[i, 2] + Ea[i, 2] + Eb[i, 2] + R[i, 2]) +# dim(comps_eligible_dose1) <- c(n_group) +# comps_eligible_dose2[] <- (S[i, 3] + Ea[i, 3] + Eb[i, 3] + R[i, 3]) +# dim(comps_eligible_dose2) <- c(n_group) + +## create indicator for which groups are allowed to be vaccinated based on: +## 1) if they are an adult or a child (via is_child) +## 2) if they are in a prioritised group, dependent on being adult or child (prioritisation_strategy_X) +## 3) if the target for this group has not yet been met (target_met_X_t) + +dim(give_dose1_children) <- c(n_group) +give_dose1_children[] <- + (is_child[i] * ceiling(prioritisation_strategy_children[ + i, prioritisation_step_1st_dose_children]) * (1 - target_met_children_t[i])) + +dim(give_dose1_adults) <- c(n_group) +give_dose1_adults[] <- + ((1 - is_child[i]) * ceiling(prioritisation_strategy_adults[ + i, prioritisation_step_1st_dose_adults]) * (1 - target_met_adults_t[i, 3])) +dim(give_dose2_adults) <- c(n_group) +give_dose2_adults[] <- + ((1 - is_child[i]) * ceiling(prioritisation_strategy_adults[ + i, prioritisation_step_2nd_dose_adults]) * (1 - target_met_adults_t[i, 4])) ## allocate the doses to the unvaccinated by age group, prioritisation strategy ## and across S, E, R -## hacky fix for now + +children_dose1_denom[] <- (S[i, 2] + Ea[i, 2] + Eb[i, 2] + R[i, 2]) * give_dose1_children[i] +dim(children_dose1_denom) <- c(n_group) + +adults_dose1_denom[] <- (S[i, 2] + Ea[i, 2] + Eb[i, 2] + R[i, 2]) * give_dose1_adults[i] +dim(adults_dose1_denom) <- c(n_group) + +adults_dose2_denom[] <- (S[i, 3] + Ea[i, 3] + Eb[i, 3] + R[i, 3]) * give_dose2_adults[i] +dim(adults_dose2_denom) <- c(n_group) ### allocate to S ## children 1st doses n_vaccination_t_S_children[] <- 0 n_vaccination_t_S_children[] <- - if (sum(n_eligible_for_dose1_children[]) == 0) 0 else - min(floor((daily_doses_children_t[2] * S[i, 2] * - prioritisation_strategy_children[ - i, prioritisation_step_1st_dose_children]) / - sum(n_eligible_for_dose1_children[])), + if (sum(children_dose1_denom[]) == 0) 0 else + min(floor(((daily_doses_children_t[2] * S[i, 2]) / + sum(children_dose1_denom[])) * give_dose1_children[i]), S[i, 2]) +##### current concern about denominator here - I think all of the compartments will be summed over before the ones that are able to be vaccinated have been taken account of, so we won't be giving out all of the compartments + + n_vaccination_t_S_adults[] <- 0 ## adults 1st doses n_vaccination_t_S_adults[] <- - if (sum(n_eligible_for_dose1_adults[]) == 0) 0 else - min(floor((daily_doses_adults_t[2] * S[i, 2] * - prioritisation_strategy_adults[ - i, prioritisation_step_1st_dose_adults]) / - sum(n_eligible_for_dose1_adults[])), + if (sum(adults_dose1_denom[]) == 0) 0 else + min(floor(((daily_doses_adults_t[2] * S[i, 2]) / + sum(adults_dose1_denom[])) * give_dose1_adults[i]), S[i, 2]) ## combine total first doses @@ -216,39 +252,34 @@ n_vaccination_t_S[, ] <- 0 n_vaccination_t_S[, 2] <- n_vaccination_t_S_children[i] + n_vaccination_t_S_adults[i] -## for the boundary case do an extra check that we haven't gone over the number -## of people in each compartment -n_vaccination_t_S[3, 2] <- min(n_vaccination_t_S[3, 2], S[3, 2]) +## this should no longer be necessary +# ## for the boundary case do an extra check that we haven't gone over the number +# ## of people in each compartment +# n_vaccination_t_S[3, 2] <- min(n_vaccination_t_S[3, 2], S[3, 2]) ## allocate 2nd doses (adults only for now) -n_vaccination_t_S[, 3] <- if (sum(n_eligible_for_dose2_adults[]) == 0) 0 else - min(floor((daily_doses_adults_t[3] * S[i, 3] * - prioritisation_strategy_adults[ - i, prioritisation_step_2nd_dose_adults]) / - sum(n_eligible_for_dose2_adults[])), +n_vaccination_t_S[, 3] <- + if (sum(adults_dose2_denom[]) == 0) 0 else + min(floor(((daily_doses_adults_t[3] * S[i, 3]) / + sum(adults_dose2_denom[])) * give_dose2_adults[i]), S[i, 3]) - ### allocate to Ea ## children 1st doses n_vaccination_t_Ea_children[] <- 0 n_vaccination_t_Ea_children[] <- - if (sum(n_eligible_for_dose1_children[]) == 0) 0 else - min(floor((daily_doses_children_t[2] * Ea[i, 2] * - prioritisation_strategy_children[ - i, prioritisation_step_1st_dose_children]) / - sum(n_eligible_for_dose1_children[])), + if (sum(children_dose1_denom[]) == 0) 0 else + min(floor(((daily_doses_children_t[2] * Ea[i, 2]) / + sum(children_dose1_denom[])) * give_dose1_children[i]), Ea[i, 2]) ## adults 1st doses n_vaccination_t_Ea_adults[] <- 0 n_vaccination_t_Ea_adults[] <- - if (sum(n_eligible_for_dose1_adults[]) == 0) 0 else - min(floor((daily_doses_adults_t[2] * Ea[i, 2] * - prioritisation_strategy_adults[ - i, prioritisation_step_1st_dose_adults]) / - sum(n_eligible_for_dose1_adults[])), + if (sum(adults_dose1_denom[]) == 0) 0 else + min(floor(((daily_doses_adults_t[2] * Ea[i, 2]) / + sum(adults_dose1_denom[])) * give_dose1_adults[i]), Ea[i, 2]) ## combine total first doses @@ -256,17 +287,16 @@ n_vaccination_t_Ea[, ] <- 0 n_vaccination_t_Ea[, 2] <- n_vaccination_t_Ea_children[i] + n_vaccination_t_Ea_adults[i] -## for the boundary case do an extra check that we haven't gone over the number -## of people in each compartment -n_vaccination_t_Ea[3, 2] <- min(n_vaccination_t_Ea[3, 2], Ea[3, 2]) +# ## for the boundary case do an extra check that we haven't gone over the number +# ## of people in each compartment +# n_vaccination_t_Ea[3, 2] <- min(n_vaccination_t_Ea[3, 2], Ea[3, 2]) ## adults 2nd doses -n_vaccination_t_Ea[, 3] <- if (sum(n_eligible_for_dose2_adults[]) == 0) 0 else - min(floor((daily_doses_adults_t[3] * Ea[i, 3] * - prioritisation_strategy_adults[ - i, prioritisation_step_2nd_dose_adults]) / - sum(n_eligible_for_dose2_adults[])), - Ea[i, 3]) +n_vaccination_t_Ea[, 3] <- + if (sum(adults_dose2_denom[]) == 0) 0 else + min(floor(((daily_doses_adults_t[3] * Ea[i, 3]) / + sum(adults_dose2_denom[])) * give_dose2_adults[i]), + Ea[i, 3]) ### allocate to Eb @@ -274,21 +304,17 @@ n_vaccination_t_Ea[, 3] <- if (sum(n_eligible_for_dose2_adults[]) == 0) 0 else ## children 1st doses n_vaccination_t_Eb_children[] <- 0 n_vaccination_t_Eb_children[] <- - if (sum(n_eligible_for_dose1_children[]) == 0) 0 else - min(floor((daily_doses_children_t[2] * Eb[i, 2] * - prioritisation_strategy_children[ - i, prioritisation_step_1st_dose_children]) / - sum(n_eligible_for_dose1_children[])), + if (sum(children_dose1_denom[]) == 0) 0 else + min(floor(((daily_doses_children_t[2] * Eb[i, 2]) / + sum(children_dose1_denom[])) * give_dose1_children[i]), Eb[i, 2]) ## adults 1st doses n_vaccination_t_Eb_adults[] <- 0 n_vaccination_t_Eb_adults[] <- - if (sum(n_eligible_for_dose1_adults[]) == 0) 0 else - min(floor((daily_doses_adults_t[2] * Eb[i, 2] * - prioritisation_strategy_adults[ - i, prioritisation_step_1st_dose_adults]) / - sum(n_eligible_for_dose1_adults[])), + if (sum(adults_dose1_denom[]) == 0) 0 else + min(floor(((daily_doses_adults_t[2] * Eb[i, 2]) / + sum(adults_dose1_denom[])) * give_dose1_adults[i]), Eb[i, 2]) ## combine total first doses @@ -296,17 +322,16 @@ n_vaccination_t_Eb[, ] <- 0 n_vaccination_t_Eb[, 2] <- n_vaccination_t_Eb_children[i] + n_vaccination_t_Eb_adults[i] -## for the boundary case do an extra check that we haven't gone over the number -## of people in each compartment -n_vaccination_t_Eb[3, 2] <- min(n_vaccination_t_Eb[3, 2], Eb[3, 2]) +# ## for the boundary case do an extra check that we haven't gone over the number +# ## of people in each compartment +# n_vaccination_t_Eb[3, 2] <- min(n_vaccination_t_Eb[3, 2], Eb[3, 2]) ## adults 2nd doses -n_vaccination_t_Eb[, 3] <- if (sum(n_eligible_for_dose2_adults[]) == 0) 0 else - min(floor((daily_doses_adults_t[3] * Eb[i, 3] * - prioritisation_strategy_adults[ - i, prioritisation_step_2nd_dose_adults]) / - sum(n_eligible_for_dose2_adults[])), - Eb[i, 3]) +n_vaccination_t_Eb[, 3] <- + if (sum(adults_dose2_denom[]) == 0) 0 else + min(floor(((daily_doses_adults_t[3] * Eb[i, 3]) / + sum(adults_dose2_denom[])) * give_dose2_adults[i]), + Eb[i, 3]) ### allocate to R @@ -314,21 +339,17 @@ n_vaccination_t_Eb[, 3] <- if (sum(n_eligible_for_dose2_adults[]) == 0) 0 else ## children 1st doses n_vaccination_t_R_children[] <- 0 n_vaccination_t_R_children[] <- - if (sum(n_eligible_for_dose1_children[]) == 0) 0 else - min(floor((daily_doses_children_t[2] * R[i, 2] * - prioritisation_strategy_children[ - i, prioritisation_step_1st_dose_children]) / - sum(n_eligible_for_dose1_children[])), + if (sum(children_dose1_denom[]) == 0) 0 else + min(floor(((daily_doses_children_t[2] * R[i, 2]) / + sum(children_dose1_denom[])) * give_dose1_children[i]), R[i, 2]) ## adults 1st doses n_vaccination_t_R_adults[] <- 0 n_vaccination_t_R_adults[] <- - if (sum(n_eligible_for_dose1_adults[]) == 0) 0 else - min(floor((daily_doses_adults_t[2] * R[i, 2] * - prioritisation_strategy_adults[ - i, prioritisation_step_1st_dose_adults]) / - sum(n_eligible_for_dose1_adults[])), + if (sum(adults_dose1_denom[]) == 0) 0 else + min(floor(((daily_doses_adults_t[2] * R[i, 2]) / + sum(adults_dose1_denom[])) * give_dose1_adults[i]), R[i, 2]) ## combine total first doses @@ -336,17 +357,16 @@ n_vaccination_t_R[, ] <- 0 n_vaccination_t_R[, 2] <- n_vaccination_t_R_children[i] + n_vaccination_t_R_adults[i] -## for the boundary case do an extra check that we haven't gone over the number -## of people in each compartment -n_vaccination_t_R[3, 2] <- min(n_vaccination_t_R[3, 2], R[3, 2]) +# ## for the boundary case do an extra check that we haven't gone over the number +# ## of people in each compartment +# n_vaccination_t_R[3, 2] <- min(n_vaccination_t_R[3, 2], R[3, 2]) ## adults 2nd doses -n_vaccination_t_R[, 3] <- if (sum(n_eligible_for_dose2_adults[]) == 0) 0 else - min(floor((daily_doses_adults_t[3] * R[i, 3] * - prioritisation_strategy_adults[ - i, prioritisation_step_2nd_dose_adults]) / - sum(n_eligible_for_dose2_adults[])), - R[i, 3]) +n_vaccination_t_R[, 3] <- + if (sum(adults_dose2_denom[]) == 0) 0 else + min(floor(((daily_doses_adults_t[3] * R[i, 3]) / + sum(adults_dose2_denom[])) * give_dose2_adults[i]), + R[i, 3]) ## net vaccination change for relevant classes (S, Ea, Eb, R) From c1987adf8d15d974f8ea16ddcb007ce0b9a72817 Mon Sep 17 00:00:00 2001 From: ruthmccabe Date: Fri, 10 Jan 2025 18:25:22 +0000 Subject: [PATCH 11/26] regenerate model for testing purposes --- inst/dust/model-targeted-vax.cpp | 363 ++++++++++++++++--------------- src/model-targeted-vax.cpp | 363 ++++++++++++++++--------------- 2 files changed, 380 insertions(+), 346 deletions(-) diff --git a/inst/dust/model-targeted-vax.cpp b/inst/dust/model-targeted-vax.cpp index 4284ced..919c97e 100644 --- a/inst/dust/model-targeted-vax.cpp +++ b/inst/dust/model-targeted-vax.cpp @@ -71,9 +71,12 @@ class model_targeted_vax { dust2::array::dimensions<1> coverage_target_1st_dose_children; dust2::array::dimensions<1> coverage_target_1st_dose_adults; dust2::array::dimensions<1> coverage_target_2nd_dose_adults; - dust2::array::dimensions<1> n_eligible_for_dose1_children; - dust2::array::dimensions<1> n_eligible_for_dose1_adults; - dust2::array::dimensions<1> n_eligible_for_dose2_adults; + dust2::array::dimensions<1> give_dose1_children; + dust2::array::dimensions<1> give_dose1_adults; + dust2::array::dimensions<1> give_dose2_adults; + dust2::array::dimensions<1> children_dose1_denom; + dust2::array::dimensions<1> adults_dose1_denom; + dust2::array::dimensions<1> adults_dose2_denom; dust2::array::dimensions<2> N; dust2::array::dimensions<2> S; dust2::array::dimensions<2> S0; @@ -197,14 +200,22 @@ class model_targeted_vax { std::vector coverage_target_1st_dose_children; std::vector coverage_target_1st_dose_adults; std::vector coverage_target_2nd_dose_adults; - std::vector n_eligible_for_dose1_children; - std::vector n_eligible_for_dose1_adults; - std::vector n_eligible_for_dose2_adults; std::vector I_infectious; std::vector delta_R; std::vector delta_D; std::vector daily_doses_children_t; std::vector daily_doses_adults_t; + std::vector give_dose1_children; + std::vector give_dose1_adults; + std::vector give_dose2_adults; + std::vector new_D; + std::vector prop_infectious; + std::vector lambda_hc; + std::vector children_dose1_denom; + std::vector adults_dose1_denom; + std::vector adults_dose2_denom; + std::vector s_ij_gen_pop; + std::vector s_ij_sex; std::vector n_vaccination_t_S_children; std::vector n_vaccination_t_S_adults; std::vector n_vaccination_t_Ea_children; @@ -213,49 +224,44 @@ class model_targeted_vax { std::vector n_vaccination_t_Eb_adults; std::vector n_vaccination_t_R_children; std::vector n_vaccination_t_R_adults; - std::vector new_D; - std::vector prop_infectious; - std::vector lambda_hc; + std::vector lambda_hh; + std::vector lambda_s; std::vector n_vaccination_t_S; std::vector n_vaccination_t_Ea; std::vector n_vaccination_t_Eb; std::vector n_vaccination_t_R; - std::vector s_ij_gen_pop; - std::vector s_ij_sex; + std::vector lambda; std::vector delta_S_n_vaccination; std::vector delta_Ea_n_vaccination; std::vector delta_Eb_n_vaccination; std::vector delta_R_n_vaccination; std::vector n_vaccination_t; - std::vector lambda_hh; - std::vector lambda_s; - std::vector new_R; - std::vector lambda; - std::vector n_EaEb; - std::vector n_EbI; std::vector p_SE; std::vector p_hh; std::vector p_s; std::vector p_hc; + std::vector new_R; + std::vector n_SEa; + std::vector n_EaEb; + std::vector n_EbI; + std::vector new_S; + std::vector n_SEa_hh; std::vector n_EbId; + std::vector delta_Ea; std::vector delta_Eb; + std::vector new_Ea; std::vector new_Eb; - std::vector n_SEa; + std::vector n_SEa_s; std::vector n_EbIr; std::vector delta_Id; - std::vector new_S; std::vector new_Id; - std::vector n_SEa_hh; - std::vector delta_Ea; + std::vector new_E; + std::vector n_SEa_hc; std::vector delta_Ir; - std::vector new_Ea; std::vector new_Ir; - std::vector n_SEa_s; - std::vector new_E; + std::vector n_SEa_z; std::vector new_I; std::vector new_N; - std::vector n_SEa_hc; - std::vector n_SEa_z; }; struct data_type { real_type cases; @@ -319,9 +325,12 @@ class model_targeted_vax { dim.coverage_target_1st_dose_children.set({static_cast(n_group)}); dim.coverage_target_1st_dose_adults.set({static_cast(n_group)}); dim.coverage_target_2nd_dose_adults.set({static_cast(n_group)}); - dim.n_eligible_for_dose1_children.set({static_cast(n_group)}); - dim.n_eligible_for_dose1_adults.set({static_cast(n_group)}); - dim.n_eligible_for_dose2_adults.set({static_cast(n_group)}); + dim.give_dose1_children.set({static_cast(n_group)}); + dim.give_dose1_adults.set({static_cast(n_group)}); + dim.give_dose2_adults.set({static_cast(n_group)}); + dim.children_dose1_denom.set({static_cast(n_group)}); + dim.adults_dose1_denom.set({static_cast(n_group)}); + dim.adults_dose2_denom.set({static_cast(n_group)}); dim.N.set({static_cast(n_group), static_cast(n_vax)}); dim.S.set({static_cast(n_group), static_cast(n_vax)}); dim.S0.set({static_cast(n_group), static_cast(n_vax)}); @@ -558,14 +567,22 @@ class model_targeted_vax { std::vector coverage_target_1st_dose_children(shared.dim.coverage_target_1st_dose_children.size); std::vector coverage_target_1st_dose_adults(shared.dim.coverage_target_1st_dose_adults.size); std::vector coverage_target_2nd_dose_adults(shared.dim.coverage_target_2nd_dose_adults.size); - std::vector n_eligible_for_dose1_children(shared.dim.n_eligible_for_dose1_children.size); - std::vector n_eligible_for_dose1_adults(shared.dim.n_eligible_for_dose1_adults.size); - std::vector n_eligible_for_dose2_adults(shared.dim.n_eligible_for_dose2_adults.size); std::vector I_infectious(shared.dim.I_infectious.size); std::vector delta_R(shared.dim.delta_R.size); std::vector delta_D(shared.dim.delta_D.size); std::vector daily_doses_children_t(shared.dim.daily_doses_children_t.size); std::vector daily_doses_adults_t(shared.dim.daily_doses_adults_t.size); + std::vector give_dose1_children(shared.dim.give_dose1_children.size); + std::vector give_dose1_adults(shared.dim.give_dose1_adults.size); + std::vector give_dose2_adults(shared.dim.give_dose2_adults.size); + std::vector new_D(shared.dim.D.size); + std::vector prop_infectious(shared.dim.prop_infectious.size); + std::vector lambda_hc(shared.dim.lambda_hc.size); + std::vector children_dose1_denom(shared.dim.children_dose1_denom.size); + std::vector adults_dose1_denom(shared.dim.adults_dose1_denom.size); + std::vector adults_dose2_denom(shared.dim.adults_dose2_denom.size); + std::vector s_ij_gen_pop(shared.dim.s_ij_gen_pop.size); + std::vector s_ij_sex(shared.dim.s_ij_sex.size); std::vector n_vaccination_t_S_children(shared.dim.n_vaccination_t_S_children.size); std::vector n_vaccination_t_S_adults(shared.dim.n_vaccination_t_S_adults.size); std::vector n_vaccination_t_Ea_children(shared.dim.n_vaccination_t_Ea_children.size); @@ -574,50 +591,45 @@ class model_targeted_vax { std::vector n_vaccination_t_Eb_adults(shared.dim.n_vaccination_t_Eb_adults.size); std::vector n_vaccination_t_R_children(shared.dim.n_vaccination_t_R_children.size); std::vector n_vaccination_t_R_adults(shared.dim.n_vaccination_t_R_adults.size); - std::vector new_D(shared.dim.D.size); - std::vector prop_infectious(shared.dim.prop_infectious.size); - std::vector lambda_hc(shared.dim.lambda_hc.size); + std::vector lambda_hh(shared.dim.lambda_hh.size); + std::vector lambda_s(shared.dim.lambda_s.size); std::vector n_vaccination_t_S(shared.dim.n_vaccination_t_S.size); std::vector n_vaccination_t_Ea(shared.dim.n_vaccination_t_Ea.size); std::vector n_vaccination_t_Eb(shared.dim.n_vaccination_t_Eb.size); std::vector n_vaccination_t_R(shared.dim.n_vaccination_t_R.size); - std::vector s_ij_gen_pop(shared.dim.s_ij_gen_pop.size); - std::vector s_ij_sex(shared.dim.s_ij_sex.size); + std::vector lambda(shared.dim.lambda.size); std::vector delta_S_n_vaccination(shared.dim.delta_S_n_vaccination.size); std::vector delta_Ea_n_vaccination(shared.dim.delta_Ea_n_vaccination.size); std::vector delta_Eb_n_vaccination(shared.dim.delta_Eb_n_vaccination.size); std::vector delta_R_n_vaccination(shared.dim.delta_R_n_vaccination.size); std::vector n_vaccination_t(shared.dim.n_vaccination_t.size); - std::vector lambda_hh(shared.dim.lambda_hh.size); - std::vector lambda_s(shared.dim.lambda_s.size); - std::vector new_R(shared.dim.R.size); - std::vector lambda(shared.dim.lambda.size); - std::vector n_EaEb(shared.dim.n_EaEb.size); - std::vector n_EbI(shared.dim.n_EbI.size); std::vector p_SE(shared.dim.p_SE.size); std::vector p_hh(shared.dim.p_hh.size); std::vector p_s(shared.dim.p_s.size); std::vector p_hc(shared.dim.p_hc.size); + std::vector new_R(shared.dim.R.size); + std::vector n_SEa(shared.dim.n_SEa.size); + std::vector n_EaEb(shared.dim.n_EaEb.size); + std::vector n_EbI(shared.dim.n_EbI.size); + std::vector new_S(shared.dim.S.size); + std::vector n_SEa_hh(shared.dim.n_SEa_hh.size); std::vector n_EbId(shared.dim.n_EbId.size); + std::vector delta_Ea(shared.dim.delta_Ea.size); std::vector delta_Eb(shared.dim.delta_Eb.size); + std::vector new_Ea(shared.dim.Ea.size); std::vector new_Eb(shared.dim.Eb.size); - std::vector n_SEa(shared.dim.n_SEa.size); + std::vector n_SEa_s(shared.dim.n_SEa_s.size); std::vector n_EbIr(shared.dim.n_EbIr.size); std::vector delta_Id(shared.dim.delta_Id.size); - std::vector new_S(shared.dim.S.size); std::vector new_Id(shared.dim.Id.size); - std::vector n_SEa_hh(shared.dim.n_SEa_hh.size); - std::vector delta_Ea(shared.dim.delta_Ea.size); + std::vector new_E(shared.dim.E.size); + std::vector n_SEa_hc(shared.dim.n_SEa_hc.size); std::vector delta_Ir(shared.dim.delta_Ir.size); - std::vector new_Ea(shared.dim.Ea.size); std::vector new_Ir(shared.dim.Ir.size); - std::vector n_SEa_s(shared.dim.n_SEa_s.size); - std::vector new_E(shared.dim.E.size); + std::vector n_SEa_z(shared.dim.n_SEa_z.size); std::vector new_I(shared.dim.I.size); std::vector new_N(shared.dim.N.size); - std::vector n_SEa_hc(shared.dim.n_SEa_hc.size); - std::vector n_SEa_z(shared.dim.n_SEa_z.size); - return internal_state{n_IrR, n_IdD, target_met_children_t, target_met_adults_t, coverage_target_1st_dose_children, coverage_target_1st_dose_adults, coverage_target_2nd_dose_adults, n_eligible_for_dose1_children, n_eligible_for_dose1_adults, n_eligible_for_dose2_adults, I_infectious, delta_R, delta_D, daily_doses_children_t, daily_doses_adults_t, n_vaccination_t_S_children, n_vaccination_t_S_adults, n_vaccination_t_Ea_children, n_vaccination_t_Ea_adults, n_vaccination_t_Eb_children, n_vaccination_t_Eb_adults, n_vaccination_t_R_children, n_vaccination_t_R_adults, new_D, prop_infectious, lambda_hc, n_vaccination_t_S, n_vaccination_t_Ea, n_vaccination_t_Eb, n_vaccination_t_R, s_ij_gen_pop, s_ij_sex, delta_S_n_vaccination, delta_Ea_n_vaccination, delta_Eb_n_vaccination, delta_R_n_vaccination, n_vaccination_t, lambda_hh, lambda_s, new_R, lambda, n_EaEb, n_EbI, p_SE, p_hh, p_s, p_hc, n_EbId, delta_Eb, new_Eb, n_SEa, n_EbIr, delta_Id, new_S, new_Id, n_SEa_hh, delta_Ea, delta_Ir, new_Ea, new_Ir, n_SEa_s, new_E, new_I, new_N, n_SEa_hc, n_SEa_z}; + return internal_state{n_IrR, n_IdD, target_met_children_t, target_met_adults_t, coverage_target_1st_dose_children, coverage_target_1st_dose_adults, coverage_target_2nd_dose_adults, I_infectious, delta_R, delta_D, daily_doses_children_t, daily_doses_adults_t, give_dose1_children, give_dose1_adults, give_dose2_adults, new_D, prop_infectious, lambda_hc, children_dose1_denom, adults_dose1_denom, adults_dose2_denom, s_ij_gen_pop, s_ij_sex, n_vaccination_t_S_children, n_vaccination_t_S_adults, n_vaccination_t_Ea_children, n_vaccination_t_Ea_adults, n_vaccination_t_Eb_children, n_vaccination_t_Eb_adults, n_vaccination_t_R_children, n_vaccination_t_R_adults, lambda_hh, lambda_s, n_vaccination_t_S, n_vaccination_t_Ea, n_vaccination_t_Eb, n_vaccination_t_R, lambda, delta_S_n_vaccination, delta_Ea_n_vaccination, delta_Eb_n_vaccination, delta_R_n_vaccination, n_vaccination_t, p_SE, p_hh, p_s, p_hc, new_R, n_SEa, n_EaEb, n_EbI, new_S, n_SEa_hh, n_EbId, delta_Ea, delta_Eb, new_Ea, new_Eb, n_SEa_s, n_EbIr, delta_Id, new_Id, new_E, n_SEa_hc, delta_Ir, new_Ir, n_SEa_z, new_I, new_N}; } static data_type build_data(cpp11::list data, const shared_state& shared) { auto cases = dust2::r::read_real(data, "cases", NA_REAL); @@ -970,15 +982,6 @@ class model_targeted_vax { for (size_t i = 1; i <= shared.dim.coverage_target_2nd_dose_adults.size; ++i) { internal.coverage_target_2nd_dose_adults[i - 1] = monty::math::ceil(shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]); } - for (size_t i = 1; i <= shared.dim.n_eligible_for_dose1_children.size; ++i) { - internal.n_eligible_for_dose1_children[i - 1] = (S[i - 1 + shared.dim.S.mult[1]] + Ea[i - 1 + shared.dim.Ea.mult[1]] + Eb[i - 1 + shared.dim.Eb.mult[1]] + R[i - 1 + shared.dim.R.mult[1]]) * shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]]; - } - for (size_t i = 1; i <= shared.dim.n_eligible_for_dose1_adults.size; ++i) { - internal.n_eligible_for_dose1_adults[i - 1] = (S[i - 1 + shared.dim.S.mult[1]] + Ea[i - 1 + shared.dim.Ea.mult[1]] + Eb[i - 1 + shared.dim.Eb.mult[1]] + R[i - 1 + shared.dim.R.mult[1]]) * shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_1st_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]; - } - for (size_t i = 1; i <= shared.dim.n_eligible_for_dose2_adults.size; ++i) { - internal.n_eligible_for_dose2_adults[i - 1] = (S[i - 1 + 2 * shared.dim.S.mult[1]] + Ea[i - 1 + 2 * shared.dim.Ea.mult[1]] + Eb[i - 1 + 2 * shared.dim.Eb.mult[1]] + R[i - 1 + 2 * shared.dim.R.mult[1]]) * shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]; - } const real_type new_deaths_00_04 = dust2::array::sum(internal.n_IdD.data(), shared.dim.n_IdD, {0, 0}, {0, shared.dim.n_IdD.dim[1] - 1}); const real_type new_deaths_SW_12_14 = monty::random::binomial(rng_state, dust2::array::sum(internal.n_IdD.data(), shared.dim.n_IdD, {16, 16}, {0, shared.dim.n_IdD.dim[1] - 1}), static_cast(0.5)); const real_type new_deaths_CSW = dust2::array::sum(internal.n_IdD.data(), shared.dim.n_IdD, {16, 16}, {0, shared.dim.n_IdD.dim[1] - 1}); @@ -1005,68 +1008,107 @@ class model_targeted_vax { const real_type prioritisation_step_1st_dose_children_proposal = (dust2::array::sum(internal.target_met_children_t.data(), shared.dim.target_met_children_t) == dust2::array::sum(internal.coverage_target_1st_dose_children.data(), shared.dim.coverage_target_1st_dose_children) ? prioritisation_step_1st_dose_children + 1 : prioritisation_step_1st_dose_children); const real_type prioritisation_step_1st_dose_adults_proposal = (dust2::array::sum(internal.target_met_adults_t.data(), shared.dim.target_met_adults_t, {0, shared.dim.target_met_adults_t.dim[0] - 1}, {2, 2}) == dust2::array::sum(internal.coverage_target_1st_dose_adults.data(), shared.dim.coverage_target_1st_dose_adults) ? prioritisation_step_1st_dose_adults + 1 : prioritisation_step_1st_dose_adults); const real_type prioritisation_step_2nd_dose_adults_proposal = (dust2::array::sum(internal.target_met_adults_t.data(), shared.dim.target_met_adults_t, {0, shared.dim.target_met_adults_t.dim[0] - 1}, {3, 3}) == dust2::array::sum(internal.coverage_target_2nd_dose_adults.data(), shared.dim.coverage_target_2nd_dose_adults) ? prioritisation_step_2nd_dose_adults + 1 : prioritisation_step_2nd_dose_adults); + for (size_t i = 1; i <= shared.dim.give_dose1_children.size; ++i) { + internal.give_dose1_children[i - 1] = (shared.is_child[i - 1] * monty::math::ceil(shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]]) * (1 - internal.target_met_children_t[i - 1])); + } + for (size_t i = 1; i <= shared.dim.give_dose1_adults.size; ++i) { + internal.give_dose1_adults[i - 1] = ((1 - shared.is_child[i - 1]) * monty::math::ceil(shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_1st_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) * (1 - internal.target_met_adults_t[i - 1 + 2 * shared.dim.target_met_adults_t.mult[1]])); + } + for (size_t i = 1; i <= shared.dim.give_dose2_adults.size; ++i) { + internal.give_dose2_adults[i - 1] = ((1 - shared.is_child[i - 1]) * monty::math::ceil(shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) * (1 - internal.target_met_adults_t[i - 1 + 3 * shared.dim.target_met_adults_t.mult[1]])); + } + for (size_t i = 1; i <= shared.dim.D.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.D.dim[1]; ++j) { + internal.new_D[i - 1 + (j - 1) * shared.dim.D.mult[1]] = D[i - 1 + (j - 1) * shared.dim.D.mult[1]] + internal.delta_D[i - 1 + (j - 1) * shared.dim.delta_D.mult[1]]; + } + } + const real_type new_deaths_SW_15_17 = dust2::array::sum(internal.n_IdD.data(), shared.dim.n_IdD, {16, 16}, {0, shared.dim.n_IdD.dim[1] - 1}) - new_deaths_SW_12_14; + const real_type new_deaths_05_14 = dust2::array::sum(internal.n_IdD.data(), shared.dim.n_IdD, {1, 2}, {0, shared.dim.n_IdD.dim[1] - 1}) + new_deaths_SW_12_14; + const real_type new_deaths_SW = new_deaths_CSW + new_deaths_ASW; + for (size_t i = 1; i <= shared.dim.prop_infectious.size; ++i) { + internal.prop_infectious[i - 1] = (dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1}) == 0 ? 0 : dust2::array::sum(internal.I_infectious.data(), shared.dim.I_infectious, {i - 1, i - 1}, {0, shared.dim.I_infectious.dim[1] - 1}) / dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1})); + } + for (size_t i = 1; i <= shared.dim.lambda_hc.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.lambda_hc.dim[1]; ++j) { + internal.lambda_hc[i - 1 + (j - 1) * shared.dim.lambda_hc.mult[1]] = (i == 20 ? shared.beta_hcw * dust2::array::sum(internal.I_infectious.data(), shared.dim.I_infectious) / dust2::array::sum(N, shared.dim.N) * (1 - shared.ve_I[i - 1 + (j - 1) * shared.dim.ve_I.mult[1]]) : 0); + } + } + for (size_t i = 1; i <= shared.dim.children_dose1_denom.size; ++i) { + internal.children_dose1_denom[i - 1] = (S[i - 1 + shared.dim.S.mult[1]] + Ea[i - 1 + shared.dim.Ea.mult[1]] + Eb[i - 1 + shared.dim.Eb.mult[1]] + R[i - 1 + shared.dim.R.mult[1]]) * internal.give_dose1_children[i - 1]; + } + for (size_t i = 1; i <= shared.dim.adults_dose1_denom.size; ++i) { + internal.adults_dose1_denom[i - 1] = (S[i - 1 + shared.dim.S.mult[1]] + Ea[i - 1 + shared.dim.Ea.mult[1]] + Eb[i - 1 + shared.dim.Eb.mult[1]] + R[i - 1 + shared.dim.R.mult[1]]) * internal.give_dose1_adults[i - 1]; + } + for (size_t i = 1; i <= shared.dim.adults_dose2_denom.size; ++i) { + internal.adults_dose2_denom[i - 1] = (S[i - 1 + 2 * shared.dim.S.mult[1]] + Ea[i - 1 + 2 * shared.dim.Ea.mult[1]] + Eb[i - 1 + 2 * shared.dim.Eb.mult[1]] + R[i - 1 + 2 * shared.dim.R.mult[1]]) * internal.give_dose2_adults[i - 1]; + } + const real_type new_deaths_15_plus = dust2::array::sum(internal.n_IdD.data(), shared.dim.n_IdD, {3, 15}, {0, shared.dim.n_IdD.dim[1] - 1}) + new_deaths_SW_15_17 + dust2::array::sum(internal.n_IdD.data(), shared.dim.n_IdD, {17, 19}, {0, shared.dim.n_IdD.dim[1] - 1}); + for (size_t i = 1; i <= shared.dim.s_ij_gen_pop.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.s_ij_gen_pop.dim[1]; ++j) { + internal.s_ij_gen_pop[i - 1 + (j - 1) * shared.dim.s_ij_gen_pop.mult[1]] = shared.m_gen_pop[i - 1 + (j - 1) * shared.dim.m_gen_pop.mult[1]] * internal.prop_infectious[j - 1]; + } + } + for (size_t i = 1; i <= shared.dim.s_ij_sex.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.s_ij_sex.dim[1]; ++j) { + internal.s_ij_sex[i - 1 + (j - 1) * shared.dim.s_ij_sex.mult[1]] = shared.m_sex[i - 1 + (j - 1) * shared.dim.m_sex.mult[1]] * internal.prop_infectious[j - 1]; + } + } for (size_t i = 1; i <= shared.dim.n_vaccination_t_S_children.size; ++i) { internal.n_vaccination_t_S_children[i - 1] = 0; } for (size_t i = 1; i <= shared.dim.n_vaccination_t_S_children.size; ++i) { - internal.n_vaccination_t_S_children[i - 1] = (dust2::array::sum(internal.n_eligible_for_dose1_children.data(), shared.dim.n_eligible_for_dose1_children) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_children_t[1] * S[i - 1 + shared.dim.S.mult[1]] * shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose1_children.data(), shared.dim.n_eligible_for_dose1_children)), S[i - 1 + shared.dim.S.mult[1]])); + internal.n_vaccination_t_S_children[i - 1] = (dust2::array::sum(internal.children_dose1_denom.data(), shared.dim.children_dose1_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_children_t[1] * S[i - 1 + shared.dim.S.mult[1]]) / dust2::array::sum(internal.children_dose1_denom.data(), shared.dim.children_dose1_denom)) * internal.give_dose1_children[i - 1]), S[i - 1 + shared.dim.S.mult[1]])); } for (size_t i = 1; i <= shared.dim.n_vaccination_t_S_adults.size; ++i) { internal.n_vaccination_t_S_adults[i - 1] = 0; } for (size_t i = 1; i <= shared.dim.n_vaccination_t_S_adults.size; ++i) { - internal.n_vaccination_t_S_adults[i - 1] = (dust2::array::sum(internal.n_eligible_for_dose1_adults.data(), shared.dim.n_eligible_for_dose1_adults) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_adults_t[1] * S[i - 1 + shared.dim.S.mult[1]] * shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_1st_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose1_adults.data(), shared.dim.n_eligible_for_dose1_adults)), S[i - 1 + shared.dim.S.mult[1]])); + internal.n_vaccination_t_S_adults[i - 1] = (dust2::array::sum(internal.adults_dose1_denom.data(), shared.dim.adults_dose1_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_adults_t[1] * S[i - 1 + shared.dim.S.mult[1]]) / dust2::array::sum(internal.adults_dose1_denom.data(), shared.dim.adults_dose1_denom)) * internal.give_dose1_adults[i - 1]), S[i - 1 + shared.dim.S.mult[1]])); } for (size_t i = 1; i <= shared.dim.n_vaccination_t_Ea_children.size; ++i) { internal.n_vaccination_t_Ea_children[i - 1] = 0; } for (size_t i = 1; i <= shared.dim.n_vaccination_t_Ea_children.size; ++i) { - internal.n_vaccination_t_Ea_children[i - 1] = (dust2::array::sum(internal.n_eligible_for_dose1_children.data(), shared.dim.n_eligible_for_dose1_children) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_children_t[1] * Ea[i - 1 + shared.dim.Ea.mult[1]] * shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose1_children.data(), shared.dim.n_eligible_for_dose1_children)), Ea[i - 1 + shared.dim.Ea.mult[1]])); + internal.n_vaccination_t_Ea_children[i - 1] = (dust2::array::sum(internal.children_dose1_denom.data(), shared.dim.children_dose1_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_children_t[1] * Ea[i - 1 + shared.dim.Ea.mult[1]]) / dust2::array::sum(internal.children_dose1_denom.data(), shared.dim.children_dose1_denom)) * internal.give_dose1_children[i - 1]), Ea[i - 1 + shared.dim.Ea.mult[1]])); } for (size_t i = 1; i <= shared.dim.n_vaccination_t_Ea_adults.size; ++i) { internal.n_vaccination_t_Ea_adults[i - 1] = 0; } for (size_t i = 1; i <= shared.dim.n_vaccination_t_Ea_adults.size; ++i) { - internal.n_vaccination_t_Ea_adults[i - 1] = (dust2::array::sum(internal.n_eligible_for_dose1_adults.data(), shared.dim.n_eligible_for_dose1_adults) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_adults_t[1] * Ea[i - 1 + shared.dim.Ea.mult[1]] * shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_1st_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose1_adults.data(), shared.dim.n_eligible_for_dose1_adults)), Ea[i - 1 + shared.dim.Ea.mult[1]])); + internal.n_vaccination_t_Ea_adults[i - 1] = (dust2::array::sum(internal.adults_dose1_denom.data(), shared.dim.adults_dose1_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_adults_t[1] * Ea[i - 1 + shared.dim.Ea.mult[1]]) / dust2::array::sum(internal.adults_dose1_denom.data(), shared.dim.adults_dose1_denom)) * internal.give_dose1_adults[i - 1]), Ea[i - 1 + shared.dim.Ea.mult[1]])); } for (size_t i = 1; i <= shared.dim.n_vaccination_t_Eb_children.size; ++i) { internal.n_vaccination_t_Eb_children[i - 1] = 0; } for (size_t i = 1; i <= shared.dim.n_vaccination_t_Eb_children.size; ++i) { - internal.n_vaccination_t_Eb_children[i - 1] = (dust2::array::sum(internal.n_eligible_for_dose1_children.data(), shared.dim.n_eligible_for_dose1_children) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_children_t[1] * Eb[i - 1 + shared.dim.Eb.mult[1]] * shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose1_children.data(), shared.dim.n_eligible_for_dose1_children)), Eb[i - 1 + shared.dim.Eb.mult[1]])); + internal.n_vaccination_t_Eb_children[i - 1] = (dust2::array::sum(internal.children_dose1_denom.data(), shared.dim.children_dose1_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_children_t[1] * Eb[i - 1 + shared.dim.Eb.mult[1]]) / dust2::array::sum(internal.children_dose1_denom.data(), shared.dim.children_dose1_denom)) * internal.give_dose1_children[i - 1]), Eb[i - 1 + shared.dim.Eb.mult[1]])); } for (size_t i = 1; i <= shared.dim.n_vaccination_t_Eb_adults.size; ++i) { internal.n_vaccination_t_Eb_adults[i - 1] = 0; } for (size_t i = 1; i <= shared.dim.n_vaccination_t_Eb_adults.size; ++i) { - internal.n_vaccination_t_Eb_adults[i - 1] = (dust2::array::sum(internal.n_eligible_for_dose1_adults.data(), shared.dim.n_eligible_for_dose1_adults) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_adults_t[1] * Eb[i - 1 + shared.dim.Eb.mult[1]] * shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_1st_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose1_adults.data(), shared.dim.n_eligible_for_dose1_adults)), Eb[i - 1 + shared.dim.Eb.mult[1]])); + internal.n_vaccination_t_Eb_adults[i - 1] = (dust2::array::sum(internal.adults_dose1_denom.data(), shared.dim.adults_dose1_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_adults_t[1] * Eb[i - 1 + shared.dim.Eb.mult[1]]) / dust2::array::sum(internal.adults_dose1_denom.data(), shared.dim.adults_dose1_denom)) * internal.give_dose1_adults[i - 1]), Eb[i - 1 + shared.dim.Eb.mult[1]])); } for (size_t i = 1; i <= shared.dim.n_vaccination_t_R_children.size; ++i) { internal.n_vaccination_t_R_children[i - 1] = 0; } for (size_t i = 1; i <= shared.dim.n_vaccination_t_R_children.size; ++i) { - internal.n_vaccination_t_R_children[i - 1] = (dust2::array::sum(internal.n_eligible_for_dose1_children.data(), shared.dim.n_eligible_for_dose1_children) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_children_t[1] * R[i - 1 + shared.dim.R.mult[1]] * shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose1_children.data(), shared.dim.n_eligible_for_dose1_children)), R[i - 1 + shared.dim.R.mult[1]])); + internal.n_vaccination_t_R_children[i - 1] = (dust2::array::sum(internal.children_dose1_denom.data(), shared.dim.children_dose1_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_children_t[1] * R[i - 1 + shared.dim.R.mult[1]]) / dust2::array::sum(internal.children_dose1_denom.data(), shared.dim.children_dose1_denom)) * internal.give_dose1_children[i - 1]), R[i - 1 + shared.dim.R.mult[1]])); } for (size_t i = 1; i <= shared.dim.n_vaccination_t_R_adults.size; ++i) { internal.n_vaccination_t_R_adults[i - 1] = 0; } for (size_t i = 1; i <= shared.dim.n_vaccination_t_R_adults.size; ++i) { - internal.n_vaccination_t_R_adults[i - 1] = (dust2::array::sum(internal.n_eligible_for_dose1_adults.data(), shared.dim.n_eligible_for_dose1_adults) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_adults_t[1] * R[i - 1 + shared.dim.R.mult[1]] * shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_1st_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose1_adults.data(), shared.dim.n_eligible_for_dose1_adults)), R[i - 1 + shared.dim.R.mult[1]])); + internal.n_vaccination_t_R_adults[i - 1] = (dust2::array::sum(internal.adults_dose1_denom.data(), shared.dim.adults_dose1_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_adults_t[1] * R[i - 1 + shared.dim.R.mult[1]]) / dust2::array::sum(internal.adults_dose1_denom.data(), shared.dim.adults_dose1_denom)) * internal.give_dose1_adults[i - 1]), R[i - 1 + shared.dim.R.mult[1]])); } - for (size_t i = 1; i <= shared.dim.D.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.D.dim[1]; ++j) { - internal.new_D[i - 1 + (j - 1) * shared.dim.D.mult[1]] = D[i - 1 + (j - 1) * shared.dim.D.mult[1]] + internal.delta_D[i - 1 + (j - 1) * shared.dim.delta_D.mult[1]]; + for (size_t i = 1; i <= shared.dim.lambda_hh.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.lambda_hh.dim[1]; ++j) { + internal.lambda_hh[i - 1 + (j - 1) * shared.dim.lambda_hh.mult[1]] = shared.beta_h * dust2::array::sum(internal.s_ij_gen_pop.data(), shared.dim.s_ij_gen_pop, {i - 1, i - 1}, {0, shared.dim.s_ij_gen_pop.dim[1] - 1}) * (1 - shared.ve_I[i - 1 + (j - 1) * shared.dim.ve_I.mult[1]]); } } - const real_type new_deaths_SW_15_17 = dust2::array::sum(internal.n_IdD.data(), shared.dim.n_IdD, {16, 16}, {0, shared.dim.n_IdD.dim[1] - 1}) - new_deaths_SW_12_14; - const real_type new_deaths_05_14 = dust2::array::sum(internal.n_IdD.data(), shared.dim.n_IdD, {1, 2}, {0, shared.dim.n_IdD.dim[1] - 1}) + new_deaths_SW_12_14; - const real_type new_deaths_SW = new_deaths_CSW + new_deaths_ASW; - for (size_t i = 1; i <= shared.dim.prop_infectious.size; ++i) { - internal.prop_infectious[i - 1] = (dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1}) == 0 ? 0 : dust2::array::sum(internal.I_infectious.data(), shared.dim.I_infectious, {i - 1, i - 1}, {0, shared.dim.I_infectious.dim[1] - 1}) / dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1})); - } - for (size_t i = 1; i <= shared.dim.lambda_hc.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.lambda_hc.dim[1]; ++j) { - internal.lambda_hc[i - 1 + (j - 1) * shared.dim.lambda_hc.mult[1]] = (i == 20 ? shared.beta_hcw * dust2::array::sum(internal.I_infectious.data(), shared.dim.I_infectious) / dust2::array::sum(N, shared.dim.N) * (1 - shared.ve_I[i - 1 + (j - 1) * shared.dim.ve_I.mult[1]]) : 0); + for (size_t i = 1; i <= shared.dim.lambda_s.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.lambda_s.dim[1]; ++j) { + internal.lambda_s[i - 1 + (j - 1) * shared.dim.lambda_s.mult[1]] = shared.beta_s * dust2::array::sum(internal.s_ij_sex.data(), shared.dim.s_ij_sex, {i - 1, i - 1}, {0, shared.dim.s_ij_sex.dim[1] - 1}) * (1 - shared.ve_I[i - 1 + (j - 1) * shared.dim.ve_I.mult[1]]); } } for (size_t i = 1; i <= shared.dim.n_vaccination_t_S.dim[0]; ++i) { @@ -1077,9 +1119,8 @@ class model_targeted_vax { for (size_t i = 1; i <= shared.dim.n_vaccination_t_S.dim[0]; ++i) { internal.n_vaccination_t_S[i - 1 + shared.dim.n_vaccination_t_S.mult[1]] = internal.n_vaccination_t_S_children[i - 1] + internal.n_vaccination_t_S_adults[i - 1]; } - internal.n_vaccination_t_S[2 + shared.dim.n_vaccination_t_S.mult[1]] = monty::math::min(internal.n_vaccination_t_S[2 + shared.dim.n_vaccination_t_S.mult[1]], S[2 + shared.dim.S.mult[1]]); for (size_t i = 1; i <= shared.dim.n_vaccination_t_S.dim[0]; ++i) { - internal.n_vaccination_t_S[i - 1 + 2 * shared.dim.n_vaccination_t_S.mult[1]] = (dust2::array::sum(internal.n_eligible_for_dose2_adults.data(), shared.dim.n_eligible_for_dose2_adults) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_adults_t[2] * S[i - 1 + 2 * shared.dim.S.mult[1]] * shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose2_adults.data(), shared.dim.n_eligible_for_dose2_adults)), S[i - 1 + 2 * shared.dim.S.mult[1]])); + internal.n_vaccination_t_S[i - 1 + 2 * shared.dim.n_vaccination_t_S.mult[1]] = (dust2::array::sum(internal.adults_dose2_denom.data(), shared.dim.adults_dose2_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_adults_t[2] * S[i - 1 + 2 * shared.dim.S.mult[1]]) / dust2::array::sum(internal.adults_dose2_denom.data(), shared.dim.adults_dose2_denom)) * internal.give_dose2_adults[i - 1]), S[i - 1 + 2 * shared.dim.S.mult[1]])); } for (size_t i = 1; i <= shared.dim.n_vaccination_t_Ea.dim[0]; ++i) { for (size_t j = 1; j <= shared.dim.n_vaccination_t_Ea.dim[1]; ++j) { @@ -1089,9 +1130,8 @@ class model_targeted_vax { for (size_t i = 1; i <= shared.dim.n_vaccination_t_Ea.dim[0]; ++i) { internal.n_vaccination_t_Ea[i - 1 + shared.dim.n_vaccination_t_Ea.mult[1]] = internal.n_vaccination_t_Ea_children[i - 1] + internal.n_vaccination_t_Ea_adults[i - 1]; } - internal.n_vaccination_t_Ea[2 + shared.dim.n_vaccination_t_Ea.mult[1]] = monty::math::min(internal.n_vaccination_t_Ea[2 + shared.dim.n_vaccination_t_Ea.mult[1]], Ea[2 + shared.dim.Ea.mult[1]]); for (size_t i = 1; i <= shared.dim.n_vaccination_t_Ea.dim[0]; ++i) { - internal.n_vaccination_t_Ea[i - 1 + 2 * shared.dim.n_vaccination_t_Ea.mult[1]] = (dust2::array::sum(internal.n_eligible_for_dose2_adults.data(), shared.dim.n_eligible_for_dose2_adults) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_adults_t[2] * Ea[i - 1 + 2 * shared.dim.Ea.mult[1]] * shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose2_adults.data(), shared.dim.n_eligible_for_dose2_adults)), Ea[i - 1 + 2 * shared.dim.Ea.mult[1]])); + internal.n_vaccination_t_Ea[i - 1 + 2 * shared.dim.n_vaccination_t_Ea.mult[1]] = (dust2::array::sum(internal.adults_dose2_denom.data(), shared.dim.adults_dose2_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_adults_t[2] * Ea[i - 1 + 2 * shared.dim.Ea.mult[1]]) / dust2::array::sum(internal.adults_dose2_denom.data(), shared.dim.adults_dose2_denom)) * internal.give_dose2_adults[i - 1]), Ea[i - 1 + 2 * shared.dim.Ea.mult[1]])); } for (size_t i = 1; i <= shared.dim.n_vaccination_t_Eb.dim[0]; ++i) { for (size_t j = 1; j <= shared.dim.n_vaccination_t_Eb.dim[1]; ++j) { @@ -1101,9 +1141,8 @@ class model_targeted_vax { for (size_t i = 1; i <= shared.dim.n_vaccination_t_Eb.dim[0]; ++i) { internal.n_vaccination_t_Eb[i - 1 + shared.dim.n_vaccination_t_Eb.mult[1]] = internal.n_vaccination_t_Eb_children[i - 1] + internal.n_vaccination_t_Eb_adults[i - 1]; } - internal.n_vaccination_t_Eb[2 + shared.dim.n_vaccination_t_Eb.mult[1]] = monty::math::min(internal.n_vaccination_t_Eb[2 + shared.dim.n_vaccination_t_Eb.mult[1]], Eb[2 + shared.dim.Eb.mult[1]]); for (size_t i = 1; i <= shared.dim.n_vaccination_t_Eb.dim[0]; ++i) { - internal.n_vaccination_t_Eb[i - 1 + 2 * shared.dim.n_vaccination_t_Eb.mult[1]] = (dust2::array::sum(internal.n_eligible_for_dose2_adults.data(), shared.dim.n_eligible_for_dose2_adults) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_adults_t[2] * Eb[i - 1 + 2 * shared.dim.Eb.mult[1]] * shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose2_adults.data(), shared.dim.n_eligible_for_dose2_adults)), Eb[i - 1 + 2 * shared.dim.Eb.mult[1]])); + internal.n_vaccination_t_Eb[i - 1 + 2 * shared.dim.n_vaccination_t_Eb.mult[1]] = (dust2::array::sum(internal.adults_dose2_denom.data(), shared.dim.adults_dose2_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_adults_t[2] * Eb[i - 1 + 2 * shared.dim.Eb.mult[1]]) / dust2::array::sum(internal.adults_dose2_denom.data(), shared.dim.adults_dose2_denom)) * internal.give_dose2_adults[i - 1]), Eb[i - 1 + 2 * shared.dim.Eb.mult[1]])); } for (size_t i = 1; i <= shared.dim.n_vaccination_t_R.dim[0]; ++i) { for (size_t j = 1; j <= shared.dim.n_vaccination_t_R.dim[1]; ++j) { @@ -1113,19 +1152,12 @@ class model_targeted_vax { for (size_t i = 1; i <= shared.dim.n_vaccination_t_R.dim[0]; ++i) { internal.n_vaccination_t_R[i - 1 + shared.dim.n_vaccination_t_R.mult[1]] = internal.n_vaccination_t_R_children[i - 1] + internal.n_vaccination_t_R_adults[i - 1]; } - internal.n_vaccination_t_R[2 + shared.dim.n_vaccination_t_R.mult[1]] = monty::math::min(internal.n_vaccination_t_R[2 + shared.dim.n_vaccination_t_R.mult[1]], R[2 + shared.dim.R.mult[1]]); for (size_t i = 1; i <= shared.dim.n_vaccination_t_R.dim[0]; ++i) { - internal.n_vaccination_t_R[i - 1 + 2 * shared.dim.n_vaccination_t_R.mult[1]] = (dust2::array::sum(internal.n_eligible_for_dose2_adults.data(), shared.dim.n_eligible_for_dose2_adults) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_adults_t[2] * R[i - 1 + 2 * shared.dim.R.mult[1]] * shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose2_adults.data(), shared.dim.n_eligible_for_dose2_adults)), R[i - 1 + 2 * shared.dim.R.mult[1]])); - } - const real_type new_deaths_15_plus = dust2::array::sum(internal.n_IdD.data(), shared.dim.n_IdD, {3, 15}, {0, shared.dim.n_IdD.dim[1] - 1}) + new_deaths_SW_15_17 + dust2::array::sum(internal.n_IdD.data(), shared.dim.n_IdD, {17, 19}, {0, shared.dim.n_IdD.dim[1] - 1}); - for (size_t i = 1; i <= shared.dim.s_ij_gen_pop.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.s_ij_gen_pop.dim[1]; ++j) { - internal.s_ij_gen_pop[i - 1 + (j - 1) * shared.dim.s_ij_gen_pop.mult[1]] = shared.m_gen_pop[i - 1 + (j - 1) * shared.dim.m_gen_pop.mult[1]] * internal.prop_infectious[j - 1]; - } + internal.n_vaccination_t_R[i - 1 + 2 * shared.dim.n_vaccination_t_R.mult[1]] = (dust2::array::sum(internal.adults_dose2_denom.data(), shared.dim.adults_dose2_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_adults_t[2] * R[i - 1 + 2 * shared.dim.R.mult[1]]) / dust2::array::sum(internal.adults_dose2_denom.data(), shared.dim.adults_dose2_denom)) * internal.give_dose2_adults[i - 1]), R[i - 1 + 2 * shared.dim.R.mult[1]])); } - for (size_t i = 1; i <= shared.dim.s_ij_sex.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.s_ij_sex.dim[1]; ++j) { - internal.s_ij_sex[i - 1 + (j - 1) * shared.dim.s_ij_sex.mult[1]] = shared.m_sex[i - 1 + (j - 1) * shared.dim.m_sex.mult[1]] * internal.prop_infectious[j - 1]; + for (size_t i = 1; i <= shared.dim.lambda.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.lambda.dim[1]; ++j) { + internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] = internal.lambda_hh[i - 1 + (j - 1) * shared.dim.lambda_hh.mult[1]] + internal.lambda_s[i - 1 + (j - 1) * shared.dim.lambda_s.mult[1]] + internal.lambda_hc[i - 1 + (j - 1) * shared.dim.lambda_hc.mult[1]] + shared.lambda_z[i - 1 + (j - 1) * shared.dim.lambda_z.mult[1]]; } } for (size_t i = 1; i <= shared.dim.delta_S_n_vaccination.dim[0]; ++i) { @@ -1153,14 +1185,24 @@ class model_targeted_vax { internal.n_vaccination_t[i - 1 + (j - 1) * shared.dim.n_vaccination_t.mult[1]] = internal.n_vaccination_t_S[i - 1 + (j - 1) * shared.dim.n_vaccination_t_S.mult[1]] + internal.n_vaccination_t_Ea[i - 1 + (j - 1) * shared.dim.n_vaccination_t_Ea.mult[1]] + internal.n_vaccination_t_Eb[i - 1 + (j - 1) * shared.dim.n_vaccination_t_Eb.mult[1]] + internal.n_vaccination_t_R[i - 1 + (j - 1) * shared.dim.n_vaccination_t_R.mult[1]]; } } - for (size_t i = 1; i <= shared.dim.lambda_hh.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.lambda_hh.dim[1]; ++j) { - internal.lambda_hh[i - 1 + (j - 1) * shared.dim.lambda_hh.mult[1]] = shared.beta_h * dust2::array::sum(internal.s_ij_gen_pop.data(), shared.dim.s_ij_gen_pop, {i - 1, i - 1}, {0, shared.dim.s_ij_gen_pop.dim[1] - 1}) * (1 - shared.ve_I[i - 1 + (j - 1) * shared.dim.ve_I.mult[1]]); + for (size_t i = 1; i <= shared.dim.p_SE.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.p_SE.dim[1]; ++j) { + internal.p_SE[i - 1 + (j - 1) * shared.dim.p_SE.mult[1]] = 1 - monty::math::exp(-internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] * dt); } } - for (size_t i = 1; i <= shared.dim.lambda_s.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.lambda_s.dim[1]; ++j) { - internal.lambda_s[i - 1 + (j - 1) * shared.dim.lambda_s.mult[1]] = shared.beta_s * dust2::array::sum(internal.s_ij_sex.data(), shared.dim.s_ij_sex, {i - 1, i - 1}, {0, shared.dim.s_ij_sex.dim[1] - 1}) * (1 - shared.ve_I[i - 1 + (j - 1) * shared.dim.ve_I.mult[1]]); + for (size_t i = 1; i <= shared.dim.p_hh.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.p_hh.dim[1]; ++j) { + internal.p_hh[i - 1 + (j - 1) * shared.dim.p_hh.mult[1]] = (internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] > 0 ? internal.lambda_hh[i - 1 + (j - 1) * shared.dim.lambda_hh.mult[1]] / internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] : 0); + } + } + for (size_t i = 1; i <= shared.dim.p_s.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.p_s.dim[1]; ++j) { + internal.p_s[i - 1 + (j - 1) * shared.dim.p_s.mult[1]] = (internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] > 0 ? internal.lambda_s[i - 1 + (j - 1) * shared.dim.lambda_s.mult[1]] / internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] : 0); + } + } + for (size_t i = 1; i <= shared.dim.p_hc.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.p_hc.dim[1]; ++j) { + internal.p_hc[i - 1 + (j - 1) * shared.dim.p_hc.mult[1]] = (internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] > 0 ? internal.lambda_hc[i - 1 + (j - 1) * shared.dim.lambda_hc.mult[1]] / internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] : 0); } } for (size_t i = 1; i <= shared.dim.R.dim[0]; ++i) { @@ -1182,9 +1224,9 @@ class model_targeted_vax { const real_type new_dose2_ASW = internal.n_vaccination_t[17 + 2 * shared.dim.n_vaccination_t.mult[1]]; const real_type new_dose2_PBS = internal.n_vaccination_t[18 + 2 * shared.dim.n_vaccination_t.mult[1]]; const real_type new_dose2_HCW = internal.n_vaccination_t[19 + 2 * shared.dim.n_vaccination_t.mult[1]]; - for (size_t i = 1; i <= shared.dim.lambda.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.lambda.dim[1]; ++j) { - internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] = internal.lambda_hh[i - 1 + (j - 1) * shared.dim.lambda_hh.mult[1]] + internal.lambda_s[i - 1 + (j - 1) * shared.dim.lambda_s.mult[1]] + internal.lambda_hc[i - 1 + (j - 1) * shared.dim.lambda_hc.mult[1]] + shared.lambda_z[i - 1 + (j - 1) * shared.dim.lambda_z.mult[1]]; + for (size_t i = 1; i <= shared.dim.n_SEa.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.n_SEa.dim[1]; ++j) { + internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]] = monty::random::binomial(rng_state, S[i - 1 + (j - 1) * shared.dim.S.mult[1]] + internal.delta_S_n_vaccination[i - 1 + (j - 1) * shared.dim.delta_S_n_vaccination.mult[1]], internal.p_SE[i - 1 + (j - 1) * shared.dim.p_SE.mult[1]]); } } for (size_t i = 1; i <= shared.dim.n_EaEb.dim[0]; ++i) { @@ -1197,30 +1239,26 @@ class model_targeted_vax { internal.n_EbI[i - 1 + (j - 1) * shared.dim.n_EbI.mult[1]] = monty::random::binomial(rng_state, Eb[i - 1 + (j - 1) * shared.dim.Eb.mult[1]] + internal.delta_Eb_n_vaccination[i - 1 + (j - 1) * shared.dim.delta_Eb_n_vaccination.mult[1]], p_EI); } } + for (size_t i = 1; i <= shared.dim.S.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.S.dim[1]; ++j) { + internal.new_S[i - 1 + (j - 1) * shared.dim.S.mult[1]] = S[i - 1 + (j - 1) * shared.dim.S.mult[1]] + internal.delta_S_n_vaccination[i - 1 + (j - 1) * shared.dim.delta_S_n_vaccination.mult[1]] - internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]]; + } + } + const real_type new_cases_00_04 = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {0, 0}, {0, shared.dim.n_SEa.dim[1] - 1}); + const real_type new_cases_SW_12_14 = monty::random::binomial(rng_state, dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {16, 16}, {0, shared.dim.n_SEa.dim[1] - 1}), static_cast(0.5)); + const real_type new_cases_CSW = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {16, 16}, {0, shared.dim.n_SEa.dim[1] - 1}); + const real_type new_cases_ASW = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {17, 17}, {0, shared.dim.n_SEa.dim[1] - 1}); + const real_type new_cases_PBS = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {18, 18}, {0, shared.dim.n_SEa.dim[1] - 1}); + const real_type new_cases_HCW = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {19, 19}, {0, shared.dim.n_SEa.dim[1] - 1}); const real_type new_dose1_SW_15_17 = internal.n_vaccination_t[16 + shared.dim.n_vaccination_t.mult[1]] - new_dose1_SW_12_14; const real_type new_dose1_05_14 = dust2::array::sum(internal.n_vaccination_t.data(), shared.dim.n_vaccination_t, {1, 2}, {1, 1}) + new_dose1_SW_12_14; const real_type new_dose1_SW = new_dose1_CSW + new_dose1_ASW; const real_type new_dose2_SW_15_17 = internal.n_vaccination_t[16 + 2 * shared.dim.n_vaccination_t.mult[1]] - new_dose2_SW_12_14; const real_type new_dose2_05_14 = dust2::array::sum(internal.n_vaccination_t.data(), shared.dim.n_vaccination_t, {1, 2}, {2, 2}) + new_dose2_SW_12_14; const real_type new_dose2_SW = new_dose2_CSW + new_dose2_ASW; - for (size_t i = 1; i <= shared.dim.p_SE.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.p_SE.dim[1]; ++j) { - internal.p_SE[i - 1 + (j - 1) * shared.dim.p_SE.mult[1]] = 1 - monty::math::exp(-internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] * dt); - } - } - for (size_t i = 1; i <= shared.dim.p_hh.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.p_hh.dim[1]; ++j) { - internal.p_hh[i - 1 + (j - 1) * shared.dim.p_hh.mult[1]] = (internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] > 0 ? internal.lambda_hh[i - 1 + (j - 1) * shared.dim.lambda_hh.mult[1]] / internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] : 0); - } - } - for (size_t i = 1; i <= shared.dim.p_s.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.p_s.dim[1]; ++j) { - internal.p_s[i - 1 + (j - 1) * shared.dim.p_s.mult[1]] = (internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] > 0 ? internal.lambda_s[i - 1 + (j - 1) * shared.dim.lambda_s.mult[1]] / internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] : 0); - } - } - for (size_t i = 1; i <= shared.dim.p_hc.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.p_hc.dim[1]; ++j) { - internal.p_hc[i - 1 + (j - 1) * shared.dim.p_hc.mult[1]] = (internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] > 0 ? internal.lambda_hc[i - 1 + (j - 1) * shared.dim.lambda_hc.mult[1]] / internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] : 0); + for (size_t i = 1; i <= shared.dim.n_SEa_hh.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.n_SEa_hh.dim[1]; ++j) { + internal.n_SEa_hh[i - 1 + (j - 1) * shared.dim.n_SEa_hh.mult[1]] = monty::random::binomial(rng_state, internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]], internal.p_hh[i - 1 + (j - 1) * shared.dim.p_hh.mult[1]]); } } for (size_t i = 1; i <= shared.dim.n_EbId.dim[0]; ++i) { @@ -1228,21 +1266,34 @@ class model_targeted_vax { internal.n_EbId[i - 1 + (j - 1) * shared.dim.n_EbId.mult[1]] = monty::random::binomial(rng_state, internal.n_EbI[i - 1 + (j - 1) * shared.dim.n_EbI.mult[1]], shared.CFR[i - 1 + (j - 1) * shared.dim.CFR.mult[1]]); } } + for (size_t i = 1; i <= shared.dim.delta_Ea.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.delta_Ea.dim[1]; ++j) { + internal.delta_Ea[i - 1 + (j - 1) * shared.dim.delta_Ea.mult[1]] = internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]] - internal.n_EaEb[i - 1 + (j - 1) * shared.dim.n_EaEb.mult[1]]; + } + } for (size_t i = 1; i <= shared.dim.delta_Eb.dim[0]; ++i) { for (size_t j = 1; j <= shared.dim.delta_Eb.dim[1]; ++j) { internal.delta_Eb[i - 1 + (j - 1) * shared.dim.delta_Eb.mult[1]] = internal.n_EaEb[i - 1 + (j - 1) * shared.dim.n_EaEb.mult[1]] - internal.n_EbI[i - 1 + (j - 1) * shared.dim.n_EbI.mult[1]]; } } + for (size_t i = 1; i <= shared.dim.Ea.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.Ea.dim[1]; ++j) { + internal.new_Ea[i - 1 + (j - 1) * shared.dim.Ea.mult[1]] = Ea[i - 1 + (j - 1) * shared.dim.Ea.mult[1]] + internal.delta_Ea_n_vaccination[i - 1 + (j - 1) * shared.dim.delta_Ea_n_vaccination.mult[1]] + internal.delta_Ea[i - 1 + (j - 1) * shared.dim.delta_Ea.mult[1]]; + } + } for (size_t i = 1; i <= shared.dim.Eb.dim[0]; ++i) { for (size_t j = 1; j <= shared.dim.Eb.dim[1]; ++j) { internal.new_Eb[i - 1 + (j - 1) * shared.dim.Eb.mult[1]] = Eb[i - 1 + (j - 1) * shared.dim.Eb.mult[1]] + internal.delta_Eb_n_vaccination[i - 1 + (j - 1) * shared.dim.delta_Eb_n_vaccination.mult[1]] + internal.delta_Eb[i - 1 + (j - 1) * shared.dim.delta_Eb.mult[1]]; } } + const real_type new_cases_SW_15_17 = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {16, 16}, {0, shared.dim.n_SEa.dim[1] - 1}) - new_cases_SW_12_14; + const real_type new_cases_05_14 = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {1, 2}, {0, shared.dim.n_SEa.dim[1] - 1}) + new_cases_SW_12_14; + const real_type new_cases_SW = new_cases_CSW + new_cases_ASW; const real_type new_dose1_15_plus = dust2::array::sum(internal.n_vaccination_t.data(), shared.dim.n_vaccination_t, {3, 15}, {1, 1}) + new_dose1_SW_15_17 + dust2::array::sum(internal.n_vaccination_t.data(), shared.dim.n_vaccination_t, {17, 19}, {1, 1}); const real_type new_dose2_15_plus = dust2::array::sum(internal.n_vaccination_t.data(), shared.dim.n_vaccination_t, {3, 15}, {2, 2}) + new_dose2_SW_15_17 + dust2::array::sum(internal.n_vaccination_t.data(), shared.dim.n_vaccination_t, {17, 19}, {2, 2}); - for (size_t i = 1; i <= shared.dim.n_SEa.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.n_SEa.dim[1]; ++j) { - internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]] = monty::random::binomial(rng_state, S[i - 1 + (j - 1) * shared.dim.S.mult[1]] + internal.delta_S_n_vaccination[i - 1 + (j - 1) * shared.dim.delta_S_n_vaccination.mult[1]], internal.p_SE[i - 1 + (j - 1) * shared.dim.p_SE.mult[1]]); + for (size_t i = 1; i <= shared.dim.n_SEa_s.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.n_SEa_s.dim[1]; ++j) { + internal.n_SEa_s[i - 1 + (j - 1) * shared.dim.n_SEa_s.mult[1]] = monty::random::binomial(rng_state, internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]] - internal.n_SEa_hh[i - 1 + (j - 1) * shared.dim.n_SEa_hh.mult[1]], internal.p_s[i - 1 + (j - 1) * shared.dim.p_s.mult[1]]); } } for (size_t i = 1; i <= shared.dim.n_EbIr.dim[0]; ++i) { @@ -1255,30 +1306,20 @@ class model_targeted_vax { internal.delta_Id[i - 1 + (j - 1) * shared.dim.delta_Id.mult[1]] = internal.n_EbId[i - 1 + (j - 1) * shared.dim.n_EbId.mult[1]] - internal.n_IdD[i - 1 + (j - 1) * shared.dim.n_IdD.mult[1]]; } } - for (size_t i = 1; i <= shared.dim.S.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.S.dim[1]; ++j) { - internal.new_S[i - 1 + (j - 1) * shared.dim.S.mult[1]] = S[i - 1 + (j - 1) * shared.dim.S.mult[1]] + internal.delta_S_n_vaccination[i - 1 + (j - 1) * shared.dim.delta_S_n_vaccination.mult[1]] - internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]]; - } - } for (size_t i = 1; i <= shared.dim.Id.dim[0]; ++i) { for (size_t j = 1; j <= shared.dim.Id.dim[1]; ++j) { internal.new_Id[i - 1 + (j - 1) * shared.dim.Id.mult[1]] = Id[i - 1 + (j - 1) * shared.dim.Id.mult[1]] + internal.delta_Id[i - 1 + (j - 1) * shared.dim.delta_Id.mult[1]]; } } - const real_type new_cases_00_04 = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {0, 0}, {0, shared.dim.n_SEa.dim[1] - 1}); - const real_type new_cases_SW_12_14 = monty::random::binomial(rng_state, dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {16, 16}, {0, shared.dim.n_SEa.dim[1] - 1}), static_cast(0.5)); - const real_type new_cases_CSW = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {16, 16}, {0, shared.dim.n_SEa.dim[1] - 1}); - const real_type new_cases_ASW = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {17, 17}, {0, shared.dim.n_SEa.dim[1] - 1}); - const real_type new_cases_PBS = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {18, 18}, {0, shared.dim.n_SEa.dim[1] - 1}); - const real_type new_cases_HCW = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {19, 19}, {0, shared.dim.n_SEa.dim[1] - 1}); - for (size_t i = 1; i <= shared.dim.n_SEa_hh.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.n_SEa_hh.dim[1]; ++j) { - internal.n_SEa_hh[i - 1 + (j - 1) * shared.dim.n_SEa_hh.mult[1]] = monty::random::binomial(rng_state, internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]], internal.p_hh[i - 1 + (j - 1) * shared.dim.p_hh.mult[1]]); + for (size_t i = 1; i <= shared.dim.E.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.E.dim[1]; ++j) { + internal.new_E[i - 1 + (j - 1) * shared.dim.E.mult[1]] = internal.new_Ea[i - 1 + (j - 1) * shared.dim.Ea.mult[1]] + internal.new_Eb[i - 1 + (j - 1) * shared.dim.Eb.mult[1]]; } } - for (size_t i = 1; i <= shared.dim.delta_Ea.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.delta_Ea.dim[1]; ++j) { - internal.delta_Ea[i - 1 + (j - 1) * shared.dim.delta_Ea.mult[1]] = internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]] - internal.n_EaEb[i - 1 + (j - 1) * shared.dim.n_EaEb.mult[1]]; + const real_type new_cases_15_plus = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {3, 15}, {0, shared.dim.n_SEa.dim[1] - 1}) + new_cases_SW_15_17 + dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {17, 19}, {0, shared.dim.n_SEa.dim[1] - 1}); + for (size_t i = 1; i <= shared.dim.n_SEa_hc.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.n_SEa_hc.dim[1]; ++j) { + internal.n_SEa_hc[i - 1 + (j - 1) * shared.dim.n_SEa_hc.mult[1]] = monty::random::binomial(rng_state, internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]] - internal.n_SEa_hh[i - 1 + (j - 1) * shared.dim.n_SEa_hh.mult[1]] - internal.n_SEa_s[i - 1 + (j - 1) * shared.dim.n_SEa_s.mult[1]], internal.p_hc[i - 1 + (j - 1) * shared.dim.p_hc.mult[1]]); } } for (size_t i = 1; i <= shared.dim.delta_Ir.dim[0]; ++i) { @@ -1286,27 +1327,14 @@ class model_targeted_vax { internal.delta_Ir[i - 1 + (j - 1) * shared.dim.delta_Ir.mult[1]] = internal.n_EbIr[i - 1 + (j - 1) * shared.dim.n_EbIr.mult[1]] - internal.n_IrR[i - 1 + (j - 1) * shared.dim.n_IrR.mult[1]]; } } - for (size_t i = 1; i <= shared.dim.Ea.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.Ea.dim[1]; ++j) { - internal.new_Ea[i - 1 + (j - 1) * shared.dim.Ea.mult[1]] = Ea[i - 1 + (j - 1) * shared.dim.Ea.mult[1]] + internal.delta_Ea_n_vaccination[i - 1 + (j - 1) * shared.dim.delta_Ea_n_vaccination.mult[1]] + internal.delta_Ea[i - 1 + (j - 1) * shared.dim.delta_Ea.mult[1]]; - } - } for (size_t i = 1; i <= shared.dim.Ir.dim[0]; ++i) { for (size_t j = 1; j <= shared.dim.Ir.dim[1]; ++j) { internal.new_Ir[i - 1 + (j - 1) * shared.dim.Ir.mult[1]] = Ir[i - 1 + (j - 1) * shared.dim.Ir.mult[1]] + internal.delta_Ir[i - 1 + (j - 1) * shared.dim.delta_Ir.mult[1]]; } } - const real_type new_cases_SW_15_17 = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {16, 16}, {0, shared.dim.n_SEa.dim[1] - 1}) - new_cases_SW_12_14; - const real_type new_cases_05_14 = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {1, 2}, {0, shared.dim.n_SEa.dim[1] - 1}) + new_cases_SW_12_14; - const real_type new_cases_SW = new_cases_CSW + new_cases_ASW; - for (size_t i = 1; i <= shared.dim.n_SEa_s.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.n_SEa_s.dim[1]; ++j) { - internal.n_SEa_s[i - 1 + (j - 1) * shared.dim.n_SEa_s.mult[1]] = monty::random::binomial(rng_state, internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]] - internal.n_SEa_hh[i - 1 + (j - 1) * shared.dim.n_SEa_hh.mult[1]], internal.p_s[i - 1 + (j - 1) * shared.dim.p_s.mult[1]]); - } - } - for (size_t i = 1; i <= shared.dim.E.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.E.dim[1]; ++j) { - internal.new_E[i - 1 + (j - 1) * shared.dim.E.mult[1]] = internal.new_Ea[i - 1 + (j - 1) * shared.dim.Ea.mult[1]] + internal.new_Eb[i - 1 + (j - 1) * shared.dim.Eb.mult[1]]; + for (size_t i = 1; i <= shared.dim.n_SEa_z.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.n_SEa_z.dim[1]; ++j) { + internal.n_SEa_z[i - 1 + (j - 1) * shared.dim.n_SEa_z.mult[1]] = internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]] - internal.n_SEa_hh[i - 1 + (j - 1) * shared.dim.n_SEa_hh.mult[1]] - internal.n_SEa_s[i - 1 + (j - 1) * shared.dim.n_SEa_s.mult[1]] - internal.n_SEa_hc[i - 1 + (j - 1) * shared.dim.n_SEa_hc.mult[1]]; } } for (size_t i = 1; i <= shared.dim.I.dim[0]; ++i) { @@ -1319,17 +1347,6 @@ class model_targeted_vax { internal.new_N[i - 1 + (j - 1) * shared.dim.N.mult[1]] = internal.new_S[i - 1 + (j - 1) * shared.dim.S.mult[1]] + internal.new_Ea[i - 1 + (j - 1) * shared.dim.Ea.mult[1]] + internal.new_Eb[i - 1 + (j - 1) * shared.dim.Eb.mult[1]] + internal.new_Ir[i - 1 + (j - 1) * shared.dim.Ir.mult[1]] + internal.new_Id[i - 1 + (j - 1) * shared.dim.Id.mult[1]] + internal.new_R[i - 1 + (j - 1) * shared.dim.R.mult[1]] + internal.new_D[i - 1 + (j - 1) * shared.dim.D.mult[1]]; } } - const real_type new_cases_15_plus = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {3, 15}, {0, shared.dim.n_SEa.dim[1] - 1}) + new_cases_SW_15_17 + dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {17, 19}, {0, shared.dim.n_SEa.dim[1] - 1}); - for (size_t i = 1; i <= shared.dim.n_SEa_hc.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.n_SEa_hc.dim[1]; ++j) { - internal.n_SEa_hc[i - 1 + (j - 1) * shared.dim.n_SEa_hc.mult[1]] = monty::random::binomial(rng_state, internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]] - internal.n_SEa_hh[i - 1 + (j - 1) * shared.dim.n_SEa_hh.mult[1]] - internal.n_SEa_s[i - 1 + (j - 1) * shared.dim.n_SEa_s.mult[1]], internal.p_hc[i - 1 + (j - 1) * shared.dim.p_hc.mult[1]]); - } - } - for (size_t i = 1; i <= shared.dim.n_SEa_z.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.n_SEa_z.dim[1]; ++j) { - internal.n_SEa_z[i - 1 + (j - 1) * shared.dim.n_SEa_z.mult[1]] = internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]] - internal.n_SEa_hh[i - 1 + (j - 1) * shared.dim.n_SEa_hh.mult[1]] - internal.n_SEa_s[i - 1 + (j - 1) * shared.dim.n_SEa_s.mult[1]] - internal.n_SEa_hc[i - 1 + (j - 1) * shared.dim.n_SEa_hc.mult[1]]; - } - } state_next[0] = (prioritisation_step_1st_dose_children_proposal > shared.N_prioritisation_steps_children ? shared.N_prioritisation_steps_children : prioritisation_step_1st_dose_children_proposal); state_next[1] = (prioritisation_step_1st_dose_adults_proposal > shared.N_prioritisation_steps_adults ? shared.N_prioritisation_steps_adults : prioritisation_step_1st_dose_adults_proposal); state_next[2] = (prioritisation_step_2nd_dose_adults_proposal > shared.N_prioritisation_steps_adults ? shared.N_prioritisation_steps_adults : prioritisation_step_2nd_dose_adults_proposal); diff --git a/src/model-targeted-vax.cpp b/src/model-targeted-vax.cpp index 5668e49..30ff1cc 100644 --- a/src/model-targeted-vax.cpp +++ b/src/model-targeted-vax.cpp @@ -73,9 +73,12 @@ class model_targeted_vax { dust2::array::dimensions<1> coverage_target_1st_dose_children; dust2::array::dimensions<1> coverage_target_1st_dose_adults; dust2::array::dimensions<1> coverage_target_2nd_dose_adults; - dust2::array::dimensions<1> n_eligible_for_dose1_children; - dust2::array::dimensions<1> n_eligible_for_dose1_adults; - dust2::array::dimensions<1> n_eligible_for_dose2_adults; + dust2::array::dimensions<1> give_dose1_children; + dust2::array::dimensions<1> give_dose1_adults; + dust2::array::dimensions<1> give_dose2_adults; + dust2::array::dimensions<1> children_dose1_denom; + dust2::array::dimensions<1> adults_dose1_denom; + dust2::array::dimensions<1> adults_dose2_denom; dust2::array::dimensions<2> N; dust2::array::dimensions<2> S; dust2::array::dimensions<2> S0; @@ -199,14 +202,22 @@ class model_targeted_vax { std::vector coverage_target_1st_dose_children; std::vector coverage_target_1st_dose_adults; std::vector coverage_target_2nd_dose_adults; - std::vector n_eligible_for_dose1_children; - std::vector n_eligible_for_dose1_adults; - std::vector n_eligible_for_dose2_adults; std::vector I_infectious; std::vector delta_R; std::vector delta_D; std::vector daily_doses_children_t; std::vector daily_doses_adults_t; + std::vector give_dose1_children; + std::vector give_dose1_adults; + std::vector give_dose2_adults; + std::vector new_D; + std::vector prop_infectious; + std::vector lambda_hc; + std::vector children_dose1_denom; + std::vector adults_dose1_denom; + std::vector adults_dose2_denom; + std::vector s_ij_gen_pop; + std::vector s_ij_sex; std::vector n_vaccination_t_S_children; std::vector n_vaccination_t_S_adults; std::vector n_vaccination_t_Ea_children; @@ -215,49 +226,44 @@ class model_targeted_vax { std::vector n_vaccination_t_Eb_adults; std::vector n_vaccination_t_R_children; std::vector n_vaccination_t_R_adults; - std::vector new_D; - std::vector prop_infectious; - std::vector lambda_hc; + std::vector lambda_hh; + std::vector lambda_s; std::vector n_vaccination_t_S; std::vector n_vaccination_t_Ea; std::vector n_vaccination_t_Eb; std::vector n_vaccination_t_R; - std::vector s_ij_gen_pop; - std::vector s_ij_sex; + std::vector lambda; std::vector delta_S_n_vaccination; std::vector delta_Ea_n_vaccination; std::vector delta_Eb_n_vaccination; std::vector delta_R_n_vaccination; std::vector n_vaccination_t; - std::vector lambda_hh; - std::vector lambda_s; - std::vector new_R; - std::vector lambda; - std::vector n_EaEb; - std::vector n_EbI; std::vector p_SE; std::vector p_hh; std::vector p_s; std::vector p_hc; + std::vector new_R; + std::vector n_SEa; + std::vector n_EaEb; + std::vector n_EbI; + std::vector new_S; + std::vector n_SEa_hh; std::vector n_EbId; + std::vector delta_Ea; std::vector delta_Eb; + std::vector new_Ea; std::vector new_Eb; - std::vector n_SEa; + std::vector n_SEa_s; std::vector n_EbIr; std::vector delta_Id; - std::vector new_S; std::vector new_Id; - std::vector n_SEa_hh; - std::vector delta_Ea; + std::vector new_E; + std::vector n_SEa_hc; std::vector delta_Ir; - std::vector new_Ea; std::vector new_Ir; - std::vector n_SEa_s; - std::vector new_E; + std::vector n_SEa_z; std::vector new_I; std::vector new_N; - std::vector n_SEa_hc; - std::vector n_SEa_z; }; struct data_type { real_type cases; @@ -321,9 +327,12 @@ class model_targeted_vax { dim.coverage_target_1st_dose_children.set({static_cast(n_group)}); dim.coverage_target_1st_dose_adults.set({static_cast(n_group)}); dim.coverage_target_2nd_dose_adults.set({static_cast(n_group)}); - dim.n_eligible_for_dose1_children.set({static_cast(n_group)}); - dim.n_eligible_for_dose1_adults.set({static_cast(n_group)}); - dim.n_eligible_for_dose2_adults.set({static_cast(n_group)}); + dim.give_dose1_children.set({static_cast(n_group)}); + dim.give_dose1_adults.set({static_cast(n_group)}); + dim.give_dose2_adults.set({static_cast(n_group)}); + dim.children_dose1_denom.set({static_cast(n_group)}); + dim.adults_dose1_denom.set({static_cast(n_group)}); + dim.adults_dose2_denom.set({static_cast(n_group)}); dim.N.set({static_cast(n_group), static_cast(n_vax)}); dim.S.set({static_cast(n_group), static_cast(n_vax)}); dim.S0.set({static_cast(n_group), static_cast(n_vax)}); @@ -560,14 +569,22 @@ class model_targeted_vax { std::vector coverage_target_1st_dose_children(shared.dim.coverage_target_1st_dose_children.size); std::vector coverage_target_1st_dose_adults(shared.dim.coverage_target_1st_dose_adults.size); std::vector coverage_target_2nd_dose_adults(shared.dim.coverage_target_2nd_dose_adults.size); - std::vector n_eligible_for_dose1_children(shared.dim.n_eligible_for_dose1_children.size); - std::vector n_eligible_for_dose1_adults(shared.dim.n_eligible_for_dose1_adults.size); - std::vector n_eligible_for_dose2_adults(shared.dim.n_eligible_for_dose2_adults.size); std::vector I_infectious(shared.dim.I_infectious.size); std::vector delta_R(shared.dim.delta_R.size); std::vector delta_D(shared.dim.delta_D.size); std::vector daily_doses_children_t(shared.dim.daily_doses_children_t.size); std::vector daily_doses_adults_t(shared.dim.daily_doses_adults_t.size); + std::vector give_dose1_children(shared.dim.give_dose1_children.size); + std::vector give_dose1_adults(shared.dim.give_dose1_adults.size); + std::vector give_dose2_adults(shared.dim.give_dose2_adults.size); + std::vector new_D(shared.dim.D.size); + std::vector prop_infectious(shared.dim.prop_infectious.size); + std::vector lambda_hc(shared.dim.lambda_hc.size); + std::vector children_dose1_denom(shared.dim.children_dose1_denom.size); + std::vector adults_dose1_denom(shared.dim.adults_dose1_denom.size); + std::vector adults_dose2_denom(shared.dim.adults_dose2_denom.size); + std::vector s_ij_gen_pop(shared.dim.s_ij_gen_pop.size); + std::vector s_ij_sex(shared.dim.s_ij_sex.size); std::vector n_vaccination_t_S_children(shared.dim.n_vaccination_t_S_children.size); std::vector n_vaccination_t_S_adults(shared.dim.n_vaccination_t_S_adults.size); std::vector n_vaccination_t_Ea_children(shared.dim.n_vaccination_t_Ea_children.size); @@ -576,50 +593,45 @@ class model_targeted_vax { std::vector n_vaccination_t_Eb_adults(shared.dim.n_vaccination_t_Eb_adults.size); std::vector n_vaccination_t_R_children(shared.dim.n_vaccination_t_R_children.size); std::vector n_vaccination_t_R_adults(shared.dim.n_vaccination_t_R_adults.size); - std::vector new_D(shared.dim.D.size); - std::vector prop_infectious(shared.dim.prop_infectious.size); - std::vector lambda_hc(shared.dim.lambda_hc.size); + std::vector lambda_hh(shared.dim.lambda_hh.size); + std::vector lambda_s(shared.dim.lambda_s.size); std::vector n_vaccination_t_S(shared.dim.n_vaccination_t_S.size); std::vector n_vaccination_t_Ea(shared.dim.n_vaccination_t_Ea.size); std::vector n_vaccination_t_Eb(shared.dim.n_vaccination_t_Eb.size); std::vector n_vaccination_t_R(shared.dim.n_vaccination_t_R.size); - std::vector s_ij_gen_pop(shared.dim.s_ij_gen_pop.size); - std::vector s_ij_sex(shared.dim.s_ij_sex.size); + std::vector lambda(shared.dim.lambda.size); std::vector delta_S_n_vaccination(shared.dim.delta_S_n_vaccination.size); std::vector delta_Ea_n_vaccination(shared.dim.delta_Ea_n_vaccination.size); std::vector delta_Eb_n_vaccination(shared.dim.delta_Eb_n_vaccination.size); std::vector delta_R_n_vaccination(shared.dim.delta_R_n_vaccination.size); std::vector n_vaccination_t(shared.dim.n_vaccination_t.size); - std::vector lambda_hh(shared.dim.lambda_hh.size); - std::vector lambda_s(shared.dim.lambda_s.size); - std::vector new_R(shared.dim.R.size); - std::vector lambda(shared.dim.lambda.size); - std::vector n_EaEb(shared.dim.n_EaEb.size); - std::vector n_EbI(shared.dim.n_EbI.size); std::vector p_SE(shared.dim.p_SE.size); std::vector p_hh(shared.dim.p_hh.size); std::vector p_s(shared.dim.p_s.size); std::vector p_hc(shared.dim.p_hc.size); + std::vector new_R(shared.dim.R.size); + std::vector n_SEa(shared.dim.n_SEa.size); + std::vector n_EaEb(shared.dim.n_EaEb.size); + std::vector n_EbI(shared.dim.n_EbI.size); + std::vector new_S(shared.dim.S.size); + std::vector n_SEa_hh(shared.dim.n_SEa_hh.size); std::vector n_EbId(shared.dim.n_EbId.size); + std::vector delta_Ea(shared.dim.delta_Ea.size); std::vector delta_Eb(shared.dim.delta_Eb.size); + std::vector new_Ea(shared.dim.Ea.size); std::vector new_Eb(shared.dim.Eb.size); - std::vector n_SEa(shared.dim.n_SEa.size); + std::vector n_SEa_s(shared.dim.n_SEa_s.size); std::vector n_EbIr(shared.dim.n_EbIr.size); std::vector delta_Id(shared.dim.delta_Id.size); - std::vector new_S(shared.dim.S.size); std::vector new_Id(shared.dim.Id.size); - std::vector n_SEa_hh(shared.dim.n_SEa_hh.size); - std::vector delta_Ea(shared.dim.delta_Ea.size); + std::vector new_E(shared.dim.E.size); + std::vector n_SEa_hc(shared.dim.n_SEa_hc.size); std::vector delta_Ir(shared.dim.delta_Ir.size); - std::vector new_Ea(shared.dim.Ea.size); std::vector new_Ir(shared.dim.Ir.size); - std::vector n_SEa_s(shared.dim.n_SEa_s.size); - std::vector new_E(shared.dim.E.size); + std::vector n_SEa_z(shared.dim.n_SEa_z.size); std::vector new_I(shared.dim.I.size); std::vector new_N(shared.dim.N.size); - std::vector n_SEa_hc(shared.dim.n_SEa_hc.size); - std::vector n_SEa_z(shared.dim.n_SEa_z.size); - return internal_state{n_IrR, n_IdD, target_met_children_t, target_met_adults_t, coverage_target_1st_dose_children, coverage_target_1st_dose_adults, coverage_target_2nd_dose_adults, n_eligible_for_dose1_children, n_eligible_for_dose1_adults, n_eligible_for_dose2_adults, I_infectious, delta_R, delta_D, daily_doses_children_t, daily_doses_adults_t, n_vaccination_t_S_children, n_vaccination_t_S_adults, n_vaccination_t_Ea_children, n_vaccination_t_Ea_adults, n_vaccination_t_Eb_children, n_vaccination_t_Eb_adults, n_vaccination_t_R_children, n_vaccination_t_R_adults, new_D, prop_infectious, lambda_hc, n_vaccination_t_S, n_vaccination_t_Ea, n_vaccination_t_Eb, n_vaccination_t_R, s_ij_gen_pop, s_ij_sex, delta_S_n_vaccination, delta_Ea_n_vaccination, delta_Eb_n_vaccination, delta_R_n_vaccination, n_vaccination_t, lambda_hh, lambda_s, new_R, lambda, n_EaEb, n_EbI, p_SE, p_hh, p_s, p_hc, n_EbId, delta_Eb, new_Eb, n_SEa, n_EbIr, delta_Id, new_S, new_Id, n_SEa_hh, delta_Ea, delta_Ir, new_Ea, new_Ir, n_SEa_s, new_E, new_I, new_N, n_SEa_hc, n_SEa_z}; + return internal_state{n_IrR, n_IdD, target_met_children_t, target_met_adults_t, coverage_target_1st_dose_children, coverage_target_1st_dose_adults, coverage_target_2nd_dose_adults, I_infectious, delta_R, delta_D, daily_doses_children_t, daily_doses_adults_t, give_dose1_children, give_dose1_adults, give_dose2_adults, new_D, prop_infectious, lambda_hc, children_dose1_denom, adults_dose1_denom, adults_dose2_denom, s_ij_gen_pop, s_ij_sex, n_vaccination_t_S_children, n_vaccination_t_S_adults, n_vaccination_t_Ea_children, n_vaccination_t_Ea_adults, n_vaccination_t_Eb_children, n_vaccination_t_Eb_adults, n_vaccination_t_R_children, n_vaccination_t_R_adults, lambda_hh, lambda_s, n_vaccination_t_S, n_vaccination_t_Ea, n_vaccination_t_Eb, n_vaccination_t_R, lambda, delta_S_n_vaccination, delta_Ea_n_vaccination, delta_Eb_n_vaccination, delta_R_n_vaccination, n_vaccination_t, p_SE, p_hh, p_s, p_hc, new_R, n_SEa, n_EaEb, n_EbI, new_S, n_SEa_hh, n_EbId, delta_Ea, delta_Eb, new_Ea, new_Eb, n_SEa_s, n_EbIr, delta_Id, new_Id, new_E, n_SEa_hc, delta_Ir, new_Ir, n_SEa_z, new_I, new_N}; } static data_type build_data(cpp11::list data, const shared_state& shared) { auto cases = dust2::r::read_real(data, "cases", NA_REAL); @@ -972,15 +984,6 @@ class model_targeted_vax { for (size_t i = 1; i <= shared.dim.coverage_target_2nd_dose_adults.size; ++i) { internal.coverage_target_2nd_dose_adults[i - 1] = monty::math::ceil(shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]); } - for (size_t i = 1; i <= shared.dim.n_eligible_for_dose1_children.size; ++i) { - internal.n_eligible_for_dose1_children[i - 1] = (S[i - 1 + shared.dim.S.mult[1]] + Ea[i - 1 + shared.dim.Ea.mult[1]] + Eb[i - 1 + shared.dim.Eb.mult[1]] + R[i - 1 + shared.dim.R.mult[1]]) * shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]]; - } - for (size_t i = 1; i <= shared.dim.n_eligible_for_dose1_adults.size; ++i) { - internal.n_eligible_for_dose1_adults[i - 1] = (S[i - 1 + shared.dim.S.mult[1]] + Ea[i - 1 + shared.dim.Ea.mult[1]] + Eb[i - 1 + shared.dim.Eb.mult[1]] + R[i - 1 + shared.dim.R.mult[1]]) * shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_1st_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]; - } - for (size_t i = 1; i <= shared.dim.n_eligible_for_dose2_adults.size; ++i) { - internal.n_eligible_for_dose2_adults[i - 1] = (S[i - 1 + 2 * shared.dim.S.mult[1]] + Ea[i - 1 + 2 * shared.dim.Ea.mult[1]] + Eb[i - 1 + 2 * shared.dim.Eb.mult[1]] + R[i - 1 + 2 * shared.dim.R.mult[1]]) * shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]; - } const real_type new_deaths_00_04 = dust2::array::sum(internal.n_IdD.data(), shared.dim.n_IdD, {0, 0}, {0, shared.dim.n_IdD.dim[1] - 1}); const real_type new_deaths_SW_12_14 = monty::random::binomial(rng_state, dust2::array::sum(internal.n_IdD.data(), shared.dim.n_IdD, {16, 16}, {0, shared.dim.n_IdD.dim[1] - 1}), static_cast(0.5)); const real_type new_deaths_CSW = dust2::array::sum(internal.n_IdD.data(), shared.dim.n_IdD, {16, 16}, {0, shared.dim.n_IdD.dim[1] - 1}); @@ -1007,68 +1010,107 @@ class model_targeted_vax { const real_type prioritisation_step_1st_dose_children_proposal = (dust2::array::sum(internal.target_met_children_t.data(), shared.dim.target_met_children_t) == dust2::array::sum(internal.coverage_target_1st_dose_children.data(), shared.dim.coverage_target_1st_dose_children) ? prioritisation_step_1st_dose_children + 1 : prioritisation_step_1st_dose_children); const real_type prioritisation_step_1st_dose_adults_proposal = (dust2::array::sum(internal.target_met_adults_t.data(), shared.dim.target_met_adults_t, {0, shared.dim.target_met_adults_t.dim[0] - 1}, {2, 2}) == dust2::array::sum(internal.coverage_target_1st_dose_adults.data(), shared.dim.coverage_target_1st_dose_adults) ? prioritisation_step_1st_dose_adults + 1 : prioritisation_step_1st_dose_adults); const real_type prioritisation_step_2nd_dose_adults_proposal = (dust2::array::sum(internal.target_met_adults_t.data(), shared.dim.target_met_adults_t, {0, shared.dim.target_met_adults_t.dim[0] - 1}, {3, 3}) == dust2::array::sum(internal.coverage_target_2nd_dose_adults.data(), shared.dim.coverage_target_2nd_dose_adults) ? prioritisation_step_2nd_dose_adults + 1 : prioritisation_step_2nd_dose_adults); + for (size_t i = 1; i <= shared.dim.give_dose1_children.size; ++i) { + internal.give_dose1_children[i - 1] = (shared.is_child[i - 1] * monty::math::ceil(shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]]) * (1 - internal.target_met_children_t[i - 1])); + } + for (size_t i = 1; i <= shared.dim.give_dose1_adults.size; ++i) { + internal.give_dose1_adults[i - 1] = ((1 - shared.is_child[i - 1]) * monty::math::ceil(shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_1st_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) * (1 - internal.target_met_adults_t[i - 1 + 2 * shared.dim.target_met_adults_t.mult[1]])); + } + for (size_t i = 1; i <= shared.dim.give_dose2_adults.size; ++i) { + internal.give_dose2_adults[i - 1] = ((1 - shared.is_child[i - 1]) * monty::math::ceil(shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) * (1 - internal.target_met_adults_t[i - 1 + 3 * shared.dim.target_met_adults_t.mult[1]])); + } + for (size_t i = 1; i <= shared.dim.D.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.D.dim[1]; ++j) { + internal.new_D[i - 1 + (j - 1) * shared.dim.D.mult[1]] = D[i - 1 + (j - 1) * shared.dim.D.mult[1]] + internal.delta_D[i - 1 + (j - 1) * shared.dim.delta_D.mult[1]]; + } + } + const real_type new_deaths_SW_15_17 = dust2::array::sum(internal.n_IdD.data(), shared.dim.n_IdD, {16, 16}, {0, shared.dim.n_IdD.dim[1] - 1}) - new_deaths_SW_12_14; + const real_type new_deaths_05_14 = dust2::array::sum(internal.n_IdD.data(), shared.dim.n_IdD, {1, 2}, {0, shared.dim.n_IdD.dim[1] - 1}) + new_deaths_SW_12_14; + const real_type new_deaths_SW = new_deaths_CSW + new_deaths_ASW; + for (size_t i = 1; i <= shared.dim.prop_infectious.size; ++i) { + internal.prop_infectious[i - 1] = (dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1}) == 0 ? 0 : dust2::array::sum(internal.I_infectious.data(), shared.dim.I_infectious, {i - 1, i - 1}, {0, shared.dim.I_infectious.dim[1] - 1}) / dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1})); + } + for (size_t i = 1; i <= shared.dim.lambda_hc.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.lambda_hc.dim[1]; ++j) { + internal.lambda_hc[i - 1 + (j - 1) * shared.dim.lambda_hc.mult[1]] = (i == 20 ? shared.beta_hcw * dust2::array::sum(internal.I_infectious.data(), shared.dim.I_infectious) / dust2::array::sum(N, shared.dim.N) * (1 - shared.ve_I[i - 1 + (j - 1) * shared.dim.ve_I.mult[1]]) : 0); + } + } + for (size_t i = 1; i <= shared.dim.children_dose1_denom.size; ++i) { + internal.children_dose1_denom[i - 1] = (S[i - 1 + shared.dim.S.mult[1]] + Ea[i - 1 + shared.dim.Ea.mult[1]] + Eb[i - 1 + shared.dim.Eb.mult[1]] + R[i - 1 + shared.dim.R.mult[1]]) * internal.give_dose1_children[i - 1]; + } + for (size_t i = 1; i <= shared.dim.adults_dose1_denom.size; ++i) { + internal.adults_dose1_denom[i - 1] = (S[i - 1 + shared.dim.S.mult[1]] + Ea[i - 1 + shared.dim.Ea.mult[1]] + Eb[i - 1 + shared.dim.Eb.mult[1]] + R[i - 1 + shared.dim.R.mult[1]]) * internal.give_dose1_adults[i - 1]; + } + for (size_t i = 1; i <= shared.dim.adults_dose2_denom.size; ++i) { + internal.adults_dose2_denom[i - 1] = (S[i - 1 + 2 * shared.dim.S.mult[1]] + Ea[i - 1 + 2 * shared.dim.Ea.mult[1]] + Eb[i - 1 + 2 * shared.dim.Eb.mult[1]] + R[i - 1 + 2 * shared.dim.R.mult[1]]) * internal.give_dose2_adults[i - 1]; + } + const real_type new_deaths_15_plus = dust2::array::sum(internal.n_IdD.data(), shared.dim.n_IdD, {3, 15}, {0, shared.dim.n_IdD.dim[1] - 1}) + new_deaths_SW_15_17 + dust2::array::sum(internal.n_IdD.data(), shared.dim.n_IdD, {17, 19}, {0, shared.dim.n_IdD.dim[1] - 1}); + for (size_t i = 1; i <= shared.dim.s_ij_gen_pop.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.s_ij_gen_pop.dim[1]; ++j) { + internal.s_ij_gen_pop[i - 1 + (j - 1) * shared.dim.s_ij_gen_pop.mult[1]] = shared.m_gen_pop[i - 1 + (j - 1) * shared.dim.m_gen_pop.mult[1]] * internal.prop_infectious[j - 1]; + } + } + for (size_t i = 1; i <= shared.dim.s_ij_sex.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.s_ij_sex.dim[1]; ++j) { + internal.s_ij_sex[i - 1 + (j - 1) * shared.dim.s_ij_sex.mult[1]] = shared.m_sex[i - 1 + (j - 1) * shared.dim.m_sex.mult[1]] * internal.prop_infectious[j - 1]; + } + } for (size_t i = 1; i <= shared.dim.n_vaccination_t_S_children.size; ++i) { internal.n_vaccination_t_S_children[i - 1] = 0; } for (size_t i = 1; i <= shared.dim.n_vaccination_t_S_children.size; ++i) { - internal.n_vaccination_t_S_children[i - 1] = (dust2::array::sum(internal.n_eligible_for_dose1_children.data(), shared.dim.n_eligible_for_dose1_children) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_children_t[1] * S[i - 1 + shared.dim.S.mult[1]] * shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose1_children.data(), shared.dim.n_eligible_for_dose1_children)), S[i - 1 + shared.dim.S.mult[1]])); + internal.n_vaccination_t_S_children[i - 1] = (dust2::array::sum(internal.children_dose1_denom.data(), shared.dim.children_dose1_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_children_t[1] * S[i - 1 + shared.dim.S.mult[1]]) / dust2::array::sum(internal.children_dose1_denom.data(), shared.dim.children_dose1_denom)) * internal.give_dose1_children[i - 1]), S[i - 1 + shared.dim.S.mult[1]])); } for (size_t i = 1; i <= shared.dim.n_vaccination_t_S_adults.size; ++i) { internal.n_vaccination_t_S_adults[i - 1] = 0; } for (size_t i = 1; i <= shared.dim.n_vaccination_t_S_adults.size; ++i) { - internal.n_vaccination_t_S_adults[i - 1] = (dust2::array::sum(internal.n_eligible_for_dose1_adults.data(), shared.dim.n_eligible_for_dose1_adults) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_adults_t[1] * S[i - 1 + shared.dim.S.mult[1]] * shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_1st_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose1_adults.data(), shared.dim.n_eligible_for_dose1_adults)), S[i - 1 + shared.dim.S.mult[1]])); + internal.n_vaccination_t_S_adults[i - 1] = (dust2::array::sum(internal.adults_dose1_denom.data(), shared.dim.adults_dose1_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_adults_t[1] * S[i - 1 + shared.dim.S.mult[1]]) / dust2::array::sum(internal.adults_dose1_denom.data(), shared.dim.adults_dose1_denom)) * internal.give_dose1_adults[i - 1]), S[i - 1 + shared.dim.S.mult[1]])); } for (size_t i = 1; i <= shared.dim.n_vaccination_t_Ea_children.size; ++i) { internal.n_vaccination_t_Ea_children[i - 1] = 0; } for (size_t i = 1; i <= shared.dim.n_vaccination_t_Ea_children.size; ++i) { - internal.n_vaccination_t_Ea_children[i - 1] = (dust2::array::sum(internal.n_eligible_for_dose1_children.data(), shared.dim.n_eligible_for_dose1_children) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_children_t[1] * Ea[i - 1 + shared.dim.Ea.mult[1]] * shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose1_children.data(), shared.dim.n_eligible_for_dose1_children)), Ea[i - 1 + shared.dim.Ea.mult[1]])); + internal.n_vaccination_t_Ea_children[i - 1] = (dust2::array::sum(internal.children_dose1_denom.data(), shared.dim.children_dose1_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_children_t[1] * Ea[i - 1 + shared.dim.Ea.mult[1]]) / dust2::array::sum(internal.children_dose1_denom.data(), shared.dim.children_dose1_denom)) * internal.give_dose1_children[i - 1]), Ea[i - 1 + shared.dim.Ea.mult[1]])); } for (size_t i = 1; i <= shared.dim.n_vaccination_t_Ea_adults.size; ++i) { internal.n_vaccination_t_Ea_adults[i - 1] = 0; } for (size_t i = 1; i <= shared.dim.n_vaccination_t_Ea_adults.size; ++i) { - internal.n_vaccination_t_Ea_adults[i - 1] = (dust2::array::sum(internal.n_eligible_for_dose1_adults.data(), shared.dim.n_eligible_for_dose1_adults) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_adults_t[1] * Ea[i - 1 + shared.dim.Ea.mult[1]] * shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_1st_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose1_adults.data(), shared.dim.n_eligible_for_dose1_adults)), Ea[i - 1 + shared.dim.Ea.mult[1]])); + internal.n_vaccination_t_Ea_adults[i - 1] = (dust2::array::sum(internal.adults_dose1_denom.data(), shared.dim.adults_dose1_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_adults_t[1] * Ea[i - 1 + shared.dim.Ea.mult[1]]) / dust2::array::sum(internal.adults_dose1_denom.data(), shared.dim.adults_dose1_denom)) * internal.give_dose1_adults[i - 1]), Ea[i - 1 + shared.dim.Ea.mult[1]])); } for (size_t i = 1; i <= shared.dim.n_vaccination_t_Eb_children.size; ++i) { internal.n_vaccination_t_Eb_children[i - 1] = 0; } for (size_t i = 1; i <= shared.dim.n_vaccination_t_Eb_children.size; ++i) { - internal.n_vaccination_t_Eb_children[i - 1] = (dust2::array::sum(internal.n_eligible_for_dose1_children.data(), shared.dim.n_eligible_for_dose1_children) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_children_t[1] * Eb[i - 1 + shared.dim.Eb.mult[1]] * shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose1_children.data(), shared.dim.n_eligible_for_dose1_children)), Eb[i - 1 + shared.dim.Eb.mult[1]])); + internal.n_vaccination_t_Eb_children[i - 1] = (dust2::array::sum(internal.children_dose1_denom.data(), shared.dim.children_dose1_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_children_t[1] * Eb[i - 1 + shared.dim.Eb.mult[1]]) / dust2::array::sum(internal.children_dose1_denom.data(), shared.dim.children_dose1_denom)) * internal.give_dose1_children[i - 1]), Eb[i - 1 + shared.dim.Eb.mult[1]])); } for (size_t i = 1; i <= shared.dim.n_vaccination_t_Eb_adults.size; ++i) { internal.n_vaccination_t_Eb_adults[i - 1] = 0; } for (size_t i = 1; i <= shared.dim.n_vaccination_t_Eb_adults.size; ++i) { - internal.n_vaccination_t_Eb_adults[i - 1] = (dust2::array::sum(internal.n_eligible_for_dose1_adults.data(), shared.dim.n_eligible_for_dose1_adults) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_adults_t[1] * Eb[i - 1 + shared.dim.Eb.mult[1]] * shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_1st_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose1_adults.data(), shared.dim.n_eligible_for_dose1_adults)), Eb[i - 1 + shared.dim.Eb.mult[1]])); + internal.n_vaccination_t_Eb_adults[i - 1] = (dust2::array::sum(internal.adults_dose1_denom.data(), shared.dim.adults_dose1_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_adults_t[1] * Eb[i - 1 + shared.dim.Eb.mult[1]]) / dust2::array::sum(internal.adults_dose1_denom.data(), shared.dim.adults_dose1_denom)) * internal.give_dose1_adults[i - 1]), Eb[i - 1 + shared.dim.Eb.mult[1]])); } for (size_t i = 1; i <= shared.dim.n_vaccination_t_R_children.size; ++i) { internal.n_vaccination_t_R_children[i - 1] = 0; } for (size_t i = 1; i <= shared.dim.n_vaccination_t_R_children.size; ++i) { - internal.n_vaccination_t_R_children[i - 1] = (dust2::array::sum(internal.n_eligible_for_dose1_children.data(), shared.dim.n_eligible_for_dose1_children) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_children_t[1] * R[i - 1 + shared.dim.R.mult[1]] * shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose1_children.data(), shared.dim.n_eligible_for_dose1_children)), R[i - 1 + shared.dim.R.mult[1]])); + internal.n_vaccination_t_R_children[i - 1] = (dust2::array::sum(internal.children_dose1_denom.data(), shared.dim.children_dose1_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_children_t[1] * R[i - 1 + shared.dim.R.mult[1]]) / dust2::array::sum(internal.children_dose1_denom.data(), shared.dim.children_dose1_denom)) * internal.give_dose1_children[i - 1]), R[i - 1 + shared.dim.R.mult[1]])); } for (size_t i = 1; i <= shared.dim.n_vaccination_t_R_adults.size; ++i) { internal.n_vaccination_t_R_adults[i - 1] = 0; } for (size_t i = 1; i <= shared.dim.n_vaccination_t_R_adults.size; ++i) { - internal.n_vaccination_t_R_adults[i - 1] = (dust2::array::sum(internal.n_eligible_for_dose1_adults.data(), shared.dim.n_eligible_for_dose1_adults) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_adults_t[1] * R[i - 1 + shared.dim.R.mult[1]] * shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_1st_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose1_adults.data(), shared.dim.n_eligible_for_dose1_adults)), R[i - 1 + shared.dim.R.mult[1]])); + internal.n_vaccination_t_R_adults[i - 1] = (dust2::array::sum(internal.adults_dose1_denom.data(), shared.dim.adults_dose1_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_adults_t[1] * R[i - 1 + shared.dim.R.mult[1]]) / dust2::array::sum(internal.adults_dose1_denom.data(), shared.dim.adults_dose1_denom)) * internal.give_dose1_adults[i - 1]), R[i - 1 + shared.dim.R.mult[1]])); } - for (size_t i = 1; i <= shared.dim.D.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.D.dim[1]; ++j) { - internal.new_D[i - 1 + (j - 1) * shared.dim.D.mult[1]] = D[i - 1 + (j - 1) * shared.dim.D.mult[1]] + internal.delta_D[i - 1 + (j - 1) * shared.dim.delta_D.mult[1]]; + for (size_t i = 1; i <= shared.dim.lambda_hh.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.lambda_hh.dim[1]; ++j) { + internal.lambda_hh[i - 1 + (j - 1) * shared.dim.lambda_hh.mult[1]] = shared.beta_h * dust2::array::sum(internal.s_ij_gen_pop.data(), shared.dim.s_ij_gen_pop, {i - 1, i - 1}, {0, shared.dim.s_ij_gen_pop.dim[1] - 1}) * (1 - shared.ve_I[i - 1 + (j - 1) * shared.dim.ve_I.mult[1]]); } } - const real_type new_deaths_SW_15_17 = dust2::array::sum(internal.n_IdD.data(), shared.dim.n_IdD, {16, 16}, {0, shared.dim.n_IdD.dim[1] - 1}) - new_deaths_SW_12_14; - const real_type new_deaths_05_14 = dust2::array::sum(internal.n_IdD.data(), shared.dim.n_IdD, {1, 2}, {0, shared.dim.n_IdD.dim[1] - 1}) + new_deaths_SW_12_14; - const real_type new_deaths_SW = new_deaths_CSW + new_deaths_ASW; - for (size_t i = 1; i <= shared.dim.prop_infectious.size; ++i) { - internal.prop_infectious[i - 1] = (dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1}) == 0 ? 0 : dust2::array::sum(internal.I_infectious.data(), shared.dim.I_infectious, {i - 1, i - 1}, {0, shared.dim.I_infectious.dim[1] - 1}) / dust2::array::sum(N, shared.dim.N, {i - 1, i - 1}, {0, shared.dim.N.dim[1] - 1})); - } - for (size_t i = 1; i <= shared.dim.lambda_hc.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.lambda_hc.dim[1]; ++j) { - internal.lambda_hc[i - 1 + (j - 1) * shared.dim.lambda_hc.mult[1]] = (i == 20 ? shared.beta_hcw * dust2::array::sum(internal.I_infectious.data(), shared.dim.I_infectious) / dust2::array::sum(N, shared.dim.N) * (1 - shared.ve_I[i - 1 + (j - 1) * shared.dim.ve_I.mult[1]]) : 0); + for (size_t i = 1; i <= shared.dim.lambda_s.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.lambda_s.dim[1]; ++j) { + internal.lambda_s[i - 1 + (j - 1) * shared.dim.lambda_s.mult[1]] = shared.beta_s * dust2::array::sum(internal.s_ij_sex.data(), shared.dim.s_ij_sex, {i - 1, i - 1}, {0, shared.dim.s_ij_sex.dim[1] - 1}) * (1 - shared.ve_I[i - 1 + (j - 1) * shared.dim.ve_I.mult[1]]); } } for (size_t i = 1; i <= shared.dim.n_vaccination_t_S.dim[0]; ++i) { @@ -1079,9 +1121,8 @@ class model_targeted_vax { for (size_t i = 1; i <= shared.dim.n_vaccination_t_S.dim[0]; ++i) { internal.n_vaccination_t_S[i - 1 + shared.dim.n_vaccination_t_S.mult[1]] = internal.n_vaccination_t_S_children[i - 1] + internal.n_vaccination_t_S_adults[i - 1]; } - internal.n_vaccination_t_S[2 + shared.dim.n_vaccination_t_S.mult[1]] = monty::math::min(internal.n_vaccination_t_S[2 + shared.dim.n_vaccination_t_S.mult[1]], S[2 + shared.dim.S.mult[1]]); for (size_t i = 1; i <= shared.dim.n_vaccination_t_S.dim[0]; ++i) { - internal.n_vaccination_t_S[i - 1 + 2 * shared.dim.n_vaccination_t_S.mult[1]] = (dust2::array::sum(internal.n_eligible_for_dose2_adults.data(), shared.dim.n_eligible_for_dose2_adults) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_adults_t[2] * S[i - 1 + 2 * shared.dim.S.mult[1]] * shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose2_adults.data(), shared.dim.n_eligible_for_dose2_adults)), S[i - 1 + 2 * shared.dim.S.mult[1]])); + internal.n_vaccination_t_S[i - 1 + 2 * shared.dim.n_vaccination_t_S.mult[1]] = (dust2::array::sum(internal.adults_dose2_denom.data(), shared.dim.adults_dose2_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_adults_t[2] * S[i - 1 + 2 * shared.dim.S.mult[1]]) / dust2::array::sum(internal.adults_dose2_denom.data(), shared.dim.adults_dose2_denom)) * internal.give_dose2_adults[i - 1]), S[i - 1 + 2 * shared.dim.S.mult[1]])); } for (size_t i = 1; i <= shared.dim.n_vaccination_t_Ea.dim[0]; ++i) { for (size_t j = 1; j <= shared.dim.n_vaccination_t_Ea.dim[1]; ++j) { @@ -1091,9 +1132,8 @@ class model_targeted_vax { for (size_t i = 1; i <= shared.dim.n_vaccination_t_Ea.dim[0]; ++i) { internal.n_vaccination_t_Ea[i - 1 + shared.dim.n_vaccination_t_Ea.mult[1]] = internal.n_vaccination_t_Ea_children[i - 1] + internal.n_vaccination_t_Ea_adults[i - 1]; } - internal.n_vaccination_t_Ea[2 + shared.dim.n_vaccination_t_Ea.mult[1]] = monty::math::min(internal.n_vaccination_t_Ea[2 + shared.dim.n_vaccination_t_Ea.mult[1]], Ea[2 + shared.dim.Ea.mult[1]]); for (size_t i = 1; i <= shared.dim.n_vaccination_t_Ea.dim[0]; ++i) { - internal.n_vaccination_t_Ea[i - 1 + 2 * shared.dim.n_vaccination_t_Ea.mult[1]] = (dust2::array::sum(internal.n_eligible_for_dose2_adults.data(), shared.dim.n_eligible_for_dose2_adults) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_adults_t[2] * Ea[i - 1 + 2 * shared.dim.Ea.mult[1]] * shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose2_adults.data(), shared.dim.n_eligible_for_dose2_adults)), Ea[i - 1 + 2 * shared.dim.Ea.mult[1]])); + internal.n_vaccination_t_Ea[i - 1 + 2 * shared.dim.n_vaccination_t_Ea.mult[1]] = (dust2::array::sum(internal.adults_dose2_denom.data(), shared.dim.adults_dose2_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_adults_t[2] * Ea[i - 1 + 2 * shared.dim.Ea.mult[1]]) / dust2::array::sum(internal.adults_dose2_denom.data(), shared.dim.adults_dose2_denom)) * internal.give_dose2_adults[i - 1]), Ea[i - 1 + 2 * shared.dim.Ea.mult[1]])); } for (size_t i = 1; i <= shared.dim.n_vaccination_t_Eb.dim[0]; ++i) { for (size_t j = 1; j <= shared.dim.n_vaccination_t_Eb.dim[1]; ++j) { @@ -1103,9 +1143,8 @@ class model_targeted_vax { for (size_t i = 1; i <= shared.dim.n_vaccination_t_Eb.dim[0]; ++i) { internal.n_vaccination_t_Eb[i - 1 + shared.dim.n_vaccination_t_Eb.mult[1]] = internal.n_vaccination_t_Eb_children[i - 1] + internal.n_vaccination_t_Eb_adults[i - 1]; } - internal.n_vaccination_t_Eb[2 + shared.dim.n_vaccination_t_Eb.mult[1]] = monty::math::min(internal.n_vaccination_t_Eb[2 + shared.dim.n_vaccination_t_Eb.mult[1]], Eb[2 + shared.dim.Eb.mult[1]]); for (size_t i = 1; i <= shared.dim.n_vaccination_t_Eb.dim[0]; ++i) { - internal.n_vaccination_t_Eb[i - 1 + 2 * shared.dim.n_vaccination_t_Eb.mult[1]] = (dust2::array::sum(internal.n_eligible_for_dose2_adults.data(), shared.dim.n_eligible_for_dose2_adults) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_adults_t[2] * Eb[i - 1 + 2 * shared.dim.Eb.mult[1]] * shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose2_adults.data(), shared.dim.n_eligible_for_dose2_adults)), Eb[i - 1 + 2 * shared.dim.Eb.mult[1]])); + internal.n_vaccination_t_Eb[i - 1 + 2 * shared.dim.n_vaccination_t_Eb.mult[1]] = (dust2::array::sum(internal.adults_dose2_denom.data(), shared.dim.adults_dose2_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_adults_t[2] * Eb[i - 1 + 2 * shared.dim.Eb.mult[1]]) / dust2::array::sum(internal.adults_dose2_denom.data(), shared.dim.adults_dose2_denom)) * internal.give_dose2_adults[i - 1]), Eb[i - 1 + 2 * shared.dim.Eb.mult[1]])); } for (size_t i = 1; i <= shared.dim.n_vaccination_t_R.dim[0]; ++i) { for (size_t j = 1; j <= shared.dim.n_vaccination_t_R.dim[1]; ++j) { @@ -1115,19 +1154,12 @@ class model_targeted_vax { for (size_t i = 1; i <= shared.dim.n_vaccination_t_R.dim[0]; ++i) { internal.n_vaccination_t_R[i - 1 + shared.dim.n_vaccination_t_R.mult[1]] = internal.n_vaccination_t_R_children[i - 1] + internal.n_vaccination_t_R_adults[i - 1]; } - internal.n_vaccination_t_R[2 + shared.dim.n_vaccination_t_R.mult[1]] = monty::math::min(internal.n_vaccination_t_R[2 + shared.dim.n_vaccination_t_R.mult[1]], R[2 + shared.dim.R.mult[1]]); for (size_t i = 1; i <= shared.dim.n_vaccination_t_R.dim[0]; ++i) { - internal.n_vaccination_t_R[i - 1 + 2 * shared.dim.n_vaccination_t_R.mult[1]] = (dust2::array::sum(internal.n_eligible_for_dose2_adults.data(), shared.dim.n_eligible_for_dose2_adults) == 0 ? 0 : monty::math::min(monty::math::floor((internal.daily_doses_adults_t[2] * R[i - 1 + 2 * shared.dim.R.mult[1]] * shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) / dust2::array::sum(internal.n_eligible_for_dose2_adults.data(), shared.dim.n_eligible_for_dose2_adults)), R[i - 1 + 2 * shared.dim.R.mult[1]])); - } - const real_type new_deaths_15_plus = dust2::array::sum(internal.n_IdD.data(), shared.dim.n_IdD, {3, 15}, {0, shared.dim.n_IdD.dim[1] - 1}) + new_deaths_SW_15_17 + dust2::array::sum(internal.n_IdD.data(), shared.dim.n_IdD, {17, 19}, {0, shared.dim.n_IdD.dim[1] - 1}); - for (size_t i = 1; i <= shared.dim.s_ij_gen_pop.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.s_ij_gen_pop.dim[1]; ++j) { - internal.s_ij_gen_pop[i - 1 + (j - 1) * shared.dim.s_ij_gen_pop.mult[1]] = shared.m_gen_pop[i - 1 + (j - 1) * shared.dim.m_gen_pop.mult[1]] * internal.prop_infectious[j - 1]; - } + internal.n_vaccination_t_R[i - 1 + 2 * shared.dim.n_vaccination_t_R.mult[1]] = (dust2::array::sum(internal.adults_dose2_denom.data(), shared.dim.adults_dose2_denom) == 0 ? 0 : monty::math::min(monty::math::floor(((internal.daily_doses_adults_t[2] * R[i - 1 + 2 * shared.dim.R.mult[1]]) / dust2::array::sum(internal.adults_dose2_denom.data(), shared.dim.adults_dose2_denom)) * internal.give_dose2_adults[i - 1]), R[i - 1 + 2 * shared.dim.R.mult[1]])); } - for (size_t i = 1; i <= shared.dim.s_ij_sex.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.s_ij_sex.dim[1]; ++j) { - internal.s_ij_sex[i - 1 + (j - 1) * shared.dim.s_ij_sex.mult[1]] = shared.m_sex[i - 1 + (j - 1) * shared.dim.m_sex.mult[1]] * internal.prop_infectious[j - 1]; + for (size_t i = 1; i <= shared.dim.lambda.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.lambda.dim[1]; ++j) { + internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] = internal.lambda_hh[i - 1 + (j - 1) * shared.dim.lambda_hh.mult[1]] + internal.lambda_s[i - 1 + (j - 1) * shared.dim.lambda_s.mult[1]] + internal.lambda_hc[i - 1 + (j - 1) * shared.dim.lambda_hc.mult[1]] + shared.lambda_z[i - 1 + (j - 1) * shared.dim.lambda_z.mult[1]]; } } for (size_t i = 1; i <= shared.dim.delta_S_n_vaccination.dim[0]; ++i) { @@ -1155,14 +1187,24 @@ class model_targeted_vax { internal.n_vaccination_t[i - 1 + (j - 1) * shared.dim.n_vaccination_t.mult[1]] = internal.n_vaccination_t_S[i - 1 + (j - 1) * shared.dim.n_vaccination_t_S.mult[1]] + internal.n_vaccination_t_Ea[i - 1 + (j - 1) * shared.dim.n_vaccination_t_Ea.mult[1]] + internal.n_vaccination_t_Eb[i - 1 + (j - 1) * shared.dim.n_vaccination_t_Eb.mult[1]] + internal.n_vaccination_t_R[i - 1 + (j - 1) * shared.dim.n_vaccination_t_R.mult[1]]; } } - for (size_t i = 1; i <= shared.dim.lambda_hh.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.lambda_hh.dim[1]; ++j) { - internal.lambda_hh[i - 1 + (j - 1) * shared.dim.lambda_hh.mult[1]] = shared.beta_h * dust2::array::sum(internal.s_ij_gen_pop.data(), shared.dim.s_ij_gen_pop, {i - 1, i - 1}, {0, shared.dim.s_ij_gen_pop.dim[1] - 1}) * (1 - shared.ve_I[i - 1 + (j - 1) * shared.dim.ve_I.mult[1]]); + for (size_t i = 1; i <= shared.dim.p_SE.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.p_SE.dim[1]; ++j) { + internal.p_SE[i - 1 + (j - 1) * shared.dim.p_SE.mult[1]] = 1 - monty::math::exp(-internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] * dt); } } - for (size_t i = 1; i <= shared.dim.lambda_s.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.lambda_s.dim[1]; ++j) { - internal.lambda_s[i - 1 + (j - 1) * shared.dim.lambda_s.mult[1]] = shared.beta_s * dust2::array::sum(internal.s_ij_sex.data(), shared.dim.s_ij_sex, {i - 1, i - 1}, {0, shared.dim.s_ij_sex.dim[1] - 1}) * (1 - shared.ve_I[i - 1 + (j - 1) * shared.dim.ve_I.mult[1]]); + for (size_t i = 1; i <= shared.dim.p_hh.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.p_hh.dim[1]; ++j) { + internal.p_hh[i - 1 + (j - 1) * shared.dim.p_hh.mult[1]] = (internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] > 0 ? internal.lambda_hh[i - 1 + (j - 1) * shared.dim.lambda_hh.mult[1]] / internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] : 0); + } + } + for (size_t i = 1; i <= shared.dim.p_s.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.p_s.dim[1]; ++j) { + internal.p_s[i - 1 + (j - 1) * shared.dim.p_s.mult[1]] = (internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] > 0 ? internal.lambda_s[i - 1 + (j - 1) * shared.dim.lambda_s.mult[1]] / internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] : 0); + } + } + for (size_t i = 1; i <= shared.dim.p_hc.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.p_hc.dim[1]; ++j) { + internal.p_hc[i - 1 + (j - 1) * shared.dim.p_hc.mult[1]] = (internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] > 0 ? internal.lambda_hc[i - 1 + (j - 1) * shared.dim.lambda_hc.mult[1]] / internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] : 0); } } for (size_t i = 1; i <= shared.dim.R.dim[0]; ++i) { @@ -1184,9 +1226,9 @@ class model_targeted_vax { const real_type new_dose2_ASW = internal.n_vaccination_t[17 + 2 * shared.dim.n_vaccination_t.mult[1]]; const real_type new_dose2_PBS = internal.n_vaccination_t[18 + 2 * shared.dim.n_vaccination_t.mult[1]]; const real_type new_dose2_HCW = internal.n_vaccination_t[19 + 2 * shared.dim.n_vaccination_t.mult[1]]; - for (size_t i = 1; i <= shared.dim.lambda.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.lambda.dim[1]; ++j) { - internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] = internal.lambda_hh[i - 1 + (j - 1) * shared.dim.lambda_hh.mult[1]] + internal.lambda_s[i - 1 + (j - 1) * shared.dim.lambda_s.mult[1]] + internal.lambda_hc[i - 1 + (j - 1) * shared.dim.lambda_hc.mult[1]] + shared.lambda_z[i - 1 + (j - 1) * shared.dim.lambda_z.mult[1]]; + for (size_t i = 1; i <= shared.dim.n_SEa.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.n_SEa.dim[1]; ++j) { + internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]] = monty::random::binomial(rng_state, S[i - 1 + (j - 1) * shared.dim.S.mult[1]] + internal.delta_S_n_vaccination[i - 1 + (j - 1) * shared.dim.delta_S_n_vaccination.mult[1]], internal.p_SE[i - 1 + (j - 1) * shared.dim.p_SE.mult[1]]); } } for (size_t i = 1; i <= shared.dim.n_EaEb.dim[0]; ++i) { @@ -1199,30 +1241,26 @@ class model_targeted_vax { internal.n_EbI[i - 1 + (j - 1) * shared.dim.n_EbI.mult[1]] = monty::random::binomial(rng_state, Eb[i - 1 + (j - 1) * shared.dim.Eb.mult[1]] + internal.delta_Eb_n_vaccination[i - 1 + (j - 1) * shared.dim.delta_Eb_n_vaccination.mult[1]], p_EI); } } + for (size_t i = 1; i <= shared.dim.S.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.S.dim[1]; ++j) { + internal.new_S[i - 1 + (j - 1) * shared.dim.S.mult[1]] = S[i - 1 + (j - 1) * shared.dim.S.mult[1]] + internal.delta_S_n_vaccination[i - 1 + (j - 1) * shared.dim.delta_S_n_vaccination.mult[1]] - internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]]; + } + } + const real_type new_cases_00_04 = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {0, 0}, {0, shared.dim.n_SEa.dim[1] - 1}); + const real_type new_cases_SW_12_14 = monty::random::binomial(rng_state, dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {16, 16}, {0, shared.dim.n_SEa.dim[1] - 1}), static_cast(0.5)); + const real_type new_cases_CSW = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {16, 16}, {0, shared.dim.n_SEa.dim[1] - 1}); + const real_type new_cases_ASW = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {17, 17}, {0, shared.dim.n_SEa.dim[1] - 1}); + const real_type new_cases_PBS = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {18, 18}, {0, shared.dim.n_SEa.dim[1] - 1}); + const real_type new_cases_HCW = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {19, 19}, {0, shared.dim.n_SEa.dim[1] - 1}); const real_type new_dose1_SW_15_17 = internal.n_vaccination_t[16 + shared.dim.n_vaccination_t.mult[1]] - new_dose1_SW_12_14; const real_type new_dose1_05_14 = dust2::array::sum(internal.n_vaccination_t.data(), shared.dim.n_vaccination_t, {1, 2}, {1, 1}) + new_dose1_SW_12_14; const real_type new_dose1_SW = new_dose1_CSW + new_dose1_ASW; const real_type new_dose2_SW_15_17 = internal.n_vaccination_t[16 + 2 * shared.dim.n_vaccination_t.mult[1]] - new_dose2_SW_12_14; const real_type new_dose2_05_14 = dust2::array::sum(internal.n_vaccination_t.data(), shared.dim.n_vaccination_t, {1, 2}, {2, 2}) + new_dose2_SW_12_14; const real_type new_dose2_SW = new_dose2_CSW + new_dose2_ASW; - for (size_t i = 1; i <= shared.dim.p_SE.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.p_SE.dim[1]; ++j) { - internal.p_SE[i - 1 + (j - 1) * shared.dim.p_SE.mult[1]] = 1 - monty::math::exp(-internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] * dt); - } - } - for (size_t i = 1; i <= shared.dim.p_hh.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.p_hh.dim[1]; ++j) { - internal.p_hh[i - 1 + (j - 1) * shared.dim.p_hh.mult[1]] = (internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] > 0 ? internal.lambda_hh[i - 1 + (j - 1) * shared.dim.lambda_hh.mult[1]] / internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] : 0); - } - } - for (size_t i = 1; i <= shared.dim.p_s.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.p_s.dim[1]; ++j) { - internal.p_s[i - 1 + (j - 1) * shared.dim.p_s.mult[1]] = (internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] > 0 ? internal.lambda_s[i - 1 + (j - 1) * shared.dim.lambda_s.mult[1]] / internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] : 0); - } - } - for (size_t i = 1; i <= shared.dim.p_hc.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.p_hc.dim[1]; ++j) { - internal.p_hc[i - 1 + (j - 1) * shared.dim.p_hc.mult[1]] = (internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] > 0 ? internal.lambda_hc[i - 1 + (j - 1) * shared.dim.lambda_hc.mult[1]] / internal.lambda[i - 1 + (j - 1) * shared.dim.lambda.mult[1]] : 0); + for (size_t i = 1; i <= shared.dim.n_SEa_hh.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.n_SEa_hh.dim[1]; ++j) { + internal.n_SEa_hh[i - 1 + (j - 1) * shared.dim.n_SEa_hh.mult[1]] = monty::random::binomial(rng_state, internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]], internal.p_hh[i - 1 + (j - 1) * shared.dim.p_hh.mult[1]]); } } for (size_t i = 1; i <= shared.dim.n_EbId.dim[0]; ++i) { @@ -1230,21 +1268,34 @@ class model_targeted_vax { internal.n_EbId[i - 1 + (j - 1) * shared.dim.n_EbId.mult[1]] = monty::random::binomial(rng_state, internal.n_EbI[i - 1 + (j - 1) * shared.dim.n_EbI.mult[1]], shared.CFR[i - 1 + (j - 1) * shared.dim.CFR.mult[1]]); } } + for (size_t i = 1; i <= shared.dim.delta_Ea.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.delta_Ea.dim[1]; ++j) { + internal.delta_Ea[i - 1 + (j - 1) * shared.dim.delta_Ea.mult[1]] = internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]] - internal.n_EaEb[i - 1 + (j - 1) * shared.dim.n_EaEb.mult[1]]; + } + } for (size_t i = 1; i <= shared.dim.delta_Eb.dim[0]; ++i) { for (size_t j = 1; j <= shared.dim.delta_Eb.dim[1]; ++j) { internal.delta_Eb[i - 1 + (j - 1) * shared.dim.delta_Eb.mult[1]] = internal.n_EaEb[i - 1 + (j - 1) * shared.dim.n_EaEb.mult[1]] - internal.n_EbI[i - 1 + (j - 1) * shared.dim.n_EbI.mult[1]]; } } + for (size_t i = 1; i <= shared.dim.Ea.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.Ea.dim[1]; ++j) { + internal.new_Ea[i - 1 + (j - 1) * shared.dim.Ea.mult[1]] = Ea[i - 1 + (j - 1) * shared.dim.Ea.mult[1]] + internal.delta_Ea_n_vaccination[i - 1 + (j - 1) * shared.dim.delta_Ea_n_vaccination.mult[1]] + internal.delta_Ea[i - 1 + (j - 1) * shared.dim.delta_Ea.mult[1]]; + } + } for (size_t i = 1; i <= shared.dim.Eb.dim[0]; ++i) { for (size_t j = 1; j <= shared.dim.Eb.dim[1]; ++j) { internal.new_Eb[i - 1 + (j - 1) * shared.dim.Eb.mult[1]] = Eb[i - 1 + (j - 1) * shared.dim.Eb.mult[1]] + internal.delta_Eb_n_vaccination[i - 1 + (j - 1) * shared.dim.delta_Eb_n_vaccination.mult[1]] + internal.delta_Eb[i - 1 + (j - 1) * shared.dim.delta_Eb.mult[1]]; } } + const real_type new_cases_SW_15_17 = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {16, 16}, {0, shared.dim.n_SEa.dim[1] - 1}) - new_cases_SW_12_14; + const real_type new_cases_05_14 = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {1, 2}, {0, shared.dim.n_SEa.dim[1] - 1}) + new_cases_SW_12_14; + const real_type new_cases_SW = new_cases_CSW + new_cases_ASW; const real_type new_dose1_15_plus = dust2::array::sum(internal.n_vaccination_t.data(), shared.dim.n_vaccination_t, {3, 15}, {1, 1}) + new_dose1_SW_15_17 + dust2::array::sum(internal.n_vaccination_t.data(), shared.dim.n_vaccination_t, {17, 19}, {1, 1}); const real_type new_dose2_15_plus = dust2::array::sum(internal.n_vaccination_t.data(), shared.dim.n_vaccination_t, {3, 15}, {2, 2}) + new_dose2_SW_15_17 + dust2::array::sum(internal.n_vaccination_t.data(), shared.dim.n_vaccination_t, {17, 19}, {2, 2}); - for (size_t i = 1; i <= shared.dim.n_SEa.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.n_SEa.dim[1]; ++j) { - internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]] = monty::random::binomial(rng_state, S[i - 1 + (j - 1) * shared.dim.S.mult[1]] + internal.delta_S_n_vaccination[i - 1 + (j - 1) * shared.dim.delta_S_n_vaccination.mult[1]], internal.p_SE[i - 1 + (j - 1) * shared.dim.p_SE.mult[1]]); + for (size_t i = 1; i <= shared.dim.n_SEa_s.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.n_SEa_s.dim[1]; ++j) { + internal.n_SEa_s[i - 1 + (j - 1) * shared.dim.n_SEa_s.mult[1]] = monty::random::binomial(rng_state, internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]] - internal.n_SEa_hh[i - 1 + (j - 1) * shared.dim.n_SEa_hh.mult[1]], internal.p_s[i - 1 + (j - 1) * shared.dim.p_s.mult[1]]); } } for (size_t i = 1; i <= shared.dim.n_EbIr.dim[0]; ++i) { @@ -1257,30 +1308,20 @@ class model_targeted_vax { internal.delta_Id[i - 1 + (j - 1) * shared.dim.delta_Id.mult[1]] = internal.n_EbId[i - 1 + (j - 1) * shared.dim.n_EbId.mult[1]] - internal.n_IdD[i - 1 + (j - 1) * shared.dim.n_IdD.mult[1]]; } } - for (size_t i = 1; i <= shared.dim.S.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.S.dim[1]; ++j) { - internal.new_S[i - 1 + (j - 1) * shared.dim.S.mult[1]] = S[i - 1 + (j - 1) * shared.dim.S.mult[1]] + internal.delta_S_n_vaccination[i - 1 + (j - 1) * shared.dim.delta_S_n_vaccination.mult[1]] - internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]]; - } - } for (size_t i = 1; i <= shared.dim.Id.dim[0]; ++i) { for (size_t j = 1; j <= shared.dim.Id.dim[1]; ++j) { internal.new_Id[i - 1 + (j - 1) * shared.dim.Id.mult[1]] = Id[i - 1 + (j - 1) * shared.dim.Id.mult[1]] + internal.delta_Id[i - 1 + (j - 1) * shared.dim.delta_Id.mult[1]]; } } - const real_type new_cases_00_04 = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {0, 0}, {0, shared.dim.n_SEa.dim[1] - 1}); - const real_type new_cases_SW_12_14 = monty::random::binomial(rng_state, dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {16, 16}, {0, shared.dim.n_SEa.dim[1] - 1}), static_cast(0.5)); - const real_type new_cases_CSW = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {16, 16}, {0, shared.dim.n_SEa.dim[1] - 1}); - const real_type new_cases_ASW = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {17, 17}, {0, shared.dim.n_SEa.dim[1] - 1}); - const real_type new_cases_PBS = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {18, 18}, {0, shared.dim.n_SEa.dim[1] - 1}); - const real_type new_cases_HCW = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {19, 19}, {0, shared.dim.n_SEa.dim[1] - 1}); - for (size_t i = 1; i <= shared.dim.n_SEa_hh.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.n_SEa_hh.dim[1]; ++j) { - internal.n_SEa_hh[i - 1 + (j - 1) * shared.dim.n_SEa_hh.mult[1]] = monty::random::binomial(rng_state, internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]], internal.p_hh[i - 1 + (j - 1) * shared.dim.p_hh.mult[1]]); + for (size_t i = 1; i <= shared.dim.E.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.E.dim[1]; ++j) { + internal.new_E[i - 1 + (j - 1) * shared.dim.E.mult[1]] = internal.new_Ea[i - 1 + (j - 1) * shared.dim.Ea.mult[1]] + internal.new_Eb[i - 1 + (j - 1) * shared.dim.Eb.mult[1]]; } } - for (size_t i = 1; i <= shared.dim.delta_Ea.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.delta_Ea.dim[1]; ++j) { - internal.delta_Ea[i - 1 + (j - 1) * shared.dim.delta_Ea.mult[1]] = internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]] - internal.n_EaEb[i - 1 + (j - 1) * shared.dim.n_EaEb.mult[1]]; + const real_type new_cases_15_plus = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {3, 15}, {0, shared.dim.n_SEa.dim[1] - 1}) + new_cases_SW_15_17 + dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {17, 19}, {0, shared.dim.n_SEa.dim[1] - 1}); + for (size_t i = 1; i <= shared.dim.n_SEa_hc.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.n_SEa_hc.dim[1]; ++j) { + internal.n_SEa_hc[i - 1 + (j - 1) * shared.dim.n_SEa_hc.mult[1]] = monty::random::binomial(rng_state, internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]] - internal.n_SEa_hh[i - 1 + (j - 1) * shared.dim.n_SEa_hh.mult[1]] - internal.n_SEa_s[i - 1 + (j - 1) * shared.dim.n_SEa_s.mult[1]], internal.p_hc[i - 1 + (j - 1) * shared.dim.p_hc.mult[1]]); } } for (size_t i = 1; i <= shared.dim.delta_Ir.dim[0]; ++i) { @@ -1288,27 +1329,14 @@ class model_targeted_vax { internal.delta_Ir[i - 1 + (j - 1) * shared.dim.delta_Ir.mult[1]] = internal.n_EbIr[i - 1 + (j - 1) * shared.dim.n_EbIr.mult[1]] - internal.n_IrR[i - 1 + (j - 1) * shared.dim.n_IrR.mult[1]]; } } - for (size_t i = 1; i <= shared.dim.Ea.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.Ea.dim[1]; ++j) { - internal.new_Ea[i - 1 + (j - 1) * shared.dim.Ea.mult[1]] = Ea[i - 1 + (j - 1) * shared.dim.Ea.mult[1]] + internal.delta_Ea_n_vaccination[i - 1 + (j - 1) * shared.dim.delta_Ea_n_vaccination.mult[1]] + internal.delta_Ea[i - 1 + (j - 1) * shared.dim.delta_Ea.mult[1]]; - } - } for (size_t i = 1; i <= shared.dim.Ir.dim[0]; ++i) { for (size_t j = 1; j <= shared.dim.Ir.dim[1]; ++j) { internal.new_Ir[i - 1 + (j - 1) * shared.dim.Ir.mult[1]] = Ir[i - 1 + (j - 1) * shared.dim.Ir.mult[1]] + internal.delta_Ir[i - 1 + (j - 1) * shared.dim.delta_Ir.mult[1]]; } } - const real_type new_cases_SW_15_17 = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {16, 16}, {0, shared.dim.n_SEa.dim[1] - 1}) - new_cases_SW_12_14; - const real_type new_cases_05_14 = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {1, 2}, {0, shared.dim.n_SEa.dim[1] - 1}) + new_cases_SW_12_14; - const real_type new_cases_SW = new_cases_CSW + new_cases_ASW; - for (size_t i = 1; i <= shared.dim.n_SEa_s.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.n_SEa_s.dim[1]; ++j) { - internal.n_SEa_s[i - 1 + (j - 1) * shared.dim.n_SEa_s.mult[1]] = monty::random::binomial(rng_state, internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]] - internal.n_SEa_hh[i - 1 + (j - 1) * shared.dim.n_SEa_hh.mult[1]], internal.p_s[i - 1 + (j - 1) * shared.dim.p_s.mult[1]]); - } - } - for (size_t i = 1; i <= shared.dim.E.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.E.dim[1]; ++j) { - internal.new_E[i - 1 + (j - 1) * shared.dim.E.mult[1]] = internal.new_Ea[i - 1 + (j - 1) * shared.dim.Ea.mult[1]] + internal.new_Eb[i - 1 + (j - 1) * shared.dim.Eb.mult[1]]; + for (size_t i = 1; i <= shared.dim.n_SEa_z.dim[0]; ++i) { + for (size_t j = 1; j <= shared.dim.n_SEa_z.dim[1]; ++j) { + internal.n_SEa_z[i - 1 + (j - 1) * shared.dim.n_SEa_z.mult[1]] = internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]] - internal.n_SEa_hh[i - 1 + (j - 1) * shared.dim.n_SEa_hh.mult[1]] - internal.n_SEa_s[i - 1 + (j - 1) * shared.dim.n_SEa_s.mult[1]] - internal.n_SEa_hc[i - 1 + (j - 1) * shared.dim.n_SEa_hc.mult[1]]; } } for (size_t i = 1; i <= shared.dim.I.dim[0]; ++i) { @@ -1321,17 +1349,6 @@ class model_targeted_vax { internal.new_N[i - 1 + (j - 1) * shared.dim.N.mult[1]] = internal.new_S[i - 1 + (j - 1) * shared.dim.S.mult[1]] + internal.new_Ea[i - 1 + (j - 1) * shared.dim.Ea.mult[1]] + internal.new_Eb[i - 1 + (j - 1) * shared.dim.Eb.mult[1]] + internal.new_Ir[i - 1 + (j - 1) * shared.dim.Ir.mult[1]] + internal.new_Id[i - 1 + (j - 1) * shared.dim.Id.mult[1]] + internal.new_R[i - 1 + (j - 1) * shared.dim.R.mult[1]] + internal.new_D[i - 1 + (j - 1) * shared.dim.D.mult[1]]; } } - const real_type new_cases_15_plus = dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {3, 15}, {0, shared.dim.n_SEa.dim[1] - 1}) + new_cases_SW_15_17 + dust2::array::sum(internal.n_SEa.data(), shared.dim.n_SEa, {17, 19}, {0, shared.dim.n_SEa.dim[1] - 1}); - for (size_t i = 1; i <= shared.dim.n_SEa_hc.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.n_SEa_hc.dim[1]; ++j) { - internal.n_SEa_hc[i - 1 + (j - 1) * shared.dim.n_SEa_hc.mult[1]] = monty::random::binomial(rng_state, internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]] - internal.n_SEa_hh[i - 1 + (j - 1) * shared.dim.n_SEa_hh.mult[1]] - internal.n_SEa_s[i - 1 + (j - 1) * shared.dim.n_SEa_s.mult[1]], internal.p_hc[i - 1 + (j - 1) * shared.dim.p_hc.mult[1]]); - } - } - for (size_t i = 1; i <= shared.dim.n_SEa_z.dim[0]; ++i) { - for (size_t j = 1; j <= shared.dim.n_SEa_z.dim[1]; ++j) { - internal.n_SEa_z[i - 1 + (j - 1) * shared.dim.n_SEa_z.mult[1]] = internal.n_SEa[i - 1 + (j - 1) * shared.dim.n_SEa.mult[1]] - internal.n_SEa_hh[i - 1 + (j - 1) * shared.dim.n_SEa_hh.mult[1]] - internal.n_SEa_s[i - 1 + (j - 1) * shared.dim.n_SEa_s.mult[1]] - internal.n_SEa_hc[i - 1 + (j - 1) * shared.dim.n_SEa_hc.mult[1]]; - } - } state_next[0] = (prioritisation_step_1st_dose_children_proposal > shared.N_prioritisation_steps_children ? shared.N_prioritisation_steps_children : prioritisation_step_1st_dose_children_proposal); state_next[1] = (prioritisation_step_1st_dose_adults_proposal > shared.N_prioritisation_steps_adults ? shared.N_prioritisation_steps_adults : prioritisation_step_1st_dose_adults_proposal); state_next[2] = (prioritisation_step_2nd_dose_adults_proposal > shared.N_prioritisation_steps_adults ? shared.N_prioritisation_steps_adults : prioritisation_step_2nd_dose_adults_proposal); From 1365a61e7a40388d895772326fd478e7ac907d37 Mon Sep 17 00:00:00 2001 From: ruthmccabe Date: Fri, 10 Jan 2025 22:02:47 +0000 Subject: [PATCH 12/26] test name updated (not relevant to this PR just found it as checking tests) --- tests/testthat/test-model-targeted-vax.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/testthat/test-model-targeted-vax.R b/tests/testthat/test-model-targeted-vax.R index aa1ad90..6c56480 100644 --- a/tests/testthat/test-model-targeted-vax.R +++ b/tests/testthat/test-model-targeted-vax.R @@ -109,7 +109,7 @@ test_that("relevant states sum correctly", { }) -test_that("when beta_h = beta_z = beta_s = 0 there are no new infections", { +test_that("when beta_h = beta_z = beta_s = beta_hcw = 0 there are no new infections", { pars <- reference_pars_targeted_vax() pars$beta_h <- 0 pars$beta_s <- 0 From 6bc1881744bfa5920eb80bc49e2f14aee8787036 Mon Sep 17 00:00:00 2001 From: ruthmccabe Date: Fri, 10 Jan 2025 22:10:46 +0000 Subject: [PATCH 13/26] as before just tidied test name (not PR specific) --- tests/testthat/test-model-targeted-vax.R | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/tests/testthat/test-model-targeted-vax.R b/tests/testthat/test-model-targeted-vax.R index 6c56480..bf2efa8 100644 --- a/tests/testthat/test-model-targeted-vax.R +++ b/tests/testthat/test-model-targeted-vax.R @@ -193,7 +193,7 @@ test_that("when CFR = 1 everybody dies", { }) -test_that("when beta_h = 0 and beta_s = 0 there are only zoonotic infections", { +test_that("when beta_h = 0, beta_s = 0, beta_hcw = 0 there are only zoonotic infections", { pars <- reference_pars_targeted_vax() pars$beta_h <- 0 pars$beta_s <- 0 @@ -219,7 +219,7 @@ test_that("when beta_h = 0 and beta_s = 0 there are only zoonotic infections", { }) -test_that("when beta_h = 0 and beta_z = 0 infections only from sexual contact", { +test_that("when beta_h = 0, beta_z = 0, beta_hcw = 0 infections only from sexual contact", { pars <- reference_pars_targeted_vax() pars$beta_h <- 0 pars$beta_z[] <- 0 From 73c0ebde82cc6d4407587ffe7eb4b5286a0de4ff Mon Sep 17 00:00:00 2001 From: ruthmccabe Date: Fri, 10 Jan 2025 22:21:12 +0000 Subject: [PATCH 14/26] tidy model notes --- inst/odin/model-targeted-vax.R | 49 +--------------------------------- 1 file changed, 1 insertion(+), 48 deletions(-) diff --git a/inst/odin/model-targeted-vax.R b/inst/odin/model-targeted-vax.R index dc71e14..49d0380 100644 --- a/inst/odin/model-targeted-vax.R +++ b/inst/odin/model-targeted-vax.R @@ -54,8 +54,6 @@ dim(prioritisation_strategy_adults) <- c(n_group, N_prioritisation_steps_adults) ## vaccination targets: for each group, what is the level of coverage we want to ## see before we move onto the next set of groups to target? -## when we input this do the calc outside by the population size as part of -## pre-processing to save faffing around in here. ## here, j=1,...,n_vax corresponds to the coverage wanted in j, so for j=0 ## and j=1 this stays at 0 @@ -167,33 +165,6 @@ initial(prioritisation_step_2nd_dose_adults) <- 1 ## now that we know the step we are on, we know who is eligible per age group to ## be vaccinated -#### Notes from meeting 8 Jan -### if target_met == 1 then vax is 0; -## if n_eligb == 0 then 0; -## else (daily_doses * S / n_elig) * ceiling(prioiritisation) (remove prioritisation) -## give_dose_X X = {dose_1_children,dose_1_adult,dose_2_adult} length n_group -## is target met or elgiible or prioiritsation step - single indicator - -## old eligibility -# ## isolate the relevant classes -# n_eligible_for_dose1_children[] <- (S[i, 2] + Ea[i, 2] + Eb[i, 2] + R[i, 2]) -# dim(n_eligible_for_dose1_children) <- c(n_group) -# ## who can get a 1st dose -# ## now we have to look in the preceding j class (e.g. who in unvaccinated j = 1 -# ## can get a 1st dose and move to j = 2) as these are the people who will be -# ## eligible -# n_eligible_for_dose1_adults[] <- (S[i, 2] + Ea[i, 2] + Eb[i, 2] + R[i, 2]) -# dim(n_eligible_for_dose1_adults) <- c(n_group) -# ## who can get a 2nd dose -# n_eligible_for_dose2_adults[] <- (S[i, 3] + Ea[i, 3] + Eb[i, 3] + R[i, 3]) -# dim(n_eligible_for_dose2_adults) <- c(n_group) - -# ## isolate the relevant compartments who can be vaccinated (relevant groups will be highlighted further down) -# comps_eligible_dose1[] <- (S[i, 2] + Ea[i, 2] + Eb[i, 2] + R[i, 2]) -# dim(comps_eligible_dose1) <- c(n_group) -# comps_eligible_dose2[] <- (S[i, 3] + Ea[i, 3] + Eb[i, 3] + R[i, 3]) -# dim(comps_eligible_dose2) <- c(n_group) - ## create indicator for which groups are allowed to be vaccinated based on: ## 1) if they are an adult or a child (via is_child) ## 2) if they are in a prioritised group, dependent on being adult or child (prioritisation_strategy_X) @@ -235,9 +206,6 @@ n_vaccination_t_S_children[] <- sum(children_dose1_denom[])) * give_dose1_children[i]), S[i, 2]) -##### current concern about denominator here - I think all of the compartments will be summed over before the ones that are able to be vaccinated have been taken account of, so we won't be giving out all of the compartments - - n_vaccination_t_S_adults[] <- 0 ## adults 1st doses @@ -252,10 +220,7 @@ n_vaccination_t_S[, ] <- 0 n_vaccination_t_S[, 2] <- n_vaccination_t_S_children[i] + n_vaccination_t_S_adults[i] -## this should no longer be necessary -# ## for the boundary case do an extra check that we haven't gone over the number -# ## of people in each compartment -# n_vaccination_t_S[3, 2] <- min(n_vaccination_t_S[3, 2], S[3, 2]) +## no longer doing boundary checks as we won't have any groups which are a mix of children and adults ## allocate 2nd doses (adults only for now) n_vaccination_t_S[, 3] <- @@ -287,10 +252,6 @@ n_vaccination_t_Ea[, ] <- 0 n_vaccination_t_Ea[, 2] <- n_vaccination_t_Ea_children[i] + n_vaccination_t_Ea_adults[i] -# ## for the boundary case do an extra check that we haven't gone over the number -# ## of people in each compartment -# n_vaccination_t_Ea[3, 2] <- min(n_vaccination_t_Ea[3, 2], Ea[3, 2]) - ## adults 2nd doses n_vaccination_t_Ea[, 3] <- if (sum(adults_dose2_denom[]) == 0) 0 else @@ -322,10 +283,6 @@ n_vaccination_t_Eb[, ] <- 0 n_vaccination_t_Eb[, 2] <- n_vaccination_t_Eb_children[i] + n_vaccination_t_Eb_adults[i] -# ## for the boundary case do an extra check that we haven't gone over the number -# ## of people in each compartment -# n_vaccination_t_Eb[3, 2] <- min(n_vaccination_t_Eb[3, 2], Eb[3, 2]) - ## adults 2nd doses n_vaccination_t_Eb[, 3] <- if (sum(adults_dose2_denom[]) == 0) 0 else @@ -357,10 +314,6 @@ n_vaccination_t_R[, ] <- 0 n_vaccination_t_R[, 2] <- n_vaccination_t_R_children[i] + n_vaccination_t_R_adults[i] -# ## for the boundary case do an extra check that we haven't gone over the number -# ## of people in each compartment -# n_vaccination_t_R[3, 2] <- min(n_vaccination_t_R[3, 2], R[3, 2]) - ## adults 2nd doses n_vaccination_t_R[, 3] <- if (sum(adults_dose2_denom[]) == 0) 0 else From 0e2a76b24b7ef50bf02aaf340cfd7e5f87bcc024 Mon Sep 17 00:00:00 2001 From: ruthmccabe Date: Fri, 10 Jan 2025 22:22:00 +0000 Subject: [PATCH 15/26] bump version --- DESCRIPTION | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/DESCRIPTION b/DESCRIPTION index df8e5cc..bc451e9 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,6 +1,6 @@ Package: mpoxseir Title: Stochastic compartmental model of mpox transmission -Version: 0.2.6 +Version: 0.2.9 Authors@R: c(person("Lilith", "Whittles", role = c("aut", "cre"), email = "l.whittles@imperial.ac.uk"), person("Ruth", "McCabe", role = c("aut")), From 1d87379c393ca868f39228b8601c61818a4bc98b Mon Sep 17 00:00:00 2001 From: ruthmccabe Date: Fri, 10 Jan 2025 22:23:38 +0000 Subject: [PATCH 16/26] update news --- NEWS.md | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/NEWS.md b/NEWS.md index 746dfe2..016677b 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,3 +1,8 @@ +# mpoxseir 0.2.9 + +* children_ind_raw replaced with is_child, adults_ind_raw replaced with (1 - ind_child) +* vaccination allocation determined based on 3 criteria of eligibility, prioritisation and target met + # mpoxseir 0.2.6 From 104e29876f562c18abb8af7e1b09d092ac56307b Mon Sep 17 00:00:00 2001 From: ruthmccabe Date: Fri, 10 Jan 2025 22:45:46 +0000 Subject: [PATCH 17/26] tidying up code from trying to tidy FOI --- inst/odin/model-targeted-vax.R | 16 +++------------- 1 file changed, 3 insertions(+), 13 deletions(-) diff --git a/inst/odin/model-targeted-vax.R b/inst/odin/model-targeted-vax.R index 49d0380..cb343b6 100644 --- a/inst/odin/model-targeted-vax.R +++ b/inst/odin/model-targeted-vax.R @@ -616,22 +616,12 @@ lambda_hh[, ] <- beta_h * sum(s_ij_gen_pop[i, ]) * (1 - ve_I[i, j]) lambda_s[, ] <- beta_s * sum(s_ij_sex[i, ]) * (1 - ve_I[i, j]) # additional foi in HCW only (i = 20) homogeneous from infected as assumed equally # likely to attend hospital -lambda_hc[, ] <- +lambda_hc[, ] <- if (i == 20) beta_hcw * sum(I_infectious) / sum(N) * (1 - ve_I[i, j]) else 0 lambda_z[, ] <- beta_z[i] * (1 - ve_I[i, j]) -lambda[, ] <- lambda_hh[i, j] + lambda_s[i, j] + lambda_hc[i, j] + lambda_z[i, j] - -#### did not like this -# lambda_hh[, ] <- beta_h * sum(s_ij_gen_pop[i, ]) -# lambda_s[, ] <- beta_s * sum(s_ij_sex[i, ]) -# # additional foi in HCW only (i = 20) homogeneous from infected as assumed equally -# # likely to attend hospital -# lambda_hc[, ] <- -# if (i == 20) beta_hcw * sum(I_infectious[, ]) else 0 -# lambda_z[, ] <- beta_z[i] -# -# lambda[, ] <- (lambda_hh[i, j] + lambda_s[i, j] + lambda_hc[i, j] + lambda_z[i, j]) * (1 - ve_I[i, j]) +lambda[, ] <- lambda_hh[i, j] + lambda_s[i, j] + lambda_hc[i, j] + lambda_z[i, j] + ## Draws from binomial distributions for numbers changing between compartments # accounting for vaccination: From 180f99b7916c40513fdb8e60320c9e243d686ca7 Mon Sep 17 00:00:00 2001 From: ruthmccabe Date: Fri, 10 Jan 2025 22:51:08 +0000 Subject: [PATCH 18/26] re-update description --- DESCRIPTION | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/DESCRIPTION b/DESCRIPTION index c053a23..d5fc276 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,6 +1,6 @@ Package: mpoxseir Title: Stochastic compartmental model of mpox transmission -Version: 0.2.7 +Version: 0.2.9 Authors@R: c(person("Lilith", "Whittles", role = c("aut", "cre"), email = "l.whittles@imperial.ac.uk"), From 788b21628315d1ce55caa1e7ddc9b7a06a612749 Mon Sep 17 00:00:00 2001 From: ruthmccabe Date: Fri, 10 Jan 2025 22:55:29 +0000 Subject: [PATCH 19/26] remove accidental extra line --- DESCRIPTION | 1 - 1 file changed, 1 deletion(-) diff --git a/DESCRIPTION b/DESCRIPTION index d5fc276..bc451e9 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,7 +1,6 @@ Package: mpoxseir Title: Stochastic compartmental model of mpox transmission Version: 0.2.9 - Authors@R: c(person("Lilith", "Whittles", role = c("aut", "cre"), email = "l.whittles@imperial.ac.uk"), person("Ruth", "McCabe", role = c("aut")), From 43a869ac0f1daa8e4da9f1e62f5450e530231961 Mon Sep 17 00:00:00 2001 From: ruthmccabe Date: Mon, 13 Jan 2025 14:28:22 +0000 Subject: [PATCH 20/26] update packages & regenerate --- R/dust.R | 2 +- src/model-targeted-vax.cpp | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/R/dust.R b/R/dust.R index b9f5c13..069b6f1 100644 --- a/R/dust.R +++ b/R/dust.R @@ -1,4 +1,4 @@ -## Generated by dust2 (version 0.3.11) - do not edit +## Generated by dust2 (version 0.3.15) - do not edit model_targeted_vax <- structure( function() get("model_targeted_vax"), class = "dust_system_generator", diff --git a/src/model-targeted-vax.cpp b/src/model-targeted-vax.cpp index 30ff1cc..fa0a9b6 100644 --- a/src/model-targeted-vax.cpp +++ b/src/model-targeted-vax.cpp @@ -1,4 +1,4 @@ -// Generated by dust2 (version 0.3.11) - do not edit +// Generated by dust2 (version 0.3.15) - do not edit // Generated by odin2 (version 0.3.17) - do not edit #include From 65c840a769f685a2f8d1228c9724624cddf0dad8 Mon Sep 17 00:00:00 2001 From: ruthmccabe <54892632+ruthmccabe@users.noreply.github.com> Date: Wed, 15 Jan 2025 09:59:32 +0000 Subject: [PATCH 21/26] Update NEWS.md Co-authored-by: edknock <47318334+edknock@users.noreply.github.com> --- NEWS.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/NEWS.md b/NEWS.md index d268596..0c798ee 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,6 +1,6 @@ # mpoxseir 0.2.9 -* children_ind_raw replaced with is_child, adults_ind_raw replaced with (1 - ind_child) +* children_ind_raw replaced with is_child, adults_ind_raw replaced with (1 - is_child) * vaccination allocation determined based on 3 criteria of eligibility, prioritisation and target met # mpoxseir 0.2.8 From 3433d3e78053e3726f90481dce03f29b1de6010a Mon Sep 17 00:00:00 2001 From: ruthmccabe <54892632+ruthmccabe@users.noreply.github.com> Date: Wed, 15 Jan 2025 10:03:12 +0000 Subject: [PATCH 22/26] Tidy code as per Lilith's suggestion Co-authored-by: Lilith Whittles --- inst/odin/model-targeted-vax.R | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/inst/odin/model-targeted-vax.R b/inst/odin/model-targeted-vax.R index 76defdf..3467c6d 100644 --- a/inst/odin/model-targeted-vax.R +++ b/inst/odin/model-targeted-vax.R @@ -171,8 +171,7 @@ initial(prioritisation_step_2nd_dose_adults) <- 1 dim(give_dose1_children) <- c(n_group) give_dose1_children[] <- - (is_child[i] * ceiling(prioritisation_strategy_children[ - i, prioritisation_step_1st_dose_children]) * (1 - target_met_children_t[i])) + is_child[i] * coverage_target_1st_dose_children[i] * (1 - target_met_children_t[i])) dim(give_dose1_adults) <- c(n_group) give_dose1_adults[] <- From 85c47b9992d69595689776332f467f1c9e60555b Mon Sep 17 00:00:00 2001 From: ruthmccabe <54892632+ruthmccabe@users.noreply.github.com> Date: Wed, 15 Jan 2025 10:03:22 +0000 Subject: [PATCH 23/26] Tidy code as per Lilith's suggestion Co-authored-by: Lilith Whittles --- inst/odin/model-targeted-vax.R | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/inst/odin/model-targeted-vax.R b/inst/odin/model-targeted-vax.R index 3467c6d..c00ab26 100644 --- a/inst/odin/model-targeted-vax.R +++ b/inst/odin/model-targeted-vax.R @@ -175,8 +175,7 @@ give_dose1_children[] <- dim(give_dose1_adults) <- c(n_group) give_dose1_adults[] <- - ((1 - is_child[i]) * ceiling(prioritisation_strategy_adults[ - i, prioritisation_step_1st_dose_adults]) * (1 - target_met_adults_t[i, 3])) + (1 - is_child[i]) * coverage_target_1st_dose_adults[i] * (1 - target_met_adults_t[i, 3]) dim(give_dose2_adults) <- c(n_group) give_dose2_adults[] <- ((1 - is_child[i]) * ceiling(prioritisation_strategy_adults[ From 10976e783b689ec078474ee67c6e112f5520d5de Mon Sep 17 00:00:00 2001 From: ruthmccabe <54892632+ruthmccabe@users.noreply.github.com> Date: Wed, 15 Jan 2025 10:03:31 +0000 Subject: [PATCH 24/26] Tidy code as per Lilith's suggestion Co-authored-by: Lilith Whittles --- inst/odin/model-targeted-vax.R | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/inst/odin/model-targeted-vax.R b/inst/odin/model-targeted-vax.R index c00ab26..2712024 100644 --- a/inst/odin/model-targeted-vax.R +++ b/inst/odin/model-targeted-vax.R @@ -178,8 +178,7 @@ give_dose1_adults[] <- (1 - is_child[i]) * coverage_target_1st_dose_adults[i] * (1 - target_met_adults_t[i, 3]) dim(give_dose2_adults) <- c(n_group) give_dose2_adults[] <- - ((1 - is_child[i]) * ceiling(prioritisation_strategy_adults[ - i, prioritisation_step_2nd_dose_adults]) * (1 - target_met_adults_t[i, 4])) + (1 - is_child[i]) * coverage_target_2nd_dose_adults[i] * (1 - target_met_adults_t[i, 4]) ## allocate the doses to the unvaccinated by age group, prioritisation strategy From 61d8f2d792cbab8c3d0531b128425eb7e53ac6a4 Mon Sep 17 00:00:00 2001 From: ruthmccabe Date: Wed, 15 Jan 2025 10:24:08 +0000 Subject: [PATCH 25/26] extra bracket removed --- inst/odin/model-targeted-vax.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/inst/odin/model-targeted-vax.R b/inst/odin/model-targeted-vax.R index 2712024..2026b6c 100644 --- a/inst/odin/model-targeted-vax.R +++ b/inst/odin/model-targeted-vax.R @@ -171,7 +171,7 @@ initial(prioritisation_step_2nd_dose_adults) <- 1 dim(give_dose1_children) <- c(n_group) give_dose1_children[] <- - is_child[i] * coverage_target_1st_dose_children[i] * (1 - target_met_children_t[i])) + is_child[i] * coverage_target_1st_dose_children[i] * (1 - target_met_children_t[i]) dim(give_dose1_adults) <- c(n_group) give_dose1_adults[] <- From 2e55d3370b62088f82d8dd4910cbf0fcb6429de7 Mon Sep 17 00:00:00 2001 From: ruthmccabe Date: Wed, 15 Jan 2025 10:26:08 +0000 Subject: [PATCH 26/26] regen model post review comments --- inst/dust/model-targeted-vax.cpp | 6 +++--- src/model-targeted-vax.cpp | 6 +++--- 2 files changed, 6 insertions(+), 6 deletions(-) diff --git a/inst/dust/model-targeted-vax.cpp b/inst/dust/model-targeted-vax.cpp index 919c97e..b6ea736 100644 --- a/inst/dust/model-targeted-vax.cpp +++ b/inst/dust/model-targeted-vax.cpp @@ -1009,13 +1009,13 @@ class model_targeted_vax { const real_type prioritisation_step_1st_dose_adults_proposal = (dust2::array::sum(internal.target_met_adults_t.data(), shared.dim.target_met_adults_t, {0, shared.dim.target_met_adults_t.dim[0] - 1}, {2, 2}) == dust2::array::sum(internal.coverage_target_1st_dose_adults.data(), shared.dim.coverage_target_1st_dose_adults) ? prioritisation_step_1st_dose_adults + 1 : prioritisation_step_1st_dose_adults); const real_type prioritisation_step_2nd_dose_adults_proposal = (dust2::array::sum(internal.target_met_adults_t.data(), shared.dim.target_met_adults_t, {0, shared.dim.target_met_adults_t.dim[0] - 1}, {3, 3}) == dust2::array::sum(internal.coverage_target_2nd_dose_adults.data(), shared.dim.coverage_target_2nd_dose_adults) ? prioritisation_step_2nd_dose_adults + 1 : prioritisation_step_2nd_dose_adults); for (size_t i = 1; i <= shared.dim.give_dose1_children.size; ++i) { - internal.give_dose1_children[i - 1] = (shared.is_child[i - 1] * monty::math::ceil(shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]]) * (1 - internal.target_met_children_t[i - 1])); + internal.give_dose1_children[i - 1] = shared.is_child[i - 1] * internal.coverage_target_1st_dose_children[i - 1] * (1 - internal.target_met_children_t[i - 1]); } for (size_t i = 1; i <= shared.dim.give_dose1_adults.size; ++i) { - internal.give_dose1_adults[i - 1] = ((1 - shared.is_child[i - 1]) * monty::math::ceil(shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_1st_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) * (1 - internal.target_met_adults_t[i - 1 + 2 * shared.dim.target_met_adults_t.mult[1]])); + internal.give_dose1_adults[i - 1] = (1 - shared.is_child[i - 1]) * internal.coverage_target_1st_dose_adults[i - 1] * (1 - internal.target_met_adults_t[i - 1 + 2 * shared.dim.target_met_adults_t.mult[1]]); } for (size_t i = 1; i <= shared.dim.give_dose2_adults.size; ++i) { - internal.give_dose2_adults[i - 1] = ((1 - shared.is_child[i - 1]) * monty::math::ceil(shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) * (1 - internal.target_met_adults_t[i - 1 + 3 * shared.dim.target_met_adults_t.mult[1]])); + internal.give_dose2_adults[i - 1] = (1 - shared.is_child[i - 1]) * internal.coverage_target_2nd_dose_adults[i - 1] * (1 - internal.target_met_adults_t[i - 1 + 3 * shared.dim.target_met_adults_t.mult[1]]); } for (size_t i = 1; i <= shared.dim.D.dim[0]; ++i) { for (size_t j = 1; j <= shared.dim.D.dim[1]; ++j) { diff --git a/src/model-targeted-vax.cpp b/src/model-targeted-vax.cpp index fa0a9b6..02ac908 100644 --- a/src/model-targeted-vax.cpp +++ b/src/model-targeted-vax.cpp @@ -1011,13 +1011,13 @@ class model_targeted_vax { const real_type prioritisation_step_1st_dose_adults_proposal = (dust2::array::sum(internal.target_met_adults_t.data(), shared.dim.target_met_adults_t, {0, shared.dim.target_met_adults_t.dim[0] - 1}, {2, 2}) == dust2::array::sum(internal.coverage_target_1st_dose_adults.data(), shared.dim.coverage_target_1st_dose_adults) ? prioritisation_step_1st_dose_adults + 1 : prioritisation_step_1st_dose_adults); const real_type prioritisation_step_2nd_dose_adults_proposal = (dust2::array::sum(internal.target_met_adults_t.data(), shared.dim.target_met_adults_t, {0, shared.dim.target_met_adults_t.dim[0] - 1}, {3, 3}) == dust2::array::sum(internal.coverage_target_2nd_dose_adults.data(), shared.dim.coverage_target_2nd_dose_adults) ? prioritisation_step_2nd_dose_adults + 1 : prioritisation_step_2nd_dose_adults); for (size_t i = 1; i <= shared.dim.give_dose1_children.size; ++i) { - internal.give_dose1_children[i - 1] = (shared.is_child[i - 1] * monty::math::ceil(shared.prioritisation_strategy_children[i - 1 + (prioritisation_step_1st_dose_children - 1) * shared.dim.prioritisation_strategy_children.mult[1]]) * (1 - internal.target_met_children_t[i - 1])); + internal.give_dose1_children[i - 1] = shared.is_child[i - 1] * internal.coverage_target_1st_dose_children[i - 1] * (1 - internal.target_met_children_t[i - 1]); } for (size_t i = 1; i <= shared.dim.give_dose1_adults.size; ++i) { - internal.give_dose1_adults[i - 1] = ((1 - shared.is_child[i - 1]) * monty::math::ceil(shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_1st_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) * (1 - internal.target_met_adults_t[i - 1 + 2 * shared.dim.target_met_adults_t.mult[1]])); + internal.give_dose1_adults[i - 1] = (1 - shared.is_child[i - 1]) * internal.coverage_target_1st_dose_adults[i - 1] * (1 - internal.target_met_adults_t[i - 1 + 2 * shared.dim.target_met_adults_t.mult[1]]); } for (size_t i = 1; i <= shared.dim.give_dose2_adults.size; ++i) { - internal.give_dose2_adults[i - 1] = ((1 - shared.is_child[i - 1]) * monty::math::ceil(shared.prioritisation_strategy_adults[i - 1 + (prioritisation_step_2nd_dose_adults - 1) * shared.dim.prioritisation_strategy_adults.mult[1]]) * (1 - internal.target_met_adults_t[i - 1 + 3 * shared.dim.target_met_adults_t.mult[1]])); + internal.give_dose2_adults[i - 1] = (1 - shared.is_child[i - 1]) * internal.coverage_target_2nd_dose_adults[i - 1] * (1 - internal.target_met_adults_t[i - 1 + 3 * shared.dim.target_met_adults_t.mult[1]]); } for (size_t i = 1; i <= shared.dim.D.dim[0]; ++i) { for (size_t j = 1; j <= shared.dim.D.dim[1]; ++j) {