R/urrplot.R

In dbasu-umass/bate: Computes Bias-Adjusted Treatment Effect

#' Region plot to demarcate URR and NURR for the bounded box
#'
#' @param parameters A vector of parameters (real numbers) that is generated by estimating the short, intermediate and auxiliary regressions.
#' @param deltalow The lower limit for delta. 
#' @param deltahigh The upper limit for delta.
#' @param Rlow The lower limit for Rmax.
#' @param Rhigh The upper limit for Rmax. 
#' @param e The step size of the grid in the x and y directions.
#'
#' @importFrom ggplot2 "geom_vline"
#' @importFrom ggplot2 "scale_fill_manual"
#' @importFrom ggplot2 "geom_tile"
#' @importFrom dplyr "filter"
#' @importFrom dplyr "mutate"
#' @importFrom latex2exp "TeX"
#' @importFrom ggplot2 "aes"
#' @importFrom ggplot2 "ggplot"
#' @importFrom ggplot2 "guide_legend"
#' @importFrom ggplot2 "stat_contour_filled"
#' @importFrom ggplot2 "theme_minimal"
#' @importFrom ggplot2 "labs"
#' @importFrom ggplot2 "guides"
#' 
#' 
#' @return A plot object created by ggplot
#' @export
#'
#'@examples 
#' 
#' ## Load data set
#' data("NLSY_IQ")
#'  
#' ## Set age and race as factor variables
#' NLSY_IQ$age <- factor(NLSY_IQ$age)
#' NLSY_IQ$race <- factor(NLSY_IQ$race)
#'    
#' ## Collect parameters from the short, intermediate and auxiliary regressions
#' parameters <- collect_par(
#' data = NLSY_IQ, outcome = "iq_std", 
#' treatment = "BF_months", 
#' control = c("age","sex","income","motherAge","motherEDU","mom_married","race"),
#' other_regressors = c("sex","age"))
#' 
#' ## Set limits for the bounded box
#' Rlow <- parameters$Rtilde
#' Rhigh <- 0.61
#' deltalow <- 0.01
#' deltahigh <- 0.99
#' e <- 0.01
#' 
#' ## Create region plot for bounded box
#' p1 <- urrplot(parameters, deltalow, deltahigh, Rlow, Rhigh, e=e)
#' 
#' ## See plot
#' print(p1)
#'  
urrplot <- function(parameters, deltalow, deltahigh, Rlow, Rhigh, e){
  
  # Check if delta parameters are 1
  if((deltalow==1) | (deltahigh==1)) stop("Values for delta cannot be exactly equal to 1")
  
  # Create grid
  consol <- expand.grid(delta = seq(deltalow,deltahigh,by=e),
                        Rmax = seq(Rlow,Rhigh,by=e)) 
  
  # Check for whether there are any values for where delta = 1
  if(nrow(filter(consol,delta==1))>0){
    while(nrow(filter(consol,delta==1))>0){
      e <- 0.99*e
      consol <- expand.grid(delta = seq(deltalow,deltahigh,by=e),
                            Rmax = seq(Rlow,Rhigh,by=e))
    }
    warning(paste0("The grid includes delta = 1. e has been modified to ",e))
  }
  
  # Determine whether each coordinate has one (D=1) or three (D=0) real roots
  consol <- consol %>% 
    mutate(D = mydisc(parameters,delta,Rmax))
  
  p <- ggplot(consol, aes(x = delta, y = Rmax)) +
    geom_tile(aes(fill=factor(D))) +
    scale_fill_manual("Region",
                      values=c("0" = "firebrick4",
                               "1" = "aquamarine2"),
                      labels=c("0" = "NURR",
                               "1" = "URR")) +
    # White line at delta=1 (if delta=1 falls within chosen box)
    {if(deltalow <1 & 1 < deltahigh)
      geom_vline(xintercept = 1, color="white")} +
    theme_minimal() +
    labs(
      x=TeX("$\\delta$"),
      y=TeX("$R_{max}$"),
      title="URR and NURR regions",
      subtitle=paste0(sum(consol$D),
                      " (",
                      nrow(consol) - sum(consol$D),
                      ") points in URR (NURR)")
    )
  # Return plot
  return(p)
}