From 933e87cf6dd4d7bc4313d078473694d07b3a149a Mon Sep 17 00:00:00 2001 From: Rob Ashton Date: Fri, 15 Mar 2024 10:56:17 +0000 Subject: [PATCH 1/3] Migrate to github actions, qualify package names --- .Rbuildignore | 8 +- .github/.gitignore | 1 + .github/workflows/R-CMD-check.yaml | 58 +++ .travis.yml | 8 - DESCRIPTION | 2 +- NAMESPACE | 1 - NEWS.md | 2 +- R/create_hts_param.R | 29 +- R/extract-pjnz.R | 10 +- R/first90-package.R | 10 + R/inputs.R | 111 ++--- R/likelihood.R | 114 ++--- R/plot_functions.R | 286 ++++++------- R/plot_outputs.R | 666 ++++++++++++++--------------- R/retest.R | 52 ++- R/simul.R | 17 +- R/time_functions.R | 10 +- README.Rmd | 8 +- README.md | 3 +- data-raw/initial-theta0-values.R | 14 +- first90release.Rproj | 20 + man/get_pjn_country.Rd | 2 +- man/get_pjn_region.Rd | 2 +- man/prepare_inputs.Rd | 3 +- tests/testthat/Rplots.pdf | Bin 0 -> 11830 bytes 25 files changed, 759 insertions(+), 678 deletions(-) create mode 100644 .github/.gitignore create mode 100644 .github/workflows/R-CMD-check.yaml delete mode 100644 .travis.yml create mode 100644 first90release.Rproj create mode 100644 tests/testthat/Rplots.pdf diff --git a/.Rbuildignore b/.Rbuildignore index 6b2aa7c..c5d67ae 100644 --- a/.Rbuildignore +++ b/.Rbuildignore @@ -1,5 +1,11 @@ ^README\.Rmd$ +^README\.md$ +^\.idea$ ^data-raw$ ^LICENSE\.md$ ^dev$ -^.travis.yml$ \ No newline at end of file +^.*\.Rproj$ +^\.Rproj\.user$ +^.github$ +^.ProgramDataGuidance\.md$ +^.SurveyDataGuidance\.md$ diff --git a/.github/.gitignore b/.github/.gitignore new file mode 100644 index 0000000..2d19fc7 --- /dev/null +++ b/.github/.gitignore @@ -0,0 +1 @@ +*.html diff --git a/.github/workflows/R-CMD-check.yaml b/.github/workflows/R-CMD-check.yaml new file mode 100644 index 0000000..87c3af0 --- /dev/null +++ b/.github/workflows/R-CMD-check.yaml @@ -0,0 +1,58 @@ +# Workflow derived from https://github.com/r-lib/actions/tree/v2/examples +# Need help debugging build failures? Start at https://github.com/r-lib/actions#where-to-find-help +on: + push: + branches: [main, master] + pull_request: + branches: [main, master] + +name: R-CMD-check + +jobs: + R-CMD-check: + runs-on: ${{ matrix.config.os }} + + name: ${{ matrix.config.os }} (${{ matrix.config.r }}) + + strategy: + fail-fast: false + matrix: + config: + - {os: macos-latest, r: 'release'} + - {os: windows-latest, r: 'release'} + - {os: ubuntu-latest, r: 'devel', http-user-agent: 'release'} + - {os: ubuntu-latest, r: 'release'} + - {os: ubuntu-latest, r: 'oldrel-1'} + + env: + GITHUB_PAT: ${{ secrets.GITHUB_TOKEN }} + R_KEEP_PKG_SOURCE: yes + + steps: + - uses: actions/checkout@v3 + + - uses: r-lib/actions/setup-pandoc@v2 + + - uses: r-lib/actions/setup-r@v2 + with: + r-version: ${{ matrix.config.r }} + http-user-agent: ${{ matrix.config.http-user-agent }} + use-public-rspm: true + + - uses: r-lib/actions/setup-r-dependencies@v2 + with: + extra-packages: any::rcmdcheck + needs: check + + - uses: actions/checkout@master + with: + name: mrc-ide/shiny90_sample_files + ## ${{ github.token }} is scoped to the current repository, so we + ## need to provide our own PAT + token: ${{ secrets.SAMPLE_FILES_GH_PAT }} + path: tests/testthat/sample_files + + - uses: r-lib/actions/check-r-package@v2 + with: + upload-snapshots: true + error-on: '"error"' diff --git a/.travis.yml b/.travis.yml deleted file mode 100644 index 359af69..0000000 --- a/.travis.yml +++ /dev/null @@ -1,8 +0,0 @@ -# R for travis: see documentation at https://docs.travis-ci.com/user/languages/r - -language: R -warnings_are_errors: false -before_install: - - git clone https://$github_pat@github.com/mrc-ide/shiny90_sample_files tests/testthat/sample_files -cache: - packages: true diff --git a/DESCRIPTION b/DESCRIPTION index 05c7ffc..1fb7a01 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -25,7 +25,7 @@ Encoding: UTF-8 LazyData: true ByteCompile: true Roxygen: list(markdown = TRUE) -RoxygenNote: 7.2.3 +RoxygenNote: 7.3.1 Depends: R (>= 2.10) Imports: diff --git a/NAMESPACE b/NAMESPACE index 7b1af33..f905a74 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -1,6 +1,5 @@ # Generated by roxygen2: do not edit by hand -export() export(add_ss_indices) export(artcov15to49) export(combine_inputs) diff --git a/NEWS.md b/NEWS.md index 4f6cd3c..4352059 100644 --- a/NEWS.md +++ b/NEWS.md @@ -122,7 +122,7 @@ Updates for 2022 UNAIDS estimates: # first90 1.4.3 -* Patch: in function `add_ss_indices()`, add argument `type.convert(..., as.is = TRUE)` +* Patch: in function `add_ss_indices()`, add argument `utils::type.convert(..., as.is = TRUE)` to suppress R 4.0 warning. * Bugfix: remove duplicate declaration of `double incrate_g[NG];` in `EPP_DIRECTINCID` incidence option. This did not affect any results because this option is not used in Shiny90 application; diff --git a/R/create_hts_param.R b/R/create_hts_param.R index 5ea4a1a..6f26b90 100644 --- a/R/create_hts_param.R +++ b/R/create_hts_param.R @@ -60,16 +60,16 @@ create_hts_param <- function(theta, fp) { n_k1 <- length(knots) - 1 # we remove 1995 n_k2 <- n_k1 + length(knots_rr_dxunt) rate_f <- exp(theta[1:n_k1]) - rr_dxunt <- plogis(theta[(n_k1 + 1):n_k2]) * 8 # Range is 0-8 - rr_m <- plogis(theta[(n_k2 + 1):(n_k2 + 2)]) * 1.1 - rr_test <- 0.95 + plogis(theta[(n_k2 + 3):(n_k2 + 4)]) * 7.05 # Range is 1-8 - rr_plhiv <- 0.05 + plogis(theta[n_k2 + 5]) * (1.95 - 0.05) - rr_dxart <- plogis(theta[n_k2 + 6]) - rr_25m <- 0.1 + plogis(theta[n_k2 + 7]) * (6 - 0.1) - rr_35m <- 0.1 + plogis(theta[n_k2 + 8]) * (6 - 0.1) - rr_25f <- 0.1 + plogis(theta[n_k2 + 9]) * (6 - 0.1) - rr_35f <- 0.1 + plogis(theta[n_k2 + 10]) * (6 - 0.1) - pr_oidx <- 0.25 + (plogis(theta[n_k2 + 11]) * (1.75 - 0.25)) + rr_dxunt <- stats::plogis(theta[(n_k1 + 1):n_k2]) * 8 # Range is 0-8 + rr_m <- stats::plogis(theta[(n_k2 + 1):(n_k2 + 2)]) * 1.1 + rr_test <- 0.95 + stats::plogis(theta[(n_k2 + 3):(n_k2 + 4)]) * 7.05 # Range is 1-8 + rr_plhiv <- 0.05 + stats::plogis(theta[n_k2 + 5]) * (1.95 - 0.05) + rr_dxart <- stats::plogis(theta[n_k2 + 6]) + rr_25m <- 0.1 + stats::plogis(theta[n_k2 + 7]) * (6 - 0.1) + rr_35m <- 0.1 + stats::plogis(theta[n_k2 + 8]) * (6 - 0.1) + rr_25f <- 0.1 + stats::plogis(theta[n_k2 + 9]) * (6 - 0.1) + rr_35f <- 0.1 + stats::plogis(theta[n_k2 + 10]) * (6 - 0.1) + pr_oidx <- 0.25 + (stats::plogis(theta[n_k2 + 11]) * (1.75 - 0.25)) # Age function for males and females agefn_m <- c(rep(1, 3), rep(rr_25m, 2), rep(rr_35m, 3), rr_35m * 0.8112) # last age group is 50+ @@ -77,15 +77,15 @@ create_hts_param <- function(theta, fp) { # Time trends in testing rates fp$t_hts_start <- as.integer(1995 - fp$ss$proj_start + 1L) - base_rate_f <- approx(knots, c(0, rate_f), seq_len(fp$ss$PROJ_YEARS), rule = 2)$y + base_rate_f <- stats::approx(knots, c(0, rate_f), seq_len(fp$ss$PROJ_YEARS), rule = 2)$y knots_rr_m <- c(2005, 2012) - fp$ss$proj_start + 1L - rr_m <- approx(knots_rr_m, rr_m, seq_len(fp$ss$PROJ_YEARS), rule = 2)$y + rr_m <- stats::approx(knots_rr_m, rr_m, seq_len(fp$ss$PROJ_YEARS), rule = 2)$y base_rate_m <- base_rate_f * rr_m # Retest with time knots_rr_test <- c(2005, 2010, 2015) - fp$ss$proj_start + 1L - rate_rr_test <- exp(approx(knots_rr_test, log(c(1, rr_test)), seq_len(fp$ss$PROJ_YEARS), rule = 2)$y) + rate_rr_test <- exp(stats::approx(knots_rr_test, log(c(1, rr_test)), seq_len(fp$ss$PROJ_YEARS), rule = 2)$y) ## Testing rate for HIV negative population: ## 1: never tested @@ -147,7 +147,7 @@ create_hts_param <- function(theta, fp) { ## Relative testing among aware # Re-testing rate among PLHIV aware (not on ART) - rate_dxunt <- approx(knots_rr_dxunt, rr_dxunt, seq_len(fp$ss$PROJ_YEARS), rule = 2)$y + rate_dxunt <- stats::approx(knots_rr_dxunt, rr_dxunt, seq_len(fp$ss$PROJ_YEARS), rule = 2)$y diagn_rate_dxunt <- sweep(diagn_rate[,,,3,, drop = FALSE], 5, rate_dxunt, "*") diagn_rate[,,,3,] <- diagn_rate_dxunt @@ -162,4 +162,3 @@ create_hts_param <- function(theta, fp) { fp } - diff --git a/R/extract-pjnz.R b/R/extract-pjnz.R index b92c4d9..1f11849 100644 --- a/R/extract-pjnz.R +++ b/R/extract-pjnz.R @@ -134,9 +134,9 @@ extract_pjnz <- function(pjnz = NULL, dp_file= NULL, pjn_file = NULL){ v$art15plus_num <- v$art15plus_num * adult_artadj_factor } - v$art15plus_eligthresh <- setNames(as.numeric(dpsub(dp, "", 2, col_idx)), proj_years) + v$art15plus_eligthresh <- stats::setNames(as.numeric(dpsub(dp, "", 2, col_idx)), proj_years) - artelig_specpop <- setNames(dpsub(dp, "", 3:9, 2:6), + artelig_specpop <- stats::setNames(dpsub(dp, "", 3:9, 2:6), c("description", "pop", "elig", "percent", "year")) artelig_specpop$pop <- c("PW", "TBHIV", "DC", "FSW", "MSM", "IDU", "OTHER") artelig_specpop$elig <- as.logical(as.integer(artelig_specpop$elig)) @@ -146,8 +146,8 @@ extract_pjnz <- function(pjnz = NULL, dp_file= NULL, pjn_file = NULL){ rownames(artelig_specpop) <- artelig_specpop$pop v$artelig_specpop <- artelig_specpop - v$median_cd4init <- setNames(as.numeric(dpsub(dp, "", 2, col_idx)), proj_years) - v$art_dropout <- setNames(as.numeric(dpsub(dp, "", 2, col_idx)), proj_years) + v$median_cd4init <- stats::setNames(as.numeric(dpsub(dp, "", 2, col_idx)), proj_years) + v$art_dropout <- stats::setNames(as.numeric(dpsub(dp, "", 2, col_idx)), proj_years) v$art_alloc_method <- get_dp_art_alloc_method(dp) v$art_prop_alloc <- get_dp_art_prop_alloc(dp) @@ -763,7 +763,7 @@ get_pjn_iso3 <- function(pjn){ #' @export file_in_zip <- function(zfile, ext){ - ff <- grep(ext, unzip(zfile, list = TRUE)$Name, value = TRUE) + ff <- grep(ext, utils::unzip(zfile, list = TRUE)$Name, value = TRUE) unz(zfile, ff) } diff --git a/R/first90-package.R b/R/first90-package.R index b7c13e5..bc134b9 100644 --- a/R/first90-package.R +++ b/R/first90-package.R @@ -5,3 +5,13 @@ NULL #' @importFrom Rcpp sourceCpp ## usethis namespace: end NULL + +## Suppress R CMD check NOTEs about "no visible binding for +## global variable". These NOTEs aren't that important as we're using +## but the thousands of lines of NOTE they print does obscure real issues. +utils::globalVariables(c( + "DT", "PROJ_YEARS", "ag.idx", "agegr", "aglast.idx", "country", "f.idx", + "h.ag.span", "h.age15plus.idx", "h.fert.idx", "hAG", "hDS", "hTS", + "hiv_steps_per_year", "hivn.idx", "hivp.idx", "hivstatus", "logit", "m.idx", + "nearPD", "outcome", "p.age15plus.idx", "p.age15to49.idx", "p.fert.idx", + "pAG", "pDS", "quantile", "sex", "tot", "totpos", "year")) diff --git a/R/inputs.R b/R/inputs.R index 0ffea3b..323b9dc 100644 --- a/R/inputs.R +++ b/R/inputs.R @@ -1,5 +1,5 @@ #' Process programmatic data on number of tests -#' +#' #' @export select_prgmdata <- function(prgm_dat, cnt, age_group) { # any country with HTS data? if so, we select them @@ -24,19 +24,19 @@ select_prgmdata <- function(prgm_dat, cnt, age_group) { anc = prg.i$anc[prg.i$sex == 'female'], ancpos = prg.i$ancpos[prg.i$sex == 'female']) value_verif <- prg_dat_both$totpos[prg_dat_both$year == yr_sex[i]] - if (length(value_verif) > 0 & is.na(tot.i$totpos)) { + if (length(value_verif) > 0 & is.na(tot.i$totpos)) { tot.i$totpos <- prg_dat_both$totpos[prg_dat_both$year == yr_sex[i]] } prg_dat_sex_both <- rbind(prg_dat_sex_both, tot.i) } prg_dat <- rbind(prg_dat_both[!(prg_dat_both$year %in% yr_sex), ], prg_dat_sex_both, prg_dat_sex) - } } - + } } + if (!any(prgm_dat$country == cnt)) { ## -- UPDATE HERE -- ## * year vector needs to be extended to output results to current year - - prg_dat <- data.frame(country = cnt, + + prg_dat <- data.frame(country = cnt, year = 2010:2022, agegr = '15-99', sex = 'both', tot = NA, totpos = NA, @@ -48,65 +48,65 @@ select_prgmdata <- function(prgm_dat, cnt, age_group) { } #' Process survey data on hiv testing behaviors -#' +#' #' @export select_hts <- function(survey_hts, cnt, age_group) { require(first90) # We first select stratified by HIV status for females - dat_f <- survey_hts[survey_hts$country == cnt & (survey_hts$hivstatus != "all") & - (survey_hts$agegr %in% age_group) & survey_hts$sex == "female" & + dat_f <- survey_hts[survey_hts$country == cnt & (survey_hts$hivstatus != "all") & + (survey_hts$agegr %in% age_group) & survey_hts$sex == "female" & survey_hts$outcome == "evertest", ] # We first select stratified by HIV status for males - dat_m <- survey_hts[survey_hts$country == cnt & (survey_hts$hivstatus != "all") & - (survey_hts$agegr %in% age_group) & survey_hts$sex == "male" & + dat_m <- survey_hts[survey_hts$country == cnt & (survey_hts$hivstatus != "all") & + (survey_hts$agegr %in% age_group) & survey_hts$sex == "male" & survey_hts$outcome == "evertest", ] if (cnt == 'Equatorial Guinea') { # special case of GE b/c data for HIV entered manually (not age disaggregated) - dat_f <- survey_hts[survey_hts$country == cnt & survey_hts$hivstatus == "all" & + dat_f <- survey_hts[survey_hts$country == cnt & survey_hts$hivstatus == "all" & (survey_hts$agegr %in% age_group) & survey_hts$sex == "female" & survey_hts$outcome == "evertest", ] - - dat_m <- survey_hts[survey_hts$country == cnt & survey_hts$hivstatus == "all" & + + dat_m <- survey_hts[survey_hts$country == cnt & survey_hts$hivstatus == "all" & (survey_hts$agegr %in% age_group) & survey_hts$sex == "male" & survey_hts$outcome == "evertest", ] - - dat_f_hiv <- survey_hts[survey_hts$country == cnt & survey_hts$hivstatus == "positive" & + + dat_f_hiv <- survey_hts[survey_hts$country == cnt & survey_hts$hivstatus == "positive" & survey_hts$agegr == '15-49' & survey_hts$sex == "female" & survey_hts$outcome == "evertest", ] - - dat_m_hiv <- survey_hts[survey_hts$country == cnt & survey_hts$hivstatus == "positive" & + + dat_m_hiv <- survey_hts[survey_hts$country == cnt & survey_hts$hivstatus == "positive" & survey_hts$agegr == '15-49' & survey_hts$sex == "male" & survey_hts$outcome == "evertest", ] - + dat_f <- rbind(dat_f, dat_f_hiv) dat_m <- rbind(dat_m, dat_m_hiv) } - + # We find surveys that collected information on HIV status but for different age groups all_svy_hiv <- unique(survey_hts$surveyid[survey_hts$country == cnt & (survey_hts$hivstatus != "all") & survey_hts$outcome == "evertest"]) other_svy_hiv <- all_svy_hiv[!(all_svy_hiv %in% unique(rbind(dat_f, dat_m)$surveyid))] - + dat_both <- NULL if (length(other_svy_hiv) > 0) { - dat_f2 <- survey_hts[(survey_hts$surveyid %in% other_svy_hiv) & (survey_hts$hivstatus != "all") & - survey_hts$agegr == '15-49' & survey_hts$sex == "female" & + dat_f2 <- survey_hts[(survey_hts$surveyid %in% other_svy_hiv) & (survey_hts$hivstatus != "all") & + survey_hts$agegr == '15-49' & survey_hts$sex == "female" & survey_hts$outcome == "evertest", ] - dat_m2 <- survey_hts[(survey_hts$surveyid %in% other_svy_hiv) & (survey_hts$hivstatus != "all") & - survey_hts$agegr == '15-49' & survey_hts$sex == "male" & + dat_m2 <- survey_hts[(survey_hts$surveyid %in% other_svy_hiv) & (survey_hts$hivstatus != "all") & + survey_hts$agegr == '15-49' & survey_hts$sex == "male" & survey_hts$outcome == "evertest", ] if (dim(dat_f2)[1] > 0 | dim(dat_m2)[1] > 0) { dat_f <- rbind(dat_f, dat_f2) - dat_m <- rbind(dat_m, dat_m2) + dat_m <- rbind(dat_m, dat_m2) } else { - dat_both <- survey_hts[(survey_hts$surveyid %in% other_svy_hiv) & (survey_hts$hivstatus != "all") & - survey_hts$agegr == '15-49' & survey_hts$sex == "both" & + dat_both <- survey_hts[(survey_hts$surveyid %in% other_svy_hiv) & (survey_hts$hivstatus != "all") & + survey_hts$agegr == '15-49' & survey_hts$sex == "both" & survey_hts$outcome == "evertest", ] } } - + # We find surveys that did not collect information on HIV status all_survey <- unique(survey_hts$surveyid[survey_hts$country == cnt & survey_hts$outcome == "evertest"]) other_survey <- all_survey[!(all_survey %in% unique(rbind(dat_f, dat_m, dat_both)$surveyid))] - + dat_fall <- NULL dat_mall <- NULL @@ -115,7 +115,7 @@ select_hts <- function(survey_hts, cnt, age_group) { age_strat_f <- unique(survey_hts$agegr[survey_hts$country == cnt & survey_hts$hivstatus == "all" & survey_hts$surveyid == other_survey[i] & survey_hts$sex == "female" & survey_hts$outcome == "evertest"]) age_strat_m <- unique(survey_hts$agegr[survey_hts$country == cnt & survey_hts$hivstatus == "all" & survey_hts$surveyid == other_survey[i] & - survey_hts$sex == "male" & survey_hts$outcome == "evertest"]) + survey_hts$sex == "male" & survey_hts$outcome == "evertest"]) ## Female - If we have the correct age disagregation, we stick with it if (all(age_group %in% age_strat_f)) { dat_f_i <- survey_hts[survey_hts$country == cnt & survey_hts$hivstatus == "all" & survey_hts$surveyid == other_survey[i] & @@ -129,7 +129,7 @@ select_hts <- function(survey_hts, cnt, age_group) { } else { dat_f_i <- survey_hts[survey_hts$country == cnt & survey_hts$hivstatus == "all" & survey_hts$surveyid == other_survey[i] & (survey_hts$agegr %in% c('15-99')) & survey_hts$sex == "female" & survey_hts$outcome == "evertest",] - } + } ## Female - If we have the correct age disagregation, we stick with it if (all(age_group %in% age_strat_m)) { dat_m_i <- survey_hts[survey_hts$country == cnt & survey_hts$hivstatus == "all" & survey_hts$surveyid == other_survey[i] & @@ -154,7 +154,7 @@ select_hts <- function(survey_hts, cnt, age_group) { } #' Function to present the survey data on shinny interface -#' +#' #' @export svy_hts_interface <- function(survey_hts, cnt, age_group = c('15-24','25-34','35-49')) { # Data used in the fitting @@ -172,13 +172,13 @@ svy_hts_interface <- function(survey_hts, cnt, age_group = c('15-24','25-34','35 sex == 'both' & (outcome %in% c('aware', 'art_both', 'art_self') | hivstatus == 'positive' & outcome == 'evertest')) - + svy_plot <- rbind(svy_plot1549, svy_plotage, svy_plot90s) svy_plot <- svy_plot[order(svy_plot$year, svy_plot$sex, svy_plot$hivstatus, svy_plot$agegr), ] svy_plot$color <- "Plotting" - + svy <- rbind(svy_dat, svy_plot) - + return(svy) } @@ -201,9 +201,9 @@ add_ss_indices <- function(dat, ss) { ## Convert 15+ to 15-99 for parsing below agegr_tmp <- sub("15\\+", "15-99", df$agegr) - - df$agestart <- type.convert(sub("([0-9]+)-([0-9]+)", "\\1", agegr_tmp), as.is = TRUE) - ageend <- type.convert(sub("([0-9]+)-([0-9]+)", "\\2", agegr_tmp), as.is = TRUE)+1L + + df$agestart <- utils::type.convert(sub("([0-9]+)-([0-9]+)", "\\1", agegr_tmp), as.is = TRUE) + ageend <- utils::type.convert(sub("([0-9]+)-([0-9]+)", "\\2", agegr_tmp), as.is = TRUE)+1L df$aidx <- df$agestart - ss$AGE_START + 1L df$agspan <- ageend - df$agestart @@ -215,14 +215,14 @@ add_ss_indices <- function(dat, ss) { ## if max age is 65 or greater, assume age 50+ end_haidx[ageend >= 65] <- length(hag_breaks) + 1L - + df$hagspan <- end_haidx - df$haidx - + ## testing history is homogenous among oldest age group df$haidx[df$agestart >= max(hag_breaks) & ageend >= max(hag_breaks)] <- length(hag_breaks) df$hagspan[df$agestart >= max(hag_breaks) & ageend >= max(hag_breaks)] <- 1L - + df$sidx <- match(df$sex, c("both", "male", "female")) - 1L df$hvidx <- match(df$hivstatus, c("all", "negative", "positive")) - 1L df$yidx <- df$year - ss$proj_start + 1L @@ -242,32 +242,32 @@ add_ss_indices <- function(dat, ss) { prepare_inputs_from_extracts <- function(pjnz_in){ pjnz_aggr <- combine_inputs(pjnz_in) - + demp <- list(basepop = pjnz_aggr$totpop, Sx = pjnz_aggr$Sx, asfr = pjnz_aggr$asfr, srb = pjnz_aggr$srb, netmigr = pjnz_aggr$netmigr) - + projp <- pjnz_aggr projp$age14totpop <- projp$totpop["14",,] dimnames(projp$fert_rat)[[1]] <- gsub("(.*)-.*", "\\1", dimnames(projp$fert_rat)[[1]]) - + fp <- create_fp(projp, demp) ## Set pop adjust fp$popadjust <- pjnz_aggr$popadjust - + ## Calculate incidence input fp$infections <- pjnz_aggr$infections[fp$ss$AGE_START + fp$ss$p.age15plus.idx, , ] fp$eppmod <- "directinfections_hts" - + ## initialise to no testing fp$hts_rate <- array(0.0, c(fp$ss$hAG, fp$ss$NG, fp$ss$pDS, fp$ss$PROJ_YEARS)) fp$diagn_rate <- array(0.0, c(fp$ss$hDS, fp$ss$hAG, fp$ss$NG, 4, fp$ss$PROJ_YEARS)) - + fp$t_hts_start <- fp$ss$PROJ_YEARS+1L - + fp } @@ -285,9 +285,10 @@ prepare_inputs_from_extracts <- function(pjnz_in){ #' #' @examples #' -#' pjnzlist <- list.files("~/Documents/Data/Spectrum files/2018 final/SSA/", "CotedIvoire.*PJNZ$", full.names=TRUE, ignore.case=TRUE) +#' pjnzlist <- list.files("~/Documents/Data/Spectrum files/2018 final/SSA/", +#' "CotedIvoire.*PJNZ$", full.names=TRUE, ignore.case=TRUE) #' pjnzlist <- "~/Documents/Data/Spectrum files/2018 final/SSA/Malawi_2018_version_8.PJNZ" -#' +#' #' @export prepare_inputs <- function(pjnzlist){ @@ -300,13 +301,13 @@ add_fp <- function(v, fp){ fp[names(v)] <- v fp } - + create_paedsurv_inputs <- function(projp, ss, artcd4elig_idx) { proj_start <- ss$proj_start proj_end <- proj_start + ss$PROJ_YEARS - 1L - + val <- list() ## HIV prevalence and ART coverage among age 15 entrants @@ -337,7 +338,7 @@ create_paedsurv_inputs <- function(projp, ss, artcd4elig_idx) { ## sum to 1 in each sex and year. for(g in 1:ss$NG) for(i in 2:ss$PROJ_YEARS){ - + if((hiv14[g,i-1] - art14[g,i-1]) > 0) val$paedsurv_cd4dist[,g,i] <- hiv_noart14[,g,i-1] %*% cd4convert / (hiv14[g,i-1] - art14[g,i-1]) if(art14[g,i-1]){ @@ -359,7 +360,7 @@ create_paedsurv_inputs <- function(projp, ss, artcd4elig_idx) { .fp_targetpop <- function(totpop, ss){ v <- list() - + proj_start <- ss$proj_start proj_end <- proj_start + ss$PROJ_YEARS - 1L diff --git a/R/likelihood.R b/R/likelihood.R index 35809a1..735f9dd 100644 --- a/R/likelihood.R +++ b/R/likelihood.R @@ -1,9 +1,9 @@ #' @export ll_hts <- function(theta, fp, likdat) { - + fp <- create_hts_param(theta, fp) mod <- simmod(fp) - + val1 <- ll_evertest(mod, fp, likdat$evertest) if (!is.null(likdat$hts)) { val2 <- ll_prgdat(mod, fp, likdat$hts) @@ -11,11 +11,11 @@ ll_hts <- function(theta, fp, likdat) { if (!is.null(likdat$hts_pos)) { val3 <- ll_prgdat(mod, fp, likdat$hts_pos) } else { val3 <- 0 } - + val_prior <- lprior_hts(theta, mod, fp) val <- val1 + val2 + val3 + val_prior - + return(val) } @@ -34,11 +34,11 @@ ll_evertest <- function(mod, fp, dat) { } if (any(is.na(mu)) || any(mu < 0) || any(mu > 1)) { - llk <- log(0) + llk <- log(0) } else { - llk <- sum(dbinom(x = dat$nsuccess, size = dat$neff, prob = mu, log = TRUE)) + llk <- sum(stats::dbinom(x = dat$nsuccess, size = dat$neff, prob = mu, log = TRUE)) } - + return(llk) } @@ -47,16 +47,17 @@ ll_evertest <- function(mod, fp, dat) { ll_prgdat <- function(mod, fp, dat) { if (as.character(dat$hivstatus[1]) == 'all') { mu <- total_tests(mod, df = add_ss_indices(dat, fp$ss)) - llk <- sum(dnorm(x = dat$tot, mean = mu, sd = dat$l_est_se, log = TRUE)) - } + llk <- sum(stats::dnorm(x = dat$tot, mean = mu, sd = dat$l_est_se, log = TRUE)) + } if (as.character(dat$hivstatus[1]) == 'positive') { mu <- total_tests(mod, df = add_ss_indices(dat, fp$ss)) - llk <- sum(dnorm(x = dat$tot, mean = mu, sd = dat$l_est_se, log = TRUE)) + llk <- sum(stats::dnorm(x = dat$tot, mean = mu, sd = dat$l_est_se, log = TRUE)) } if (!(as.character(dat$hivstatus[1]) %in% c('all', 'positive'))) { - print('Error - HIV status is incorrect'); break + print('Error - HIV status is incorrect') + return() } - + return(llk) } @@ -68,56 +69,56 @@ art_constraint_penalty <- function(mod, fp, max_year = 2022) { tot_late <- apply(attr(mod, "late_diagnoses")[,,, ind_year], 4, sum) tot_untreated_pop <- apply(attr(mod, "hivpop")[,,, ind_year], 4, sum) ratio_late_per_1000_untreated <- tot_late / tot_untreated_pop * 1000 - penalty <- sum(dnorm(x = 0, mean = ratio_late_per_1000_untreated, + penalty <- sum(stats::dnorm(x = 0, mean = ratio_late_per_1000_untreated, sd = 20, log = TRUE)) return(penalty) } # Include this in ll_hts if you want to incorporate the likelihood constraint on ART. # val_art_penalty <- art_constraint_penalty(mod, fp, max_year = 2022) # val <- val1 + val2 + val3 + val_prior + val_art_penalty - + # Function to prepare the data for input in the likelihood function. #' @export prepare_hts_likdat <- function(dat_evertest, dat_prg, fp) { - + # Testing behavior data dat_evertest$est <- ifelse(is.na(dat_evertest$est), NA, dat_evertest$est) - dat_evertest <- dat_evertest[complete.cases(dat_evertest$est),] - dat_evertest$l_est <- qlogis(dat_evertest$est) + dat_evertest <- dat_evertest[stats::complete.cases(dat_evertest$est),] + dat_evertest$l_est <- stats::qlogis(dat_evertest$est) dat_evertest$l_est_se <- dat_evertest$se / (dat_evertest$est * (1 - dat_evertest$est)) - + # if using binomial, we calculate Neff from SE - dat_evertest$neff <- (dat_evertest$est * (1 - dat_evertest$est)) / dat_evertest$se^2 + dat_evertest$neff <- (dat_evertest$est * (1 - dat_evertest$est)) / dat_evertest$se^2 dat_evertest$neff <- ifelse(dat_evertest$est < 1e-5, dat_evertest$counts, dat_evertest$neff) dat_evertest$neff <- ifelse(dat_evertest$est > 0.999, dat_evertest$counts, dat_evertest$neff) dat_evertest$neff <- ifelse(dat_evertest$neff == Inf, dat_evertest$counts, dat_evertest$neff) dat_evertest$nsuccess <- round(dat_evertest$est * dat_evertest$neff, 0) dat_evertest$neff <- ifelse(round(dat_evertest$neff, 0) < 1, 1, round(dat_evertest$neff, 0)) dat_evertest <- add_ss_indices(dat_evertest, fp$ss) - + # Verifying if input data is OK if(any(dat_prg$agegr != '15-99')) { stop(print('Error, HTS program data should be for the 15+ age group only /n other age-grouping not supported at the moment'))} - + # We remove years with NA from programmatic data dat_prg_pos <- dat_prg - dat_prg_pos <- dat_prg_pos[complete.cases(dat_prg_pos$totpos),] + dat_prg_pos <- dat_prg_pos[stats::complete.cases(dat_prg_pos$totpos),] dat_prg_pos_sex <- subset(dat_prg_pos, sex != 'both') yr_sex <- unique(dat_prg_pos_sex$year) dat_prg_pos1 <- subset(dat_prg_pos, sex == 'both' & !(year %in% yr_sex)) dat_prg_pos2 <- subset(dat_prg_pos, sex != 'both' & (year %in% yr_sex)) dat_prg_pos <- rbind(dat_prg_pos1, dat_prg_pos2) - + dat_prg_pos$tot <- dat_prg_pos$totpos - + dat_prg_sex <- subset(dat_prg, sex != 'both') yr_sex <- unique(dat_prg_sex$year) dat_prg1 <- subset(dat_prg, sex == 'both' & !(year %in% yr_sex)) dat_prg2 <- subset(dat_prg, sex != 'both' & (year %in% yr_sex)) dat_prg <- rbind(dat_prg1, dat_prg2) - dat_prg <- dat_prg[complete.cases(dat_prg$tot),] - - # Programmatic data - nb of test (total: VCT + PMTCT + other; ideally among 15+) + dat_prg <- dat_prg[stats::complete.cases(dat_prg$tot),] + + # Programmatic data - nb of test (total: VCT + PMTCT + other; ideally among 15+) if (dim(dat_prg)[1] > 0) { dat_prg$l_est_se <- ifelse(dat_prg$sex == 'both', dat_prg$tot * 0.05, dat_prg$tot * 0.05 * sqrt(1 + 1 + 2*1)) # the sqrt() is to maintain same SE as if fitting on both sex - assuming perfect correlation dat_prg$hivstatus <- 'all' @@ -130,7 +131,7 @@ prepare_hts_likdat <- function(dat_evertest, dat_prg, fp) { dat_prg_pos$hivstatus <- 'positive' dat_prg_pos <- add_ss_indices(dat_prg_pos, fp$ss) } else { dat_prg_pos <- NULL } - + return(list(evertest = dat_evertest, hts = dat_prg, hts_pos = dat_prg_pos)) } @@ -143,47 +144,46 @@ lprior_hts <- function(theta, mod, fp) { ## * Extend knots by 1 year to current year knots <- 2000:2023 - fp$ss$proj_start + 1L ## -- UPDATE ABOVE -- - + n_k1 <- length(knots) - n_k2 <- n_k1 * 2 - 10 + n_k2 <- n_k1 * 2 - 10 penalty_f <- log(fp$hts_rate[1,2,1, knots[-1]]) - log(fp$hts_rate[1,2,1, knots[-n_k1]]) penalty_rr_dxunt <- theta[c((n_k1 + 2):n_k2)] - theta[c((n_k1 + 1):(n_k2 - 1))] - penalty_rr_m <- log(plogis(theta[n_k2 + 2]) * 1.1) - log(plogis(theta[n_k2 + 1]) * 1.1) + penalty_rr_m <- log(stats::plogis(theta[n_k2 + 2]) * 1.1) - log(stats::plogis(theta[n_k2 + 1]) * 1.1) penalty_rr_test <- theta[n_k2 + 4] - theta[n_k2 + 3] - lprior <- + lprior <- ## Prior for first baseline rate for females # exp(log(0.001) + 1.96*0.25) - dnorm(x = theta[1], mean = log(0.005), sd = 0.25, log = TRUE) + - ## Relative testing among PLHIV diagnosed, untreated. 1.50 (95%CI: 0.14-6.00) # plogis(qlogis(1.5/8) - 1.96*1.31)*8 - dnorm(x = theta[n_k2], mean = qlogis(1.5/8), sd = 1.31, log = TRUE) + - ## Prior for male RR 0.6 (95%CI: 0.07-1.05)# plogis(qlogis(0.6/1.1) + 1.96*1.46) * 1.1 - sum(dnorm(x = theta[n_k2 + 1], mean = qlogis(0.6/1.1), sd = 1.46, log = TRUE)) + - ## Relative increase among previously tested. 1.93 (95%CI: 1.08-5.00) # 0.95 + plogis(qlogis((1.93 - 0.95)/7.05) - qnorm(0.975)*1.084)*7.05 - dnorm(x = theta[n_k2 + 3], mean = qlogis((1.93 - 0.95)/7.05), sd = 1.084, log = TRUE) + - ## Relative factor for re-testing among PLHIV. 1.00 (95%CI: 0.10-1.90) # 0.05 + plogis(qlogis(0.95 / 1.90) - 1.96*1.85) * (1.95 - 0.05) - dnorm(x = theta[n_k2 + 5], mean = qlogis(0.95 / 1.90), sd = 1.85, log = TRUE) + - ## Relative testing among PLHIV already on ART (95%CI: 0.01-0.90) # plogis(qlogis(0.25) - 1.96*1.68) - dnorm(x = theta[n_k2 + 6], mean = qlogis(0.25), sd = 1.68, log = TRUE) + + stats::dnorm(x = theta[1], mean = log(0.005), sd = 0.25, log = TRUE) + + ## Relative testing among PLHIV diagnosed, untreated. 1.50 (95%CI: 0.14-6.00) # stats::plogis(stats::qlogis(1.5/8) - 1.96*1.31)*8 + stats::dnorm(x = theta[n_k2], mean = stats::qlogis(1.5/8), sd = 1.31, log = TRUE) + + ## Prior for male RR 0.6 (95%CI: 0.07-1.05)# stats::plogis(stats::qlogis(0.6/1.1) + 1.96*1.46) * 1.1 + sum(stats::dnorm(x = theta[n_k2 + 1], mean = stats::qlogis(0.6/1.1), sd = 1.46, log = TRUE)) + + ## Relative increase among previously tested. 1.93 (95%CI: 1.08-5.00) # 0.95 + stats::plogis(stats::qlogis((1.93 - 0.95)/7.05) - stats::qnorm(0.975)*1.084)*7.05 + stats::dnorm(x = theta[n_k2 + 3], mean = stats::qlogis((1.93 - 0.95)/7.05), sd = 1.084, log = TRUE) + + ## Relative factor for re-testing among PLHIV. 1.00 (95%CI: 0.10-1.90) # 0.05 + stats::plogis(stats::qlogis(0.95 / 1.90) - 1.96*1.85) * (1.95 - 0.05) + stats::dnorm(x = theta[n_k2 + 5], mean = stats::qlogis(0.95 / 1.90), sd = 1.85, log = TRUE) + + ## Relative testing among PLHIV already on ART (95%CI: 0.01-0.90) # stats::plogis(stats::qlogis(0.25) - 1.96*1.68) + stats::dnorm(x = theta[n_k2 + 6], mean = stats::qlogis(0.25), sd = 1.68, log = TRUE) + ## Prior for age (95% CI is 0.14-5.0)# 0.1 + invlogit(logit(0.9/5.9) - 1.96*1.5) * (6 - 0.1) - sum(dnorm(x = theta[c((n_k2 + 7):(n_k2 + 10))], mean = qlogis(0.9/5.9), sd = 1.685, log = TRUE)) + - ## RR OI diangosed for HIV relative to ART coverage 1.0 (0.3-1.7) # 0.25 + plogis(qlogis(0.5) - 1.96*1.75) * (1.75 - 0.25) - dnorm(x = theta[n_k2 + 11], mean = qlogis(0.5), sd = 1.75, log = TRUE) - + sum(stats::dnorm(x = theta[c((n_k2 + 7):(n_k2 + 10))], mean = stats::qlogis(0.9/5.9), sd = 1.685, log = TRUE)) + + ## RR OI diangosed for HIV relative to ART coverage 1.0 (0.3-1.7) # 0.25 + stats::plogis(stats::qlogis(0.5) - 1.96*1.75) * (1.75 - 0.25) + stats::dnorm(x = theta[n_k2 + 11], mean = stats::qlogis(0.5), sd = 1.75, log = TRUE) + prior_sd_f <- 0.205 # exp(log(0.5) - 1.96*0.205); previous of 0.35 when double-counting prior_sd_rr_m <- 0.26 # exp(log(0.6) + 1.96*0.26) prior_sd_rr_dxunt <- 0.25 prior_sd_rr_test <- 0.25 - - hyperprior <- + + hyperprior <- ## Penalty for 1st degree difference (female baseline rate) - sum(dnorm(x = penalty_f, mean = 0, sd = prior_sd_f, log = TRUE)) + - ## Penalty for 1st degree difference (male RR) - sum(dnorm(x = penalty_rr_m, mean = 0, sd = prior_sd_rr_m, log = TRUE)) + - ## Penalty for 1st degree difference (RR_dxunt) - sum(dnorm(x = penalty_rr_dxunt, mean = 0, sd = prior_sd_rr_dxunt, log = TRUE)) + - ## Penalty for 1st degree difference (RR_test) - sum(dnorm(x = penalty_rr_test, mean = 0, sd = prior_sd_rr_test, log = TRUE)) + sum(stats::dnorm(x = penalty_f, mean = 0, sd = prior_sd_f, log = TRUE)) + + ## Penalty for 1st degree difference (male RR) + sum(stats::dnorm(x = penalty_rr_m, mean = 0, sd = prior_sd_rr_m, log = TRUE)) + + ## Penalty for 1st degree difference (RR_dxunt) + sum(stats::dnorm(x = penalty_rr_dxunt, mean = 0, sd = prior_sd_rr_dxunt, log = TRUE)) + + ## Penalty for 1st degree difference (RR_test) + sum(stats::dnorm(x = penalty_rr_test, mean = 0, sd = prior_sd_rr_test, log = TRUE)) return(lprior + hyperprior) } - diff --git a/R/plot_functions.R b/R/plot_functions.R index e21656e..015ea38 100644 --- a/R/plot_functions.R +++ b/R/plot_functions.R @@ -20,61 +20,61 @@ get_pjnz_summary_data <- function(fp) { ## -- UPDATE HERE -- ## * update yr_pred to current year plot_pjnz_prv <- function(pjnz_summary, yr_pred = 2022) { - pjnz_summary <- na.omit(data.frame(year = pjnz_summary[["year"]], prv = pjnz_summary[["prevalence"]]*100)) + pjnz_summary <- stats::na.omit(data.frame(year = pjnz_summary[["year"]], prv = pjnz_summary[["prevalence"]]*100)) pjnz_summary$year <- pjnz_summary$year + 0.5 - plot(pjnz_summary$prv ~ pjnz_summary$year, + plot(pjnz_summary$prv ~ pjnz_summary$year, type = 'l', lwd = 2, xlim = c(2000, yr_pred), ylim = c(0, max(pjnz_summary$prv) * 1.2), main = 'HIV Prevalence (15-49 years)', ylab = 'HIV Prevalence (%)', xlab = 'Year', col='grey40') - polygon(x = c(pjnz_summary$year, rev(pjnz_summary$year)), + graphics::polygon(x = c(pjnz_summary$year, rev(pjnz_summary$year)), y = c(rep(0, length(pjnz_summary$prv)), rev(pjnz_summary$prv)), - col = rgb(155, 0, 50, 150, max = 255), border = NA) + col = grDevices::rgb(155, 0, 50, 150, max = 255), border = NA) } #' @export ## -- UPDATE HERE -- ## * update yr_pred to current year plot_pjnz_inc <- function(pjnz_summary, yr_pred = 2022) { - pjnz_summary <- na.omit(data.frame(year = pjnz_summary[["year"]], inc = pjnz_summary[["incidence"]]*1000)) + pjnz_summary <- stats::na.omit(data.frame(year = pjnz_summary[["year"]], inc = pjnz_summary[["incidence"]]*1000)) pjnz_summary$year <- pjnz_summary$year + 0.5 pjnz_summary <- subset(pjnz_summary, year >= 2000) - plot(pjnz_summary$inc ~ pjnz_summary$year, + plot(pjnz_summary$inc ~ pjnz_summary$year, type = 'l', lwd = 2, xlim = c(2000, yr_pred), ylim = c(0, max(pjnz_summary$inc) * 1.1), main='HIV Incidence (15-49 years)', ylab = 'HIV Incidence (per 1000)', xlab='Year', col='grey40') - polygon(x = c(pjnz_summary$year, rev(pjnz_summary$year)), + graphics::polygon(x = c(pjnz_summary$year, rev(pjnz_summary$year)), y = c(rep(0, length(pjnz_summary$inc)), rev(pjnz_summary$inc)), - col = rgb(255, 155, 0, 200, max = 255), border = NA) + col = grDevices::rgb(255, 155, 0, 200, max = 255), border = NA) } #' @export ## -- UPDATE HERE -- ## * update yr_pred to current year plot_pjnz_pop <- function(pjnz_summary, yr_pred = 2022) { - pjnz_summary <- na.omit(data.frame(year = pjnz_summary[["year"]], pop = pjnz_summary[["pop"]]/1000)) + pjnz_summary <- stats::na.omit(data.frame(year = pjnz_summary[["year"]], pop = pjnz_summary[["pop"]]/1000)) pjnz_summary$year <- pjnz_summary$year + 0.5 - plot(pjnz_summary$pop ~ pjnz_summary$year, + plot(pjnz_summary$pop ~ pjnz_summary$year, type = 'l', lwd = 2, xlim = c(2000, yr_pred), ylim = c(0, max(pjnz_summary$pop) * 1.2), main = 'Total Population (15-49 years)', ylab = 'Population (in 1,000s)', xlab='Year', col='grey40') - polygon(x = c(pjnz_summary$year, rev(pjnz_summary$year)), + graphics::polygon(x = c(pjnz_summary$year, rev(pjnz_summary$year)), y = c(rep(0, length(pjnz_summary$pop)), rev(pjnz_summary$pop)), - col = rgb(90, 170, 240, 150, max = 255), border = NA) + col = grDevices::rgb(90, 170, 240, 150, max = 255), border = NA) } #' @export ## -- UPDATE HERE -- ## * update yr_pred to current year plot_pjnz_plhiv <- function(pjnz_summary, yr_pred = 2022) { - pjnz_summary <- na.omit(data.frame(year = pjnz_summary[["year"]], plhiv = pjnz_summary[["plhiv"]]/1000)) + pjnz_summary <- stats::na.omit(data.frame(year = pjnz_summary[["year"]], plhiv = pjnz_summary[["plhiv"]]/1000)) pjnz_summary$year <- pjnz_summary$year + 0.5 - plot(pjnz_summary$plhiv ~ pjnz_summary$year, + plot(pjnz_summary$plhiv ~ pjnz_summary$year, type = 'l', lwd = 2, xlim = c(2000, yr_pred), ylim = c(0, max(pjnz_summary$plhiv) * 1.2), main = 'Total Number of PLHIV (15-49 years)', ylab = 'Population of PLHIV (in 1,000s)', xlab='Year', col='grey40') - polygon(x = c(pjnz_summary$year, rev(pjnz_summary$year)), + graphics::polygon(x = c(pjnz_summary$year, rev(pjnz_summary$year)), y = c(rep(0, length(pjnz_summary$plhiv)), rev(pjnz_summary$plhiv)), - col = rgb(250, 50, 10, 150, max = 255), border = NA) + col = grDevices::rgb(250, 50, 10, 150, max = 255), border = NA) } @@ -85,7 +85,7 @@ plot_pjnz_plhiv <- function(pjnz_summary, yr_pred = 2022) { ## * update yr_pred to current year plot_pjnz <- function(fp, yr_pred = 2022) { summary <- get_pjnz_summary_data(fp) - par(mfrow=c(2,2)) + graphics::par(mfrow=c(2,2)) plot_pjnz_prv(summary, yr_pred) plot_pjnz_inc(summary, yr_pred) plot_pjnz_pop(summary, yr_pred) @@ -142,15 +142,15 @@ plot_input_tot <- function(prgdat, fp, yr_pred = 2022) { xlim = c(2005, yr_pred), lwd = 2, col = 'steelblue3', xlab = 'Year', ylab = 'Total number of tests (in 1,000s)', pch = 16, main = expression(bold(paste("Total ", N^o, " of Tests Performed (VCT & ANC)")))) - lines(I(prg_t$tot/1000) ~ prg_t$year, lwd = 2, col = 'steelblue4') + graphics::lines(I(prg_t$tot/1000) ~ prg_t$year, lwd = 2, col = 'steelblue4') p_test_pop <- round(tot_prg$tot/pop[(tot_prg$year - start)]*100, 0) p_test_pop_r <- c(min(p_test_pop, na.rm = TRUE), max(p_test_pop, na.rm = TRUE)) - text(paste('No Test/Pop = ', p_test_pop_r[1], '-', p_test_pop_r[2], '%'), + graphics::text(paste('No Test/Pop = ', p_test_pop_r[1], '-', p_test_pop_r[2], '%'), x = 2005, y = max(tot_prg$tot/1000, na.rm = TRUE)/20, pos = 4) - lines(I(prg_f$tot/1000) ~ prg_f$year, lty = 2, lwd = 2, col = 'steelblue4') - points(I(prg_f$tot/1000) ~ prg_f$year, pch = 'f') - lines(I(prg_m$tot/1000) ~ prg_m$year, lty = 2, lwd = 2, col = 'steelblue4') - points(I(prg_m$tot/1000) ~ prg_m$year, pch = 'm') + graphics::lines(I(prg_f$tot/1000) ~ prg_f$year, lty = 2, lwd = 2, col = 'steelblue4') + graphics::points(I(prg_f$tot/1000) ~ prg_f$year, pch = 'f') + graphics::lines(I(prg_m$tot/1000) ~ prg_m$year, lty = 2, lwd = 2, col = 'steelblue4') + graphics::points(I(prg_m$tot/1000) ~ prg_m$year, pch = 'm') } } @@ -179,21 +179,21 @@ plot_input_totpos <- function(prgdat, fp, yr_pred = 2022) { ylab = 'Total of HIV+ tests (in 1,000s)', main = expression(bold(paste("Total ", N^o, " of HIV+ Tests (VCT & ANC)")))) - lines(I(prg_t$totpos/1000) ~ prg_t$year, lwd = 2, col = 'violetred4') + graphics::lines(I(prg_t$totpos/1000) ~ prg_t$year, lwd = 2, col = 'violetred4') yield <- round(tot_prg$totpos/tot_prg$tot * 100, 1) y_range <- c(min(yield, na.rm = TRUE), max(yield, na.rm = TRUE)) p_pos_pop <- round(tot_prg$totpos/plhiv[(tot_prg$year - start)]*100, 0) p_pos_pop_r <- c(min(p_pos_pop, na.rm = TRUE), max(p_pos_pop, na.rm = TRUE)) - text(paste('Positivity = ', y_range[1], '-', y_range[2], '%', sep = ''), + graphics::text(paste('Positivity = ', y_range[1], '-', y_range[2], '%', sep = ''), x = 2005, y = max(tot_prg$totpos/1000, na.rm = TRUE)/6, pos = 4) - text(paste('No Positive/PLHIV = ', p_pos_pop_r[1], '-', + graphics::text(paste('No Positive/PLHIV = ', p_pos_pop_r[1], '-', p_pos_pop_r[2], '%', sep=''), x = 2005, y = max(tot_prg$totpos/1000, na.rm = TRUE)/20, pos = 4) - lines(I(prg_f$totpos/1000) ~ prg_f$year, lty = 2, lwd = 2, col = 'violetred4') - points(I(prg_f$totpos/1000) ~ prg_f$year, pch = 'f') - lines(I(prg_m$totpos/1000) ~ prg_m$year, lty = 2, lwd = 2, col = 'violetred4') - points(I(prg_m$totpos/1000) ~ prg_m$year, pch = 'm') + graphics::lines(I(prg_f$totpos/1000) ~ prg_f$year, lty = 2, lwd = 2, col = 'violetred4') + graphics::points(I(prg_f$totpos/1000) ~ prg_f$year, pch = 'f') + graphics::lines(I(prg_m$totpos/1000) ~ prg_m$year, lty = 2, lwd = 2, col = 'violetred4') + graphics::points(I(prg_m$totpos/1000) ~ prg_m$year, pch = 'm') } } } @@ -210,7 +210,7 @@ plot_input_anctot <- function(prgdat, fp, yr_pred = 2022) { lwd = 2, col = 'darkorange2', xlab = 'Year', ylab = 'Total number of ANC tests (in 1,000s)', pch=16, main = expression(bold(paste("Total ", N^o, " of HIV Tests Performed at ANC")))) - lines(I(prgdat$anc/1000) ~ prgdat$year, lwd = 2, col = 'darkorange3') + graphics::lines(I(prgdat$anc/1000) ~ prgdat$year, lwd = 2, col = 'darkorange3') } } @@ -221,12 +221,12 @@ plot_input_ancpos <- function(prgdat, fp, yr_pred = 2022) { if (sum(prgdat$ancpos, na.rm = TRUE) > 0) { prgdat <- subset(prgdat, sex != 'male') - plot(I(prgdat$ancpos/1000) ~ prgdat$year, + plot(I(prgdat$ancpos/1000) ~ prgdat$year, ylim = c(0, max(prgdat$ancpos/1000, na.rm = TRUE)*1.1), xlim = c(2005, yr_pred ), - lwd = 2, col = 'deeppink2', xlab='Year', + lwd = 2, col = 'deeppink2', xlab='Year', ylab = 'Total number of positive ANC tests (in 1,000s)', pch = 17, main = expression(bold(paste("Total ", N^o, " of HIV+ Tests at ANC")))) - lines(I(prgdat$ancpos/1000) ~ prgdat$year, lwd = 2, col = 'deeppink2') + graphics::lines(I(prgdat$ancpos/1000) ~ prgdat$year, lwd = 2, col = 'deeppink2') } } @@ -236,42 +236,42 @@ plot_input_ancpos <- function(prgdat, fp, yr_pred = 2022) { ## -- UPDATE HERE -- ## * update yr_pred to current year plot_inputdata <- function(prgm_dat, fp, yr_pred = 2022) { - par(mfrow = c(2,2)) + graphics::par(mfrow = c(2,2)) plot_input_tot(prgm_dat, fp, yr_pred) plot_input_totpos(prgm_dat, fp, yr_pred) plot_input_anctot(prgm_dat, fp, yr_pred) - plot_input_ancpos(prgm_dat, fp, yr_pred) + plot_input_ancpos(prgm_dat, fp, yr_pred) } plot_prior <- function(plotprior = TRUE, sd = c(1.625, 1.85, 1.31, 1.68)){ if (plotprior == TRUE){ x <- logit(seq(0.00001, 1, 0.001)) - par(mfrow = c(1,4)) - rr_test <- dnorm(x, mean = logit(0.7/7.2), sd = sd[1]) - rr_plhiv <- dnorm(x, mean = logit(0.5), sd = sd[2]) - rr_dxunt <- dnorm(x, mean = logit(1.5/8), sd = sd[3]) - rr_dxart <- dnorm(x, mean = logit(0.25), sd = sd[4]) - plot(rr_test ~ I(0.8 + plogis(x) * 7.2), type = 'l', xlab = 'RR re-testing', ylab = 'Density', xlim = c(0,8)) - abline(v = 0.8 + 0.7/7.2 * 7.2, lty = 2, col = 'grey50') - lci <- round(0.8 + plogis(logit(0.7/7.2) - qnorm(0.975) * sd[1]) * 7.2, 1) - uci <- round(0.8 + plogis(logit(0.7/7.2) + qnorm(0.975) * sd[1]) * 7.2, 1) - text(x = 8, y = max(rr_test), paste('Prior: ', 1.5, '(95%: ', lci, '-', uci, ')'), pos = 2) - plot(rr_plhiv ~ I(0.05 + plogis(x) * 1.9), type = 'l', xlab = 'RR PLHIV', ylab = 'Density', xlim = c(0,2)) - abline(v = 0.05 + 0.5 * 1.9, lty = 2, col = 'grey50') - lci <- round(0.05 + plogis(logit(0.5) - qnorm(0.975) * sd[2]) * 1.9, 1) - uci <- round(0.05 + plogis(logit(0.5) + qnorm(0.975) * sd[2]) * 1.9, 1) - text(x = 2, y = max(rr_plhiv), paste('Prior: ', 1.0, '(95%: ', lci, '-', uci, ')'), pos = 2) - plot(rr_dxunt ~ I(plogis(x) * 8), type = 'l', xlab = 'RR re-testing for untreated diagnosed', ylab = 'Density', xlim = c(0,8)) - abline(v = 1.5/8 * 8, lty = 2, col = 'grey50') - lci <- round(plogis(logit(1.5/8) - qnorm(0.975) * sd[3]) * 8, 1) - uci <- round(plogis(logit(1.5/8) + qnorm(0.975) * sd[3]) * 8, 1) - text(x = 8, y = max(rr_test), paste('Prior: ', 1, '(95%: ', lci, '-', uci, ')'), pos = 2) - plot(rr_dxart ~ I(plogis(x)), type = 'l', xlab = 'RR re-testing on ART', ylab = 'Density', xlim = c(0,1)) - abline(v = 0.25, lty = 2, col = 'grey50') - lci <- round(plogis(logit(0.25) - qnorm(0.975) * sd[4]), 1) - uci <- round(plogis(logit(0.25) + qnorm(0.975) * sd[4]), 1) - text(x = 8, y = max(rr_test), paste('Prior: ', 0.25, '(95%: ', lci, '-', uci, ')'), pos = 2) + graphics::par(mfrow = c(1,4)) + rr_test <- stats::dnorm(x, mean = logit(0.7/7.2), sd = sd[1]) + rr_plhiv <- stats::dnorm(x, mean = logit(0.5), sd = sd[2]) + rr_dxunt <- stats::dnorm(x, mean = logit(1.5/8), sd = sd[3]) + rr_dxart <- stats::dnorm(x, mean = logit(0.25), sd = sd[4]) + plot(rr_test ~ I(0.8 + stats::plogis(x) * 7.2), type = 'l', xlab = 'RR re-testing', ylab = 'Density', xlim = c(0,8)) + graphics::abline(v = 0.8 + 0.7/7.2 * 7.2, lty = 2, col = 'grey50') + lci <- round(0.8 + stats::plogis(logit(0.7/7.2) - stats::qnorm(0.975) * sd[1]) * 7.2, 1) + uci <- round(0.8 + stats::plogis(logit(0.7/7.2) + stats::qnorm(0.975) * sd[1]) * 7.2, 1) + graphics::text(x = 8, y = max(rr_test), paste('Prior: ', 1.5, '(95%: ', lci, '-', uci, ')'), pos = 2) + plot(rr_plhiv ~ I(0.05 + stats::plogis(x) * 1.9), type = 'l', xlab = 'RR PLHIV', ylab = 'Density', xlim = c(0,2)) + graphics::abline(v = 0.05 + 0.5 * 1.9, lty = 2, col = 'grey50') + lci <- round(0.05 + stats::plogis(logit(0.5) - stats::qnorm(0.975) * sd[2]) * 1.9, 1) + uci <- round(0.05 + stats::plogis(logit(0.5) + stats::qnorm(0.975) * sd[2]) * 1.9, 1) + graphics::text(x = 2, y = max(rr_plhiv), paste('Prior: ', 1.0, '(95%: ', lci, '-', uci, ')'), pos = 2) + plot(rr_dxunt ~ I(stats::plogis(x) * 8), type = 'l', xlab = 'RR re-testing for untreated diagnosed', ylab = 'Density', xlim = c(0,8)) + graphics::abline(v = 1.5/8 * 8, lty = 2, col = 'grey50') + lci <- round(stats::plogis(logit(1.5/8) - stats::qnorm(0.975) * sd[3]) * 8, 1) + uci <- round(stats::plogis(logit(1.5/8) + stats::qnorm(0.975) * sd[3]) * 8, 1) + graphics::text(x = 8, y = max(rr_test), paste('Prior: ', 1, '(95%: ', lci, '-', uci, ')'), pos = 2) + plot(rr_dxart ~ I(stats::plogis(x)), type = 'l', xlab = 'RR re-testing on ART', ylab = 'Density', xlim = c(0,1)) + graphics::abline(v = 0.25, lty = 2, col = 'grey50') + lci <- round(stats::plogis(logit(0.25) - stats::qnorm(0.975) * sd[4]), 1) + uci <- round(stats::plogis(logit(0.25) + stats::qnorm(0.975) * sd[4]), 1) + graphics::text(x = 8, y = max(rr_test), paste('Prior: ', 0.25, '(95%: ', lci, '-', uci, ')'), pos = 2) } } @@ -279,19 +279,19 @@ plot_oi <- function(plotprior = TRUE, par = c(logit(0.5), 0.56)) { if (plotprior == TRUE){ r <- 0.5 tdx <- seq(1995, 2017, by = 0.1) - date_oidx <- 2005 + plogis(rnorm(2500, mean = par[1], sd = par[2])) * (2015 - 2005) - plot(I(rnorm(2)) ~ 1, pch = '', xlab = 'Year', ylab = '% OI Diagnosed', + date_oidx <- 2005 + stats::plogis(stats::rnorm(2500, mean = par[1], sd = par[2])) * (2015 - 2005) + plot(I(stats::rnorm(2)) ~ 1, pch = '', xlab = 'Year', ylab = '% OI Diagnosed', ylim = c(0, 1), xlim = c(1995, 2017)) for (i in 1:length(date_oidx)) { pr_oi <- 0.95 / (1 + exp(-r * (tdx - date_oidx[i]))) - lines(pr_oi ~ tdx, col = rgb(150, 100, 200, 25, max = 255), lwd = 1) + graphics::lines(pr_oi ~ tdx, col = grDevices::rgb(150, 100, 200, 25, max = 255), lwd = 1) } pr_oimin <- 0.95 / (1 + exp(-r * (tdx - 2005))) pr_oimax <- 0.95 / (1 + exp(-r * (tdx - 2015))) - pr_median <- 0.95 / (1 + exp(-r * (tdx - median(date_oidx)))) - lines(pr_oimin ~ tdx, col = rgb(100, 100, 100, 150, max = 255), lwd = 1, lty = 2) - lines(pr_oimax ~ tdx, col = rgb(100, 100, 100, 150, max = 255), lwd = 1, lty = 2) - lines(pr_median ~ tdx, col = 'steelblue4', lwd = 2, lty = 1) + pr_median <- 0.95 / (1 + exp(-r * (tdx - stats::median(date_oidx)))) + graphics::lines(pr_oimin ~ tdx, col = grDevices::rgb(100, 100, 100, 150, max = 255), lwd = 1, lty = 2) + graphics::lines(pr_oimax ~ tdx, col = grDevices::rgb(100, 100, 100, 150, max = 255), lwd = 1, lty = 2) + graphics::lines(pr_median ~ tdx, col = 'steelblue4', lwd = 2, lty = 1) } } @@ -299,16 +299,15 @@ plot_oi <- function(plotprior = TRUE, par = c(logit(0.5), 0.56)) { optimized_par <- function(opt, param = NULL) { n_k2 <- length(opt$par) - 11 n_k1 <- n_k2 - 10 - require(Matrix) - if (is.null(opt$hessian)) { - print('You forgot to specify Hessian = TRUE; 95%CrI not calculated.') - se <- rep(NA, length(opt$par)) + if (is.null(opt$hessian)) { + print('You forgot to specify Hessian = TRUE; 95%CrI not calculated.') + se <- rep(NA, length(opt$par)) } else { # We first test is system is singular (TEMPORARY FIX) if (any(diag(opt$hessian) == 0)) { index <- which(diag(opt$hessian) == 0) opt$hessian[index, index] <- opt$hessian[index, index] - 1e-3 - } + } # From the hessian, simulate the model vcova <- solve(-opt$hessian) # Test if positive semi-definite @@ -317,95 +316,90 @@ optimized_par <- function(opt, param = NULL) { if (!all(ev >= -1e-06 * abs(ev[1L]))) { vcova <- nearPD(vcova, corr = FALSE)$mat } se <- sqrt(diag(as.matrix(vcova))) } - + param_names <- c("RRm_05", "RRm_12", "RR_Test10","RR_Test15", - "RR_PLHIV", "RR_DxUnt10", "RR_DxUnt17", + "RR_PLHIV", "RR_DxUnt10", "RR_DxUnt17", "RR_DxART", "RR_25p_m", "RR_35p_m", "RR_25p_f", "RR_35p_f", "RR OI Dx") - - description <- c('RR testing: men in 2005', + + description <- c('RR testing: men in 2005', 'RR testing: men in 2012', - 'RR re-testing 2010', - 'RR re-testing 2015', + 'RR re-testing 2010', + 'RR re-testing 2015', 'RR testing: PLHIV unaware', 'RR re-testing: PLHIV aware (not ART) 2010', ## -- UPDATE HERE -- ## * update label to current year - 'RR re-testing: PLHIV aware (not ART) 2022', - 'RR re-testing: PLHIV on ART (*RR not ART)', + 'RR re-testing: PLHIV aware (not ART) 2022', + 'RR re-testing: PLHIV on ART (*RR not ART)', 'RR among 25-34 men', 'RR among 35+ men', 'RR among 25-34 women', - 'RR among 35+ women', + 'RR among 35+ women', 'RR OI Dx (ART Cov)') - - pt <- c(round(plogis(opt$par[n_k2 + 1]) * 1.1, 2), - round(plogis(opt$par[n_k2 + 2]) * 1.1, 2), - round(0.95 + plogis(opt$par[n_k2 + 3]) * 7.05, 2), - round(0.95 + plogis(opt$par[n_k2 + 4]) * 7.05, 2), - round(0.05 + plogis(opt$par[n_k2 + 5]) * (1.95 - 0.05), 2), - round(plogis(opt$par[n_k1 + 1]) * 8, 2), - round(plogis(opt$par[n_k2 - 1]) * 8, 2), - round(plogis(opt$par[n_k2 + 6]), 2), - round(0.1 + plogis(opt$par[n_k2 + 7]) * (6 - 0.1), 2), - round(0.1 + plogis(opt$par[n_k2 + 8]) * (6 - 0.1), 2), - round(0.1 + plogis(opt$par[n_k2 + 9]) * (6 - 0.1), 2), - round(0.1 + plogis(opt$par[n_k2 + 10]) * (6 - 0.1), 2), - round(0.25 + (plogis(opt$par[n_k2 + 11])) * (1.75 - 0.25), 1)) - + + pt <- c(round(stats::plogis(opt$par[n_k2 + 1]) * 1.1, 2), + round(stats::plogis(opt$par[n_k2 + 2]) * 1.1, 2), + round(0.95 + stats::plogis(opt$par[n_k2 + 3]) * 7.05, 2), + round(0.95 + stats::plogis(opt$par[n_k2 + 4]) * 7.05, 2), + round(0.05 + stats::plogis(opt$par[n_k2 + 5]) * (1.95 - 0.05), 2), + round(stats::plogis(opt$par[n_k1 + 1]) * 8, 2), + round(stats::plogis(opt$par[n_k2 - 1]) * 8, 2), + round(stats::plogis(opt$par[n_k2 + 6]), 2), + round(0.1 + stats::plogis(opt$par[n_k2 + 7]) * (6 - 0.1), 2), + round(0.1 + stats::plogis(opt$par[n_k2 + 8]) * (6 - 0.1), 2), + round(0.1 + stats::plogis(opt$par[n_k2 + 9]) * (6 - 0.1), 2), + round(0.1 + stats::plogis(opt$par[n_k2 + 10]) * (6 - 0.1), 2), + round(0.25 + (stats::plogis(opt$par[n_k2 + 11])) * (1.75 - 0.25), 1)) + if (is.null(param)) { - lci <- c(round(plogis(opt$par[n_k2 + 1] - qnorm(0.975) * se[n_k2 + 1]) * 1.1, 2), - round(plogis(opt$par[n_k2 + 2] - qnorm(0.975) * se[n_k2 + 2]) * 1.1, 2), - round(0.95 + plogis(opt$par[n_k2 + 3] - qnorm(0.975) * se[n_k2 + 3]) * 7.05, 2), - round(0.95 + plogis(opt$par[n_k2 + 4] - qnorm(0.975) * se[n_k2 + 4]) * 7.05, 2), - round(0.05 + plogis(opt$par[n_k2 + 5] - qnorm(0.975) * se[n_k2 + 5]) * (1.95 - 0.05), 2), - round(plogis(opt$par[n_k1 + 1] - qnorm(0.975) * se[n_k1 + 1]) * 8, 2), - round(plogis(opt$par[n_k2 - 1] - qnorm(0.975) * se[n_k2 - 1]) * 8, 2), - round(plogis(opt$par[n_k2 + 6] - qnorm(0.975) * se[n_k2 + 6]), 2), - round(0.1 + plogis(opt$par[n_k2 + 7] - qnorm(0.975) * se[n_k2 + 7]) * (6 - 0.1), 2), - round(0.1 + plogis(opt$par[n_k2 + 8] - qnorm(0.975) * se[n_k2 + 8]) * (6 - 0.1), 2), - round(0.1 + plogis(opt$par[n_k2 + 9] - qnorm(0.975) * se[n_k2 + 9]) * (6 - 0.1), 2), - round(0.1 + plogis(opt$par[n_k2 + 10] - qnorm(0.975) * se[n_k2 + 10]) * (6 - 0.1), 2), - round(0.25 + (plogis(opt$par[n_k2 + 11] - qnorm(0.975) * se[n_k2 + 11])) * (1.75 - 0.25), 1)) - - uci <- c(round(plogis(opt$par[n_k2 + 1] + qnorm(0.975) * se[n_k2 + 1]) * 1.1, 2), - round(plogis(opt$par[n_k2 + 2] + qnorm(0.975) * se[n_k2 + 2]) * 1.1, 2), - round(0.95 + plogis(opt$par[n_k2 + 3] + qnorm(0.975) * se[n_k2 + 3]) * 7.05, 2), - round(0.95 + plogis(opt$par[n_k2 + 4] + qnorm(0.975) * se[n_k2 + 4]) * 7.05, 2), - round(0.05 + plogis(opt$par[n_k2 + 5] + qnorm(0.975) * se[n_k2 + 5]) * (1.95 - 0.05), 2), - round(plogis(opt$par[n_k1 + 1] + qnorm(0.975) * se[n_k1 + 1]) * 8, 2), - round(plogis(opt$par[n_k2 - 1] + qnorm(0.975) * se[n_k2 - 1]) * 8, 2), - round(plogis(opt$par[n_k2 + 6] + qnorm(0.975) * se[n_k2 + 6]), 2), - round(0.1 + plogis(opt$par[n_k2 + 7] + qnorm(0.975) * se[n_k2 + 7]) * (6 - 0.1), 2), - round(0.1 + plogis(opt$par[n_k2 + 8] + qnorm(0.975) * se[n_k2 + 8]) * (6 - 0.1), 2), - round(0.1 + plogis(opt$par[n_k2 + 9] + qnorm(0.975) * se[n_k2 + 9]) * (6 - 0.1), 2), - round(0.1 + plogis(opt$par[n_k2 + 10] + qnorm(0.975) * se[n_k2 + 10]) * (6 - 0.1), 2), - round(0.25 + (plogis(opt$par[n_k2 + 11] + qnorm(0.975) * se[n_k2 + 11])) * (1.75 - 0.25), 1)) + lci <- c(round(stats::plogis(opt$par[n_k2 + 1] - stats::qnorm(0.975) * se[n_k2 + 1]) * 1.1, 2), + round(stats::plogis(opt$par[n_k2 + 2] - stats::qnorm(0.975) * se[n_k2 + 2]) * 1.1, 2), + round(0.95 + stats::plogis(opt$par[n_k2 + 3] - stats::qnorm(0.975) * se[n_k2 + 3]) * 7.05, 2), + round(0.95 + stats::plogis(opt$par[n_k2 + 4] - stats::qnorm(0.975) * se[n_k2 + 4]) * 7.05, 2), + round(0.05 + stats::plogis(opt$par[n_k2 + 5] - stats::qnorm(0.975) * se[n_k2 + 5]) * (1.95 - 0.05), 2), + round(stats::plogis(opt$par[n_k1 + 1] - stats::qnorm(0.975) * se[n_k1 + 1]) * 8, 2), + round(stats::plogis(opt$par[n_k2 - 1] - stats::qnorm(0.975) * se[n_k2 - 1]) * 8, 2), + round(stats::plogis(opt$par[n_k2 + 6] - stats::qnorm(0.975) * se[n_k2 + 6]), 2), + round(0.1 + stats::plogis(opt$par[n_k2 + 7] - stats::qnorm(0.975) * se[n_k2 + 7]) * (6 - 0.1), 2), + round(0.1 + stats::plogis(opt$par[n_k2 + 8] - stats::qnorm(0.975) * se[n_k2 + 8]) * (6 - 0.1), 2), + round(0.1 + stats::plogis(opt$par[n_k2 + 9] - stats::qnorm(0.975) * se[n_k2 + 9]) * (6 - 0.1), 2), + round(0.1 + stats::plogis(opt$par[n_k2 + 10] - stats::qnorm(0.975) * se[n_k2 + 10]) * (6 - 0.1), 2), + round(0.25 + (stats::plogis(opt$par[n_k2 + 11] - stats::qnorm(0.975) * se[n_k2 + 11])) * (1.75 - 0.25), 1)) + + uci <- c(round(stats::plogis(opt$par[n_k2 + 1] + stats::qnorm(0.975) * se[n_k2 + 1]) * 1.1, 2), + round(stats::plogis(opt$par[n_k2 + 2] + stats::qnorm(0.975) * se[n_k2 + 2]) * 1.1, 2), + round(0.95 + stats::plogis(opt$par[n_k2 + 3] + stats::qnorm(0.975) * se[n_k2 + 3]) * 7.05, 2), + round(0.95 + stats::plogis(opt$par[n_k2 + 4] + stats::qnorm(0.975) * se[n_k2 + 4]) * 7.05, 2), + round(0.05 + stats::plogis(opt$par[n_k2 + 5] + stats::qnorm(0.975) * se[n_k2 + 5]) * (1.95 - 0.05), 2), + round(stats::plogis(opt$par[n_k1 + 1] + stats::qnorm(0.975) * se[n_k1 + 1]) * 8, 2), + round(stats::plogis(opt$par[n_k2 - 1] + stats::qnorm(0.975) * se[n_k2 - 1]) * 8, 2), + round(stats::plogis(opt$par[n_k2 + 6] + stats::qnorm(0.975) * se[n_k2 + 6]), 2), + round(0.1 + stats::plogis(opt$par[n_k2 + 7] + stats::qnorm(0.975) * se[n_k2 + 7]) * (6 - 0.1), 2), + round(0.1 + stats::plogis(opt$par[n_k2 + 8] + stats::qnorm(0.975) * se[n_k2 + 8]) * (6 - 0.1), 2), + round(0.1 + stats::plogis(opt$par[n_k2 + 9] + stats::qnorm(0.975) * se[n_k2 + 9]) * (6 - 0.1), 2), + round(0.1 + stats::plogis(opt$par[n_k2 + 10] + stats::qnorm(0.975) * se[n_k2 + 10]) * (6 - 0.1), 2), + round(0.25 + (stats::plogis(opt$par[n_k2 + 11] + stats::qnorm(0.975) * se[n_k2 + 11])) * (1.75 - 0.25), 1)) } if (!is.null(param)){ - est <- rbind(plogis(param[, n_k2 + 1]) * 1.1, - plogis(param[, n_k2 + 2]) * 1.1, - 0.95 + plogis(param[, n_k2 + 3]) * 7.05, - 0.95 + plogis(param[, n_k2 + 4]) * 7.05, - 0.05 + plogis(param[, n_k2 + 5]) * (1.95 - 0.05), - plogis(param[, n_k1 + 1]) * 8, - plogis(param[, n_k2 - 1]) * 8, - plogis(param[, n_k2 + 6]), - 0.1 + plogis(param[, n_k2 + 7]) * (6 - 0.1), - 0.1 + plogis(param[, n_k2 + 8]) * (6 - 0.1), - 0.1 + plogis(param[, n_k2 + 9]) * (6 - 0.1), - 0.1 + plogis(param[, n_k2 + 10]) * (6 - 0.1), - 0.25 + (plogis(param[, n_k2 + 11])) * (1.75 - 0.25)) - lci <- round(apply(est, 1, quantile, probs = 0.025), 2) - uci <- round(apply(est, 1, quantile, probs = 0.975), 2) - + est <- rbind(stats::plogis(param[, n_k2 + 1]) * 1.1, + stats::plogis(param[, n_k2 + 2]) * 1.1, + 0.95 + stats::plogis(param[, n_k2 + 3]) * 7.05, + 0.95 + stats::plogis(param[, n_k2 + 4]) * 7.05, + 0.05 + stats::plogis(param[, n_k2 + 5]) * (1.95 - 0.05), + stats::plogis(param[, n_k1 + 1]) * 8, + stats::plogis(param[, n_k2 - 1]) * 8, + stats::plogis(param[, n_k2 + 6]), + 0.1 + stats::plogis(param[, n_k2 + 7]) * (6 - 0.1), + 0.1 + stats::plogis(param[, n_k2 + 8]) * (6 - 0.1), + 0.1 + stats::plogis(param[, n_k2 + 9]) * (6 - 0.1), + 0.1 + stats::plogis(param[, n_k2 + 10]) * (6 - 0.1), + 0.25 + (stats::plogis(param[, n_k2 + 11])) * (1.75 - 0.25)) + lci <- round(apply(est, 1, stats::quantile, probs = 0.025), 2) + uci <- round(apply(est, 1, stats::quantile, probs = 0.975), 2) + } RR_opt <- data.frame(Parameter_Name = description, Estimate = pt, LCI = lci, UCI = uci) - + return(RR_opt) } - - - - - diff --git a/R/plot_outputs.R b/R/plot_outputs.R index 43c1627..bc3331b 100644 --- a/R/plot_outputs.R +++ b/R/plot_outputs.R @@ -3,7 +3,7 @@ get_out_evertest <- function(mod, fp, agegr = c("15-24", '25-34','35-49', '15-49', '15-64', '15+'), sex = c("both", "female", "male"), hivstatus = c("all", "negative", "positive") ) { - + end_date <- fp$ss$proj_start + fp$ss$PROJ_YEARS - 1L out_evertest <- expand.grid(year = 2000:end_date, outcome = "evertest", @@ -12,13 +12,13 @@ get_out_evertest <- function(mod, fp, hivstatus = hivstatus) out_evertest$value <- evertest(mod, fp, add_ss_indices(out_evertest, fp$ss)) - + out_evertest } #' @export get_out_nbtest <- function(mod, fp) { - + end_date <- fp$ss$proj_start + fp$ss$PROJ_YEARS - 1L out_nbtest <- expand.grid(year = 2000:end_date, outcome = "numbertests", @@ -32,7 +32,7 @@ get_out_nbtest <- function(mod, fp) { #' @export get_out_nbtest_sex <- function(mod, fp) { - + end_date <- fp$ss$proj_start + fp$ss$PROJ_YEARS - 1L out_nbtest <- expand.grid(year = 2000:end_date, outcome = "numbertests", @@ -54,7 +54,7 @@ get_out_nbtest_pos <- function(mod, fp) { hivstatus = 'positive') out_nbtest_pos$hivstatus <- as.character(out_nbtest_pos$hivstatus) out_nbtest_pos$value <- total_tests(mod, add_ss_indices(out_nbtest_pos, fp$ss)) - + out_nbtest_pos } @@ -68,7 +68,7 @@ get_out_nbtest_pos_sex <- function(mod, fp) { hivstatus = 'positive') out_nbtest_pos$hivstatus <- as.character(out_nbtest_pos$hivstatus) out_nbtest_pos$value <- total_tests(mod, add_ss_indices(out_nbtest_pos, fp$ss)) - + out_nbtest_pos } @@ -82,7 +82,7 @@ get_out_aware <- function(mod, fp, agegr = c('15-49', '15-64', "15+"), sex = sex, hivstatus = 'positive') out_aware$value <- diagnosed(mod, fp, add_ss_indices(out_aware, fp$ss)) - + out_aware } @@ -96,7 +96,7 @@ get_out_nbaware <- function(mod, fp, agegr = '15-49', sex = sex, hivstatus = 'positive') out_nbaware$value <- nb_aware(mod, fp, add_ss_indices(out_nbaware, fp$ss)) - + out_nbaware } @@ -107,7 +107,7 @@ get_out_art <- function(mod, fp, out_artcov <- data.frame(year = 1970:end_date, outcome = "artcov", agegr = "15-49", sex = gender, hivstatus = "positive", value = artcov15to49(mod, sex = gender)) out_artcov <- subset(out_artcov) - + out_artcov } @@ -122,18 +122,18 @@ get_out_pregprev <- function(mod, fp) { undiagn_w <- colSums((attr(mod, "hivpop") - attr(mod, "diagnpop"))[, fp$ss$h.fert.idx, fp$ss$f.idx, ]) diagn_w <- colSums(attr(mod, "diagnpop")[, fp$ss$h.fert.idx, fp$ss$f.idx, ]) artpop_w <- colSums(attr(mod, "artpop")[, , fp$ss$h.fert.idx, fp$ss$f.idx, ],,2) - + tot_w <- (hivneg_w + undiagn_w + diagn_w + artpop_w) - + ## Number of births by HIV, awareness, and ART status hivneg_births <- colSums(births_ha * hivneg_w / tot_w) undiagn_births <- colSums(births_ha * undiagn_w / tot_w) diagn_births <- colSums(births_ha * diagn_w / tot_w) artpop_births <- colSums(births_ha * artpop_w / tot_w) - + births <- colSums(births_ha) hivp_births <- births - hivneg_births - + data.frame(year = fp$ss$proj_start - 1 + seq_len(fp$ss$PROJ_YEARS), prev = hivp_births / births, aware = (diagn_births + artpop_births) / hivp_births, @@ -143,7 +143,7 @@ get_out_pregprev <- function(mod, fp) { # ---- Individuals functions to plot model outputs ---- #' @export -## +## ## -- UPDATE HERE -- ## * update yr_pred to current year plot_out_nbtest <- function(mod, fp, likdat, cnt, simul = NULL, yr_pred = 2022, @@ -151,30 +151,30 @@ plot_out_nbtest <- function(mod, fp, likdat, cnt, simul = NULL, yr_pred = 2022, # if fitting with HTS program data stratified by sex, we add both sex back ld <- likdat$hts if (!is.null(ld)) { - likdat$hts <- aggregate(ld$tot, by = list(ld$year), FUN = sum) - colnames(likdat$hts) <- c('year','tot') + likdat$hts <- stats::aggregate(ld$tot, by = list(ld$year), FUN = sum) + colnames(likdat$hts) <- c('year','tot') } # redact <- c('Namibia','Uganda','Zambia','Zimbabwe') redact <- c('XXX') - + cnt_to_plot <- ifelse(cnt == "Cote d'Ivoire", "Côte d'Ivoire", - ifelse(cnt == "Swaziland", "eSwatini", cnt)) + ifelse(cnt == "Swaziland", "eSwatini", cnt)) start <- fp$ss$proj_start mod <- simmod(fp) - plhiv <- apply(attr(mod, "hivpop")[,1:8,,], 4, FUN = sum) + + plhiv <- apply(attr(mod, "hivpop")[,1:8,,], 4, FUN = sum) + apply(attr(mod, "artpop")[,,1:8,,], 5, FUN = sum) pop <- apply(mod[1:35,,,], 4, FUN = sum) out_nbtest <- subset(get_out_nbtest(mod, fp), year <= yr_pred) out_nbtest$year <- out_nbtest$year + 0.5 - + if (!is.null(ld) & any(ld$sex != 'both')) { ld_sex <- subset(ld, sex != 'both') ld_sex_pred <- total_tests(mod, df = add_ss_indices(ld_sex, fp$ss)) } - + # Decide if we plot CI or not if (!is.null(simul)) { ci <- getCI(simul$number.test) @@ -186,76 +186,76 @@ plot_out_nbtest <- function(mod, fp, likdat, cnt, simul = NULL, yr_pred = 2022, } else { plot.ci <- FALSE } - + if (plot_title == TRUE) { main_title <- expression(bold(paste("Total ", N^o, " of Tests"))) } else { main_title <- "" } - - col_ci <- rgb(150, 150, 150, 125, max = 255) - + + col_ci <- grDevices::rgb(150, 150, 150, 125, max = 255) + # Total tested if (plot.ci == T){ plot(I(out_nbtest$value / 1000) ~ out_nbtest$year, pch = '', ylim = c(0, max(out_nbtest$value, na.rm = TRUE) / 1000 * 1.25), - cex = 1, ylab = expression(paste(N^o, " of Tests (in 1,000)")), - xlab = 'Year', xlim = c(2000, yr_pred), + cex = 1, ylab = expression(paste(N^o, " of Tests (in 1,000)")), + xlab = 'Year', xlim = c(2000, yr_pred), main = main_title) - polygon(x = c(ci$year, rev(ci$year)), + graphics::polygon(x = c(ci$year, rev(ci$year)), y = c(ci$upper / 1000, rev(ci$lower / 1000)), col = col_ci, border = NA) - lines(I(out_nbtest$value / 1000) ~ out_nbtest$year, lwd = 1, col = 'seagreen3') - text(x = 2000, y = max(out_nbtest$value, na.rm = TRUE) / 1000 * 1.15, cnt_to_plot, cex = 1.25, pos = 4) - + graphics::lines(I(out_nbtest$value / 1000) ~ out_nbtest$year, lwd = 1, col = 'seagreen3') + graphics::text(x = 2000, y = max(out_nbtest$value, na.rm = TRUE) / 1000 * 1.15, cnt_to_plot, cex = 1.25, pos = 4) + if(length(likdat$hts$tot) > 0) { if (cnt %in% redact) { pchpt <- '' } else { pchpt <- 16 } - points(I(likdat$hts$tot/1000) ~ I(likdat$hts$year + 0.5), pch = pchpt) + graphics::points(I(likdat$hts$tot/1000) ~ I(likdat$hts$year + 0.5), pch = pchpt) p_test_pop <- round(likdat$hts$tot / pop[likdat$hts$year - start] * 100, 0) p_test_pop_r <- c(min(p_test_pop, na.rm = TRUE), max(p_test_pop, na.rm = TRUE)) - text(paste('No Test/Pop = ', p_test_pop_r[1], '-', p_test_pop_r[2], '%'), + graphics::text(paste('No Test/Pop = ', p_test_pop_r[1], '-', p_test_pop_r[2], '%'), x = yr_pred, y = max(likdat$hts$tot / 1000, na.rm = TRUE)/25, pos = 2) } } else { plot(I(out_nbtest$value/1000) ~ out_nbtest$year, pch = '', ylim = c(0, max(out_nbtest$value, na.rm = TRUE) / 1000 + max(out_nbtest$value, na.rm = TRUE)/1000*0.25), - cex = 1, ylab = expression(paste(N^o, " of Tests (in 1,000)")), - xlab = 'Year', xlim = c(2000, yr_pred), + cex = 1, ylab = expression(paste(N^o, " of Tests (in 1,000)")), + xlab = 'Year', xlim = c(2000, yr_pred), main = main_title) - lines(I(out_nbtest$value/1000) ~ out_nbtest$year, lwd = 1, col = 'seagreen3') - text(x = 2000, y = max(out_nbtest$value, na.rm = TRUE) / 1000, cnt_to_plot, cex = 1.5, pos = 2) + graphics::lines(I(out_nbtest$value/1000) ~ out_nbtest$year, lwd = 1, col = 'seagreen3') + graphics::text(x = 2000, y = max(out_nbtest$value, na.rm = TRUE) / 1000, cnt_to_plot, cex = 1.5, pos = 2) if(length(likdat$hts$tot) > 0) { if (cnt %in% redact) { pchpt <- '' } else { pchpt <- 16 } - points(I(likdat$hts$tot / 1000) ~ I(likdat$hts$year + 0.5), pch = pchpt) + graphics::points(I(likdat$hts$tot / 1000) ~ I(likdat$hts$year + 0.5), pch = pchpt) p_test_pop <- round(likdat$hts$tot/pop[(likdat$hts$year - start)] * 100, 0) p_test_pop_r <- c(min(p_test_pop, na.rm = TRUE), max(p_test_pop, na.rm = TRUE)) - text(paste('No Test/Pop = ', p_test_pop_r[1], '-', p_test_pop_r[2], '%'), + graphics::text(paste('No Test/Pop = ', p_test_pop_r[1], '-', p_test_pop_r[2], '%'), x = yr_pred, y = max(likdat$hts$tot / 1000, na.rm = TRUE) / 20, pos = 2) } } } -## +## ## -- UPDATE HERE -- ## * update yr_pred to current year plot_out_nbtest_sex <- function(mod, fp, likdat, cnt, simul = NULL, yr_pred = 2022) { - + # redact <- c('Namibia','Uganda','Zambia','Zimbabwe') redact <- c('XXX') - + cnt_to_plot <- ifelse(cnt == "Cote d'Ivoire", "Côte d'Ivoire", - ifelse(cnt == "Swaziland", "eSwatini", cnt)) + ifelse(cnt == "Swaziland", "eSwatini", cnt)) start <- fp$ss$proj_start mod <- simmod(fp) - plhiv <- apply(attr(mod, "hivpop")[,1:8,,], 4, FUN=sum) + + plhiv <- apply(attr(mod, "hivpop")[,1:8,,], 4, FUN=sum) + apply(attr(mod, "artpop")[,,1:8,,], 5, FUN=sum) pop <- apply(mod[1:35,,,], 4, FUN=sum) - + out_nbtest <- subset(get_out_nbtest_sex(mod, fp), year <= yr_pred) out_nbtest$year <- out_nbtest$year + 0.5 out_nbtest_f <- subset(out_nbtest, sex == 'female') out_nbtest_m <- subset(out_nbtest, sex == 'male') max_ylim <- max(c(out_nbtest_f$value, out_nbtest_m$value), na.rm = TRUE) - + lik_sex <- likdat$hts if (!is.null(lik_sex)){ lik_f <- subset(lik_sex, sex == 'female') @@ -274,51 +274,51 @@ plot_out_nbtest_sex <- function(mod, fp, likdat, cnt, simul = NULL, yr_pred = 20 } else { plot.ci <- FALSE } - - col_ci <- rgb(150, 150, 150, 125, max = 255) - + + col_ci <- grDevices::rgb(150, 150, 150, 125, max = 255) + # Total tested if (plot.ci == T){ plot(I(out_nbtest_f$value / 1000) ~ out_nbtest_f$year, pch = '', ylim = c(0, max_ylim / 1000 * 1.25), - cex = 1, ylab = expression(paste(N^o, " of Tests (in 1,000)")), - xlab = 'Year', xlim = c(2000, yr_pred), + cex = 1, ylab = expression(paste(N^o, " of Tests (in 1,000)")), + xlab = 'Year', xlim = c(2000, yr_pred), main = expression(bold(paste("Total ", N^o, " of Tests")))) - polygon(x = c(ci_f$year, rev(ci_f$year)), + graphics::polygon(x = c(ci_f$year, rev(ci_f$year)), y = c(ci_f$upper / 1000, rev(ci_f$lower / 1000)), col = col_ci, border = NA) - polygon(x = c(ci_m$year, rev(ci_m$year)), + graphics::polygon(x = c(ci_m$year, rev(ci_m$year)), y = c(ci_m$upper / 1000, rev(ci_m$lower / 1000)), col = col_ci, border = NA) - lines(I(out_nbtest_f$value / 1000) ~ out_nbtest_f$year, lwd = 1, col = 'seagreen3') - lines(I(out_nbtest_m$value / 1000) ~ out_nbtest_m$year, lwd = 1, col = 'seagreen3') - text(x = 2000, y = max(out_nbtest$value, na.rm = TRUE) / 1000 * 1.15, cnt_to_plot, cex = 1.25, pos = 4) - + graphics::lines(I(out_nbtest_f$value / 1000) ~ out_nbtest_f$year, lwd = 1, col = 'seagreen3') + graphics::lines(I(out_nbtest_m$value / 1000) ~ out_nbtest_m$year, lwd = 1, col = 'seagreen3') + graphics::text(x = 2000, y = max(out_nbtest$value, na.rm = TRUE) / 1000 * 1.15, cnt_to_plot, cex = 1.25, pos = 4) + if(length(lik_f$tot) > 0 | length(lik_m$tot) > 0) { if (cnt %in% redact) { pchpt <- '' } else { pchpt <- 16 } - points(I(lik_f$tot / 1000) ~ I(lik_f$year + 0.5), pch = pchpt) - points(I(lik_m$tot / 1000) ~ I(lik_m$year + 0.5), pch = pchpt + 1) + graphics::points(I(lik_f$tot / 1000) ~ I(lik_f$year + 0.5), pch = pchpt) + graphics::points(I(lik_m$tot / 1000) ~ I(lik_m$year + 0.5), pch = pchpt + 1) } } else { plot(I(out_nbtest_f$value / 1000) ~ out_nbtest_f$year, pch = '', ylim = c(0, max_ylim / 1000 * 1.25), - cex = 1, ylab = expression(paste(N^o, " of Tests (in 1,000)")), - xlab = 'Year', xlim = c(2000, yr_pred), + cex = 1, ylab = expression(paste(N^o, " of Tests (in 1,000)")), + xlab = 'Year', xlim = c(2000, yr_pred), main = expression(bold(paste("Total ", N^o, " of Tests")))) - lines(I(out_nbtest_f$value / 1000) ~ out_nbtest_f$year, lwd = 1, col = 'seagreen3') - lines(I(out_nbtest_m$value / 1000) ~ out_nbtest_m$year, lwd = 1, col = 'seagreen3') - text(x = 2000, y = max(out_nbtest$value, na.rm = TRUE) / 1000, cnt_to_plot, cex = 1.5, pos = 2) - + graphics::lines(I(out_nbtest_f$value / 1000) ~ out_nbtest_f$year, lwd = 1, col = 'seagreen3') + graphics::lines(I(out_nbtest_m$value / 1000) ~ out_nbtest_m$year, lwd = 1, col = 'seagreen3') + graphics::text(x = 2000, y = max(out_nbtest$value, na.rm = TRUE) / 1000, cnt_to_plot, cex = 1.5, pos = 2) + if(length(likdat$hts$tot) > 0) { if (cnt %in% redact) { pchpt <- '' } else { pchpt <- 16 } - points(I(lik_f$tot / 1000) ~ I(lik_f$year + 0.5), pch = pchpt) - points(I(lik_m$tot / 1000) ~ I(lik_m$year + 0.5), pch = pchpt + 1) } + graphics::points(I(lik_f$tot / 1000) ~ I(lik_f$year + 0.5), pch = pchpt) + graphics::points(I(lik_m$tot / 1000) ~ I(lik_m$year + 0.5), pch = pchpt + 1) } } } #' @export -## +## ## -- UPDATE HERE -- ## * update yr_pred to current year plot_out_nbpostest <- function(mod, fp, likdat, cnt, simul = NULL, yr_pred = 2022, @@ -326,29 +326,29 @@ plot_out_nbpostest <- function(mod, fp, likdat, cnt, simul = NULL, yr_pred = 202 # if fitting with HTS program data stratified by sex, we add both sex back ld <- likdat$hts if (!is.null(ld)) { - likdat$hts <- aggregate(ld$tot, by = list(ld$year), FUN = sum) + likdat$hts <- stats::aggregate(ld$tot, by = list(ld$year), FUN = sum) colnames(likdat$hts) <- c('year','tot') - likdat$hts$totpos <- aggregate(ld$totpos, by = list(ld$year), FUN = sum)$x + likdat$hts$totpos <- stats::aggregate(ld$totpos, by = list(ld$year), FUN = sum)$x } ldp <- likdat$hts_pos if (!is.null(ldp)) { - likdat$hts_pos <- aggregate(ldp$tot, by = list(ldp$year), FUN = sum) + likdat$hts_pos <- stats::aggregate(ldp$tot, by = list(ldp$year), FUN = sum) colnames(likdat$hts_pos) <- c('year','tot') } - + # redact <- c('Namibia','Uganda','Zambia','Zimbabwe') redact <- c('XXX') - + start <- fp$ss$proj_start mod <- simmod(fp) - plhiv <- apply(attr(mod, "hivpop")[, 1:8,,], 4, FUN = sum) + + plhiv <- apply(attr(mod, "hivpop")[, 1:8,,], 4, FUN = sum) + apply(attr(mod, "artpop")[,, 1:8,,], 5, FUN = sum) pop <- apply(mod[1:35,,,], 4, FUN = sum) out_nbtest_pos <- get_out_nbtest_pos(mod, fp) out_nbtest_pos <- subset(out_nbtest_pos, year <= yr_pred) - out_nbtest_pos$year <- out_nbtest_pos$year + 0.5 - + out_nbtest_pos$year <- out_nbtest_pos$year + 0.5 + # Decide if we plot CI or not if (!is.null(simul)){ ci <- getCI(simul$number.test.pos) @@ -366,77 +366,77 @@ plot_out_nbpostest <- function(mod, fp, likdat, cnt, simul = NULL, yr_pred = 202 } else { main_title <- "" } - - col_ci <- rgb(150, 150, 150, 125, max = 255) - + + col_ci <- grDevices::rgb(150, 150, 150, 125, max = 255) + # Postive tests if (plot.ci == T){ plot(I(out_nbtest_pos$value/1000) ~ out_nbtest_pos$year, pch = '', ylim = c(0, max(out_nbtest_pos$value, na.rm = TRUE) / 1000 * 1.25), - cex = 1, ylab = expression(paste(N^o, " of Positive Tests (in 1,000)")), - xlab='Year', xlim = c(2000, yr_pred), + cex = 1, ylab = expression(paste(N^o, " of Positive Tests (in 1,000)")), + xlab='Year', xlim = c(2000, yr_pred), main = main_title) - polygon(x = c(ci$year, rev(ci$year)), + graphics::polygon(x = c(ci$year, rev(ci$year)), y = c(ci$upper / 1000, rev(ci$lower / 1000)), col = col_ci, border = NA) - lines(I(out_nbtest_pos$value / 1000) ~ out_nbtest_pos$year, lwd = 1, col = 'orangered2') + graphics::lines(I(out_nbtest_pos$value / 1000) ~ out_nbtest_pos$year, lwd = 1, col = 'orangered2') if (!is.null(likdat$hts_pos$tot)) { if (cnt %in% redact) { pchpt <- '' } else { pchpt <- 16 } - points(I(likdat$hts_pos$tot / 1000) ~ I(likdat$hts_pos$year + 0.5), pch = pchpt) + graphics::points(I(likdat$hts_pos$tot / 1000) ~ I(likdat$hts_pos$year + 0.5), pch = pchpt) yield <- round(likdat$hts$totpos / likdat$hts$tot * 100, 1) y_range <- c(min(yield, na.rm = TRUE), max(yield, na.rm = TRUE)) p_pos_pop <- round(likdat$hts_pos$tot / plhiv[(likdat$hts_pos$year - start)] * 100, 0) p_pos_pop_r <- c(min(p_pos_pop, na.rm = TRUE), max(p_pos_pop, na.rm = TRUE)) - text(paste('Positivity = ', y_range[1], '-', y_range[2], '%', sep=''), + graphics::text(paste('Positivity = ', y_range[1], '-', y_range[2], '%', sep=''), x = yr_pred, y = max(likdat$hts_pos$tot / 1000, na.rm = TRUE) / 6, pos = 2) - text(paste('No Positive/PLHIV = ', p_pos_pop_r[1], '-', - p_pos_pop_r[2], '%', sep=''), x = yr_pred, + graphics::text(paste('No Positive/PLHIV = ', p_pos_pop_r[1], '-', + p_pos_pop_r[2], '%', sep=''), x = yr_pred, y = max(likdat$hts_pos$tot / 1000, na.rm = TRUE)/20, pos = 2) } } else { plot(I(out_nbtest_pos$value / 1000) ~ out_nbtest_pos$year, pch = '', ylim = c(0, max(out_nbtest_pos$value, na.rm = TRUE) / 1000 * 1.25), - cex = 1, ylab = expression(paste(N^o, " of Positive Tests (in 1,000)")), - xlab = 'Year', xlim = c(2000, yr_pred), + cex = 1, ylab = expression(paste(N^o, " of Positive Tests (in 1,000)")), + xlab = 'Year', xlim = c(2000, yr_pred), main = main_title) - lines(I(out_nbtest_pos$value / 1000) ~ out_nbtest_pos$year, lwd = 1, col = 'orangered2') + graphics::lines(I(out_nbtest_pos$value / 1000) ~ out_nbtest_pos$year, lwd = 1, col = 'orangered2') if (!is.null(likdat$hts_pos$tot)) { if (cnt %in% redact) { pchpt <- '' } else { pchpt <- 16 } - points(I(likdat$hts_pos$tot / 1000) ~ I(likdat$hts_pos$year + 0.5), pch = pchpt) + graphics::points(I(likdat$hts_pos$tot / 1000) ~ I(likdat$hts_pos$year + 0.5), pch = pchpt) yield <- round(likdat$hts$totpos / likdat$hts$tot * 100, 1) y_range <- c(min(yield, na.rm = TRUE), max(yield, na.rm = TRUE)) p_pos_pop <- round(likdat$hts_pos$tot / plhiv[(likdat$hts_pos$year - start)] * 100, 0) p_pos_pop_r <- c(min(p_pos_pop, na.rm = TRUE), max(p_pos_pop, na.rm = TRUE)) - text(paste('Positivity = ', y_range[1], '-', y_range[2], '%', sep=''), + graphics::text(paste('Positivity = ', y_range[1], '-', y_range[2], '%', sep=''), x = yr_pred, y = max(likdat$hts_pos$tot / 1000, na.rm = TRUE) / 6, pos = 2) - text(paste('No Positive/PLHIV = ', p_pos_pop_r[1], '-', - p_pos_pop_r[2], '%', sep = ''), x = yr_pred, + graphics::text(paste('No Positive/PLHIV = ', p_pos_pop_r[1], '-', + p_pos_pop_r[2], '%', sep = ''), x = yr_pred, y = max(likdat$hts_pos$tot / 1000, na.rm = TRUE) / 20, pos = 2) } } } #' @export -## +## ## -- UPDATE HERE -- ## * update yr_pred to current year plot_out_nbpostest_sex <- function(mod, fp, likdat, cnt, simul = NULL, yr_pred = 2022) { - + # redact <- c('Namibia','Uganda','Zambia','Zimbabwe') redact <- c('XXX') - + start <- fp$ss$proj_start mod <- simmod(fp) - plhiv <- apply(attr(mod, "hivpop")[,1:8,,], 4, FUN = sum) + + plhiv <- apply(attr(mod, "hivpop")[,1:8,,], 4, FUN = sum) + apply(attr(mod, "artpop")[,, 1:8,,], 5, FUN = sum) pop <- apply(mod[1:35,,,], 4, FUN = sum) - + out_nbtest_pos <- subset(get_out_nbtest_pos_sex(mod, fp), year <= yr_pred) out_nbtest_pos$year <- out_nbtest_pos$year + 0.5 out_nbtest_pos_f <- subset(out_nbtest_pos, sex == 'female') out_nbtest_pos_m <- subset(out_nbtest_pos, sex == 'male') max_ylim <- max(c(out_nbtest_pos_f$value, out_nbtest_pos_m$value), na.rm = TRUE) - + lik_sex <- likdat$hts_pos if (!is.null(lik_sex)){ lik_f <- subset(lik_sex, sex == 'female') @@ -455,57 +455,57 @@ plot_out_nbpostest_sex <- function(mod, fp, likdat, cnt, simul = NULL, yr_pred = } else { plot.ci <- FALSE } - - col_ci <- rgb(150, 150, 150, 125, max = 255) - + + col_ci <- grDevices::rgb(150, 150, 150, 125, max = 255) + # Postive tests if (plot.ci == T){ plot(I(out_nbtest_pos_f$value / 1000) ~ out_nbtest_pos_f$year, pch = '', ylim = c(0, max_ylim / 1000 * 1.25), - cex = 1, ylab = expression(paste(N^o, " of Positive Tests (in 1,000)")), - xlab = 'Year', xlim = c(2000, yr_pred), + cex = 1, ylab = expression(paste(N^o, " of Positive Tests (in 1,000)")), + xlab = 'Year', xlim = c(2000, yr_pred), main = expression(bold(paste("Total ", N^o, " of Positive Tests")))) - polygon(x = c(ci_f$year, rev(ci_f$year)), + graphics::polygon(x = c(ci_f$year, rev(ci_f$year)), y = c(ci_f$upper / 1000, rev(ci_f$lower / 1000)), col = col_ci, border = NA) - polygon(x = c(ci_m$year, rev(ci_m$year)), + graphics::polygon(x = c(ci_m$year, rev(ci_m$year)), y = c(ci_m$upper / 1000, rev(ci_m$lower / 1000)), col = col_ci, border = NA) - lines(I(out_nbtest_pos_f$value / 1000) ~ out_nbtest_pos_f$year, lwd = 1, col = 'orangered2') - lines(I(out_nbtest_pos_m$value / 1000) ~ out_nbtest_pos_m$year, lwd = 1, col = 'orangered2') - + graphics::lines(I(out_nbtest_pos_f$value / 1000) ~ out_nbtest_pos_f$year, lwd = 1, col = 'orangered2') + graphics::lines(I(out_nbtest_pos_m$value / 1000) ~ out_nbtest_pos_m$year, lwd = 1, col = 'orangered2') + if (length(lik_f) > 0 || length(lik_m) > 0) { if (cnt %in% redact) { pchpt <- '' } else { pchpt <- 16 } - points(I(lik_f$totpos / 1000) ~ I(lik_f$year + 0.5), pch = pchpt) - points(I(lik_m$totpos / 1000) ~ I(lik_m$year + 0.5), pch = pchpt + 1) + graphics::points(I(lik_f$totpos / 1000) ~ I(lik_f$year + 0.5), pch = pchpt) + graphics::points(I(lik_m$totpos / 1000) ~ I(lik_m$year + 0.5), pch = pchpt + 1) } } else { plot(I(out_nbtest_pos_f$value / 1000) ~ out_nbtest_pos_f$year, pch = '', ylim = c(0, max_ylim / 1000 * 1.25), - cex = 1, ylab = expression(paste(N^o, " of Positive Tests (in 1,000)")), - xlab = 'Year', xlim = c(2000, yr_pred), + cex = 1, ylab = expression(paste(N^o, " of Positive Tests (in 1,000)")), + xlab = 'Year', xlim = c(2000, yr_pred), main = expression(bold(paste("Total ", N^o, " of Positive Tests")))) - lines(I(out_nbtest_pos_f$value / 1000) ~ out_nbtest_pos_f$year, lwd = 1, col = 'orangered2') - lines(I(out_nbtest_pos_m$value / 1000) ~ out_nbtest_pos_m$year, lwd = 1, col = 'orangered2') + graphics::lines(I(out_nbtest_pos_f$value / 1000) ~ out_nbtest_pos_f$year, lwd = 1, col = 'orangered2') + graphics::lines(I(out_nbtest_pos_m$value / 1000) ~ out_nbtest_pos_m$year, lwd = 1, col = 'orangered2') if (!is.null(likdat$hts_pos$tot)) { if (cnt %in% redact) { pchpt <- '' } else { pchpt <- 16 } - points(I(lik_f$totpos / 1000) ~ I(lik_f$year + 0.5), pch = pchpt) - points(I(lik_m$totpos / 1000) ~ I(lik_m$year + 0.5), pch = pchpt + 1) + graphics::points(I(lik_f$totpos / 1000) ~ I(lik_f$year + 0.5), pch = pchpt) + graphics::points(I(lik_m$totpos / 1000) ~ I(lik_m$year + 0.5), pch = pchpt + 1) } } } #' @export -## +## ## -- UPDATE HERE -- ## * update yr_pred to current year -plot_out_evertestneg <- function(mod, fp, likdat, cnt, survey_hts, out_evertest, +plot_out_evertestneg <- function(mod, fp, likdat, cnt, survey_hts, out_evertest, simul = NULL, plot_legend = TRUE, yr_pred = 2022) { - + out_evertest <- subset(out_evertest, year <= yr_pred) out_evertest$year <- out_evertest$year + 0.5 survey_hts$year <- survey_hts$year + 0.5 - + # Decide if we plot CI or not if (!is.null(simul)){ ci <- getCI(simul$ever.test) @@ -517,14 +517,14 @@ plot_out_evertestneg <- function(mod, fp, likdat, cnt, survey_hts, out_evertest, } if (is.null(simul)) plot.ci <- FALSE - col_ci <- rgb(150, 150, 150, 125, max = 255) - + col_ci <- grDevices::rgb(150, 150, 150, 125, max = 255) + # Plots of Ever Test by Sex (among negatives) out_f <- subset(out_evertest, agegr == '15-49' & outcome == 'evertest' & sex =='female' & hivstatus == 'negative') out_m <- subset(out_evertest, agegr == '15-49' & outcome == 'evertest' & sex =='male' & hivstatus == 'negative') dat_f <- subset(survey_hts, country == cnt & outcome == 'evertest' & agegr == '15-49' & sex == 'female' & hivstatus == 'negative') dat_m <- subset(survey_hts, country == cnt & outcome == 'evertest' & agegr == '15-49' & sex == 'male' & hivstatus == 'negative') - + if (plot.ci == T){ ci_f <- subset(ci, agegr == '15-49' & outcome == 'evertest' & sex =='female' & hivstatus == 'negative') ci_m <- subset(ci, agegr == '15-49' & outcome == 'evertest' & sex =='male' & hivstatus == 'negative') @@ -532,51 +532,51 @@ plot_out_evertestneg <- function(mod, fp, likdat, cnt, survey_hts, out_evertest, plot(I(out_f$value * 100) ~ out_f$year, pch = "", ylim = c(0, 100), col = "maroon", main = "Negative Ever Tested", xlim = c(2000, yr_pred), ylab = "Proportion Ever Tested Among Susceptibles", xlab = "Year", lwd = 1) - polygon( + graphics::polygon( x = c(ci_f$year, rev(ci_f$year)), y = c(I(ci_f$upper * 100), rev(I(ci_f$lower * 100))), col = col_ci, border = NA) - polygon( + graphics::polygon( x = c(ci_m$year, rev(ci_m$year)), y = c(I(ci_m$upper * 100), rev(I(ci_m$lower * 100))), col = col_ci, border = NA) - lines(I(out_f$value * 100) ~ out_f$year, col = "maroon", lwd = 1) - lines(I(out_m$value * 100) ~ out_m$year, col = "navy", lwd = 1) - points(I(dat_f$est * 100) ~ dat_f$year, pch = 15, col = "maroon") - points(I(dat_m$est * 100) ~ dat_m$year, pch = 16, col = "navy") - segments( + graphics::lines(I(out_f$value * 100) ~ out_f$year, col = "maroon", lwd = 1) + graphics::lines(I(out_m$value * 100) ~ out_m$year, col = "navy", lwd = 1) + graphics::points(I(dat_f$est * 100) ~ dat_f$year, pch = 15, col = "maroon") + graphics::points(I(dat_m$est * 100) ~ dat_m$year, pch = 16, col = "navy") + graphics::segments( x0 = dat_f$year, y0 = I(dat_f$ci_l * 100), x1 = dat_f$year, y1 = I(dat_f$ci_u * 100), lwd = 1, col = "maroon") - segments( + graphics::segments( x0 = dat_m$year, y0 = I(dat_m$ci_l * 100), x1 = dat_m$year, y1 = I(dat_m$ci_u * 100), lwd = 1, col = "navy") - if (plot_legend) { - legend(x = 2000, y = 100, - legend = c('Women (15-49 years)','Men (15-49 years)'), + if (plot_legend) { + graphics::legend(x = 2000, y = 100, + legend = c('Women (15-49 years)','Men (15-49 years)'), col = c('maroon','navy'), bty = 'n', lwd = 1, pch = c(15,16)) } } else { plot(I(out_f$value * 100) ~ out_f$year, type = "l", ylim = c(0, 100), col = "maroon", main = "Negative Ever Tested", xlab = "Year", lwd = 1, xlim = c(2000, yr_pred), ylab = "Proportion Ever Tested Among Susceptibles") - lines(I(out_m$value * 100) ~ out_m$year, col = "navy", lwd = 1) - points(I(dat_f$est * 100) ~ dat_f$year, pch = 15, col = "maroon") - points(I(dat_m$est * 100) ~ dat_m$year, pch = 16, col = "navy") - segments( + graphics::lines(I(out_m$value * 100) ~ out_m$year, col = "navy", lwd = 1) + graphics::points(I(dat_f$est * 100) ~ dat_f$year, pch = 15, col = "maroon") + graphics::points(I(dat_m$est * 100) ~ dat_m$year, pch = 16, col = "navy") + graphics::segments( x0 = dat_f$year, y0 = I(dat_f$ci_l * 100), x1 = dat_f$year, y1 = I(dat_f$ci_u * 100), lwd = 1, col = "maroon") - segments( + graphics::segments( x0 = dat_m$year, y0 = I(dat_m$ci_l * 100), x1 = dat_m$year, y1 = I(dat_m$ci_u * 100), lwd = 1, col = "navy") - if (plot_legend) { - legend(x = 2000, y = 100, - legend = c('Women (15-49 years)','Men (15-49 years)'), + if (plot_legend) { + graphics::legend(x = 2000, y = 100, + legend = c('Women (15-49 years)','Men (15-49 years)'), col = c('maroon','navy'), bty = 'n', lwd = 2, pch = c(15, 16)) } } } #' @export -## +## ## -- UPDATE HERE -- ## * update yr_pred to current year plot_out_evertestpos <- function(mod, fp, likdat, cnt, survey_hts, out_evertest, @@ -586,7 +586,7 @@ plot_out_evertestpos <- function(mod, fp, likdat, cnt, survey_hts, out_evertest, out_evertest <- subset(out_evertest, year <= yr_pred) out_evertest$year <- out_evertest$year + 0.5 survey_hts$year <- survey_hts$year + 0.5 - + # Decide if we plot CI or not if (!is.null(simul)){ ci <- getCI(simul$ever.test) @@ -597,71 +597,71 @@ plot_out_evertestpos <- function(mod, fp, likdat, cnt, survey_hts, out_evertest, plot.ci <- TRUE } if (is.null(simul)) plot.ci <- FALSE - + if (plot_title == TRUE) { main_title <- "PLHIV Ever Tested" } else { main_title <- "" } - col_ci <- rgb(150, 150, 150, 125, max = 255) - + col_ci <- grDevices::rgb(150, 150, 150, 125, max = 255) + # Plots of Ever Test by Sex (among PLHIV) out_f <- subset(out_evertest, agegr == '15-49' & outcome == 'evertest' & sex =='female' & hivstatus == 'positive') out_m <- subset(out_evertest, agegr == '15-49' & outcome == 'evertest' & sex =='male' & hivstatus == 'positive') dat_f <- subset(survey_hts, country == cnt & agegr == '15-49' & sex == 'female' & hivstatus == 'positive' & outcome == 'evertest') dat_m <- subset(survey_hts, country == cnt & agegr == '15-49' & sex == 'male' & hivstatus == 'positive' & outcome == 'evertest') - + if (plot.ci == T){ ci_f <- subset(ci, agegr == '15-49' & outcome == 'evertest' & sex =='female' & hivstatus == 'positive') ci_m <- subset(ci, agegr == '15-49' & outcome == 'evertest' & sex =='male' & hivstatus == 'positive') plot(I(out_f$value * 100) ~ out_f$year, pch = '', ylim = c(0, 100), col = 'maroon', main = main_title, xlim=c(2000, yr_pred), ylab = 'Proportion PLHIV Ever Tested', xlab = 'Year', lwd = 1) - polygon(x = c(ci_f$year, rev(ci_f$year)), + graphics::polygon(x = c(ci_f$year, rev(ci_f$year)), y = c(I(ci_f$upper * 100), rev(I(ci_f$lower*100))), col = col_ci, border = NA) - polygon(x = c(ci_m$year, rev(ci_m$year)), + graphics::polygon(x = c(ci_m$year, rev(ci_m$year)), y = c(I(ci_m$upper * 100), rev(I(ci_m$lower * 100))), col = col_ci, border = NA) - lines(I(out_f$value * 100) ~ out_f$year, col = 'maroon', lwd = 1) - lines(I(out_m$value * 100) ~ out_m$year, col = 'navy', lwd = 1) - points(I(dat_f$est * 100) ~ I(dat_f$year - 0.1), pch = 15, col = 'maroon') - points(I(dat_m$est * 100) ~ I(dat_m$year + 0.1), pch = 16, col = 'navy') - segments(x0 = dat_f$year - 0.1, y0 = I(dat_f$ci_l * 100), + graphics::lines(I(out_f$value * 100) ~ out_f$year, col = 'maroon', lwd = 1) + graphics::lines(I(out_m$value * 100) ~ out_m$year, col = 'navy', lwd = 1) + graphics::points(I(dat_f$est * 100) ~ I(dat_f$year - 0.1), pch = 15, col = 'maroon') + graphics::points(I(dat_m$est * 100) ~ I(dat_m$year + 0.1), pch = 16, col = 'navy') + graphics::segments(x0 = dat_f$year - 0.1, y0 = I(dat_f$ci_l * 100), x1 = dat_f$year - 0.1, y1 = I(dat_f$ci_u * 100), lwd = 1, col = 'maroon') - segments(x0 = dat_m$year + 0.1, y0 = I(dat_m$ci_l * 100), + graphics::segments(x0 = dat_m$year + 0.1, y0 = I(dat_m$ci_l * 100), x1 = dat_m$year + 0.1, y1 = I(dat_m$ci_u * 100), lwd = 1, col = 'navy') - if (plot_legend) { - legend(x = 2000, y = 100, - legend = c('Women (15-49 years)','Men (15-49 years)'), + if (plot_legend) { + graphics::legend(x = 2000, y = 100, + legend = c('Women (15-49 years)','Men (15-49 years)'), col = c('maroon','navy'), bty = 'n', lwd = 1, pch = c(15, 16)) } } else { plot(I(out_f$value * 100) ~ out_f$year, type='l', ylim = c(0,100), col = 'maroon', main = main_title, xlim = c(2000, yr_pred), ylab = 'Proportion PLHIV Ever Tested', xlab = 'Year', lwd = 1) - lines(I(out_m$value * 100) ~ out_m$year, col = 'navy', lwd = 1) - points(I(dat_f$est * 100) ~ I(dat_f$year - 0.1), pch = 15, col = 'maroon') - points(I(dat_m$est * 100) ~ I(dat_m$year + 0.1), pch = 16, col = 'navy') - segments(x0 = dat_f$year - 0.1, y0 = I(dat_f$ci_l * 100), + graphics::lines(I(out_m$value * 100) ~ out_m$year, col = 'navy', lwd = 1) + graphics::points(I(dat_f$est * 100) ~ I(dat_f$year - 0.1), pch = 15, col = 'maroon') + graphics::points(I(dat_m$est * 100) ~ I(dat_m$year + 0.1), pch = 16, col = 'navy') + graphics::segments(x0 = dat_f$year - 0.1, y0 = I(dat_f$ci_l * 100), x1 = dat_f$year - 0.1, y1 = I(dat_f$ci_u * 100), lwd = 1, col = 'maroon') - segments(x0 = dat_m$year + 0.1, y0 = I(dat_m$ci_l * 100), + graphics::segments(x0 = dat_m$year + 0.1, y0 = I(dat_m$ci_l * 100), x1 = dat_m$year + 0.1, y1 = I(dat_m$ci_u * 100), lwd = 1, col = 'navy') - if (plot_legend) { - legend(x = 2000, y = 100, - legend = c('Women (15-49 years)', 'Men (15-49 years)'), + if (plot_legend) { + graphics::legend(x = 2000, y = 100, + legend = c('Women (15-49 years)', 'Men (15-49 years)'), col = c('maroon', 'navy'), bty = 'n', lwd = 2, pch = c(15, 16)) } } } #' @export -## +## ## -- UPDATE HERE -- ## * update yr_pred to current year -plot_out_evertest <- function(mod, fp, likdat, cnt, survey_hts, out_evertest, +plot_out_evertest <- function(mod, fp, likdat, cnt, survey_hts, out_evertest, simul = NULL, plot_legend = TRUE, yr_pred = 2022, plot_title = TRUE) { out_evertest <- subset(out_evertest, year <= yr_pred) out_evertest$year <- out_evertest$year + 0.5 survey_hts$year <- survey_hts$year + 0.5 - + # Decide if we plot CI or not if (!is.null(simul)){ ci <- getCI(simul$ever.test) @@ -672,7 +672,7 @@ plot_out_evertest <- function(mod, fp, likdat, cnt, survey_hts, out_evertest, plot.ci <- TRUE } if (is.null(simul)) plot.ci <- FALSE - + if (plot_title == TRUE) { main_title <- "Population Ever Tested" } else { main_title <- "" } @@ -681,57 +681,57 @@ plot_out_evertest <- function(mod, fp, likdat, cnt, survey_hts, out_evertest, out_m <- subset(out_evertest, agegr == '15-49' & outcome == 'evertest' & sex =='male' & hivstatus == 'all') dat_f <- subset(survey_hts, country == cnt & agegr == '15-49' & sex == 'female' & hivstatus == 'all' & outcome == 'evertest') dat_m <- subset(survey_hts, country == cnt & agegr == '15-49' & sex == 'male' & hivstatus == 'all' & outcome == 'evertest') - - col_ci <- rgb(150, 150, 150, 125, max = 255) - + + col_ci <- grDevices::rgb(150, 150, 150, 125, max = 255) + if (plot.ci == T){ ci_f <- subset(ci, agegr == '15-49' & outcome == 'evertest' & sex =='female' & hivstatus == 'all') ci_m <- subset(ci, agegr == '15-49' & outcome == 'evertest' & sex =='male' & hivstatus == 'all') plot(I(out_f$value * 100) ~ out_f$year, pch='', ylim = c(0,100), col = 'maroon', main = main_title, xlim=c(2000, yr_pred), ylab = 'Proportion Ever Tested', xlab = 'Year', lwd = 1) - polygon(x = c(ci_f$year, rev(ci_f$year)), + graphics::polygon(x = c(ci_f$year, rev(ci_f$year)), y = c(I(ci_f$upper * 100), rev(I(ci_f$lower*100))), col = col_ci, border = NA) - lines(I(out_f$value * 100) ~ out_m$year, col='maroon', lwd=1) - polygon(x = c(ci_m$year, rev(ci_m$year)), + graphics::lines(I(out_f$value * 100) ~ out_m$year, col='maroon', lwd=1) + graphics::polygon(x = c(ci_m$year, rev(ci_m$year)), y = c(I(ci_m$upper * 100), rev(I(ci_m$lower*100))), col = col_ci, border = NA) - lines(I(out_m$value * 100) ~ out_m$year, col = 'navy', lwd = 1) - points(I(dat_f$est * 100) ~ I(dat_f$year + 0.1), pch = 15, col='maroon') - points(I(dat_m$est * 100) ~ I(dat_m$year - 0.1), pch = 16, col = 'navy') - segments(x0 = dat_f$year + 0.1, y0 = I(dat_f$ci_l*100), + graphics::lines(I(out_m$value * 100) ~ out_m$year, col = 'navy', lwd = 1) + graphics::points(I(dat_f$est * 100) ~ I(dat_f$year + 0.1), pch = 15, col='maroon') + graphics::points(I(dat_m$est * 100) ~ I(dat_m$year - 0.1), pch = 16, col = 'navy') + graphics::segments(x0 = dat_f$year + 0.1, y0 = I(dat_f$ci_l*100), x1 = dat_f$year + 0.1, y1 = I(dat_f$ci_u*100), lwd = 1, col = 'maroon') - segments(x0 = dat_m$year - 0.1, y0 = I(dat_m$ci_l*100), + graphics::segments(x0 = dat_m$year - 0.1, y0 = I(dat_m$ci_l*100), x1 = dat_m$year - 0.1, y1 = I(dat_m$ci_u*100), lwd = 1, col = 'navy') - if (plot_legend) { - legend(x = 2000, y = 100, - legend = c('Women (15-49 years)','Men (15-49 years)'), + if (plot_legend) { + graphics::legend(x = 2000, y = 100, + legend = c('Women (15-49 years)','Men (15-49 years)'), col = c('maroon','navy'), bty = 'n', lwd = 1, pch = c(15,16)) } } else { - plot(I(out_f$value * 100) ~ out_f$year, type = 'l', ylim = c(0, 100), + plot(I(out_f$value * 100) ~ out_f$year, type = 'l', ylim = c(0, 100), col = 'maroon', main = main_title, - xlim = c(2000, yr_pred), ylab = 'Proportion Ever Tested', + xlim = c(2000, yr_pred), ylab = 'Proportion Ever Tested', xlab = 'Year', lwd = 1) - lines(I(out_m$value*100) ~ out_m$year, col = 'navy', lwd = 1) - points(I(dat_f$est*100) ~ I(dat_f$year + 0.1), pch = 15, col = 'maroon') - points(I(dat_m$est*100) ~ I(dat_m$year - 0.1), pch = 16, col = 'navy') - segments(x0 = dat_f$year + 0.1, y0 = I(dat_f$ci_l*100), + graphics::lines(I(out_m$value*100) ~ out_m$year, col = 'navy', lwd = 1) + graphics::points(I(dat_f$est*100) ~ I(dat_f$year + 0.1), pch = 15, col = 'maroon') + graphics::points(I(dat_m$est*100) ~ I(dat_m$year - 0.1), pch = 16, col = 'navy') + graphics::segments(x0 = dat_f$year + 0.1, y0 = I(dat_f$ci_l*100), x1 = dat_f$year + 0.1, y1 = I(dat_f$ci_u*100), lwd = 1, col = 'maroon') - segments(x0 = dat_m$year - 0.1, y0 = I(dat_m$ci_l*100), + graphics::segments(x0 = dat_m$year - 0.1, y0 = I(dat_m$ci_l*100), x1 = dat_m$year - 0.1, y1 = I(dat_m$ci_u*100), lwd = 1, col = 'navy') - if (plot_legend) { - legend(x = 2000, y = 100, - legend = c('Women (15-49 years)','Men (15-49 years)'), + if (plot_legend) { + graphics::legend(x = 2000, y = 100, + legend = c('Women (15-49 years)','Men (15-49 years)'), col = c('maroon','navy'), bty = 'n', lwd = 2, pch = c(15,16)) } } } #' @export -## +## ## -- UPDATE HERE -- ## * update yr_pred to current year -plot_out_90s <- function(mod, fp, likdat, cnt, out_evertest, survey_hts, +plot_out_90s <- function(mod, fp, likdat, cnt, out_evertest, survey_hts, simul = NULL, plot_legend = TRUE, yr_pred = 2022) { phia_aware <- subset(survey_hts, country == cnt & agegr == '15-49' & @@ -740,16 +740,16 @@ plot_out_90s <- function(mod, fp, likdat, cnt, out_evertest, survey_hts, sex == 'both' & outcome == 'art_both') phia_aware$year <- phia_aware$year + 0.5 phia_art$year <- phia_art$year + 0.5 - + if (all(is.na(phia_art$est))) { phia_art <- subset(survey_hts, country == cnt & agegr == '15-49' & - sex == 'both' & outcome == 'art_self') + sex == 'both' & outcome == 'art_self') phia_art$year <- phia_art$year + 0.5 } svy_evertest <- subset(survey_hts, country == cnt & agegr == '15-49' & hivstatus == 'positive' & sex == 'both' & outcome == 'evertest') svy_evertest$year <- svy_evertest$year + 0.5 - + out_evertest <- subset(out_evertest, year <= yr_pred) out_aware <- subset(get_out_aware(mod, fp), year <= yr_pred) out_art <- subset(get_out_art(mod, fp), year >= 2000 & year <= yr_pred) @@ -771,15 +771,15 @@ plot_out_90s <- function(mod, fp, likdat, cnt, out_evertest, survey_hts, } if (is.null(simul)) plot.ci <- FALSE - col_ci <- rgb(150, 150, 150, 125, max = 255) - + col_ci <- grDevices::rgb(150, 150, 150, 125, max = 255) + # First 90 out_f <- subset(out_aware, agegr == '15-49' & outcome == 'aware' & sex =='female') out_m <- subset(out_aware, agegr == '15-49' & outcome == 'aware' & sex =='male') out_a <- subset(out_aware, agegr == '15-49' & outcome == 'aware' & sex =='both') out_ever_all <- subset(out_evertest, agegr == '15-49' & outcome == 'evertest' & sex == 'both' & hivstatus == 'positive') - + if (plot.ci == TRUE){ ci_f <- subset(ci$diagnoses, agegr == '15-49' & outcome == 'aware' & sex =='female') ci_m <- subset(ci$diagnoses, agegr == '15-49' & outcome == 'aware' & sex =='male') @@ -787,78 +787,78 @@ plot_out_90s <- function(mod, fp, likdat, cnt, out_evertest, survey_hts, ci_ever_all <- subset(ci$ever.test, agegr == '15-49' & outcome == 'evertest' & sex == 'both' & hivstatus == 'positive') - plot(I(out_a$value*100) ~ out_a$year, pch = '', ylim = c(0,100), + plot(I(out_a$value*100) ~ out_a$year, pch = '', ylim = c(0,100), col = 'darkslategrey', main = 'PLHIV Ever Tested, 1st 90, and ART Coverage', xlim = c(2000, yr_pred), ylab = 'Proportion', xlab='Year', lwd = 1) - polygon(x = c(ci_a$year, rev(ci_a$year)), + graphics::polygon(x = c(ci_a$year, rev(ci_a$year)), y = c(I(ci_a$upper*100), rev(I(ci_a$lower*100))), col = col_ci, border = NA) - lines(I(out_a$value*100) ~ out_a$year, col = 'darkslategrey', lwd = 1) - lines(I(out_art$value*100) ~ out_art$year, col = 'deepskyblue2', lwd = 1, lty = 1) - polygon(x = c(ci_ever_all$year, rev(ci_ever_all$year)), + graphics::lines(I(out_a$value*100) ~ out_a$year, col = 'darkslategrey', lwd = 1) + graphics::lines(I(out_art$value*100) ~ out_art$year, col = 'deepskyblue2', lwd = 1, lty = 1) + graphics::polygon(x = c(ci_ever_all$year, rev(ci_ever_all$year)), y = c(I(ci_ever_all$upper*100), rev(I(ci_ever_all$lower*100))), col = col_ci,border = NA) - lines(I(out_ever_all$value*100) ~ out_ever_all$year, col = 'orange3', lwd = 1, lty = 1) - legend(x = 2000, y = 100, legend = c('PLHIV Ever Tested','PLHIV Aware','ART Coverage'), + graphics::lines(I(out_ever_all$value*100) ~ out_ever_all$year, col = 'orange3', lwd = 1, lty = 1) + graphics::legend(x = 2000, y = 100, legend = c('PLHIV Ever Tested','PLHIV Aware','ART Coverage'), col = c('orange3','darkslategrey','deepskyblue2'), pch = c(25, 1, 2), lwd = 1, bty = 'n') # We plot ever tested if (dim(svy_evertest)[1] > 0) { - points(I(svy_evertest$est*100) ~ svy_evertest$year, pch = 25, bg = 'grey15') - segments(x0 = svy_evertest$year, y0 = I(svy_evertest$ci_l*100), + graphics::points(I(svy_evertest$est*100) ~ svy_evertest$year, pch = 25, bg = 'grey15') + graphics::segments(x0 = svy_evertest$year, y0 = I(svy_evertest$ci_l*100), x1 = svy_evertest$year, y1 = I(svy_evertest$ci_u*100), lwd = 1, col = 'bisque4') } # We plot PHIA if any if (dim(phia_aware)[1] > 0) { - points(I(phia_aware$est*100) ~ phia_aware$year, pch = 1) - segments(x0 = phia_aware$year, y0 = I(phia_aware$ci_l*100), + graphics::points(I(phia_aware$est*100) ~ phia_aware$year, pch = 1) + graphics::segments(x0 = phia_aware$year, y0 = I(phia_aware$ci_l*100), x1 = phia_aware$year, y1 = I(phia_aware$ci_u*100), lwd = 1, col = 'bisque4') - text('Empty symbols are \n direct survey measures \n (not included in fitting)', x = 2015, y = 12.5, cex = 0.8) + graphics::text('Empty symbols are \n direct survey measures \n (not included in fitting)', x = 2015, y = 12.5, cex = 0.8) if (any(cnt == c('Swaziland','Zambia','Zimbabwe'))) { - text('(not ARV-corrected)', x = 2015, y = 1, cex = 0.8) } + graphics::text('(not ARV-corrected)', x = 2015, y = 1, cex = 0.8) } } if (dim(phia_art)[1] > 0) { - points(I(phia_art$est*100) ~ phia_art$year, pch = 2) - segments(x0 = phia_art$year, y0 = I(phia_art$ci_l*100), + graphics::points(I(phia_art$est*100) ~ phia_art$year, pch = 2) + graphics::segments(x0 = phia_art$year, y0 = I(phia_art$ci_l*100), x1 = phia_art$year, y1 = I(phia_art$ci_u*100), lwd = 1, col = 'bisque4') } } else { - plot(I(out_a$value*100) ~ out_a$year, type = 'l', ylim = c(0,100), + plot(I(out_a$value*100) ~ out_a$year, type = 'l', ylim = c(0,100), col = 'darkslategrey', main='PLHIV Ever Tested, 1st 90, and ART Coverage', xlim = c(2000, yr_pred), ylab = 'Proportion', xlab = 'Year', lwd = 1) - lines(I(out_art$value*100) ~ out_art$year, col = 'deepskyblue2', lwd = 1, lty = 1) - lines(I(out_ever_all$value*100) ~ out_ever_all$year, col='orange3', lwd = 1, lty = 1) - legend(x = 2000, y = 100, legend = c('PLHIV Ever Tested','PLHIV Aware','ART Coverage'), + graphics::lines(I(out_art$value*100) ~ out_art$year, col = 'deepskyblue2', lwd = 1, lty = 1) + graphics::lines(I(out_ever_all$value*100) ~ out_ever_all$year, col='orange3', lwd = 1, lty = 1) + graphics::legend(x = 2000, y = 100, legend = c('PLHIV Ever Tested','PLHIV Aware','ART Coverage'), col = c('orange3','darkslategrey','deepskyblue2'), pch = c(25, 1, 2), lwd = 1, bty = 'n') # We plot ever tested if (dim(svy_evertest)[1] > 0) { - points(I(svy_evertest$est*100) ~ svy_evertest$year, pch = 25, bg = 'grey15') - segments(x0 = svy_evertest$year, y0 = I(svy_evertest$ci_l*100), + graphics::points(I(svy_evertest$est*100) ~ svy_evertest$year, pch = 25, bg = 'grey15') + graphics::segments(x0 = svy_evertest$year, y0 = I(svy_evertest$ci_l*100), x1 = svy_evertest$year, y1 = I(svy_evertest$ci_u*100), lwd = 1, col = 'bisque4') } # We plot PHIA if any if (dim(phia_aware)[1]>0) { - points(I(phia_aware$est*100) ~ phia_aware$year, pch = 1) - segments(x0 = phia_aware$year, y0 = I(phia_aware$ci_l*100), + graphics::points(I(phia_aware$est*100) ~ phia_aware$year, pch = 1) + graphics::segments(x0 = phia_aware$year, y0 = I(phia_aware$ci_l*100), x1 = phia_aware$year, y1 = I(phia_aware$ci_u*100), lwd = 1, col = 'bisque4') - text('X-Validation', x = 2015, y = 12.5, cex = 0.8) + graphics::text('X-Validation', x = 2015, y = 12.5, cex = 0.8) if (any(cnt == c('Swaziland','Zambia','Zimbabwe'))) { - text('(not ARV-corrected)', x = 2015, y = 1, cex = 0.8) } + graphics::text('(not ARV-corrected)', x = 2015, y = 1, cex = 0.8) } } if (dim(phia_art)[1]>0) { - points(I(phia_art$est*100) ~ phia_art$year, pch = 2) - segments(x0 = phia_art$year, y0 = I(phia_art$ci_l*100), + graphics::points(I(phia_art$est*100) ~ phia_art$year, pch = 2) + graphics::segments(x0 = phia_art$year, y0 = I(phia_art$ci_l*100), x1 = phia_art$year, y1 = I(phia_art$ci_u*100), lwd = 1, col = 'bisque4') } } } #' @export -## +## ## -- UPDATE HERE -- ## * update yr_pred to current year -plot_out_evertest_fbyage <- function(mod, fp, likdat, cnt, survey_hts, out_evertest, +plot_out_evertest_fbyage <- function(mod, fp, likdat, cnt, survey_hts, out_evertest, simul = NULL, plot_legend = TRUE, yr_pred = 2022) { out_evertest <- subset(out_evertest, year <= yr_pred) out_evertest$year <- out_evertest$year + 0.5 survey_hts$year <- survey_hts$year + 0.5 - + # Decide if we plot CI or not if (!is.null(simul)){ ci <- getCI(simul$ever.test) @@ -880,79 +880,79 @@ plot_out_evertest_fbyage <- function(mod, fp, likdat, cnt, survey_hts, out_evert agegr == "25-34" & sex == "female" & outcome == "evertest") datage35f <- subset(survey_hts, country == cnt & hivstatus == "all" & agegr == "35-49" & sex == "female" & outcome == "evertest") - - col_ci <- rgb(150, 150, 150, 125, max = 255) - + + col_ci <- grDevices::rgb(150, 150, 150, 125, max = 255) + if (plot.ci == T){ ci_1524f <- subset(ci, agegr == '15-24' & outcome == 'evertest' & hivstatus=='all' & sex=='female') ci_2534f <- subset(ci, agegr == '25-34' & outcome == 'evertest' & hivstatus=='all' & sex=='female') ci_3549f <- subset(ci, agegr == '35-49' & outcome == 'evertest' & hivstatus=='all' & sex=='female') - plot(I(out_1524f$value*100) ~ out_1524f$year, pch = '', ylim = c(0,100), + plot(I(out_1524f$value*100) ~ out_1524f$year, pch = '', ylim = c(0,100), xlim=c(2000, yr_pred), - ylab = 'Proportion ever tested', xlab = 'Year', + ylab = 'Proportion ever tested', xlab = 'Year', main = 'Women Ever Tested, by Age') - polygon(x = c(ci_1524f$year, rev(ci_1524f$year)), + graphics::polygon(x = c(ci_1524f$year, rev(ci_1524f$year)), y = c(I(ci_1524f$upper*100), rev(I(ci_1524f$lower*100))), col = col_ci, border = NA) - polygon(x = c(ci_2534f$year, rev(ci_2534f$year)), + graphics::polygon(x = c(ci_2534f$year, rev(ci_2534f$year)), y = c(I(ci_2534f$upper*100), rev(I(ci_2534f$lower*100))), col = col_ci, border = NA) - polygon(x = c(ci_3549f$year, rev(ci_3549f$year)), + graphics::polygon(x = c(ci_3549f$year, rev(ci_3549f$year)), y = c(I(ci_3549f$upper*100), rev(I(ci_3549f$lower*100))), col = col_ci, border = NA) - lines(I(out_1524f$value*100) ~ out_1524f$year, col = 'maroon', lwd = 1, lty = 3) - lines(I(out_2534f$value*100) ~ out_2534f$year, col = 'maroon', lwd = 1, lty = 2) - lines(I(out_3549f$value*100) ~ out_3549f$year, col = 'maroon', lwd = 1, lty = 1) - if (plot_legend) { - legend(x = 2000, y = 100, - legend = c('15-24 years', '25-34 years','35-49 years'), + graphics::lines(I(out_1524f$value*100) ~ out_1524f$year, col = 'maroon', lwd = 1, lty = 3) + graphics::lines(I(out_2534f$value*100) ~ out_2534f$year, col = 'maroon', lwd = 1, lty = 2) + graphics::lines(I(out_3549f$value*100) ~ out_3549f$year, col = 'maroon', lwd = 1, lty = 1) + if (plot_legend) { + graphics::legend(x = 2000, y = 100, + legend = c('15-24 years', '25-34 years','35-49 years'), col = c('maroon','maroon'), lwd = 2, lty = c(3,2,1), pch = c(16,15,17), pt.cex = 0.7, bty = 'n') } - points(I(datage15f$est*100) ~ I(datage15f$year-0.1), pch = 16, col = 'goldenrod3', cex = 0.8) - points(I(datage25f$est*100) ~ datage25f$year, pch=15, col = 'goldenrod3', cex = 0.8) - points(I(datage35f$est*100) ~ I(datage35f$year+0.1), pch = 17, col = 'goldenrod3', cex = 0.8) - segments(x0 = datage15f$year-0.1, y0 = I(datage15f$ci_l*100), + graphics::points(I(datage15f$est*100) ~ I(datage15f$year-0.1), pch = 16, col = 'goldenrod3', cex = 0.8) + graphics::points(I(datage25f$est*100) ~ datage25f$year, pch=15, col = 'goldenrod3', cex = 0.8) + graphics::points(I(datage35f$est*100) ~ I(datage35f$year+0.1), pch = 17, col = 'goldenrod3', cex = 0.8) + graphics::segments(x0 = datage15f$year-0.1, y0 = I(datage15f$ci_l*100), x1 = datage15f$year-0.1, y1 = I(datage15f$ci_u*100), col = 'goldenrod3', lwd=1) - segments(x0 = datage25f$year, y0 = I(datage25f$ci_l*100), + graphics::segments(x0 = datage25f$year, y0 = I(datage25f$ci_l*100), x1 = datage25f$year, y1 = I(datage25f$ci_u*100), col = 'goldenrod3', lwd = 1) - segments(x0 = datage35f$year+0.1, y0 = I(datage35f$ci_l*100), + graphics::segments(x0 = datage35f$year+0.1, y0 = I(datage35f$ci_l*100), x1 = datage35f$year+0.1, y1 = I(datage35f$ci_u*100), col = 'goldenrod3', lwd = 1) } else { plot(I(out_1524f$value*100) ~ out_1524f$year, type = 'l', col = 'maroon', ylim = c(0,100), xlim = c(2000, yr_pred), lwd = 1, lty = 3, - ylab = 'Proportion ever tested', xlab = 'Year', + ylab = 'Proportion ever tested', xlab = 'Year', main = 'Women Ever Tested, by Age') - lines(I(out_2534f$value*100) ~ out_2534f$year, col = 'maroon', lwd = 1, lty = 2) - lines(I(out_3549f$value*100) ~ out_3549f$year, col = 'maroon', lwd = 1, lty = 1) - if (plot_legend) { - legend(x = 2000, y = 100, - legend = c('15-24 years', '25-34 years','35-49 years'), + graphics::lines(I(out_2534f$value*100) ~ out_2534f$year, col = 'maroon', lwd = 1, lty = 2) + graphics::lines(I(out_3549f$value*100) ~ out_3549f$year, col = 'maroon', lwd = 1, lty = 1) + if (plot_legend) { + graphics::legend(x = 2000, y = 100, + legend = c('15-24 years', '25-34 years','35-49 years'), col = c('maroon','maroon'), lwd = 2, lty = c(3,2,1), pch = c(16,15,17), pt.cex = 0.7, bty = 'n') } - points(I(datage15f$est*100) ~ I(datage15f$year-0.1), pch = 16, col = 'goldenrod3', cex = 0.8) - points(I(datage25f$est*100) ~ datage25f$year, pch = 15, col = 'goldenrod3', cex = 0.8) - points(I(datage35f$est*100) ~ I(datage35f$year+0.1), pch = 17, col = 'goldenrod3', cex = 0.8) - segments(x0 = datage15f$year-0.1, y0 = I(datage15f$ci_l*100), + graphics::points(I(datage15f$est*100) ~ I(datage15f$year-0.1), pch = 16, col = 'goldenrod3', cex = 0.8) + graphics::points(I(datage25f$est*100) ~ datage25f$year, pch = 15, col = 'goldenrod3', cex = 0.8) + graphics::points(I(datage35f$est*100) ~ I(datage35f$year+0.1), pch = 17, col = 'goldenrod3', cex = 0.8) + graphics::segments(x0 = datage15f$year-0.1, y0 = I(datage15f$ci_l*100), x1 = datage15f$yea-0.1, y1 = I(datage15f$ci_u*100), col = 'goldenrod3', lwd = 1) - segments(x0 = datage25f$year, y0 = I(datage25f$ci_l*100), + graphics::segments(x0 = datage25f$year, y0 = I(datage25f$ci_l*100), x1 = datage25f$year, y1 = I(datage25f$ci_u*100), col = 'goldenrod3', lwd = 1) - segments(x0 = datage35f$year+0.1, y0 = I(datage35f$ci_l*100), + graphics::segments(x0 = datage35f$year+0.1, y0 = I(datage35f$ci_l*100), x1 = datage35f$year+0.1, y1 = I(datage35f$ci_u*100), col = 'goldenrod3', lwd = 1) } } #' @export -## +## ## -- UPDATE HERE -- ## * update yr_pred to current year -plot_out_evertest_mbyage <- function(mod, fp, likdat, cnt, survey_hts, - out_evertest, simul = NULL, +plot_out_evertest_mbyage <- function(mod, fp, likdat, cnt, survey_hts, + out_evertest, simul = NULL, plot_legend = TRUE, yr_pred = 2022) { out_evertest <- subset(out_evertest, year <= yr_pred) out_evertest$year <- out_evertest$year + 0.5 survey_hts$year <- survey_hts$year + 0.5 - + # Decide if we plot CI or not if (!is.null(simul)){ ci <- getCI(simul$ever.test) @@ -974,61 +974,61 @@ plot_out_evertest_mbyage <- function(mod, fp, likdat, cnt, survey_hts, datage35m <- subset(survey_hts, country == cnt & hivstatus == "all" & agegr == "35-49" & sex == "male" & outcome == "evertest") - col_ci <- rgb(150, 150, 150, 125, max = 255) - + col_ci <- grDevices::rgb(150, 150, 150, 125, max = 255) + if (plot.ci == T){ ci_1524m <- subset(ci, agegr == '15-24' & outcome == 'evertest' & hivstatus=='all' & sex=='male') ci_2534m <- subset(ci, agegr == '25-34' & outcome == 'evertest' & hivstatus=='all' & sex=='male') ci_3549m <- subset(ci, agegr == '35-49' & outcome == 'evertest' & hivstatus=='all' & sex=='male') - plot(I(out_1524m$value*100) ~ out_1524m$year, pch = '', + plot(I(out_1524m$value*100) ~ out_1524m$year, pch = '', ylim = c(0,100), xlim = c(2000, yr_pred), ylab = 'Proportion ever tested', xlab = 'Year', main = 'Men Ever Tested, by Age') - polygon(x = c(ci_1524m$year, rev(ci_1524m$year)), + graphics::polygon(x = c(ci_1524m$year, rev(ci_1524m$year)), y = c(I(ci_1524m$upper*100), rev(I(ci_1524m$lower*100))), col = col_ci, border = NA) - polygon(x = c(ci_2534m$year, rev(ci_2534m$year)), + graphics::polygon(x = c(ci_2534m$year, rev(ci_2534m$year)), y = c(I(ci_2534m$upper*100), rev(I(ci_2534m$lower*100))), col = col_ci, border = NA) - polygon(x = c(ci_3549m$year, rev(ci_3549m$year)), + graphics::polygon(x = c(ci_3549m$year, rev(ci_3549m$year)), y = c(I(ci_3549m$upper*100), rev(I(ci_3549m$lower*100))), col = col_ci, border = NA) - lines(I(out_1524m$value*100) ~ I(out_1524m$year-0.1), col = 'navy', lwd = 1, lty = 3) - lines(I(out_2534m$value*100) ~ out_2534m$year, col = 'navy', lwd = 1, lty = 2) - lines(I(out_3549m$value*100) ~ I(out_3549m$year+0.1), col = 'navy', lwd = 1, lty = 1) - if (plot_legend) { - legend(x = 2000, y = 100, + graphics::lines(I(out_1524m$value*100) ~ I(out_1524m$year-0.1), col = 'navy', lwd = 1, lty = 3) + graphics::lines(I(out_2534m$value*100) ~ out_2534m$year, col = 'navy', lwd = 1, lty = 2) + graphics::lines(I(out_3549m$value*100) ~ I(out_3549m$year+0.1), col = 'navy', lwd = 1, lty = 1) + if (plot_legend) { + graphics::legend(x = 2000, y = 100, legend = c('15-24', '25-34', '35-49'), col = c('navy','navy'), lwd = 2, lty = c(3,2,1), pch = c(16,15,17), pt.cex = 0.7, bty = 'n') } - points(I(datage15m$est*100) ~ I(datage15m$year-0.1), pch = 16, col = 'goldenrod3', cex = 0.8) - points(I(datage25m$est*100) ~ datage25m$year, pch = 15, col = 'goldenrod3', cex = 0.8) - points(I(datage35m$est*100) ~ I(datage25m$year+0.1), pch = 17, col = 'goldenrod3', cex = 0.8) - segments(x0 = datage15m$year-0.1, y0 = I(datage15m$ci_l*100), + graphics::points(I(datage15m$est*100) ~ I(datage15m$year-0.1), pch = 16, col = 'goldenrod3', cex = 0.8) + graphics::points(I(datage25m$est*100) ~ datage25m$year, pch = 15, col = 'goldenrod3', cex = 0.8) + graphics::points(I(datage35m$est*100) ~ I(datage25m$year+0.1), pch = 17, col = 'goldenrod3', cex = 0.8) + graphics::segments(x0 = datage15m$year-0.1, y0 = I(datage15m$ci_l*100), x1 = datage15m$year-0.1, y1 = I(datage15m$ci_u*100), col = 'goldenrod3', lwd = 1) - segments(x0 = datage25m$year, y0 = I(datage25m$ci_l*100), + graphics::segments(x0 = datage25m$year, y0 = I(datage25m$ci_l*100), x1=datage25m$year, y1 = I(datage25m$ci_u*100), col = 'goldenrod3', lwd = 1) - segments(x0 = datage35m$year+0.1, y0 = I(datage35m$ci_l*100), + graphics::segments(x0 = datage35m$year+0.1, y0 = I(datage35m$ci_l*100), x1 = datage35m$year+0.1, y1 = I(datage35m$ci_u*100), col = 'goldenrod3', lwd = 1) } else { - plot(I(out_1524m$value*100) ~ out_1524m$year, type = 'l', col = 'navy', + plot(I(out_1524m$value*100) ~ out_1524m$year, type = 'l', col = 'navy', ylim = c(0,100), xlim = c(2000, yr_pred), lwd = 1, lty = 3, - ylab = 'Proportion ever tested', xlab = 'Year', + ylab = 'Proportion ever tested', xlab = 'Year', main = 'Men Ever Tested, by Age') - lines(I(out_2534m$value*100) ~ out_2534m$year, col = 'navy', lwd = 1, lty = 2) - lines(I(out_3549m$value*100) ~ out_3549m$year, col = 'navy', lwd = 1, lty = 1) - if (plot_legend) { - legend(x = 2000, y = 100, + graphics::lines(I(out_2534m$value*100) ~ out_2534m$year, col = 'navy', lwd = 1, lty = 2) + graphics::lines(I(out_3549m$value*100) ~ out_3549m$year, col = 'navy', lwd = 1, lty = 1) + if (plot_legend) { + graphics::legend(x = 2000, y = 100, legend = c('15-24', '25-34', '35-49'), col = c('navy','navy'), lwd = 2, lty = c(3,2,1), pch = c(16,15,17), pt.cex = 0.7, bty = 'n') } - points(I(datage15m$est*100) ~ I(datage15m$year-0.1), pch = 16, col = 'goldenrod3', cex = 0.8) - points(I(datage25m$est*100) ~ datage25m$year, pch = 15, col = 'goldenrod3', cex = 0.8) - points(I(datage35m$est*100) ~ I(datage25m$year+0.1), pch = 17, col = 'goldenrod3', cex = 0.8) - segments(x0=datage15m$year-0.1, y0=I(datage15m$ci_l*100), + graphics::points(I(datage15m$est*100) ~ I(datage15m$year-0.1), pch = 16, col = 'goldenrod3', cex = 0.8) + graphics::points(I(datage25m$est*100) ~ datage25m$year, pch = 15, col = 'goldenrod3', cex = 0.8) + graphics::points(I(datage35m$est*100) ~ I(datage25m$year+0.1), pch = 17, col = 'goldenrod3', cex = 0.8) + graphics::segments(x0=datage15m$year-0.1, y0=I(datage15m$ci_l*100), x1=datage15m$year-0.1, y1=I(datage15m$ci_u*100), col = 'goldenrod3', lwd = 1) - segments(x0=datage25m$year, y0=I(datage25m$ci_l*100), + graphics::segments(x0=datage25m$year, y0=I(datage25m$ci_l*100), x1=datage25m$year, y1=I(datage25m$ci_u*100), col = 'goldenrod3', lwd = 1) - segments(x0=datage35m$year+0.1, y0=I(datage35m$ci_l*100), + graphics::segments(x0=datage35m$year+0.1, y0=I(datage35m$ci_l*100), x1=datage35m$year+0.1, y1=I(datage25m$ci_u*100), col = 'goldenrod3', lwd = 1) } } @@ -1039,49 +1039,49 @@ plot_out_evertest_mbyage <- function(mod, fp, likdat, cnt, survey_hts, #' @export ## -- UPDATE HERE -- ## * update yr_pred to current year -plot_out <- function(mod, fp, likdat, cnt, survey_hts, out_evertest, simul = NULL, +plot_out <- function(mod, fp, likdat, cnt, survey_hts, out_evertest, simul = NULL, plot_legend = TRUE, yr_pred = 2022) { -par(mfrow = c(3,2), mar = c(4,4,2,2)) +graphics::par(mfrow = c(3,2), mar = c(4,4,2,2)) plot_out_nbtest(mod, fp, likdat, cnt, simul, yr_pred) - plot_out_nbpostest(mod, fp, likdat, cnt, simul, yr_pred) + plot_out_nbpostest(mod, fp, likdat, cnt, simul, yr_pred) plot_out_evertestneg(mod, fp, likdat, cnt, survey_hts, out_evertest, simul, plot_legend, yr_pred) plot_out_evertestpos(mod, fp, likdat, cnt, survey_hts, out_evertest, simul, plot_legend, yr_pred) plot_out_evertest(mod, fp, likdat, cnt, survey_hts, out_evertest, simul, plot_legend, yr_pred) - plot_out_90s(mod, fp, likdat, cnt, out_evertest, survey_hts, simul, plot_legend, yr_pred) + plot_out_90s(mod, fp, likdat, cnt, out_evertest, survey_hts, simul, plot_legend, yr_pred) } #' @export -## +## ## -- UPDATE HERE -- ## * update yr_pred to current year -plot_out_strat <- function(mod, fp, likdat, cnt, survey_hts, out_evertest, simul = NULL, +plot_out_strat <- function(mod, fp, likdat, cnt, survey_hts, out_evertest, simul = NULL, plot_legend = TRUE, yr_pred = 2022) { - par(mfrow = c(1,2), mar = c(4,4,2,2)) - plot_out_evertest_mbyage(mod, fp, likdat, cnt, survey_hts, out_evertest, simul, plot_legend, yr_pred) - plot_out_evertest_fbyage(mod, fp, likdat, cnt, survey_hts, out_evertest, simul, plot_legend, yr_pred) + graphics::par(mfrow = c(1,2), mar = c(4,4,2,2)) + plot_out_evertest_mbyage(mod, fp, likdat, cnt, survey_hts, out_evertest, simul, plot_legend, yr_pred) + plot_out_evertest_fbyage(mod, fp, likdat, cnt, survey_hts, out_evertest, simul, plot_legend, yr_pred) } # ---- Output Tabular Format ---- #' @export end_of_year <- function(year, value){ - if (length(unique(year)) != length(year)) { print('non unique years'); break } + if (length(unique(year)) != length(year)) { print('non unique years'); return() } new_x <- year + 0.5 - new_value <- approx(year, value, new_x, method = 'linear', rule = 2)$y + new_value <- stats::approx(year, value, new_x, method = 'linear', rule = 2)$y return(new_value) -} +} #' @export -## +## ## -- UPDATE HERE -- ## * update year_range to include current year -tab_out_evertest <- function(mod, fp, age_grp = '15-49', gender = 'both', - hiv = 'all', year_range = c(2010, 2022), +tab_out_evertest <- function(mod, fp, age_grp = '15-49', gender = 'both', + hiv = 'all', year_range = c(2010, 2022), simul = NULL, end_year = TRUE) { interpolate_output <- end_year && fp$projection_period == "midyear" || !end_year && fp$projection_period == "calendar" - + if (length(year_range) == 1) { year_range <- c(year_range, year_range) } if (is.null(simul)) { out <- get_out_evertest(mod, fp, age_grp, gender, hiv) @@ -1104,7 +1104,7 @@ tab_out_evertest <- function(mod, fp, age_grp = '15-49', gender = 'both', } outci$lower <- round(outci$lower * 100, 1) outci$upper <- round(outci$upper * 100, 1) - outci <- subset(outci, year >= year_range[1] & year <= year_range[2]) + outci <- subset(outci, year >= year_range[1] & year <= year_range[2]) tab_evertest <- merge(out, outci) } row.names(tab_evertest) <- NULL @@ -1112,11 +1112,11 @@ tab_out_evertest <- function(mod, fp, age_grp = '15-49', gender = 'both', } #' @export -## +## ## -- UPDATE HERE -- ## * update year_range to include current year -tab_out_aware <- function(mod, fp, age_grp = '15-49', gender = 'both', - year_range = c(2010, 2022), simul = NULL, +tab_out_aware <- function(mod, fp, age_grp = '15-49', gender = 'both', + year_range = c(2010, 2022), simul = NULL, end_year = TRUE) { interpolate_output <- end_year && fp$projection_period == "midyear" || @@ -1125,7 +1125,7 @@ tab_out_aware <- function(mod, fp, age_grp = '15-49', gender = 'both', if (age_grp == "15-99") { age_grp <- "15+" } - + if (length(year_range) == 1) { year_range <- c(year_range, year_range) } @@ -1144,30 +1144,30 @@ tab_out_aware <- function(mod, fp, age_grp = '15-49', gender = 'both', } out$value <- round(out$value * 100, 1) outci <- getCI(simul$diagnoses) - outci <- subset(outci, agegr == age_grp & sex == gender) + outci <- subset(outci, agegr == age_grp & sex == gender) - if (interpolate_output) { + if (interpolate_output) { outci$lower <- end_of_year(outci$year, outci$lower) outci$upper <- end_of_year(outci$year, outci$upper) } outci$lower <- round(outci$lower * 100, 1) outci$upper <- round(outci$upper * 100, 1) - outci <- subset(outci, year >= year_range[1] & year <= year_range[2]) + outci <- subset(outci, year >= year_range[1] & year <= year_range[2]) tab_aware <- merge(out, outci) } - + row.names(tab_aware) <- NULL tab_aware } #' @export -## +## ## -- UPDATE HERE -- ## * update year_range to include current year -tab_out_nbaware <- function(mod, fp, age_grp = '15-49', - gender = 'both', year_range = c(2010, 2022), +tab_out_nbaware <- function(mod, fp, age_grp = '15-49', + gender = 'both', year_range = c(2010, 2022), end_year = TRUE) { interpolate_output <- end_year && fp$projection_period == "midyear" || @@ -1191,13 +1191,13 @@ tab_out_nbaware <- function(mod, fp, age_grp = '15-49', #' @export -## +## ## -- UPDATE HERE -- ## * update year_range to include current year -tab_out_artcov <- function(mod, fp, gender = 'both', +tab_out_artcov <- function(mod, fp, gender = 'both', year_range = c(2010, 2022)) { ## ART coverage is already end-of-year, no need to adjust - + if (length(year_range) == 1) { year_range <- c(year_range, year_range) } @@ -1232,10 +1232,10 @@ tab_out_artcov <- function(mod, fp, gender = 'both', } #' @export -## +## ## -- UPDATE HERE -- ## * update year_range to include current year -tab_out_pregprev <- function(mod, fp, year_range = c(2010, 2022), +tab_out_pregprev <- function(mod, fp, year_range = c(2010, 2022), end_year = TRUE) { if (length(year_range) == 1) { diff --git a/R/retest.R b/R/retest.R index e3b8476..ac3ce16 100644 --- a/R/retest.R +++ b/R/retest.R @@ -120,8 +120,8 @@ plot_retest_test_neg <- function(mod, fp, likdat, cnt, relative = F, add_ss_indices(out_retest, fp$ss))$tests out_retest$first_test <- out_retest$tests - out_retest$retests - col1 <- rgb(230, 200, 0, 200, max=255) - col2 <- rgb(50, 100, 155, 200, max=255) + col1 <- grDevices::rgb(230, 200, 0, 200, max=255) + col2 <- grDevices::rgb(50, 100, 155, 200, max=255) if(relative == F) { if (plot_title == TRUE) { main_title <- expression(bold(paste("Total ", N^o, " of HIV- Tests (in 1,000)"))) @@ -131,13 +131,13 @@ plot_retest_test_neg <- function(mod, fp, likdat, cnt, relative = F, main = main_title, xlim = c(2000, yr_pred), ylim = c(0, max(out_retest$tests/1000)*1.2), ylab=expression(paste(N^o, " of HIV- Tests (in 1,000)")), xlab = 'Year', lty = 1, lwd = 1) - polygon(x = c(out_retest$year, rev(out_retest$year)), + graphics::polygon(x = c(out_retest$year, rev(out_retest$year)), y = c(rep(0, length = length(out_retest$tests)), rev(out_retest$first_test/1000)), col = col1, border = NA) - polygon(x = c(out_retest$year, rev(out_retest$year)), + graphics::polygon(x = c(out_retest$year, rev(out_retest$year)), y = c(out_retest$tests/1000, rev(out_retest$first_test/1000)), col = col2, border = NA) - legend('topleft', inset = 0.02, legend = c('Repeat Testers','First-Time Testers'), + graphics::legend('topleft', inset = 0.02, legend = c('Repeat Testers','First-Time Testers'), col = c(col2, col1), bty = 'n', pch = c(15,15), pt.cex = 2, cex = 0.9) } else if(relative == T) { if (plot_title == TRUE) { main_title <- "Proportion of HIV- Tests \n Among Never and Ever Tested" @@ -146,13 +146,13 @@ plot_retest_test_neg <- function(mod, fp, likdat, cnt, relative = F, plot((out_retest$tests/out_retest$tests)*100 ~ out_retest$year, type='l', ylim=c(0,100), col=col1, main = main_title, xlim = c(2000, yr_pred), ylab='Proportion of Total HIV- Tests (%)', xlab='Year', lty=1, lwd=1) - polygon(x = c(out_retest$year, rev(out_retest$year)), + graphics::polygon(x = c(out_retest$year, rev(out_retest$year)), y = c((out_retest$retest/out_retest$tests)*100, rep(100, length = length(out_retest$tests))), col = col1, border = NA) - polygon(x = c(out_retest$year, rev(out_retest$year)), + graphics::polygon(x = c(out_retest$year, rev(out_retest$year)), y = c(rep(0, length = length(out_retest$tests)), rev((out_retest$retest/out_retest$tests)*100)), col = col2, border = NA) - legend('topleft', inset = 0.02, legend = c('% Repeat Testers','% First-Time Testers'), + graphics::legend('topleft', inset = 0.02, legend = c('% Repeat Testers','% First-Time Testers'), col = c(col2, col1), bg = 'grey90', box.col = 'grey90', pch = c(15,15), pt.cex = 2, cex = 0.8) } else { print('True for relative scale, False for absolute scale') @@ -184,9 +184,9 @@ plot_retest_test_pos <- function(mod, fp, likdat, cnt, relative = F, ylim <- c(0, max(out_retest$tests/1000, na.rm = TRUE) * 1.3) - col1 <- rgb(230, 180, 205, 200, max=255) - col2 <- rgb(150, 0, 50, 200, max=255) - col3 <- rgb(100, 0, 100, 200, max=255) + col1 <- grDevices::rgb(230, 180, 205, 200, max=255) + col2 <- grDevices::rgb(150, 0, 50, 200, max=255) + col3 <- grDevices::rgb(100, 0, 100, 200, max=255) if (plot_title == TRUE) { main_title <- "Distribution of HIV+ Tests by \n Awareness and ART Status" @@ -199,32 +199,32 @@ plot_retest_test_pos <- function(mod, fp, likdat, cnt, relative = F, main = main_title, xlim = c(2000, yr_pred), ylab = expression(paste(N^o, " of HIV+ Tests (in 1,000)")), xlab = 'Year') - polygon(x = c(out_retest$year, rev(out_retest$year)), + graphics::polygon(x = c(out_retest$year, rev(out_retest$year)), y = c(out_retest$newuna/1000, rev(out_retest$tests/1000)), col = col1, border = NA) - polygon(x = c(out_retest$year, rev(out_retest$year)), + graphics::polygon(x = c(out_retest$year, rev(out_retest$year)), y = c(out_retest$newdiag/1000, rev(out_retest$newuna/1000)), col = col2, border = NA) - polygon(x = c(out_retest$year, rev(out_retest$year)), + graphics::polygon(x = c(out_retest$year, rev(out_retest$year)), y = c(rep(0, length = length(out_retest$tests)), rev(out_retest$newdiag/1000)), col = col3, border = NA) if (plot_legend) { - legend('topleft', inset = 0.02, legend = c('Retests - PLHIV on ART','Retests - Aware Not on ART','New Diagnoses'), + graphics::legend('topleft', inset = 0.02, legend = c('Retests - PLHIV on ART','Retests - Aware Not on ART','New Diagnoses'), col=c(col1, col2, col3), bty = 'n', pch = c(15,15, 15), pt.cex = 2, cex = 0.9) } } else if(relative == T) { plot((out_retest$tests/out_retest$tests)*100 ~ out_retest$year, type='n', ylim = c(0,100), main = "Distribution of HIV+ Tests by \n Awareness and ART Status", xlim = c(2000, yr_pred), ylab = "Proportion of Total HIV+ Tests (%)", xlab = 'Year', lty=1, lwd=1) - polygon(x = c(out_retest$year, rev(out_retest$year)), + graphics::polygon(x = c(out_retest$year, rev(out_retest$year)), y = c((out_retest$newuna/out_retest$tests)*100, rev((out_retest$tests/out_retest$tests)*100)), col = col1, border = NA) - polygon(x = c(out_retest$year, rev(out_retest$year)), + graphics::polygon(x = c(out_retest$year, rev(out_retest$year)), y = c((out_retest$newdiag/out_retest$tests)*100, rev((out_retest$newuna/out_retest$tests)*100)), col = col2, border = NA) - polygon(x = c(out_retest$year, rev(out_retest$year)), + graphics::polygon(x = c(out_retest$year, rev(out_retest$year)), y = c(rep(0, length = length(out_retest$tests)), rev((out_retest$newdiag/out_retest$tests)*100)), col = col3, border = NA) - legend('bottomleft', inset = 0.02, legend = c('% Retests - PLHIV on ART','% Retests - Aware Not on ART','% New Diagnoses'), + graphics::legend('bottomleft', inset = 0.02, legend = c('% Retests - PLHIV on ART','% Retests - Aware Not on ART','% New Diagnoses'), col = c(col1, col2, col3), bg = 'grey90', box.col = 'grey90', pch = c(15,15, 15), pt.cex = 2, cex = 0.8) } else { print('TRUE for relative scale; FALSE for absolute scale') @@ -277,9 +277,9 @@ plot_prv_pos_yld <- function(mod, fp, likdat, cnt, yr_pred = 2022, ylim <- c(0, ifelse(max(c(prv, yld, ndx) * 1.3, na.rm = TRUE) > 1, 1, max(c(prv, yld, ndx) * 1.3, na.rm = TRUE))) * 100 - col1 <- rgb(255, 155, 50, 250, max = 255) - col2 <- rgb(255, 0, 130, 250, max = 255) - col3 <- rgb(100, 0, 100, 250, max = 255) + col1 <- grDevices::rgb(255, 155, 50, 250, max = 255) + col2 <- grDevices::rgb(255, 0, 130, 250, max = 255) + col3 <- grDevices::rgb(100, 0, 100, 250, max = 255) if (plot_title == TRUE) { main_title <- "HIV Prevalence (15+), Positivity, and \n Yield of New HIV Diagnoses" @@ -291,12 +291,10 @@ plot_prv_pos_yld <- function(mod, fp, likdat, cnt, yr_pred = 2022, main = main_title, xlim = c(2000, yr_pred), col = col1, lty = 2, lwd = 1.5, ylab = "Proportion (%)", xlab = 'Year') - lines(I(yld * 100) ~ out_test$year, col = col2, lwd = 1.5) - lines(I(ndx * 100) ~ out_postest$year, col = col3, lwd = 1.5) + graphics::lines(I(yld * 100) ~ out_test$year, col = col2, lwd = 1.5) + graphics::lines(I(ndx * 100) ~ out_postest$year, col = col3, lwd = 1.5) if (plot_legend) { - legend('topright', inset = 0.01, legend = c('HIV Prevalence (15+)','Positivity','Yield of New Diagnoses'), + graphics::legend('topright', inset = 0.01, legend = c('HIV Prevalence (15+)','Positivity','Yield of New Diagnoses'), col=c(col1, col2, col3), bty = 'n', lty = c(2, 1, 1), lwd = 1.5, cex = 0.7) } } - - diff --git a/R/simul.R b/R/simul.R index ae50735..1563111 100644 --- a/R/simul.R +++ b/R/simul.R @@ -19,8 +19,8 @@ simul.sample <- function(hessian, par, fp, sim = 3000, likdat = NULL, SIR = FALS ev <- eS$values if (!all(ev >= -1e-06 * abs(ev[1L]))) { vcova <- Matrix::nearPD(vcova, corr = FALSE)$mat - vcova <- matrix(vcova@x, - nrow = vcova@Dim[1], + vcova <- matrix(vcova@x, + nrow = vcova@Dim[1], ncol = vcova@Dim[2] ) } @@ -33,7 +33,10 @@ simul.sample <- function(hessian, par, fp, sim = 3000, likdat = NULL, SIR = FALS # Sampling Importance Resampling # http://www.sumsar.net/blog/2013/12/shaping_up_laplace_approximation/ if (SIR == TRUE) { - if (is.null(likdat)) { print('SIR needs to include '); break } + if (is.null(likdat)) { + print('SIR needs to include ') + return() + } # Sample parameters from multivariate t-distribution to get thicker tails par_sir <- mvtnorm::rmvt(n = nsir, delta = par, sigma = as.matrix(vcova), df = 2) # We had the mode to the resampled values (to be sure that it is included in the CI) @@ -99,7 +102,7 @@ simul.run <- function(samp, fp, progress = NULL){ nbtestpos_ss <- add_ss_indices(out_nbtest_pos, fp$ss) nbtest_ss$hivstatus <- as.character(nbtest_ss$hivstatus) nbtestpos_ss$hivstatus <- as.character(nbtestpos_ss$hivstatus) - + # Create parameters (proper scale, etc.), and simulate model for (i in 1:nrow(samp)){ fp <- create_hts_param(samp[i,], fp) @@ -123,7 +126,7 @@ simul.run <- function(samp, fp, progress = NULL){ # Simulation of the dataset #' @export simul.test <- function(opt, fp, sim = 3000, likdat = NULL, - SIR = FALSE, nsir = 50000, with_replacement = TRUE, + SIR = FALSE, nsir = 50000, with_replacement = TRUE, sir_progress = NULL, run_progress = NULL){ samp <- simul.sample(opt$hessian, opt$par, fp, sim = sim, likdat = likdat, @@ -131,13 +134,13 @@ simul.test <- function(opt, fp, sim = 3000, likdat = NULL, simul.run(samp, fp, progress = run_progress) } -# Get the confidence interval for the data, first for parameters, second for +# Get the confidence interval for the data, first for parameters, second for # the data itself getCI <- function(df) { n_col <- ncol(df) # Function for the confidence interval of the estimates CI <- apply(X = df[,(6:n_col)], MARGIN = 1, - FUN = quantile, probs = c(0.025, 0.975), na.rm = TRUE) + FUN = stats::quantile, probs = c(0.025, 0.975), na.rm = TRUE) lower <- CI[1,] upper <- CI[2,] df1 <- data.frame(df[c(1:5)], lower, upper) diff --git a/R/time_functions.R b/R/time_functions.R index fc27a96..27a2977 100644 --- a/R/time_functions.R +++ b/R/time_functions.R @@ -54,9 +54,9 @@ prb_dx_one_yr <- function(fp, year = c(2000:2022), age = "15-24", sex = "male", ind_yr <- year[n] - fp$ss$proj_start + 1 # what is the average age of incident HIV infection - avg_age15 <- 14 + weighted.mean(1:10, w = fp$infections[1:10, ind_sex, ind_yr]) - avg_age25 <- 14 + weighted.mean(11:20, w = fp$infections[11:20, ind_sex, ind_yr]) - avg_age35 <- 14 + weighted.mean(21:35, w = fp$infections[21:35, ind_sex, ind_yr]) + avg_age15 <- 14 + stats::weighted.mean(1:10, w = fp$infections[1:10, ind_sex, ind_yr]) + avg_age25 <- 14 + stats::weighted.mean(11:20, w = fp$infections[11:20, ind_sex, ind_yr]) + avg_age35 <- 14 + stats::weighted.mean(21:35, w = fp$infections[21:35, ind_sex, ind_yr]) # Testing rates by CD4 category test_rate15 <- fp$diagn_rate[, 1, ind_sex, ind_th, ind_yr] @@ -186,7 +186,7 @@ pool_prb_dx_one_yr <- function(mod, fp, year = c(2000:2022), denom$w_e[denom$agegr == time_mat$age[i] & denom$sex == time_mat$sex[i]] } } time_mat_sel <- time_mat[time_mat$age %in% age & time_mat$sex %in% sex, ] - val$prb1yr[val$year == year[n]] <- weighted.mean(time_mat_sel$prb1, w = time_mat_sel$w) + val$prb1yr[val$year == year[n]] <- stats::weighted.mean(time_mat_sel$prb1, w = time_mat_sel$w) } return(val) } @@ -217,7 +217,7 @@ simul_pool_prb_dx_one_yr <- function(samp, mod, fp, year = c(2010:2022), prb1yr = NA, prb1yr_lci = NA, prb1yr_uci = NA) for (n in 1:n_year) { - uncertainty <- quantile(val_df$prb1yr[val_df$year == year[n]], probs = c(0.5, 0.025, 0.975)) + uncertainty <- stats::quantile(val_df$prb1yr[val_df$year == year[n]], probs = c(0.5, 0.025, 0.975)) val$prb1yr[val$year == year[n]] <- uncertainty[1] val$prb1yr_lci[val$year == year[n]] <- uncertainty[2] val$prb1yr_uci[val$year == year[n]] <- uncertainty[3] diff --git a/README.Rmd b/README.Rmd index 8a0e0ba..70a77d6 100644 --- a/README.Rmd +++ b/README.Rmd @@ -18,7 +18,7 @@ library(data.table) ``` # first90 -[![Travis-CI Build Status](https://travis-ci.com/mrc-ide/first90release.svg?branch=master)](https://travis-ci.com/mrc-ide/first90release) +![R CMD check](https://github.com/mrc-ide/first90release/actions/workflows/R-CMD-check.yml/badge.svg) UNAIDS put forward the ambitious 90-90-90 target to end the AIDS epidemic by 2030. This target aims for 90% of people living with HIV (PLHIV) to be aware of their HIV-positive status, 90% of those diagnosed to receive antiretroviral therapy, and 90% of those on treatment to have a suppressed viral load by 2020 (each reaching 95% by 2030). HIV testing remains an important bottleneck in this cascade, however, and obtaining reliable epidemiological data on the proportion of PLHIV aware of their status is difficult. Such information is nevertheless crucial to effectively monitor HIV prevention efforts. Tracking progress towards achievement of this “first 90” target could be improved by combining population-based surveys and programmatic data on the number of HIV tests performed (and yield) in a coherent deterministic/statistical model. This type of integrative systems modelling is especially useful to fully consider HIV incidence, mortality, testing behaviours, as well as to coherently combine different sources of data. @@ -146,18 +146,18 @@ first90::plot_out_evertest_mbyage(mod, fp, likdat, cnt, survey_hts, out_evertest We can compare HIV tests' positivity through time, the estimated *true* yield of new HIV diagnoses, and compare those to population-level HIV prevalence. ```{r, fig.height = 5, fig.width = 5, warning = FALSE} -par(mfrow = c(1,1)) +graphics::par(mfrow = c(1,1)) first90::plot_prv_pos_yld(mod, fp, likdat, cnt, yr_pred = 2018) ``` We can also examine some ouptuts related to the distribution of HIV tests performed in those susceptibles to HIV and PLHIV by different awareness and treatment status. (First sets of graphs is on the absolute scale, second one on the relative scale.) ```{r, fig.height = 5, fig.width = 8.5, warning = FALSE} -par(mfrow = c(1,2)) +graphics::par(mfrow = c(1,2)) first90::plot_retest_test_neg(mod, fp, likdat, cnt) first90::plot_retest_test_pos(mod, fp, likdat, cnt) -par(mfrow = c(1,2)) +graphics::par(mfrow = c(1,2)) first90::plot_retest_test_neg(mod, fp, likdat, cnt, relative = TRUE) first90::plot_retest_test_pos(mod, fp, likdat, cnt, relative = TRUE) ``` diff --git a/README.md b/README.md index 50a4f3d..92eb9f5 100644 --- a/README.md +++ b/README.md @@ -4,8 +4,7 @@ Please edit that file and run `knitr::knit("README.Rmd") from the root of this d # first90 -[![Travis-CI Build -Status](https://travis-ci.com/mrc-ide/first90release.svg?branch=master)](https://travis-ci.com/mrc-ide/first90release) +![R CMD check](https://github.com/mrc-ide/first90release/actions/workflows/R-CMD-check.yml/badge.svg) UNAIDS put forward the ambitious 90-90-90 target to end the AIDS epidemic by 2030. This target aims for 90% of people living with HIV diff --git a/data-raw/initial-theta0-values.R b/data-raw/initial-theta0-values.R index f7fe1c3..6490a11 100644 --- a/data-raw/initial-theta0-values.R +++ b/data-raw/initial-theta0-values.R @@ -28,13 +28,13 @@ theta0 <- c(par_med_f15to24rate, rep(par_med_f15to24rate[n_k_2021], n_k - n_k_20 # old theta - keeping for comparisons theta1 <- c(log(seq(0.01, 0.2, length.out = n_k)), # Females 15-24 baseline testing rate 2000 to 2020 - qlogis(rep(1.5 / 8, n_k - 10)), # Factor testing among diagnosed, untreated from 2010 to 2020 - qlogis(rep(0.6 / 1.1, 2)), # RR for males in 2005, 2012 - qlogis(rep((1.93 - 0.95) / 7.05, 2)), # Factor increase among previously tested - qlogis(0.5), # RR re-test among PLHIV - qlogis(0.25), # Factor testing among diagnosed, on ART - rep(qlogis(0.9 / 5.9), 4), # Rate ratio for 25-34 - qlogis(0.5)) # Date of inflection point for logistic growth curve OI + stats::qlogis(rep(1.5 / 8, n_k - 10)), # Factor testing among diagnosed, untreated from 2010 to 2020 + stats::qlogis(rep(0.6 / 1.1, 2)), # RR for males in 2005, 2012 + stats::qlogis(rep((1.93 - 0.95) / 7.05, 2)), # Factor increase among previously tested + stats::qlogis(0.5), # RR re-test among PLHIV + stats::qlogis(0.25), # Factor testing among diagnosed, on ART + rep(stats::qlogis(0.9 / 5.9), 4), # Rate ratio for 25-34 + stats::qlogis(0.5)) # Date of inflection point for logistic growth curve OI stopifnot(length(theta0) == n_k + n_k-10 + 11) stopifnot(length(theta1) == n_k + n_k-10 + 11) diff --git a/first90release.Rproj b/first90release.Rproj new file mode 100644 index 0000000..497f8bf --- /dev/null +++ b/first90release.Rproj @@ -0,0 +1,20 @@ +Version: 1.0 + +RestoreWorkspace: Default +SaveWorkspace: Default +AlwaysSaveHistory: Default + +EnableCodeIndexing: Yes +UseSpacesForTab: Yes +NumSpacesForTab: 2 +Encoding: UTF-8 + +RnwWeave: Sweave +LaTeX: pdfLaTeX + +AutoAppendNewline: Yes +StripTrailingWhitespace: Yes + +BuildType: Package +PackageUseDevtools: Yes +PackageInstallArgs: --no-multiarch --with-keep.source diff --git a/man/get_pjn_country.Rd b/man/get_pjn_country.Rd index abcbd9f..48924a3 100644 --- a/man/get_pjn_country.Rd +++ b/man/get_pjn_country.Rd @@ -13,5 +13,5 @@ get_pjn_country(pjn) Get country name from parsed PJN } \details{ -\code{pjn} should be via \code{read.csv(pjn_file, as.is = TRUE)} +\code{pjn} should be via \code{first90_read_csv_character(pjn_file)} } diff --git a/man/get_pjn_region.Rd b/man/get_pjn_region.Rd index 28bd522..314048d 100644 --- a/man/get_pjn_region.Rd +++ b/man/get_pjn_region.Rd @@ -13,5 +13,5 @@ get_pjn_region(pjn) Get subnational region from parsed PJN } \details{ -\code{pjn} should be via \code{read.csv(pjn_file, as.is = TRUE)} +\code{pjn} should be via \code{first90_read_csv_character(pjn_file)} } diff --git a/man/prepare_inputs.Rd b/man/prepare_inputs.Rd index 62f984d..11fed31 100644 --- a/man/prepare_inputs.Rd +++ b/man/prepare_inputs.Rd @@ -21,7 +21,8 @@ The aggregation makes a number of assumptions: } \examples{ -pjnzlist <- list.files("~/Documents/Data/Spectrum files/2018 final/SSA/", "CotedIvoire.*PJNZ$", full.names=TRUE, ignore.case=TRUE) +pjnzlist <- list.files("~/Documents/Data/Spectrum files/2018 final/SSA/", + "CotedIvoire.*PJNZ$", full.names=TRUE, ignore.case=TRUE) pjnzlist <- "~/Documents/Data/Spectrum files/2018 final/SSA/Malawi_2018_version_8.PJNZ" } diff --git a/tests/testthat/Rplots.pdf b/tests/testthat/Rplots.pdf new file mode 100644 index 0000000000000000000000000000000000000000..6b077ad6346af3a6b8365a5feb4cc8d3bba99a3d GIT binary patch literal 11830 zcmeHNcUV)|wnsWh4NZDLL68y>dQ-X}9qGM=A|XH$dhb#N5drC4L7MccNS9tj5Tq-; zgH-Va=g!O>eecfP|6ji3J73mWYp;D)_AjflX~-$?fO!Q7*#Z{=7Xt?aM~xALU?3mR z$?^fAgajc-76G+D**n4IEKpD&hnyIM4)@T;3Jn<>!T+C2b!xS2&!isBkwv{@P1el))ScWVqu)dNIk&DdL!X2#Uw0kxYh-b3d{E&poQ9;1n zGAW)+lD%Y0QqL7G6l^lXR|HMyA zuoL0g7L(PfZ|PuHO4*9@=5ws%_2oKN{`WZ(1VY$6*m~Ije4&ux$M1u~oW-j8>*+0E z0}1YgM!v~d?9Yf!;wuu2&8u4k=jL z!m^;uwdO2jFoGH6EExobPpQLFC0!#50~|H<6vj`Bi$I6O;X zr^%J0s+b^Owtdssx^PkBa$kP8&2;x7an~@q{()scxTw$ZF)<%?ctxMk6tU0N&P2!2 zVN=_uD6W@DmEv`u-lj%Xx}=4xLpEgBU51Z$D%vg2}}1h zX|o-2VIN=Eo*-*t=d+T^RRfLH{mP_T-V!>MHa0AfX9};z5_Pu^rVcYkXz!jg8r%_C zL1vqNTy<(KXPdn9-23xuG(`V(-2Vku(Lb>AJ6QiC_P~692fiO1KcP=ZP~;z=PiG2F z2PX49((AjzqEpi6T=zp8g~WF>w=gSb(06!Y2TsghULD zz#;PF8Ym`r=Ic0Z#mM-Mxoh^;`GicwR06A$SAuCbN}{h&Im&o$W|*{uTyq(JHlE8g zKEP4(>5OSZbmV2pV#Y&*mv_~UEV;EJ50O4D9cIm;cNV^o+>pS>NjguxN=2}voj*tb<_C2zeli{@qTEvP;T6#XzJ1pFTegF z_(%ZX7{ODIAc&vew7UQK>jt5$bfQ7Ld@inc7?rd#JcbFMb{w0xTyQW1n=44uak?fn z&O&vE$yZ*noOS6XVYb?z#-IT{n(^Z-{ViX&}(zpsv|W3W@jWbKLrFVv~aM`+K>?=`vIm z(-eiab6o9vK2svvn&%!GJ`?mquA#Tk4d}jt9>Lrb>Q3_AgKY|LJD`8CEVGp%(BzI4ba; zu96V+|DNjqJ=Onvs{h+lfiJ_x-yF+7gWqpq;eUePCC4wWO5pEZRlbg@69t&8dAlt6 zDhJo_3)}O0JQ4xSynAv_UP?XR%&v;adH{GvxqC#r?M356nIpika2U1hX_L{#6odrU zg?BArGPUYHC-R8t3}6Y)koK18EzfAB(O1Y_K<%%7qm3>MmC2`ET0}|2`!|Qqc6=7( zsN=#P`0Q?`ArL?sDr~2{S5C%EmqOu#kqj@X?G|0CQfbW<`sNjs5~xckdVot{gNv2q zi=EPHkPd13w9)O1Y~)s(t+3M_47@ffBVIK@ehaVj?&KOyv_7YijIUW!NoxO5r4r() zy4p=wk@VoUcSNg$JEg^hve$;{k7jhk6$kD;afttVJHtJ9dR)w6MWaS#kD$OzosZwV z-_FM=&fF#|Q9fq|*gBmR^Fp{=rjM$7TmnKEN!w-VK=?Q0%X zw&qb8MvUsL(E0M-btRB)wL3)ZV|*4ITlOtIbAPwPvbUa4}My+#5V1t z*LdSuhYP~&)5m9q_qc~!`@NQLf9D`gdFS#ttvc7#XF=^-K1RxSL!ZR45i@=fa%InBngsA28_DmX2~r!Rsz6+*50} zV#g;6x-_}vR$a^&jLw+7)tL5lDTV-D3(rj)HK35Ob@N^3nxmx^>5kXtylP zSK^NuqN*cIJJfH0#P;NO@koM7B7b#Lf`6Rqe})sGzjafW96#ek=pW%k^AP|{<~M5w z2GI%ahN(?CHj>`af?1ibZdSIc#mWUV{z8KVe_d>obz&<+u3H{hm@0oO#o4# z`IScx4P-wCxF;*D(skEC+}s9t~X}!$v$(DxVEcj)S9{Qg&KRB23yvT8pjjwSVKR? zP!62LT`0qO-3Q|u^CT7o7b?(oF76$Y5)8o1S0*91A;Dw0kDJownjm~oCo?;AA>K2U zn}Vej9ma&b}G?{GJ_k$5esIxkPH6jo)V!AZaw8Uw>Wplp^O$SKM7H ziPK92pSzy&rLojI!`W@!`s?D^>l!xfexm6cXeC-oV?v^G$)X&ZESiXieD3(7Y~kD; z#(^(4idBeJ?TAG_tB2~%NL1X9rDd`fy@h}%>wH>_Pt+E*9fy7Af`^=mD)CmiziTeVAbnE<0tAcDhI#Wj5 zTTEF^p*G16BwquV-CH`2xDFNYWL&LX|4e$#pO0b^a%nK5iV z+0otTSby*7mlGS)po0#_rsR56zFeX-%yV!g zBHbTGnt>3H@OI>gP=9sbzawDzb<(UYz1eVS33TMF9@DFXc*e%S?BM>lZHEc+R|OXz zpfWAFf+$WT!gsuWZ}NN42Mof}W6uipisP%qaoKNvMa+L$fB^q(1>y(C&xjG_|A(GM znh5P#c5t2~K3`1P_1Or-iotjRAZOzud*H=k&-(mX2)q2buS8?zloDkQF3H*MxS!Qo zOjq$7CgpUitP26DIoAzbr#M=DuorPBR~4NQ@))z&vhLf#1%u60&=3cv4E=&6zJs_C z6~-+-G2DLKl3?gX2djZJbuh(}hPX5nvAf;42Od#DtsW)CS{>tqad=ZBh>lajg)?mB z16-rgEi!$Bv#h;;KZ7dj`Q2OY0`Z6l%WUnvi#?flC{Oc2~C6 z_bhqk)q)LgdXVx>R@QbaP4=YiZEq0Rkc8p8t}p`n7~EvO(pr-&yS8MXQFRhCwv3K2 z9b!J|w6=5NaU<-xayG$o2*7>T{ZwO?%4ywZcM--W*n7Q@%j>oYUO3hEi26f#5uxlv z_WErvZ!$TAYL0PK^pVSJ`X8B4IHdCN-MN_z^DkL33SXOLgHnQ7qV64%5PTtVRBM)f zYFR3!`_-TTm+OEeD58tymb1oru5pKs1S0{z$@GRC5FNK!tW}WNSy!)`^K+G}&cfzxTqi2<^@4D?0B}_sGx`@hfT?rr# zcq(F}@3LuTJ7ndJ7P2WZH2=GxaAt?dSy)0QOFRm`!^e9%HfT!9Or zlO1L6?qV4F@6N4=#%xvXTOx+%_v&+@+4e}LRKtQ`W^}4!9b9WLK6qMual}AED1V>U zF41A1u`FBJ;4yvVfRN313RA3WB3Ps#HS(u2=kyNuI|yi|%w66lv&TA;eQzk(%!b}= z;l9AQd%FSjE1dqzngZl+YYIO&euk5<@IQf*@@1gb@>`(R-M)ClBr^0+4@Xm3R~{Xw z-HK6sFZ4hdb9bAzYbHjRP*aO!*K6VA%dV${tTPdI?T&uh+v~}Mp}2(_4nujvN*b|4 ziXIJtPCH~zA~p99+sM&@niLE6(mG0F-akEXw&R2i|5%fZVSl?OD+|MnuuwBqH zlG2a#8CB{j`;4DpN@-cc7t@AnwcJ!1WEyX)TNVJVZj{IXkf9n>Rk{x?3$mWbzt&Z> zEF}xxqI|`Gh=iDi(V9FPGT0$(Hqop)baS&*QhuquaAaOZiK`>GZktsiOEdh{xk?LD zsch(i*3K;S0lpo89-7Exo$c9H^n$G_XPl0!7ook?p6LV4zPY(?r*q4vS@T$|RL4&6 zBSN%+Z2lyRxAp2Z?x%Y08m2P7)1Mor9t+z+ZS470$=?Tc;EI$c5~%8lNAZp{O%Tqk zrDTeEz3;t3&3}BnxDbpxxp{tNMR3+CLB8=r)#77}qJuna ziUru0-jlMEJ<=wNmSK{%y#);3Tqs;3F(&Nh)s;!!6gayZYA>bU9Aq74Jf9#`Z&U5WK2 zaff4*zjX;3Z5pDEF*WriOdKf^~r^dG@TonD#<4WFHS z9bFz9ri>5S;ef-?kND51bZ{8RaLBRc$~B8G&b-`ggcCSo6urBxn=fuZ`sUfmHBWY9 z`mL2K(IazALd?Q=d_{m5&POEy!i>^6dI{~b%@3q_&Xc;B#I$4_)n9bDa5g{7vNg;0 zYfEngoZZ6>Yz-zMX$}+r>VUGepMB#xEga6MF8Y!T^PQ_$#*UHgoSv%}38tSmfaJv8 zs0Tj)lz9C{d_v1|l^{Z>rwJ>$^dtewQ?xN4dUOG!#^AxebqC{>0fuf31B3$uPnC1p zD^#=Y+fn4yka4hfGNZjxO(F+(TvWQ*j!036*VB7wO^7}yo;Tc3GI?Dl_f$}uS{}Z+12i51E%0e zIYjvpblVUpDv<`{glXl-BsI(MK3?J8`Mgt;U}xoCgJ9z%YPjySwW(VJz^xwmG@0(m z*w8-SSXuprQ-(41B9M_%Jtx+2TurCt5!Ei%^?kV&Ro1uWlWKFVP4lrG@ds(*#rySz z)kt-FVwaJ)&!a1rL4~;mxFB-@GvjG#Uy}ZLjkT|eGg&(FdBfK7g&ExA-S))nu`UQ| zlM(xt_$C+ImQ@#xB*8c5bhrJ)JaLSPW?v2KG}_I>P!RYWTT^=U;M)6lHIJxe?HFjT z-QRfYU*@@f?dyWQRoZI?T)(9)rNKTg4Lu5drXl`)s6!hd48{~|Hy%V~DpKn}Zx^)OvB)yUO@SZIXxT*e_?iGEKl zAAw84)G%iwm5h{MMr`MKS;?Hbiyr}|9~B_>7^&B7jxkH9tjNkHxp)Z?%_ck%lhM09 z*R++=5%5E$KB6kJ=f<{i?R{A^uX=&}A}^5(<6K1T#CLABgn40VS!lIif2WYc_vYoE z>KOku4wfRoTH}IjGpdF05+m>CwKhof2C@oo`#f6t=#47H(=$1ka8{`Y{BT7*-P6xZ za~$$jczSXTRuhtGpC>Ch#I;WYQK81Pa$x?Uk6A{vH^9OVmE3U+h(+aTL3p$#5~LC! z35Ee~&W|4AMiq}Kz0!|A(U-`Dish{0*tIP(4cquvjW@~LMYw^a#F-mdvhwvEiUs0u zc1eX)xw5Gu@O|}VhU+MtPTsv$`v|pjec)s0;1gL$S{z1Hn@dcaLAnM=3uyImaV7oC zHg?yfWJ_>;3yNb20IrW5cVa#C(}=mO-{o0a5c{ z`e}hETwYw;TsqUXv)z@P4ZRxQTJ4a;&zQ22%13vUNf@V1B^?j;AHnA`eY{6@jh4hp z^B?rUk*0aDLQ>%Do}`dLjksjM_vUx$M7>ITnDa--lob~Oo$ksM-h$VKQMMj<_i}IJ zcNI(4MQ3>KdQdPvS{~XrfqYWg*g57xwrXw!=J{fB%xmn{t;3hty_;Br-Ueyic)L&c zrjpu6>m28;QcLzP$ob!ag9L(fPwer40IVSyn(_e&YcD z#`42gX#pX>v6!Gcy-u#^jxE}?s6wsnEo7WLfJT?mogkmU<;UdTn!f*Ddqek@t?oln zKqC;kyNxbdKs`_7{!v`+>Cfi!Q+MtAlhomZ&BDzZkI5pR z3^~|oVJZ?gE5`HQ0uzwkdlQGjDkP`vETF;TYWcvpya zFtA`4bIz5OGStTkm?@Y(6Vb2~iesAQ1`^=O*{#t#(XANYGoREOu!1o(lB1XRjwg+2 zX?jE?_0rgC)q9d_1XKt3?hUXS1VMBYkbB%lK3`Lyq@w$f5DD5|*kiUc z0@Cjrr06B>80a(P1>7inNXu?~KP)DJB5HncP<1(~qmPe7IOh5Scp`NxY+F=^XH-empLh>yqQLr?WF^I7Pm{FpE}Zjl2Cysf#0wtcU5-o-lU zj~Y%%PCpwAy71MjQ@QA-zTqEhbyp5j2hbeYcy?}xDFDQfm%>_;qF=$Hc!OQ@hP)1g zbPfwsn&Tn1i#d2Mh%)Yq9FT@GpgQZ?lb{mwo6m7wv&aR4g`BUU0T5cWxtg@ksKQE7#vfU?~TuR`_X=$Ei?WXUxcp5U{|$|zz|Nv&sTJ!dEi zTg$dtFS^iTAicK+aMe2gt*;m--mTl*e2jsTg0_ zZA6JHFAvsGd2;m~Ggat9&eLq&Y@ysDO+1bC=eNgn#yB1j9+ImnIOOyX53FFIxW#DG zqj(kcvQmc_S1MP?4@6F>aD%h0qQSA({8bp(sVXrmX)EbRV?{%0JEJWMD)e89jWeuA z5_cfVEn(mqLT@r}fH#hKtZ(r31(b7*h?%j3+C#2(ojF72o6dU96m_&sfElvwXwr`T zo!(vOv4Rh&AO4-#OcoA`XoAuZMVVU`JlT(j6hxHkls%+x01h9jwL9ko)P?CQNytg4 zGqAO=X)4;>5aJBH!Q#lWc{82$6KgXtp+~`*6d_@tG9v?L&t#ds`C2pUmi{eljvY!} zt}zbOEn`iOf`@tn+6UUkFO$pkvpwWo6n9?$OnVd~O4>CI^R;tP+N)tBRS=y?1yOBL zS^~yOi7=4;N(TH6uCBIBYqrv_p#2q4(>@T*nmF4Aw-5W%&zmhP3SStDNkD zt2n75Bh&9y-OIQa^Uj5elQElq7;9zq0CH+^>Ut`!!=t0F)3V0?Ug$knQQ7-1&pn^- zzP~s2@x8*>%nSZ8-?7oL+86FGW=Ey$lkGzsmL6;jixlsoNE~B(8Me|Mi-s>xkB12M z3tSh_v)dYF7)l!-8Q&XyH&QermV2n3TDqB?Ueqjm7VPa092GEulu^eTYSgp@WIG}ds0Cz$6dPa}7x0j++A zkz?gr?WvukQ#$;Uoo`o|C57F6yHTtC;*u1e4>g@%KN{~D&(;?&yAArnZQ)(DQ&#Sy z84s8(U$?e=xCjo;4erkFw-_0I10Jpz-^GpUB(!q0Nq3hK_Y^I06Z>*7KG|xp;WKeU zLsOvqOgXAI*Ryl~c%N=Lrhcs+VZ3gV4=1XXt8IpBxtqAhdIz6jk@}Km+)d~TDSTmD zv@AVnYdUVqT1Ha#0cQL4-d9ei_Nl?fcJHC^1&Yu2H<7b?gTj-w@7>=G&Nzk6a{5>K z$9)(2PJhMl%5=b&fEQ9-ZxRE2rSr0Iq~A-E2K5KYw@R`&E1GwIjceQ|UcS2Q9vU%okoVeMr&n+GVwvyO_FANOQ_PUGFT>dOKTW{*-l5-ZP4qdt5)Q z<@gxT!Is8$7cd+Fj^K&d7R^q>9XJ>X$gfr|tb=XwZc$|5IpH}y-hKyKR=-ito_JI8 zo}21ShD6FjhF1KZ_8V&Ls~Lr-h2!}xYFAaKGQV<4vNua&CB-MZXe*Kz5S&4(`S&5= zm3@_wm3h5sNt?z(Y=i(osO)6atHut2$-%LJcDYxkkY|fe?7oXqPd z_f>?)mP&HVrgvH5oR_GAgezN98ZSJJlH~$}jfPTQB1^lx*pA>JO|AskV3U=cwZ4 z2Km(U+L)-v>kdE5sW|joH;ST)l)x!> z7kJBYV!vLZneWHPRnYRjY}lY~pYHtl`qXqSBZnUEAc%R2=~Xil&+2YZY?64=Mw5nx zo@Q-Sx$-MzLt5iW6R1XY^_s8iy+g$qJyq`#o6Ulc9~wOy%Qki=Ms*8y^bC*SW#-%c zUawCFaC5^VeUp7f&zOW1gunXtpElGh<-Vv|ZS$M=TMc}O<=Y?I-_P@m2U_A>oaee{ z{uyeFoE4k2daj&K#asTpgf^HMKU`XjpQ7fqdhawTO4dEQZU7D;~}e z-T)Fe*`m!5_F0HLOB~f}F3V`7cHSK6K98SY*Zkl#qgcn%+|zinqI+7hT#JMcP8%R; zd~267Mw%nfo*mKq34E(CSp!PmUXv8tF>ehPx0RWI%5J`qbSn zcHf(zJEq4T=St%VDTN((xV{^o^({@hzxFtX@KK2iUKk&okCvsHJ%2s%ePj!|jhZ<4 zrbk<5W17^Y>5rV_S$F*2nCwo$U%y?k&2=Z)Z{_6L!R}UJw0fo=)=|j**)PF5}_nI=<&(R`A8$6fFQ`x9`1?+ ziVA^TT%Ax*7}N$u2o?|kSv#R8K9NX!AjlSB;f4-!Ev#HGV_hqIgq5qK4Gii51flF< z)@YKWg%!dH4g^^upqGqjNee40^h72h$l4yQ3Tcl{BY5FXC~K$<5QO-#Nh@>~0E1cl z8{N(oZfk*Xb%a^CqJSVLTPHZw0j96bYK`!Tzex`~3LIZ<;T@MTvK_Ak;%7R5NbH$%&qQB^c zh+NJJ{*@0R@~d8cq04sKU-|fj(N_A2CLq9%PTc;=Cm{H%tRO_>7yAkF3;nWZAu#wC z`w4;3x!zAYghVfs;=k$@5*0;f$bX@M`Gom@wkeoTRP^WPK$}_Q-<}tRu&{?g5$M@& zAn3lm7ZiPlfFLEfjT7*)PJB6D=H!F|{#FD2k#AZcP(Mx^7>xcQ2-(==)fEW;2b^y! A(*OVf literal 0 HcmV?d00001 From f13b6eafacf9cf822ce9da77fe34fef9ade6e911 Mon Sep 17 00:00:00 2001 From: Rob Ashton Date: Fri, 15 Mar 2024 10:57:53 +0000 Subject: [PATCH 2/3] Bump version number --- .github/workflows/R-CMD-check.yaml | 2 +- DESCRIPTION | 2 +- NEWS.md | 6 +++++- tests/testthat/Rplots.pdf | Bin 11830 -> 0 bytes 4 files changed, 7 insertions(+), 3 deletions(-) delete mode 100644 tests/testthat/Rplots.pdf diff --git a/.github/workflows/R-CMD-check.yaml b/.github/workflows/R-CMD-check.yaml index 87c3af0..0c42ded 100644 --- a/.github/workflows/R-CMD-check.yaml +++ b/.github/workflows/R-CMD-check.yaml @@ -46,7 +46,7 @@ jobs: - uses: actions/checkout@master with: - name: mrc-ide/shiny90_sample_files + repository: mrc-ide/shiny90_sample_files ## ${{ github.token }} is scoped to the current repository, so we ## need to provide our own PAT token: ${{ secrets.SAMPLE_FILES_GH_PAT }} diff --git a/DESCRIPTION b/DESCRIPTION index 1fb7a01..18a3906 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,5 +1,5 @@ Package: first90 -Version: 1.6.9 +Version: 1.6.10 Title: The first90 model Description: Implements the Shiny90 model for estimating progress towards the UNAIDS "first 90" target for HIV awareness of status in sub-Saharan Africa. Authors@R: diff --git a/NEWS.md b/NEWS.md index 4352059..39dc33b 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,3 +1,7 @@ +# first90 1.6.10 + +* Qualify package names, fix some R CMD warnings and notes + # first90 1.6.9 * Bug fix: account for end-year net migration in the ART population in the first year of ART start (implemented in v1.6.0). @@ -122,7 +126,7 @@ Updates for 2022 UNAIDS estimates: # first90 1.4.3 -* Patch: in function `add_ss_indices()`, add argument `utils::type.convert(..., as.is = TRUE)` +* Patch: in function `add_ss_indices()`, add argument `type.convert(..., as.is = TRUE)` to suppress R 4.0 warning. * Bugfix: remove duplicate declaration of `double incrate_g[NG];` in `EPP_DIRECTINCID` incidence option. This did not affect any results because this option is not used in Shiny90 application; diff --git a/tests/testthat/Rplots.pdf b/tests/testthat/Rplots.pdf deleted file mode 100644 index 6b077ad6346af3a6b8365a5feb4cc8d3bba99a3d..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 11830 zcmeHNcUV)|wnsWh4NZDLL68y>dQ-X}9qGM=A|XH$dhb#N5drC4L7MccNS9tj5Tq-; zgH-Va=g!O>eecfP|6ji3J73mWYp;D)_AjflX~-$?fO!Q7*#Z{=7Xt?aM~xALU?3mR z$?^fAgajc-76G+D**n4IEKpD&hnyIM4)@T;3Jn<>!T+C2b!xS2&!isBkwv{@P1el))ScWVqu)dNIk&DdL!X2#Uw0kxYh-b3d{E&poQ9;1n zGAW)+lD%Y0QqL7G6l^lXR|HMyA zuoL0g7L(PfZ|PuHO4*9@=5ws%_2oKN{`WZ(1VY$6*m~Ije4&ux$M1u~oW-j8>*+0E z0}1YgM!v~d?9Yf!;wuu2&8u4k=jL z!m^;uwdO2jFoGH6EExobPpQLFC0!#50~|H<6vj`Bi$I6O;X zr^%J0s+b^Owtdssx^PkBa$kP8&2;x7an~@q{()scxTw$ZF)<%?ctxMk6tU0N&P2!2 zVN=_uD6W@DmEv`u-lj%Xx}=4xLpEgBU51Z$D%vg2}}1h zX|o-2VIN=Eo*-*t=d+T^RRfLH{mP_T-V!>MHa0AfX9};z5_Pu^rVcYkXz!jg8r%_C zL1vqNTy<(KXPdn9-23xuG(`V(-2Vku(Lb>AJ6QiC_P~692fiO1KcP=ZP~;z=PiG2F z2PX49((AjzqEpi6T=zp8g~WF>w=gSb(06!Y2TsghULD zz#;PF8Ym`r=Ic0Z#mM-Mxoh^;`GicwR06A$SAuCbN}{h&Im&o$W|*{uTyq(JHlE8g zKEP4(>5OSZbmV2pV#Y&*mv_~UEV;EJ50O4D9cIm;cNV^o+>pS>NjguxN=2}voj*tb<_C2zeli{@qTEvP;T6#XzJ1pFTegF z_(%ZX7{ODIAc&vew7UQK>jt5$bfQ7Ld@inc7?rd#JcbFMb{w0xTyQW1n=44uak?fn z&O&vE$yZ*noOS6XVYb?z#-IT{n(^Z-{ViX&}(zpsv|W3W@jWbKLrFVv~aM`+K>?=`vIm z(-eiab6o9vK2svvn&%!GJ`?mquA#Tk4d}jt9>Lrb>Q3_AgKY|LJD`8CEVGp%(BzI4ba; zu96V+|DNjqJ=Onvs{h+lfiJ_x-yF+7gWqpq;eUePCC4wWO5pEZRlbg@69t&8dAlt6 zDhJo_3)}O0JQ4xSynAv_UP?XR%&v;adH{GvxqC#r?M356nIpika2U1hX_L{#6odrU zg?BArGPUYHC-R8t3}6Y)koK18EzfAB(O1Y_K<%%7qm3>MmC2`ET0}|2`!|Qqc6=7( zsN=#P`0Q?`ArL?sDr~2{S5C%EmqOu#kqj@X?G|0CQfbW<`sNjs5~xckdVot{gNv2q zi=EPHkPd13w9)O1Y~)s(t+3M_47@ffBVIK@ehaVj?&KOyv_7YijIUW!NoxO5r4r() zy4p=wk@VoUcSNg$JEg^hve$;{k7jhk6$kD;afttVJHtJ9dR)w6MWaS#kD$OzosZwV z-_FM=&fF#|Q9fq|*gBmR^Fp{=rjM$7TmnKEN!w-VK=?Q0%X zw&qb8MvUsL(E0M-btRB)wL3)ZV|*4ITlOtIbAPwPvbUa4}My+#5V1t z*LdSuhYP~&)5m9q_qc~!`@NQLf9D`gdFS#ttvc7#XF=^-K1RxSL!ZR45i@=fa%InBngsA28_DmX2~r!Rsz6+*50} zV#g;6x-_}vR$a^&jLw+7)tL5lDTV-D3(rj)HK35Ob@N^3nxmx^>5kXtylP zSK^NuqN*cIJJfH0#P;NO@koM7B7b#Lf`6Rqe})sGzjafW96#ek=pW%k^AP|{<~M5w z2GI%ahN(?CHj>`af?1ibZdSIc#mWUV{z8KVe_d>obz&<+u3H{hm@0oO#o4# z`IScx4P-wCxF;*D(skEC+}s9t~X}!$v$(DxVEcj)S9{Qg&KRB23yvT8pjjwSVKR? zP!62LT`0qO-3Q|u^CT7o7b?(oF76$Y5)8o1S0*91A;Dw0kDJownjm~oCo?;AA>K2U zn}Vej9ma&b}G?{GJ_k$5esIxkPH6jo)V!AZaw8Uw>Wplp^O$SKM7H ziPK92pSzy&rLojI!`W@!`s?D^>l!xfexm6cXeC-oV?v^G$)X&ZESiXieD3(7Y~kD; z#(^(4idBeJ?TAG_tB2~%NL1X9rDd`fy@h}%>wH>_Pt+E*9fy7Af`^=mD)CmiziTeVAbnE<0tAcDhI#Wj5 zTTEF^p*G16BwquV-CH`2xDFNYWL&LX|4e$#pO0b^a%nK5iV z+0otTSby*7mlGS)po0#_rsR56zFeX-%yV!g zBHbTGnt>3H@OI>gP=9sbzawDzb<(UYz1eVS33TMF9@DFXc*e%S?BM>lZHEc+R|OXz zpfWAFf+$WT!gsuWZ}NN42Mof}W6uipisP%qaoKNvMa+L$fB^q(1>y(C&xjG_|A(GM znh5P#c5t2~K3`1P_1Or-iotjRAZOzud*H=k&-(mX2)q2buS8?zloDkQF3H*MxS!Qo zOjq$7CgpUitP26DIoAzbr#M=DuorPBR~4NQ@))z&vhLf#1%u60&=3cv4E=&6zJs_C z6~-+-G2DLKl3?gX2djZJbuh(}hPX5nvAf;42Od#DtsW)CS{>tqad=ZBh>lajg)?mB z16-rgEi!$Bv#h;;KZ7dj`Q2OY0`Z6l%WUnvi#?flC{Oc2~C6 z_bhqk)q)LgdXVx>R@QbaP4=YiZEq0Rkc8p8t}p`n7~EvO(pr-&yS8MXQFRhCwv3K2 z9b!J|w6=5NaU<-xayG$o2*7>T{ZwO?%4ywZcM--W*n7Q@%j>oYUO3hEi26f#5uxlv z_WErvZ!$TAYL0PK^pVSJ`X8B4IHdCN-MN_z^DkL33SXOLgHnQ7qV64%5PTtVRBM)f zYFR3!`_-TTm+OEeD58tymb1oru5pKs1S0{z$@GRC5FNK!tW}WNSy!)`^K+G}&cfzxTqi2<^@4D?0B}_sGx`@hfT?rr# zcq(F}@3LuTJ7ndJ7P2WZH2=GxaAt?dSy)0QOFRm`!^e9%HfT!9Or zlO1L6?qV4F@6N4=#%xvXTOx+%_v&+@+4e}LRKtQ`W^}4!9b9WLK6qMual}AED1V>U zF41A1u`FBJ;4yvVfRN313RA3WB3Ps#HS(u2=kyNuI|yi|%w66lv&TA;eQzk(%!b}= z;l9AQd%FSjE1dqzngZl+YYIO&euk5<@IQf*@@1gb@>`(R-M)ClBr^0+4@Xm3R~{Xw z-HK6sFZ4hdb9bAzYbHjRP*aO!*K6VA%dV${tTPdI?T&uh+v~}Mp}2(_4nujvN*b|4 ziXIJtPCH~zA~p99+sM&@niLE6(mG0F-akEXw&R2i|5%fZVSl?OD+|MnuuwBqH zlG2a#8CB{j`;4DpN@-cc7t@AnwcJ!1WEyX)TNVJVZj{IXkf9n>Rk{x?3$mWbzt&Z> zEF}xxqI|`Gh=iDi(V9FPGT0$(Hqop)baS&*QhuquaAaOZiK`>GZktsiOEdh{xk?LD zsch(i*3K;S0lpo89-7Exo$c9H^n$G_XPl0!7ook?p6LV4zPY(?r*q4vS@T$|RL4&6 zBSN%+Z2lyRxAp2Z?x%Y08m2P7)1Mor9t+z+ZS470$=?Tc;EI$c5~%8lNAZp{O%Tqk zrDTeEz3;t3&3}BnxDbpxxp{tNMR3+CLB8=r)#77}qJuna ziUru0-jlMEJ<=wNmSK{%y#);3Tqs;3F(&Nh)s;!!6gayZYA>bU9Aq74Jf9#`Z&U5WK2 zaff4*zjX;3Z5pDEF*WriOdKf^~r^dG@TonD#<4WFHS z9bFz9ri>5S;ef-?kND51bZ{8RaLBRc$~B8G&b-`ggcCSo6urBxn=fuZ`sUfmHBWY9 z`mL2K(IazALd?Q=d_{m5&POEy!i>^6dI{~b%@3q_&Xc;B#I$4_)n9bDa5g{7vNg;0 zYfEngoZZ6>Yz-zMX$}+r>VUGepMB#xEga6MF8Y!T^PQ_$#*UHgoSv%}38tSmfaJv8 zs0Tj)lz9C{d_v1|l^{Z>rwJ>$^dtewQ?xN4dUOG!#^AxebqC{>0fuf31B3$uPnC1p zD^#=Y+fn4yka4hfGNZjxO(F+(TvWQ*j!036*VB7wO^7}yo;Tc3GI?Dl_f$}uS{}Z+12i51E%0e zIYjvpblVUpDv<`{glXl-BsI(MK3?J8`Mgt;U}xoCgJ9z%YPjySwW(VJz^xwmG@0(m z*w8-SSXuprQ-(41B9M_%Jtx+2TurCt5!Ei%^?kV&Ro1uWlWKFVP4lrG@ds(*#rySz z)kt-FVwaJ)&!a1rL4~;mxFB-@GvjG#Uy}ZLjkT|eGg&(FdBfK7g&ExA-S))nu`UQ| zlM(xt_$C+ImQ@#xB*8c5bhrJ)JaLSPW?v2KG}_I>P!RYWTT^=U;M)6lHIJxe?HFjT z-QRfYU*@@f?dyWQRoZI?T)(9)rNKTg4Lu5drXl`)s6!hd48{~|Hy%V~DpKn}Zx^)OvB)yUO@SZIXxT*e_?iGEKl zAAw84)G%iwm5h{MMr`MKS;?Hbiyr}|9~B_>7^&B7jxkH9tjNkHxp)Z?%_ck%lhM09 z*R++=5%5E$KB6kJ=f<{i?R{A^uX=&}A}^5(<6K1T#CLABgn40VS!lIif2WYc_vYoE z>KOku4wfRoTH}IjGpdF05+m>CwKhof2C@oo`#f6t=#47H(=$1ka8{`Y{BT7*-P6xZ za~$$jczSXTRuhtGpC>Ch#I;WYQK81Pa$x?Uk6A{vH^9OVmE3U+h(+aTL3p$#5~LC! z35Ee~&W|4AMiq}Kz0!|A(U-`Dish{0*tIP(4cquvjW@~LMYw^a#F-mdvhwvEiUs0u zc1eX)xw5Gu@O|}VhU+MtPTsv$`v|pjec)s0;1gL$S{z1Hn@dcaLAnM=3uyImaV7oC zHg?yfWJ_>;3yNb20IrW5cVa#C(}=mO-{o0a5c{ z`e}hETwYw;TsqUXv)z@P4ZRxQTJ4a;&zQ22%13vUNf@V1B^?j;AHnA`eY{6@jh4hp z^B?rUk*0aDLQ>%Do}`dLjksjM_vUx$M7>ITnDa--lob~Oo$ksM-h$VKQMMj<_i}IJ zcNI(4MQ3>KdQdPvS{~XrfqYWg*g57xwrXw!=J{fB%xmn{t;3hty_;Br-Ueyic)L&c zrjpu6>m28;QcLzP$ob!ag9L(fPwer40IVSyn(_e&YcD z#`42gX#pX>v6!Gcy-u#^jxE}?s6wsnEo7WLfJT?mogkmU<;UdTn!f*Ddqek@t?oln zKqC;kyNxbdKs`_7{!v`+>Cfi!Q+MtAlhomZ&BDzZkI5pR z3^~|oVJZ?gE5`HQ0uzwkdlQGjDkP`vETF;TYWcvpya zFtA`4bIz5OGStTkm?@Y(6Vb2~iesAQ1`^=O*{#t#(XANYGoREOu!1o(lB1XRjwg+2 zX?jE?_0rgC)q9d_1XKt3?hUXS1VMBYkbB%lK3`Lyq@w$f5DD5|*kiUc z0@Cjrr06B>80a(P1>7inNXu?~KP)DJB5HncP<1(~qmPe7IOh5Scp`NxY+F=^XH-empLh>yqQLr?WF^I7Pm{FpE}Zjl2Cysf#0wtcU5-o-lU zj~Y%%PCpwAy71MjQ@QA-zTqEhbyp5j2hbeYcy?}xDFDQfm%>_;qF=$Hc!OQ@hP)1g zbPfwsn&Tn1i#d2Mh%)Yq9FT@GpgQZ?lb{mwo6m7wv&aR4g`BUU0T5cWxtg@ksKQE7#vfU?~TuR`_X=$Ei?WXUxcp5U{|$|zz|Nv&sTJ!dEi zTg$dtFS^iTAicK+aMe2gt*;m--mTl*e2jsTg0_ zZA6JHFAvsGd2;m~Ggat9&eLq&Y@ysDO+1bC=eNgn#yB1j9+ImnIOOyX53FFIxW#DG zqj(kcvQmc_S1MP?4@6F>aD%h0qQSA({8bp(sVXrmX)EbRV?{%0JEJWMD)e89jWeuA z5_cfVEn(mqLT@r}fH#hKtZ(r31(b7*h?%j3+C#2(ojF72o6dU96m_&sfElvwXwr`T zo!(vOv4Rh&AO4-#OcoA`XoAuZMVVU`JlT(j6hxHkls%+x01h9jwL9ko)P?CQNytg4 zGqAO=X)4;>5aJBH!Q#lWc{82$6KgXtp+~`*6d_@tG9v?L&t#ds`C2pUmi{eljvY!} zt}zbOEn`iOf`@tn+6UUkFO$pkvpwWo6n9?$OnVd~O4>CI^R;tP+N)tBRS=y?1yOBL zS^~yOi7=4;N(TH6uCBIBYqrv_p#2q4(>@T*nmF4Aw-5W%&zmhP3SStDNkD zt2n75Bh&9y-OIQa^Uj5elQElq7;9zq0CH+^>Ut`!!=t0F)3V0?Ug$knQQ7-1&pn^- zzP~s2@x8*>%nSZ8-?7oL+86FGW=Ey$lkGzsmL6;jixlsoNE~B(8Me|Mi-s>xkB12M z3tSh_v)dYF7)l!-8Q&XyH&QermV2n3TDqB?Ueqjm7VPa092GEulu^eTYSgp@WIG}ds0Cz$6dPa}7x0j++A zkz?gr?WvukQ#$;Uoo`o|C57F6yHTtC;*u1e4>g@%KN{~D&(;?&yAArnZQ)(DQ&#Sy z84s8(U$?e=xCjo;4erkFw-_0I10Jpz-^GpUB(!q0Nq3hK_Y^I06Z>*7KG|xp;WKeU zLsOvqOgXAI*Ryl~c%N=Lrhcs+VZ3gV4=1XXt8IpBxtqAhdIz6jk@}Km+)d~TDSTmD zv@AVnYdUVqT1Ha#0cQL4-d9ei_Nl?fcJHC^1&Yu2H<7b?gTj-w@7>=G&Nzk6a{5>K z$9)(2PJhMl%5=b&fEQ9-ZxRE2rSr0Iq~A-E2K5KYw@R`&E1GwIjceQ|UcS2Q9vU%okoVeMr&n+GVwvyO_FANOQ_PUGFT>dOKTW{*-l5-ZP4qdt5)Q z<@gxT!Is8$7cd+Fj^K&d7R^q>9XJ>X$gfr|tb=XwZc$|5IpH}y-hKyKR=-ito_JI8 zo}21ShD6FjhF1KZ_8V&Ls~Lr-h2!}xYFAaKGQV<4vNua&CB-MZXe*Kz5S&4(`S&5= zm3@_wm3h5sNt?z(Y=i(osO)6atHut2$-%LJcDYxkkY|fe?7oXqPd z_f>?)mP&HVrgvH5oR_GAgezN98ZSJJlH~$}jfPTQB1^lx*pA>JO|AskV3U=cwZ4 z2Km(U+L)-v>kdE5sW|joH;ST)l)x!> z7kJBYV!vLZneWHPRnYRjY}lY~pYHtl`qXqSBZnUEAc%R2=~Xil&+2YZY?64=Mw5nx zo@Q-Sx$-MzLt5iW6R1XY^_s8iy+g$qJyq`#o6Ulc9~wOy%Qki=Ms*8y^bC*SW#-%c zUawCFaC5^VeUp7f&zOW1gunXtpElGh<-Vv|ZS$M=TMc}O<=Y?I-_P@m2U_A>oaee{ z{uyeFoE4k2daj&K#asTpgf^HMKU`XjpQ7fqdhawTO4dEQZU7D;~}e z-T)Fe*`m!5_F0HLOB~f}F3V`7cHSK6K98SY*Zkl#qgcn%+|zinqI+7hT#JMcP8%R; zd~267Mw%nfo*mKq34E(CSp!PmUXv8tF>ehPx0RWI%5J`qbSn zcHf(zJEq4T=St%VDTN((xV{^o^({@hzxFtX@KK2iUKk&okCvsHJ%2s%ePj!|jhZ<4 zrbk<5W17^Y>5rV_S$F*2nCwo$U%y?k&2=Z)Z{_6L!R}UJw0fo=)=|j**)PF5}_nI=<&(R`A8$6fFQ`x9`1?+ ziVA^TT%Ax*7}N$u2o?|kSv#R8K9NX!AjlSB;f4-!Ev#HGV_hqIgq5qK4Gii51flF< z)@YKWg%!dH4g^^upqGqjNee40^h72h$l4yQ3Tcl{BY5FXC~K$<5QO-#Nh@>~0E1cl z8{N(oZfk*Xb%a^CqJSVLTPHZw0j96bYK`!Tzex`~3LIZ<;T@MTvK_Ak;%7R5NbH$%&qQB^c zh+NJJ{*@0R@~d8cq04sKU-|fj(N_A2CLq9%PTc;=Cm{H%tRO_>7yAkF3;nWZAu#wC z`w4;3x!zAYghVfs;=k$@5*0;f$bX@M`Gom@wkeoTRP^WPK$}_Q-<}tRu&{?g5$M@& zAn3lm7ZiPlfFLEfjT7*)PJB6D=H!F|{#FD2k#AZcP(Mx^7>xcQ2-(==)fEW;2b^y! A(*OVf From 607d14bbce5c94c317e32dd4ff7ee7349cbf78bb Mon Sep 17 00:00:00 2001 From: Rob Ashton Date: Mon, 18 Mar 2024 10:17:35 +0000 Subject: [PATCH 3/3] Ignore plots file generated during testing --- .gitignore | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/.gitignore b/.gitignore index 4666775..ce0729a 100644 --- a/.gitignore +++ b/.gitignore @@ -7,7 +7,8 @@ .DS_Store .idea tests/testthat/sample_files +tests/testthat/Rplots.pdf *.tar.gz *.Rcheck -README.html \ No newline at end of file +README.html