R/estimate.R

In validmeasures/wsup: An Implementation of the WSUP Urban Water and Sanitation Survey in R

################################################################################
#
#' A wrapper function for \code{bootBW} in \code{bbw} package that applies the
#' blocked weighted bootstrap to an indicators data frame with multiple
#' stratifications.
#' 
#' @param steer A data frame containing information on the various indicators
#'     to be analysed
#' @param indicators A data frame containing indicators data
#' @param pop A data frame containing population data for each primary sampling
#'     unit (PSU)
#' @param statistic A character value specifying the name of the statistic to
#'     apply on the data on bootstrap. Two statistics are available for use via
#'     \code{bbw} package: \code{bootClassic} and \code{bootProbit}.
#' @param replicates A numeric value of the number of replicates in performing
#'     blocked weighted bootstrap. Default is 399. Accepted values greater than
#'     or equal to 399
#' @return A data frame of boostrapped results
#' @examples
#'   #
#'   estimate_indicators(steer = steerIndicators[steerIndicators$varList == "jmpWater4", ], 
#'                       indicators = washdata::indicatorsDataBGD, 
#'                       pop = washdata::popBGD, 
#'                       statistic = "bootClassic", 
#'                       replicates = 9)
#'
#' @export
#'
#
################################################################################

estimate_indicators <- function(steer, 
                                indicators, 
                                pop, 
                                statistic,
                                replicates = 399) {
  ## Create empty container object for bootstrap results
  summaryTable <- NULL
  ## Cycle through each enumeration zones
  for(i in unique(indicators$zone)) {
    ## Subset population data to slum areas
    slum.pop <- pop[pop$zone == i & pop$type == 1, ]
    ## Subset population data to other areas
    other.pop <- pop[pop$zone == i & pop$type == 2, ]
    ## Subset survey data to slum areas
    slum.data <- indicators[indicators$zone == i & indicators$type == 1, ]
    ## Subset survey data to other areas
    other.data <- indicators[indicators$zone == i & indicators$type == 2, ]  
    ## Get population data for all areas in current zone
    all.pop <- pop[pop$zone == i, ]
    ## Get survey data for all areas in current zone
    all.data <- indicators[indicators$zone == i, ]
    ## Cycle through indicator codes
    for(j in steer$varList) {
      ## Select indicator and statistic pair
      INDICATOR <- j
      STATISTIC <- statistic
      ## Report progress
      cat("\n\n", INDICATOR, " : Survey Area ", i, "\n\n", sep = "")
      ## Boostrap (BBW) for current indicator and statistic pair - slum
      boot.slum <- bbw::bootBW(x = slum.data, 
                               w = slum.pop, 
                               statistic = eval(parse(text = STATISTIC)), 
                               params = INDICATOR, 
                               outputColumns = INDICATOR, 
                               replicates = replicates)
      ## Boostrap (BBW) for current indicator and statistic pair - other
      boot.other <- bbw::bootBW(x = other.data, 
                                w = other.pop, 
                                statistic = eval(parse(text = STATISTIC)), 
                                params = INDICATOR, 
                                outputColumns = INDICATOR, 
                                replicates = replicates)
      ## Boostrap (BBW) for current indicator and statistic pair - total
      boot.total <- bootBW(x = all.data, 
                           w = all.pop, 
                           statistic = eval(parse(text = STATISTIC)), 
                           params = INDICATOR, 
                           outputColumns = INDICATOR, 
                           replicates = replicates)     
      ## Concatenate various boot vectors into data.frame
      boot.all <- data.frame(boot.slum, boot.other, boot.total)
      names(boot.all) <- c("Slum", "Other", "Total")
      ## Create a row of results and insert quantiles of boot results
      rowResult <- c(unique(indicators$country), unique(indicators$ccode), unique(indicators$month), unique(indicators$year), 
                     as.character(steer$varNames[steer$varList == j]), j, paste("Survey Area ", i, sep = ""))
      for(k in names(boot.all)) {
        rowResult <- c(rowResult, quantile(boot.all[[k]], probs = c(0.5, 0.025, 0.975), na.rm = TRUE))
      }
      ## Record results
      summaryTable <- rbind(summaryTable, rowResult)
    }
    ## Tidy results (row names, column types, column names, and order of rows)
    summaryTable <- data.frame(summaryTable, row.names = 1:nrow(summaryTable), stringsAsFactors = FALSE)
    ## Convert estimate results into numeric
    for(col in 8:ncol(summaryTable)) {
      summaryTable[, col] <- as.numeric(summaryTable[, col])
    }
  }
  ## Rename summaryTable
  names(summaryTable) <- c("country", "ccode", "month", "year", 
                           "indicator.name", "indicator.code", "survey.area",
                           "slum.estimate", "slum.lcl", "slum.ucl",
                           "other.estimate", "other.lcl", "other.ucl",
                           "city.estimate", "city.lcl", "city.ucl")
  ## Return output
  return(summaryTable)
}