updates from utils/update.r scripts, adding CI paramteres to ei_iter,…

DESCRIPTION

@@ -43,11 +43,11 @@ License: GPL-3
 Depends: R (>= 3.5.0), eiPack, ei, wru (>= 1.0.0)
 Imports: bayestestR, coda, data.table, doSNOW,
         dplyr, foreach, ggplot2, graphics, magrittr, mcmcse, methods,
-        overlapping, purrr, rlang, sf, stringr, tidyr,tidyselect
+        overlapping, purrr, rlang, sf, stringr, tidyr, tidyselect, viridis
 NeedsCompilation: no
 Suggests: knitr, plyr, rmarkdown, reshape2, RColorBrewer,
         RJSONIO, testthat, tigris
-RoxygenNote: 7.2.1
+RoxygenNote: 7.3.2
 Encoding: UTF-8
 VignetteBuilder: knitr
 Packaged: 2020-09-08 07:00:35 UTC; lorencollingwood

R/ei_iter.R

@@ -38,6 +38,8 @@
 #' @param n_cores The number of cores to use in parallel computation. Defaulted to NULL, in which case parallel::detectCores() - 1 is used
 #' @param verbose A boolean indicating whether to print out status messages.
 #' @param plot_path A string to specify plot save location. If NULL, plot is not saved
+#' @param CI  Numeric. Confidence interval level (default = 0.95). Specifies the
+#'   interval width for calculation with bayestestR package.
 #' @param ... Additional arguments passed directly to ei::ei()
 #'
 #' @return If eiCompare_class = TRUE, an object of class eiCompare is returned.
@@ -76,6 +78,7 @@ ei_iter <- function(
                     n_cores = NULL,
                     verbose = FALSE,
                     plot_path = NULL,
+                    CI = .95,
                     ...) {
 
   # Preparation for parallel processing if user specifies parallelization
@@ -156,7 +159,7 @@ ei_iter <- function(
   # Loop through each 2x2 ei
   ei_results <- foreach::foreach(
     i = seq_len(n_iters),
-    .inorder = FALSE,
+    .inorder = TRUE,
     .packages = c("ei", "stats", "utils", "mvtnorm"),
     .options.snow = opts
   ) %myinfix% {
@@ -200,11 +203,6 @@ ei_iter <- function(
               )
           })
           break
-          # This was meant to enable parameterization of the ei importance sample
-          # size, but its inclusion changes results dramatically.
-          # utils::capture.output({
-          #  ei_out <- suppressMessages(ei_sim(ei_out, samples))
-          # })
         },
         error = function(cond) {
           if (ii == n_erhos) {
@@ -393,7 +391,7 @@ ei_iter <- function(
       # Both CIs
       suppressMessages({
         suppressWarnings({
-          cis <- bayestestR::ci(aggs, ci = 0.95, method = "HDI")
+          cis <- bayestestR::ci(aggs, ci = CI, method = "HDI")
         })
       })
       ci_lowers <- append(ci_lowers, cis$CI_low)
@@ -436,9 +434,9 @@ ei_iter <- function(
     estimates <- data.frame(cbind(means, sds, ci_lowers, ci_uppers))
     estimates <- cbind(cands, races, estimates)
     colnames(estimates) <- c(
-      "cand", "race", "mean", "sd", "ci_95_lower", "ci_95_upper"
+      "cand", "race", "mean", "sd",
+      "ci_lower", "ci_upper"
     )
-
     output <- list(
       "type" = "iter",
       "estimates" = estimates,

R/ei_rxc.R

@@ -162,7 +162,7 @@ ei_rxc <- function(
 
     md_mcmc <- foreach::foreach(
       chain = seq_len(n_chains),
-      .inorder = FALSE,
+      .inorder = TRUE,
       .packages = c("ei"),
       .options.snow = opts
     ) %myinfix% {
@@ -261,51 +261,32 @@ ei_rxc <- function(
       chains_pr[, race_indices] <- race_pr
     }
 
-    # Get point estimates and standard errors
     estimate <- mcmcse::mcse.mat(chains_pr)
-
-    # Get standard deviation of each distribution
     sds <- apply(chains_pr, 2, stats::sd)
 
-    # The upper and lower CI estimates also have standard errors. Here these
-    # errors are conservatively used to extend the 95% confidence bound further
-
-    # Set bounds according to
     if (eiCompare_class) {
-      # eiCompare class object reports fixed CIs
       ci_lower <- 0.025
       ci_upper <- 0.975
     } else {
-      # Get upper, lower CI limits
       ci_lower <- (1 - ci_size) / 2
       ci_upper <- 1 - ci_lower
       if (verbose) {
         message(paste("Setting CI lower bound equal to", ci_lower))
         message(paste("Setting CI upper bound equal to", ci_upper))
       }
     }
-
-    # Lower CI estimate
     lower <- mcmcse::mcse.q.mat(chains_pr, q = ci_lower)
     lower_est <- lower[, 1]
     lower_se <- lower[, 2]
     lower <- lower_est - lower_se
-
-    # Upper CI estimate
     upper <- mcmcse::mcse.q.mat(chains_pr, q = ci_upper)
     upper_est <- upper[, 1]
     upper_se <- upper[, 2]
     upper <- upper_est + upper_se
-
-    # Get race and cand cols for the final table
     cand_col <- rep(cand_cols, each = length(race_cols))
     race_col <- rep(race_cols, times = length(cand_cols))
-
-    # Put names on chains_pr
     names <- paste(cand_col, race_col, sep = "_")
     colnames(chains_pr) <- names
-
-    # Create, name an output table
     results_table <- data.frame(cbind(estimate[, 1], sds, lower, upper))
     results_table <- cbind(cand_col, race_col, results_table)
     if (!eiCompare_class) {
@@ -325,10 +306,7 @@ ei_rxc <- function(
     }
 
     if (!eiCompare_class) {
-      # Match expected output
       results_table <- get_md_bayes_gen_output(results_table)
-
-      # Return results and chains if requested
       if (ret_mcmc) {
         return(list(table = results_table, chains = chains_pr))
       } else {

R/rpv_coef_plot.R

@@ -0,0 +1,62 @@
+rpv_coef_plot <- function (rpvDF = NULL, 
+                            title = "Racially Polarized Voting Analysis Estimates", 
+                            caption = "Data: eiCompare RPV estimates", 
+                            ylab = "Pct. Voting Black-Preferred Candidate", 
+                            colors = NULL, 
+                            race_order = c("Black", "White")
+) 
+
+{
+  colnames(rpvDF) <- stringr::str_to_lower(colnames(rpvDF))
+
+  rpvDF$voter_race <- factor(rpvDF$voter_race, levels = race_order)
+  len_race <- length(unique(rpvDF$voter_race))
+  if (is.null(colors)) {
+    if (len_race == 2) {
+      race_colors <- c(viridis::viridis(10)[4], viridis::viridis(10)[7])
+      names(race_colors) <- race_order
+      ggplot_color_obj <- scale_color_manual(values = race_colors)
+    }
+    else {
+      ggplot_color_obj <- viridis::scale_color_viridis(drop = FALSE, 
+                                                       discrete = TRUE, option = "turbo", alpha = 0.8)
+    }
+  }
+  if (is.null(ylab)) {
+    prty <- unique(rpvDF$party) %>% stringr::str_to_title()
+    ylab <- paste("Percent Voting for", prty, "Candidate")
+  }
+
+  coef_plot <- ggplot(rpvDF, aes(x = 0, y = 0:100)) + 
+    scale_y_continuous(breaks = seq(0, 100, by = 10), 
+                       limits = c(0, 100), 
+                       labels = sprintf("%0.1f%%", seq(0, 100, by = 10)), 
+                       expand = c(0, 0)) + geom_hline(yintercept = 50, 
+                                                      colour = "#000000", size = 0.75) + 
+    geom_pointrange(aes(y = .data$estimate, 
+                        ymin = .data$lower_bound, 
+                        ymax = .data$upper_bound, 
+                        color = .data$voter_race),
+                    position = position_jitter(width = 0.1), 
+                    linewidth = .5, fatten = 1.5, 
+                    show.legend = F) + 
+    ggplot_color_obj + 
+    facet_grid(~panellab) + 
+    labs(y = ylab, title = title, caption = caption) + theme_minimal() + 
+    theme(legend.title = element_blank(), axis.title.x = element_blank(), 
+          axis.ticks.x = element_blank(), axis.text.x = element_blank(), 
+          panel.border = element_rect(fill = NA, colour = "grey"), 
+          panel.grid.major.x = element_blank(), 
+          panel.grid.minor.x = element_blank(), 
+          panel.grid.minor.y = element_blank(), 
+          axis.text.y = element_text(size = 20, face = "bold", 
+                                     family = "serif"), 
+          axis.title.y = element_text(size = 24, 
+                                      face = "bold", family = "serif"), 
+          strip.text.x = element_text(size = 15, 
+                                      face = "bold", family = "serif"), 
+          title = element_text(size = 25, hjust = 0.5, face = "bold", family = "serif"), 
+          plot.caption = element_text(size = 12, face = "italic", 
+                                      family = "serif"))
+  return(coef_plot)
+}

R/rpv_toDF.R

@@ -0,0 +1,42 @@
+rpv_toDF <- function (rpv_results = NULL, model = NULL, jurisdiction = "", 
+                       preferred_candidate = "", party = "", election_type = "", 
+                       year = "", contest = "", candidate = "") 
+{
+  if (inherits(rpv_results, "eiCompare")) {
+    fun <- function(x) {
+      df <- summary(rpv_results)[[x]]
+      df <- df %>% dplyr::select(-"sd")
+      colnames(df) <- paste(x, colnames(df), sep = ".")
+      df <- data.frame(original_name = row.names(df), df)
+      return(df)
+    }
+    sink(tempfile())
+    smry_dfs <- lapply(names(summary(rpv_results)), fun)
+    rpv_data <- suppressMessages(Reduce(dplyr::inner_join, 
+                                        smry_dfs))
+    sink()
+  }
+  else if (inherits(rpv_results, "data.frame")) {
+    rpv_data <- data.frame(original_name = row.names(rpv_results), 
+                           rpv_results)
+  }
+  else {
+    stop("incorrect class type for argument rpv_results")
+  }
+  rownames(rpv_data) <- 1:nrow(rpv_data)
+  colnames(rpv_data) <- colnames(rpv_data) %>% stringr::str_to_lower()
+  newcols <- gsub("mean", "Estimate", colnames(rpv_data))
+  newcols <- gsub("ci_95_lower", "Lower_Bound", newcols)
+  newcols <- gsub("ci_95_upper", "Upper_Bound", newcols)
+  colnames(rpv_data) <- newcols
+  plotDF <- rpv_data %>% dplyr::mutate(Model = model, Jurisdiction = jurisdiction, 
+                                       Election_Type = election_type, Year = as.numeric(year), 
+                                       Contest = contest, Candidate = candidate, Party = party, 
+                                       Preferred_Candidate = preferred_candidate) %>% tidyr::pivot_longer(cols = grep("\\.", 
+                                                                                                                      colnames(rpv_data), value = TRUE), names_to = c("Voter_Race", 
+                                                                                                                                                                      ".value"), names_pattern = "(.*)\\.(.*)", names_repair = "unique")
+  plotDF$Voter_Race <- gsub("^pct_", "", plotDF$Voter_Race)
+  colnames(plotDF) <- gsub("_ei", "", colnames(plotDF))
+  colnames(plotDF) <- gsub("_rxc", "", colnames(plotDF))
+  return(plotDF)
+}

R/summary.eiCompare.R

@@ -5,7 +5,7 @@
 #' @return A nicely formatted dataframe for printing results
 #' @export
 summary.eiCompare <- function(object, ...) {
-  objects <- list(object, ...)
+  objects <- list(object)
   tables <- vector("list", length(objects))
 
   for (ii in seq_along(objects)) {
@@ -25,5 +25,5 @@ summary.eiCompare <- function(object, ...) {
     race <- race_groups[[ii]]
     outputs[[race]] <- data.frame(lapply(tables, function(x) x[[race]]))
   }
-  message(outputs)
+  print(outputs)
 }