From 713d8847e7316a5a455b23d60e967b6783950a5f Mon Sep 17 00:00:00 2001 From: jeffeaton Date: Mon, 18 Nov 2024 05:40:16 -0500 Subject: [PATCH 1/3] add ANC testing outputs to T4 --- DESCRIPTION | 2 +- NEWS.md | 4 ++++ R/model.R | 4 ++++ R/outputs.R | 14 ++++++++++++++ R/tmb-model.R | 3 +++ src/tmb.cpp | 44 ++++++++++++++++++++++++++++++++++++++++++++ 6 files changed, 70 insertions(+), 1 deletion(-) diff --git a/DESCRIPTION b/DESCRIPTION index 4614b402..27626c74 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,6 +1,6 @@ Package: naomi Title: Naomi Model for Subnational HIV Estimates -Version: 2.10.1 +Version: 2.10.2 Authors@R: person(given = "Jeff", family = "Eaton", diff --git a/NEWS.md b/NEWS.md index 3e00f431..d57dacd5 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,3 +1,7 @@ +# naomi 2.10.2 + +* Add ANC testing outputs to T4 projection period for including in PEPFAR datapack output. + # naomi 2.10.1 * Patch error in reading `anc_already_art` from Spectrum PJNZ file (was errantly diff --git a/R/model.R b/R/model.R index fe55d6ef..a21d8f6e 100644 --- a/R/model.R +++ b/R/model.R @@ -567,6 +567,10 @@ naomi_model_frame <- function(area_merged, spec_prev_t4 = prevalence, spec_incid_t4 = incidence, spec_artcov_t4 = art_coverage, + spec_unaware_untreated_prop_t4 = unaware_untreated_prop, + asfr_t4 = asfr, + frr_plhiv_t4 = frr_plhiv, + frr_already_art_t4 = frr_already_art ), by = c("spectrum_region_code", "sex", "age_group") ) %>% diff --git a/R/outputs.R b/R/outputs.R index 40485f1d..8e04f485 100644 --- a/R/outputs.R +++ b/R/outputs.R @@ -184,6 +184,16 @@ extract_indicators <- function(naomi_fit, naomi_mf, na.rm = FALSE) { "anc_rho_t3_out" = "anc_prevalence", "anc_alpha_t3_out" = "anc_art_coverage") + indicators_anc_t4 <- c("anc_clients_t4_out" = "anc_clients", + "anc_plhiv_t4_out" = "anc_plhiv", + "anc_already_art_t4_out" = "anc_already_art", + "anc_art_new_t4_out" = "anc_art_new", + "anc_known_pos_t4_out" = "anc_known_pos", + "anc_tested_pos_t4_out" = "anc_tested_pos", + "anc_tested_neg_t4_out" = "anc_tested_neg", + "anc_rho_t4_out" = "anc_prevalence", + "anc_alpha_t4_out" = "anc_art_coverage") + indicator_anc_est_t1 <- Map(get_est, names(indicators_anc_t1), indicators_anc_t1, naomi_mf$calendar_quarter1, list(naomi_mf$mf_anc_out)) @@ -191,11 +201,14 @@ extract_indicators <- function(naomi_fit, naomi_mf, na.rm = FALSE) { naomi_mf$calendar_quarter2, list(naomi_mf$mf_anc_out)) indicator_anc_est_t3 <- Map(get_est, names(indicators_anc_t3), indicators_anc_t3, naomi_mf$calendar_quarter3, list(naomi_mf$mf_anc_out)) + indicator_anc_est_t4 <- Map(get_est, names(indicators_anc_t4), indicators_anc_t4, + naomi_mf$calendar_quarter4, list(naomi_mf$mf_anc_out)) indicator_anc_est_t1 <- dplyr::bind_rows(indicator_anc_est_t1) indicator_anc_est_t2 <- dplyr::bind_rows(indicator_anc_est_t2) indicator_anc_est_t3 <- dplyr::bind_rows(indicator_anc_est_t3) + indicator_anc_est_t4 <- dplyr::bind_rows(indicator_anc_est_t4) mf_anc_out <- naomi_mf$mf_areas %>% dplyr::transmute(area_id, @@ -210,6 +223,7 @@ extract_indicators <- function(naomi_fit, naomi_mf, na.rm = FALSE) { indicator_est_t3, indicator_anc_est_t3, indicator_est_t4, + indicator_anc_est_t4, indicator_est_t5 ) diff --git a/R/tmb-model.R b/R/tmb-model.R index 4b74c09c..d1b2bb83 100644 --- a/R/tmb-model.R +++ b/R/tmb-model.R @@ -294,12 +294,15 @@ prepare_tmb_inputs <- function(naomi_data, log_asfr_t1_offset = log(df$asfr_t1), log_asfr_t2_offset = log(df$asfr_t2), log_asfr_t3_offset = log(df$asfr_t3), + log_asfr_t4_offset = log(df$asfr_t4), logit_anc_rho_t1_offset = log(df$frr_plhiv_t1), logit_anc_rho_t2_offset = log(df$frr_plhiv_t2), logit_anc_rho_t3_offset = log(df$frr_plhiv_t3), + logit_anc_rho_t4_offset = log(df$frr_plhiv_t4), logit_anc_alpha_t1_offset = log(df$frr_already_art_t1), logit_anc_alpha_t2_offset = log(df$frr_already_art_t2), logit_anc_alpha_t3_offset = log(df$frr_already_art_t3), + logit_anc_alpha_t4_offset = log(df$frr_already_art_t4), ## logit_rho_offset = naomi_data$mf_model$logit_rho_offset * naomi_data$mf_model$bin_rho_model, logit_alpha_offset = naomi_data$mf_model$logit_alpha_offset, diff --git a/src/tmb.cpp b/src/tmb.cpp index 8fb3173b..9e8b79c3 100644 --- a/src/tmb.cpp +++ b/src/tmb.cpp @@ -104,14 +104,17 @@ Type objective_function::operator() () DATA_VECTOR(log_asfr_t1_offset); DATA_VECTOR(log_asfr_t2_offset); DATA_VECTOR(log_asfr_t3_offset); + DATA_VECTOR(log_asfr_t4_offset); DATA_VECTOR(logit_anc_rho_t1_offset); DATA_VECTOR(logit_anc_rho_t2_offset); DATA_VECTOR(logit_anc_rho_t3_offset); + DATA_VECTOR(logit_anc_rho_t4_offset); DATA_VECTOR(logit_anc_alpha_t1_offset); DATA_VECTOR(logit_anc_alpha_t2_offset); DATA_VECTOR(logit_anc_alpha_t3_offset); + DATA_VECTOR(logit_anc_alpha_t4_offset); DATA_SPARSE_MATRIX(Z_asfr_x); DATA_SPARSE_MATRIX(Z_ancrho_x); @@ -1019,6 +1022,24 @@ Type objective_function::operator() () vector infections_t4(lambda_t4 * (population_t4 - plhiv_t4)); + // Note: currently assuming same district effects parameters from t2 for t4 + vector mu_anc_rho_t4(logit(rho_t4) + + logit_anc_rho_t2_offset + + X_ancrho * vector(beta_anc_rho + beta_anc_rho_t2) + + Z_ancrho_x * vector(ui_anc_rho_x + ui_anc_rho_xt)); + vector anc_rho_t4(invlogit(mu_anc_rho_t4)); + + vector mu_anc_alpha_t4(mu_alpha_t4 + + logit_anc_alpha_t4_offset + + X_ancalpha * vector(beta_anc_alpha + beta_anc_alpha_t2) + + Z_ancalpha_x * vector(ui_anc_alpha_x + ui_anc_alpha_xt)); + vector anc_alpha_t4(invlogit(mu_anc_alpha_t4)); + + vector anc_clients_t4(population_t4 * exp(log_asfr_t4_offset + mu_asfr)); + vector anc_plhiv_t4(anc_clients_t4 * anc_rho_t4); + vector anc_already_art_t4(anc_plhiv_t4 * anc_alpha_t4); + + vector prop_art_ij_t4((Xart_idx * prop_art_t4) * (Xart_gamma * gamma_art_t2)); // Note: using same ART attendance as T2 vector population_ij_t4(Xart_idx * population_t4); vector artnum_ij_t4(population_ij_t4 * prop_art_ij_t4); @@ -1040,6 +1061,18 @@ Type objective_function::operator() () vector infections_t4_out(A_out * infections_t4); vector lambda_t4_out(infections_t4_out / (population_t4_out - plhiv_t4_out)); + vector anc_clients_t4_out(A_anc_out * anc_clients_t4); + vector anc_plhiv_t4_out(A_anc_out * anc_plhiv_t4); + vector anc_already_art_t4_out(A_anc_out * anc_already_art_t4); + vector anc_art_new_t4_out(anc_plhiv_t4_out - anc_already_art_t4_out); + vector anc_known_pos_t4_out(anc_already_art_t4_out); + vector anc_tested_pos_t4_out(anc_plhiv_t4_out - anc_known_pos_t4_out); + vector anc_tested_neg_t4_out(anc_clients_t4_out - anc_plhiv_t4_out); + + vector anc_rho_t4_out(anc_plhiv_t4_out / anc_clients_t4_out); + vector anc_alpha_t4_out(anc_already_art_t4_out / anc_plhiv_t4_out); + + REPORT(population_t4_out); // REPORT(rho_t4_out); REPORT(plhiv_t4_out); @@ -1053,6 +1086,17 @@ Type objective_function::operator() () REPORT(lambda_t4_out); REPORT(infections_t4_out); + REPORT(anc_clients_t4_out); + REPORT(anc_plhiv_t4_out); + REPORT(anc_already_art_t4_out); + REPORT(anc_art_new_t4_out); + REPORT(anc_known_pos_t4_out); + REPORT(anc_tested_pos_t4_out); + REPORT(anc_tested_neg_t4_out); + REPORT(anc_rho_t4_out); + REPORT(anc_alpha_t4_out); + + // Projection to time 5 // Only PLHIV, ART and infections calculated. No ANC or awareness indicators From dac898bd4de39d61670b2c6e62103468b34136ff Mon Sep 17 00:00:00 2001 From: jeffeaton Date: Mon, 18 Nov 2024 05:42:48 -0500 Subject: [PATCH 2/3] add PMTCT outputs at T4 (.T) to datapack CSV --- NEWS.md | 3 +++ R/pepfar-datapack.R | 2 +- inst/datapack/datapack_indicator_mapping.csv | 7 +++++++ 3 files changed, 11 insertions(+), 1 deletion(-) diff --git a/NEWS.md b/NEWS.md index d57dacd5..393c2bfa 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,6 +1,9 @@ # naomi 2.10.2 * Add ANC testing outputs to T4 projection period for including in PEPFAR datapack output. +* Rename Datapack input CSV in the output ZIP folder for 2025 to `"pepfar_datapack_indicators_2025.csv"`. +* Add ANC testing indicators to T4 projection reprsenting the end of one year COP planning + period. New indicators are `PMTCT_STAT*.T` and `PMTCT_ART*.T` # naomi 2.10.1 diff --git a/R/pepfar-datapack.R b/R/pepfar-datapack.R index 504cd7af..d808b69a 100644 --- a/R/pepfar-datapack.R +++ b/R/pepfar-datapack.R @@ -1,4 +1,4 @@ -PEPFAR_DATAPACK_FILENAME <- "pepfar_datapack_indicators_2024.csv" +PEPFAR_DATAPACK_FILENAME <- "pepfar_datapack_indicators_2025.csv" #' Export naomi outputs to PEPFAR Data Pack format #' diff --git a/inst/datapack/datapack_indicator_mapping.csv b/inst/datapack/datapack_indicator_mapping.csv index 76fa2b5b..cc7b083d 100644 --- a/inst/datapack/datapack_indicator_mapping.csv +++ b/inst/datapack/datapack_indicator_mapping.csv @@ -11,6 +11,13 @@ anc_clients,RM8gRoxtsNw,PMTCT_STAT_SUBNAT.D.T_1,# Pregnant Women who present at anc_known_pos,tAE7ZD7p9zu,PMTCT_STAT_SUBNAT.N.Known.Pos.T_1,# Pregnant Women who are known HIV-positive at ANC1,TRUE,3 anc_tested_neg,tAE7ZD7p9zu,PMTCT_STAT_SUBNAT.N.New.Neg.T_1,# Pregnant Women who are newly tested and found HIV-negative at ANC1,TRUE,3 anc_tested_pos,tAE7ZD7p9zu,PMTCT_STAT_SUBNAT.N.New.Pos.T_1,# Pregnant Women who are newly tested and found HIV-positive at ANC1,TRUE,3 +anc_plhiv,PMTCT_ART_SUBNAT.D.T,PMTCT_ART_SUBNAT.D.T,# HIV-positIve Pregnant Women identified at ANC1,TRUE,4 +anc_already_art,PMTCT_ART_SUBNAT.N.Already.T,PMTCT_ART_SUBNAT.N.Already.T,# HIV-positive Pregnant Women already on ART at ANC1,TRUE,4 +anc_art_new,PMTCT_ART_SUBNAT.N.New.T,PMTCT_ART_SUBNAT.N.New.T,# HIV-positive Pregnant Women newly diagnosed at ANC1 and placed on ART,TRUE,4 +anc_clients,PMTCT_STAT_SUBNAT.D.T,PMTCT_STAT_SUBNAT.D.T,# Pregnant Women who present at ANC1,TRUE,4 +anc_known_pos,PMTCT_STAT_SUBNAT.N.Known.Pos.T,PMTCT_STAT_SUBNAT.N.Known.Pos.T,# Pregnant Women who are known HIV-positive at ANC1,TRUE,4 +anc_tested_neg,PMTCT_STAT_SUBNAT.N.New.Neg.T,PMTCT_STAT_SUBNAT.N.New.Neg.T,# Pregnant Women who are newly tested and found HIV-negative at ANC1,TRUE,4 +anc_tested_pos,PMTCT_STAT_SUBNAT.N.New.Pos.T,PMTCT_STAT_SUBNAT.N.New.Pos.T,# Pregnant Women who are newly tested and found HIV-positive at ANC1,TRUE,4 art_current,xghQXueYJxu,TX_CURR_SUBNAT.T_1,Total # of PLHIV on ART,TRUE,3 art_coverage,PopART.Rt.T_1,PopART.Rt.T_1,ART coverage,FALSE,3 aware_plhiv_num,nF19GOjcnoD,DIAGNOSED_SUBNAT.T_1,Number of PLHIV aware of HIV positive status,TRUE,3 From 1938b309104bbe4bea07b373bbce9b1c465ffed2 Mon Sep 17 00:00:00 2001 From: jeffeaton Date: Mon, 18 Nov 2024 06:03:46 -0500 Subject: [PATCH 3/3] update R implementation of objective funtion; test updates --- R/tmb-model-r-implementation.R | 42 +++++++++++++++++++++++++++++- tests/testthat/test-01-run-model.R | 6 ++--- 2 files changed, 44 insertions(+), 4 deletions(-) diff --git a/R/tmb-model-r-implementation.R b/R/tmb-model-r-implementation.R index 92e1b11a..e4be9c03 100644 --- a/R/tmb-model-r-implementation.R +++ b/R/tmb-model-r-implementation.R @@ -758,6 +758,26 @@ naomi_objective_function_r <- function(d, p) { infections_t4 <- lambda_t4 * (d$population_t4 - plhiv_t4) + ## Note: currently assuming same district effects parameters from t2 for t4 + mu_anc_rho_t4 <- qlogis(rho_t4) + + d$logit_anc_rho_t2_offset + + d$X_ancrho %*% (p$beta_anc_rho + p$beta_anc_rho_t2) + + d$Z_ancrho_x %*% (p$ui_anc_rho_x + p$ui_anc_rho_xt) + mu_anc_rho_t4 <- as.vector(mu_anc_rho_t4) + anc_rho_t4 <- stats::plogis(mu_anc_rho_t4) + + mu_anc_alpha_t4 <- mu_alpha_t4 + + d$logit_anc_alpha_t4_offset + + d$X_ancalpha %*% (p$beta_anc_alpha + p$beta_anc_alpha_t2) + + d$Z_ancalpha_x %*% (p$ui_anc_alpha_x + p$ui_anc_alpha_xt) + mu_anc_alpha_t4 <- as.vector(mu_anc_alpha_t4) + anc_alpha_t4 <- plogis(mu_anc_alpha_t4) + + anc_clients_t4 <- d$population_t4 * exp(d$log_asfr_t4_offset + mu_asfr) + anc_plhiv_t4 <- anc_clients_t4 * anc_rho_t4 + anc_already_art_t4 <- anc_plhiv_t4 * anc_alpha_t4 + + prop_art_ij_t4 <- as.vector(d$Xart_idx %*% prop_art_t4) * as.vector(d$Xart_gamma %*% gamma_art_t2) ## Note: using same ART attendance as T2 population_ij_t4 <- as.vector(d$Xart_idx %*% d$population_t4) artnum_ij_t4 <- population_ij_t4 * prop_art_ij_t4 @@ -774,11 +794,31 @@ naomi_objective_function_r <- function(d, p) { infections_t4_out <- as.vector(d$A_out %*% infections_t4) lambda_t4_out <- infections_t4_out / (population_t4_out - plhiv_t4_out) + anc_clients_t4_out <- as.vector(d$A_anc_out %*% anc_clients_t4) + anc_plhiv_t4_out <- as.vector(d$A_anc_out %*% anc_plhiv_t4) + anc_already_art_t4_out <- as.vector(d$A_anc_out %*% anc_already_art_t4) + anc_art_new_t4_out <- anc_plhiv_t4_out - anc_already_art_t4_out + anc_known_pos_t4_out <- anc_already_art_t4_out + anc_tested_pos_t4_out <- anc_plhiv_t4_out - anc_known_pos_t4_out + anc_tested_neg_t4_out <- anc_clients_t4_out - anc_plhiv_t4_out + + anc_rho_t4_out <- anc_plhiv_t4_out / anc_clients_t4_out + anc_alpha_t4_out <- anc_already_art_t4_out / anc_plhiv_t4_out + report_t4 <- list(population_t4_out = population_t4_out, plhiv_t4_out = plhiv_t4_out, plhiv_attend_t4_out = plhiv_attend_t4_out, lambda_t4_out = lambda_t4_out, - infections_t4_out = infections_t4_out) + infections_t4_out = infections_t4_out, + anc_clients_t4_out = anc_clients_t4_out, + anc_plhiv_t4_out = anc_plhiv_t4_out, + anc_already_art_t4_out = anc_already_art_t4_out, + anc_art_new_t4_out = anc_art_new_t4_out, + anc_known_pos_t4_out = anc_known_pos_t4_out, + anc_tested_pos_t4_out = anc_tested_pos_t4_out, + anc_tested_neg_t4_out = anc_tested_neg_t4_out, + anc_rho_t4_out = anc_rho_t4_out, + anc_alpha_t4_out = anc_alpha_t4_out) ## Projection to time 5 diff --git a/tests/testthat/test-01-run-model.R b/tests/testthat/test-01-run-model.R index 6ef10e97..3a2229ff 100644 --- a/tests/testthat/test-01-run-model.R +++ b/tests/testthat/test-01-run-model.R @@ -311,12 +311,12 @@ test_that("model run can be calibrated", { ## * 1 output time - 4 indicators [NOT INCLUDED IN PLOT OUTPUTS] ## ## ANC indicators outputs - ## 3 = number or output times + ## 4 = number or output times ## 9 = number of ANC indicators ## 9 = number of areas ## 12 = number of ANC age groups expect_equal(nrow(calibrated_output_obj$output_package$indicators), - 33 * 3 * 9 * (3 * 16 + 1 * 5 + 1 * 4) + 3 * 9 * 9 * 12) + 33 * 3 * 9 * (3 * 16 + 1 * 5 + 1 * 4) + 4 * 9 * 9 * 12) ## Plot data output: T3 and T4 indicators not included -> fewer rows plot_data_output <- read_hintr_output(calibrated_output$plot_data_path) @@ -505,7 +505,7 @@ test_that("Model can be run without .shiny90 file", { ## Check there is some data ## 11 indicators (5 fewer because missing awareness of status indicators) expect_equal(nrow(indicators_output$output_package$indicators), - 33 * 3 * 9 * (3 * 11 + 1 * 5 + 1 * 4) + 3 * 9 * 9 * 12) + 33 * 3 * 9 * (3 * 11 + 1 * 5 + 1 * 4) + 4 * 9 * 9 * 12) }) test_that("hintr_run_model can skip validation", {