Use fixed extreme case CIs in prevalence calculations

NAMESPACE

@@ -20,10 +20,8 @@ export(anthro_zscores)
 importFrom(stats,aggregate)
 importFrom(stats,as.formula)
 importFrom(stats,confint)
-importFrom(stats,glm)
 importFrom(stats,plogis)
 importFrom(stats,qt)
-importFrom(stats,quasibinomial)
 importFrom(stats,sd)
 importFrom(stats,setNames)
 importFrom(survey,degf)

R/prevalence-simple.R

@@ -101,7 +101,7 @@ compute_and_aggregate <- function(
   data
 }
 
-#' @importFrom stats glm plogis qt quasibinomial
+#' @importFrom stats plogis qt
 logit_rate_estimate <- function(x, N, empty_data_prototype) {
   x <- x[!is.na(x)]
   if (length(x) == 0) {

R/prevalence-survey.R

@@ -65,7 +65,6 @@ compute_prevalence_zscore_summaries_by.survey_design <- function(
   )
 }
 
-
 #' @export
 compute_prevalence_sample_size_by.survey_design <- function(
     data, indicator, subset_col_name) {
@@ -104,7 +103,6 @@ compute_prevalence_sample_size_by.survey_design <- function(
   )
 }
 
-
 #' @export
 compute_prevalence_estimates_for_column_by.survey_design <- function(
     data, indicator_name, subset_col_name, prev_col_name) {
@@ -131,12 +129,23 @@ compute_prevalence_estimates_for_column_by.survey_design <- function(
     na.rm.all = TRUE,
     level = 1 - prevalence_significance_level
   )[, 3L:4L]
-  data.frame(
+  res <- data.frame(
     Group = as.character(mean_est_prev[[subset_col_name]]),
     r = mean_est_prev[[prev_col_name]] * 100,
     se = survey::SE(mean_est_prev) * 100,
     ll = mean_est_ci_prev$ci_l * 100,
     ul = mean_est_ci_prev$ci_u * 100,
     stringsAsFactors = FALSE
   )
+  # For the extreme cases of `r = 0` and `r = 1` we set the CIs
+  # to [0,0] and [1,1] respectively. Mostly for the convenience
+  # of the human user who consumes the prevalence estimates and to be
+  # in line with the method of the `simple` computation.
+  boundary_0 <- res$r == 0
+  boundary_1 <- res$r == 1
+  res$ll[boundary_0] <- 0
+  res$ul[boundary_0] <- 0
+  res$ll[boundary_1] <- 1
+  res$ul[boundary_1] <- 1
+  res
 }

anthro.Rproj

@@ -1,4 +1,5 @@
 Version: 1.0
+ProjectId: 16901b22-6b7e-4e48-b8f9-c12418c40896
 
 RestoreWorkspace: No
 SaveWorkspace: No

tests/testthat/test-prevalence.R

@@ -362,6 +362,22 @@ test_that("Cluster/strata/sw information is passed correctly to survey", {
   expect_equal(
     observed$HAZ_unwpop, as.numeric(expected_total_unweighted[2])
   )
+
+  # we also enforce that the cis
+  HA_2_WH_2_r_se_0 <- res$HA_2_WH_2_r == 0 & res$HA_2_WH_2_se == 0
+  expect_true(
+    all(res$HA_2_WH_2_ll[HA_2_WH_2_r_se_0] == 0, na.rm = TRUE)
+  )
+  expect_true(
+    all(res$HA_2_WH_2_ul[HA_2_WH_2_r_se_0] == 0, na.rm = TRUE)
+  )
+  WH_3_r_se_0 <- res$WH_3_r == 0 & res$WH_3_se == 0
+  expect_true(
+    all(res$WH_3_ll[WH_3_r_se_0] == 0, na.rm = TRUE)
+  )
+  expect_true(
+    all(res$WH_3_ul[WH_3_r_se_0] == 0, na.rm = TRUE)
+  )
 })
 
 test_that("pop/unwpop are 0 if no values in that group", {