##############################################################################
# Replication Code for Figure 4.5 
# 
# "Exit from International Organizations: Costly Negotiation for Institutional Change"
# 
#  Inken von Borzyskowski and Felicity Vabulas
# 
#  Cambridge University Press 2025
#  
#  Last Updated: March 2025
############################################################################## 

# Code adapted from Avery Schmidt's 2025
#   "Damaged Relations: How Treaty Withdrawal Impacts International Cooperation," 
#   American Journal of Political Science 69(1): 223–239.



setwd("~/Library/CloudStorage/Dropbox/IGO suspensions/Book project/Replication")

# Set Up 
# tidyverse (Version 1.3.1)
if((!"tidyverse" %in% installed.packages()) | # If package is not installed 
   packageVersion("tidyverse") != "1.3.1"){ # Or the installed version is not the correct version
  # Install and load correct version of package
  package_url <- "https://cran.r-project.org/src/contrib/Archive/tidyverse/tidyverse_1.3.1.tar.gz"
  install.packages(package_url, repos=NULL, type="source")
  library("tidyverse")
} else { # Otherwise, load package
  library("tidyverse")
}

# PanelMatch (Version 1.0.0)
if((!"PanelMatch" %in% installed.packages()) | # If package is not installed 
   packageVersion("PanelMatch") != "1.0.0"){ # Or the installed version is not the correct version
  # Install and load correct version of package
  package_url <- "https://cran.r-project.org/src/contrib/Archive/PanelMatch/PanelMatch_1.0.0.tar.gz" 
  install.packages(package_url, repos=NULL, type="source")
  library("PanelMatch")
} else { # Otherwise, load package
  library("PanelMatch")
}


# gridExtra (Version 2.2.1)
if((!"gridExtra" %in% installed.packages()) | # If package is not installed 
   packageVersion("gridExtra") != "2.2.1"){ # Or the installed version is not the correct version
  # Install and load correct version of package
  package_url <- "https://cran.r-project.org/src/contrib/Archive/gridExtra/gridExtra_2.2.1.tar.gz"
  install.packages(package_url, repos=NULL, type="source")
  
  library("gridExtra")
} else { # Otherwise, load package
  library("gridExtra")
}

# grid (Version 4.2.0)
# Note: installation not necessary, because gird is a base R package; therefore,
#       its version corresponds with the version of R installed locally
library("grid")
packageVersion("grid") == "4.2.0" # R Version 4.2.0 was used for original analysis 

# Set seed
set.seed(123) 

# Load results 
file_name <- "Data/Chapter4_Figure5.RData"
load(file = file_name)


# Make functions to simplifying presentation
alpha <- 0.05 # Set confidence levels
ci_upper_bound_fxn <- function(x){quantile(x, probs = c(1-alpha/2))} # Functions to extract CI bounds from bootstrapped estimates
ci_lower_bound_fxn <- function(x){quantile(x, probs = c(alpha/2))}
pm_extractor_fxn <- function(analysis, sample, ests){ # Functions to extract relevant data from PanelEstimate objects
  est_names <- names(ests$estimates) # Extract estimate names
  est_betas <- unname(ests$estimates) # Extract point estimates
  boots <- ests$bootstrapped.estimates # Extract matrix of bootstrapped estimates
  # Place relevant quantities in a vector
  dat_temp <- c(rep(analysis, length(est_names)), # Name of analysis
                rep(sample, length(est_names)), # Name of sample
                est_betas, # Point Estimates
                apply(boots, 2, ci_lower_bound_fxn), # 95% CI lower bound
                apply(boots, 2, ci_upper_bound_fxn), # 95% CI lower bound
                est_names) # Time period of estimate
  table_out <- matrix(data = dat_temp, nrow = length(est_betas), ncol = 6, byrow = F) # Convert vector to matrix
  colnames(table_out) <- c("pooling", "name", "estimate", "lower_95", "upper_95", "t") # Assign column names
  return(table_out) # Return output 
}

# Main Results 
graph_values <- as.data.frame(matrix(NA, ncol = 9, nrow = 0)) # Make empty dataframe to fill with relevant data
colnames(graph_values) <- c("pooling", "name", "estimate", "lower_95", "upper_95", "t") # Assign column names
graph_values <- rbind(graph_values, pm_extractor_fxn("Logged", "Estimates", ests_inken_felicity_log_init)) # Extract true estimates
graph_values <- rbind(graph_values, pm_extractor_fxn("Logged", "Placebo", ests_placebo_t3_inken_felicity_log_init )) # Extract placebo estimates
for(i in 3:5){graph_values[,i] <- as.numeric(graph_values[,i])} # Coerce estimates and CIs to numeric
graph_values$t <- as.factor(graph_values$t) # Convert t to factor 
graph_values$t_num <- as.numeric(graph_values$t) # Make numeric value of t
graph_values$t_num[graph_values$name == "Placebo"] <- graph_values$t_num[graph_values$name == "Placebo"] - 3 # Correct numeric t for placebos
graph_values$t <- graph_values$t_num - 1 # Shift t_num down by one so that the true estimate at t = 0 is also t_num = 0; overwrite t

# Figure 5 ----
# removed grid with panel.grid = element_blank(),
graph_values %>% 
  ggplot(aes(x = t, y = estimate)) +  
  geom_crossbar(aes(ymin = lower_95, ymax = upper_95), 
                width=0.01, fill = "black", color = "black") +
  geom_point(aes(x = t, y = estimate), size = 3) + 
  theme_bw() +
  geom_hline(yintercept=0, linetype="dashed", color = "black") +
  geom_vline(xintercept = seq(-0.5, 10, by = 0.01), colour="grey10", linetype = "solid", alpha = 0.05) +
  coord_cartesian(xlim = c(-3.25, 6.25)) + 
  scale_x_continuous(breaks=c(-5, -4, -3, -2, -1, 0,  1, 2, 3, 4, 5, 6),
                     labels=c("-5","-4", "-3", "-2", "-1", "0", "1", "2", "3", "4", "5", "6")) + 
  scale_y_continuous(breaks = c(-0.2, -0.15, -0.1, - 0.05, 0, 0.05, 0.1, 0.15, 0.2), 
                     labels = c("-20%", "-15%", "-10%", "-5%", "0%", "5%", "10%", "15%", "20%")) + 
  ylab("Percent change in joining treaties") +
  xlab("Years relative to withdrawal") + 
  theme(legend.position= "none", 
        legend.title=element_blank(), 
        legend.text = element_text(size=14),
        strip.text = element_text(size = 14),
        plot.title = element_text(hjust = 0.5, size=16), 
        text = element_text(size = 14),
        panel.grid = element_blank(),
        axis.title = element_text(size = 14), 
        axis.text = element_text(size = 12)) 

ggsave("Figures/Chapter4_Figure5.pdf", plot = last_plot(), width = 7, height = 4, units = "in")

