---
title: "Legitimation Amid Dependence (1960-2010): A Dataset on the Determinants of Authoritarian Legitimation Strategies - DATA ANALYSIS SCRIPT"
author: "Farah Aly"
date: "2024-05-30"
output: html_document
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```


```{r libraries}
rm(list = ls()) 
p_required <- c("tidyr", "dplyr", "tidyverse", "haven", "readxl", "countrycode", "lmtest", "car", "clarify", "data.table", "stargazer", "ggbreak", "ggtext", "scales", "ggfortify") 
packages <- rownames(installed.packages())
p_to_install <- p_required[!(p_required %in% packages)]
if (length(p_to_install) > 0) {
    install.packages(p_to_install)
}
sapply(p_required, require, character.only = TRUE)
rm(p_required, p_to_install, packages)
```



```{r wd, test & simulation}

data1960 <- read.csv("DIRECTORY/legitimation_dependence_1960_2010.csv")

model_test <- lm(normalized_extent ~ netdependence_diff + gdppc_diff + int_conflict + election_yes_or_no + v2cagenmob + leaderchange + gwf_casename + year, 
                                         data = data1960)

#summary(model_test)

###confidence intervals
#confint(model_test)

###diagnostic tests
resettest(model_test, power = 2:3) ### high omitted variable bias

vif(model_test) ### no multicollinearity

bptest(model_test) ### high heteroscedasticity

shapiro.test(residuals(model_test))

dwtest(model_test) ### no autocorrelation

plot(model_test, which = 1)

autoplot(model_test) #### ggfortify 



### SIMULATION / MARGINAL EFFECT PLOT
### using the clarify package
# King, G., Tomz, M., & Wittenberg, J. (2000). Making the Most of Statistical Analyses: Improving Interpretation and Presentation. American Journal of Political Science, 44(2), 347–361. https://doi.org/10.2307/2669316
coef <- coef(model_test)
varcov <- vcov(model_test)


s1 <- clarify::sim(fit = model_test, n=1000, coefs=coef, vcov=varcov)

sim_1 <- sim_setx (s1, 
                   x= list(netdependence_diff=seq(-4,18,2)),
                   verbose = FALSE,
                   )

p <- plot(sim_1)


p <- p + labs(
  x = "<b>∆ Net Dependence</b>", 
  y = "<b>Predicted Magnitude of Shift<br>from Ideological to Pragmatic Legitimation<br><span style='font-size:8pt;'>(Scaled 0 to 100)</span></b>"
) +
theme(
  axis.title.x = element_markdown(),
  axis.title.y = element_markdown(),
  text = element_text(family = "Times New Roman")
) +
scale_fill_manual(values = "green")

p$layers[[2]]$aes_params$fill <- "green"

p

#ggsave("clarify_plot.png", plot = p, width = 10, height = 6, units = "in")


sd(data1960$normalized_extent, na.rm = T)

plot_cont <- ggplot(data = data1960, aes(x = normalized_extent)) +
  geom_histogram()
plot_cont


g2 <- plot_cont + scale_y_break(c(150, 3350), scales = 0.3, space= 1) + theme_minimal()

g2


negative_cases <- data1960 %>%
  dplyr::filter(extent_shift_ideo_to_prag < 0) %>%
  nrow()

negative_cases

positive_cases <- data1960 %>%
  dplyr::filter(extent_shift_ideo_to_prag > 0) %>%
  nrow()

positive_cases

no_change <- data1960 %>%
  dplyr::filter(extent_shift_ideo_to_prag == 0) %>%
  nrow()

no_change

```





```{r regression table main}
model0 <- lm(normalized_extent ~ netdependence_diff + gdppc_diff  + gwf_casename + year, 
                                         data = data1960)
#summary(model0)


model1 <- lm(normalized_extent ~ netdependence_diff + gdppc_diff + int_conflict + gwf_casename + year, 
                                         data = data1960)

#summary(model1)


model2 <- lm(normalized_extent ~ netdependence_diff + gdppc_diff + int_conflict + election_yes_or_no + gwf_casename + year, 
                                         data = data1960)

#summary(model2)


model3 <- lm(normalized_extent ~ netdependence_diff + gdppc_diff + int_conflict + election_yes_or_no + v2cagenmob + gwf_casename + year, 
                                         data = data1960)

#summary(model3)


model4 <- lm(normalized_extent ~ netdependence_diff + gdppc_diff + int_conflict + election_yes_or_no + v2cagenmob + naturaldeath + gwf_casename + year, 
                                         data = data1960)

#summary(model4)


model5 <- lm(normalized_extent ~ netdependence_diff + gdppc_diff + int_conflict + election_yes_or_no + v2cagenmob + leaderchange + gwf_casename + year, 
                                         data = data1960)
#summary(model5)


stargazer(
  model0, model1, model2, model3, model4, model5, title="Results", align=TRUE, type = 'text', omit = c('gwf_casename', 'year', 'Constant'),
  add.lines = list(
    c('Regime Fixed effects', '✓', '✓', '✓', '✓', '✓', '✓'),
    c('Year Fixed effects', '✓', '✓', '✓', '✓', '✓', '✓')
  ),
  column.labels = c(
    "Model 1", "Model 2", "Model 3", "Model 4", "Model 5", "Model6"
  ), dep.var.labels=c("Shift from Ideological to Pragmatic Legitimation"), order = c("netdependence_diff", "gdppc_diff", "int_conflict", "election_yes_or_no",  "v2cagenmob", "naturaldeath", "leaderchange"),
covariate.labels=c("Δ Net Dependence", 
    "Δ GDP per capita", 
    "Interstate/Internationalized Conflict", "Election-Year",
    "Mass Mobilization",
    "Natural Death of Leader",
    "Leadership Change"), single.row=TRUE,
  out = "MODELS.htm", omit.stat = c("ser", "f"),
  out.header = FALSE
)


#resettest(model5, power = 2)

#durbinWatsonTest(model5)
```

```{r robustness checks}
FBIC_robustness_data <- read.csv("FBIC_robustness_data.csv")

dependence_alt_specifications <- FBIC_robustness_data %>%
  group_by(countryb, year) %>%
  dplyr::summarise(
    economicdependence_sum = sum(economicdependence),
    securitydependence_sum = sum(securitydependence),
    dependence_nonnorm_sum = sum(dependence_nonnorm),
    dependence_xonworld = sum(dependence)
  ) %>%
  ungroup()

dependence_alt_specifications <- dependence_alt_specifications %>% 
  dplyr::mutate(countryb = ifelse(countryb %in% names(country_replacements_FBIC), 
                           country_replacements_FBIC[countryb], 
                           countryb))


merged_for_robustness <- left_join(data1960, dependence_alt_specifications, by = c("gwf_country" = "countryb", "year"))

merged_for_robustness <- merged_for_robustness %>% 
  dplyr::mutate(economicdependence_diff = economicdependence_sum - lag(economicdependence_sum),
         securitydependence_diff = securitydependence_sum - lag(securitydependence_sum),
         dependence_nonnorm_diff = dependence_nonnorm_sum - lag(dependence_nonnorm_sum),
         dependence_xonworld_diff = dependence_xonworld.y - lag(dependence_xonworld.y))

#write.csv(merged_for_robustness, "dependence_alt_specifications.csv", row.names = FALSE)


##### changing the IV //// splitting it into its components
model8 <- lm(normalized_extent ~ dependence_xonworld.y + gdppc_diff + int_conflict + election_yes_or_no + v2cagenmob + leaderchange + gwf_casename + year, 
                                         data = merged_for_robustness)

#summary(model8)

model9 <- lm(normalized_extent ~ economicdependence_diff + gdppc_diff + int_conflict + election_yes_or_no + v2cagenmob + leaderchange + gwf_casename + year, 
                                         data = merged_for_robustness)

model10 <- lm(normalized_extent ~ securitydependence_diff + gdppc_diff + int_conflict + election_yes_or_no + v2cagenmob + leaderchange + gwf_casename + year, 
                                         data = merged_for_robustness)


#### splitting the DV
model12 <- lm(legitideol_diff ~ netdependence_diff + gdppc_diff + int_conflict + election_yes_or_no + v2cagenmob + leaderchange + gwf_casename + year, 
                                         data = data1960)

model13 <- lm(pragmatic_diff ~ netdependence_diff + gdppc_diff + int_conflict + election_yes_or_no + v2cagenmob + leaderchange + gwf_casename + year, 
                                         data = data1960)


stargazer(model8, model9, model10, model12, model13, title="Robustness Checks", align=TRUE, type = 'text', omit = c('gwf_casename', 'year'),
  add.lines = list(
    c('Regime + Year FE', '✓', '✓', '✓', '✓', '✓')),
  column.labels = c(
    "Model 7", "Model 8", "Model 9", "Model 10", "Model 11"
  ), dep.var.labels=c("Shift from Ideological to Pragmatic Legitimation", "Δ Ideological Legitimation", "Δ Pragmatic Legitimation"), order = c("netdependence_diff", "dependence_xonworld.y","economicdependence_diff", "securitydependence_diff", "gdppc_diff", "int_conflict", "election_yes_or_no",  "v2cagenmob", "leaderchange"),
covariate.labels=c("Δ Net Dependence", 
                   "Δ Dependence",
                   "Δ Economic Dependence",
                   "Δ Security Dependence",
    "Δ GDP per capita", 
    "Interstate/Internationalized Conflict", "Elections",
    "Mass Mobilization",
    "Leadership Change"), single.row=TRUE,
  out = "MODELS_robust.htm", omit.stat = c("ser", "f"),
  out.header = FALSE
)

```


```{r robustness checks splitting the sample}

##### splitting the sample to Sub–Saharan Africa
data_SSA <- data1960 %>% 
  dplyr::filter(e_regionpol == 4)

model13 <- lm(normalized_extent ~ netdependence_diff + gdppc_diff + int_conflict + election_yes_or_no + v2cagenmob + leaderchange + gwf_casename + year, 
                                         data = data_SSA)

#summary(model13)

### MENA
data_MENA <- data1960 %>% 
  dplyr::filter(e_regionpol == 3)

model14 <- lm(normalized_extent ~ netdependence_diff + gdppc_diff + int_conflict + election_yes_or_no + v2cagenmob + leaderchange + gwf_casename + year, 
                                         data = data_MENA)

#summary(model14)

### South East Asia
data_SEA <- data1960 %>% 
  dplyr::filter(e_regionpol == 7)

model15 <- lm(normalized_extent ~ netdependence_diff + gdppc_diff + int_conflict + election_yes_or_no + v2cagenmob + leaderchange + gwf_casename + year, 
                                         data = data_SEA)

#summary(model15)

### Eastern Europe
data_EE <- data1960 %>% 
  dplyr::filter(e_regionpol == 1)

model16 <- lm(normalized_extent ~ netdependence_diff + gdppc_diff + int_conflict + election_yes_or_no + v2cagenmob + leaderchange + gwf_casename + year, 
                                         data = data_EE)

#summary(model16)

### LatAm
data_LatAm <- data1960 %>% 
  dplyr::filter(e_regionpol == 2)

model17 <- lm(normalized_extent ~ netdependence_diff + gdppc_diff + int_conflict + election_yes_or_no + v2cagenmob + leaderchange + gwf_casename + year, 
                                         data = data_LatAm)

#summary(model17)


### East Asia
data_EA <- data1960 %>% 
  dplyr::filter(e_regionpol == 6)

model18 <- lm(normalized_extent ~ netdependence_diff + gdppc_diff + int_conflict + election_yes_or_no + v2cagenmob + leaderchange + gwf_casename + year, 
                                         data = data_EA)

#summary(model18)


### Southern Asia
data_SA <- data1960 %>% 
  dplyr::filter(e_regionpol == 8)

model19 <- lm(normalized_extent ~ netdependence_diff + gdppc_diff + int_conflict + election_yes_or_no + v2cagenmob + leaderchange + gwf_casename + year, 
                                         data = data_SA)

#summary(model19)


stargazer(model13, model14, model15, model16, model17, model18, model19, title="Regression Results for Different Regions", align=TRUE, type = 'text', omit = c('gwf_casename', 'year'),
  add.lines = list(
    c('Regime + Year FE', '✓', '✓', '✓', '✓', '✓', '✓', '✓')),
  column.labels = c("Sub-Saharan Africa", "MENA", "South East Asia", "Eastern Europe", "Latin America", "East Asia", "Southern Asia"), dep.var.labels=c("Shift from Ideological to Pragmatic Legitimation"), order = c("netdependence_diff", "gdppc_diff", "int_conflict", "election_yes_or_no",  "v2cagenmob", "naturaldeath", "leaderchange"),
covariate.labels=c("Δ Net Dependence", 
    "Δ GDP per capita", 
    "Interstate/Internationalized Conflict", "Election-Year",
    "Mass Mobilization",
    "Natural Death of Leader",
    "Leadership Change"), single.row=TRUE,
  out = "MODELS_regions.htm", omit.stat = c("ser", "f"),
  out.header = FALSE
)

```

```{r legitimation descriptive graphs: SEPARATE}
###CHINA
china_gwf_fail_years <- data1960 %>% 
  dplyr::filter(gwf_fail == 1, cowcode == 710) %>%
  distinct(year) %>%
  pull()

china_leaderchange_years <- data1960 %>% 
  dplyr::filter(leaderchange == 1, cowcode == 710) %>%
  distinct(year) %>%
  pull()

china_years_list <- list(china_gwf_fail_years = china_gwf_fail_years, china_leaderchange_years = china_leaderchange_years)

china_years_list

china1 <- ggplot(data1960 %>% 
         dplyr::filter(cowcode == 710, year > 1959 & year < 2011) %>% 
         dplyr::select(cowcode, gwf_country, year, v2exl_legitlead, v2exl_legitperf, v2exl_legitratio, v2exl_legitideol, gwf_fail, leaderchange) %>% 
         pivot_longer(-c(cowcode, gwf_country, year, gwf_fail, leaderchange), names_to = "variable", values_to = "value"), 
       aes(x = year, y = value, color = variable)) +
  geom_line() +
  geom_vline(xintercept = c(1976, 1980, 1997, 2003), linetype = "dashed") + # Adding vertical lines
  labs(title = "China",
       x = "Year",
       y = "Value") +
  guides(color = FALSE) +
  scale_color_manual(values = c("v2exl_legitlead" = "purple", 
                                 "v2exl_legitperf" = "blue", 
                                 "v2exl_legitratio" = "red", 
                                 "v2exl_legitideol" = "green"),
                      labels = c("v2exl_legitlead" = "Person of the Leader", 
                                 "v2exl_legitperf" = "Performance Legitimacy", 
                                 "v2exl_legitratio" = "Rational Legal Legitimation", 
                                 "v2exl_legitideol" = "Ideological Legitimacy")) +
  theme_minimal() +  theme(plot.title = element_text(hjust = 0.5))

china1

###RUSSIA
russia_gwf_fail_years <- data1960 %>% 
  dplyr::filter(gwf_fail == 1, cowcode == 365) %>%
  distinct(year) %>%
  pull()

russia_leaderchange_years <- data1960 %>% 
  dplyr::filter(leaderchange == 1, cowcode == 365) %>%
  distinct(year) %>%
  pull()

list(russia_gwf_fail_years = russia_gwf_fail_years, russia_leaderchange_years = russia_leaderchange_years)


russia1 <- ggplot(data1960 %>% 
         dplyr::filter(cowcode == 365, year > 1959 & year < 2011) %>% 
         dplyr::select(cowcode, gwf_country, year, v2exl_legitlead, v2exl_legitperf, v2exl_legitratio, v2exl_legitideol, gwf_fail, leaderchange) %>% 
         pivot_longer(-c(cowcode, gwf_country, year, gwf_fail, leaderchange), names_to = "variable", values_to = "value"), 
       aes(x = year, y = value, color = variable)) +
  geom_line() +
  geom_vline(xintercept = c(1964, 1982, 1984, 1985, 1991, 1999), linetype = "dashed") + 
  labs(title = "Russia",
       x = "Year",
       y = "Value") +
  guides(color = FALSE) +
  scale_color_manual(values = c("v2exl_legitlead" = "purple", 
                                 "v2exl_legitperf" = "blue", 
                                 "v2exl_legitratio" = "red", 
                                 "v2exl_legitideol" = "green"),
                      labels = c("v2exl_legitlead" = "Person of the Leader", 
                                 "v2exl_legitperf" = "Performance Legitimacy", 
                                 "v2exl_legitratio" = "Rational Legal Legitimation", 
                                 "v2exl_legitideol" = "Ideological Legitimacy")) +
  theme_minimal()+ theme(plot.title = element_text(hjust = 0.5))

russia1

###UAE
uae_gwf_fail_years <- data1960 %>% 
  dplyr::filter(gwf_fail == 1, cowcode == 696) %>%
  distinct(year) %>%
  pull()

uae_leaderchange_years <- data1960 %>% 
  dplyr::filter(leaderchange == 1, cowcode == 696) %>%
  distinct(year) %>%
  pull()

list(uae_gwf_fail_years = uae_gwf_fail_years, uae_leaderchange_years = uae_leaderchange_years)


uae1 <- ggplot(data1960 %>% 
         dplyr::filter(cowcode == 696, year > 1959 & year < 2011) %>% 
         dplyr::select(cowcode, gwf_country, year, v2exl_legitlead, v2exl_legitperf, v2exl_legitratio, v2exl_legitideol, gwf_fail, leaderchange) %>% 
         pivot_longer(-c(cowcode, gwf_country, year, gwf_fail, leaderchange), names_to = "variable", values_to = "value"), 
       aes(x = year, y = value, color = variable)) +
  geom_line() +
  geom_vline(xintercept = c(2004), linetype = "dashed") + 
  labs(title = "UAE",
       x = "Year",
       y = "Value",
       color = "Variable") + 
  guides(color = FALSE) +
  scale_color_manual(values = c("v2exl_legitlead" = "purple", 
                                 "v2exl_legitperf" = "blue", 
                                 "v2exl_legitratio" = "red", 
                                 "v2exl_legitideol" = "green"),
                      labels = c("v2exl_legitlead" = "Person of the Leader", 
                                 "v2exl_legitperf" = "Performance Legitimacy", 
                                 "v2exl_legitratio" = "Rational Legal Legitimation", 
                                 "v2exl_legitideol" = "Ideological Legitimacy")) +
  theme_minimal() +  theme(plot.title = element_text(hjust = 0.5))


uae1


```

```{r FINAL COMBINED graph}

combined_data <- data1960 %>%
  dplyr::filter(cowcode %in% c(702, 710, 775, 811, 812, 373, 439, 461, 483, 500, 510, 541, 625, 696, 135)) %>%
  dplyr::select(cowcode, gwf_country, year, v2exl_legitlead, v2exl_legitperf, v2exl_legitratio, v2exl_legitideol, gwf_fail, leaderchange) %>%
  pivot_longer(-c(cowcode, gwf_country, year, gwf_fail, leaderchange), names_to = "variable", values_to = "value") %>%
  dplyr::mutate(gwf_country = ifelse(gwf_country == "United Arab Emirates", "UAE", gwf_country))

#combined_plot <- ggplot(combined_data, aes(x = year, y = value, color = variable)) +
 # geom_line(size = 1) +  
  #geom_vline(data = combined_data %>% filter(leaderchange == 1), 
   #          aes(xintercept = year), 
    #         linetype = "solid",  # Change linetype to solid
     #        color = "grey30", 
      #       size = 0.4) +  labs(title = "Legitimation Types Over Time",
#       x = "\nAuthoritarian Year",
 #      y = "Legitimation Score",
  #     color = "Legitimation Type") +
#  scale_x_continuous(breaks = pretty_breaks(n = 3)) + 
#  scale_color_manual(values = c("v2exl_legitperf" = "#1f77b4",  
 #                               "v2exl_legitratio" = "#ff7f0e", 
  #                              "v2exl_legitideol" = "#2ca02c", 
   #                             "v2exl_legitlead" = "#d62728"),
    #                 labels = c("v2exl_legitperf" = "Performance", 
      #                          "v2exl_legitratio" = "Rational/Legal", 
     #                           "v2exl_legitideol" = "Ideology", 
      #                          "v2exl_legitlead" = "Person of the Leader")) +
#  theme_minimal(base_size = 12) +
 # theme(plot.title = element_text(hjust = 0.5, size = 14, face = "bold"),  
  #      axis.title = element_text(size = 11, face = "italic"),  
   #     axis.text = element_text(size = 10),  
    #    axis.text.x = element_text(size = 8),  
     #   legend.position = "bottom",
      #  legend.box = "horizontal",
       # legend.title = element_text(size = 10, face = "bold"),  
      #  legend.text = element_text(size = 10),  
#        panel.grid.major = element_line(color = "grey90"),  
 #       panel.grid.minor = element_blank(),  
  #      panel.background = element_rect(fill = "grey98", color = NA),  
   #     strip.background = element_rect(fill = "grey90", color = "grey50"),  
    #    strip.text = element_text(size = 10, face = "bold")) +  
#  guides(color = guide_legend(title.position = "top", title.hjust = 0.5)) +  
 # facet_wrap(~ gwf_country, scales = "free_y", ncol = 5)  

  
combined_plot <- ggplot(combined_data, aes(x = year, y = value, color = variable)) +
  geom_line(linewidth = 1) +  
  geom_vline(data = combined_data %>% filter(leaderchange == 1), 
             aes(xintercept = year, linetype = 'Leadership Change'),  
             colour = 'grey30', 
             size = 0.5) +  
  labs(title = "Legitimation Types Over Time",
       x = "  Authoritarian Year",
       y = "Legitimation Score",
       color = "Legitimation Type",
       linetype = "") +
  scale_x_continuous(breaks = pretty_breaks(n = 3)) +
  scale_y_continuous(breaks = pretty_breaks(n = 3)) +
  scale_color_manual(values = c("v2exl_legitperf" = "#1f77b4",  
                                "v2exl_legitratio" = "#ff7f0e", 
                                "v2exl_legitideol" = "#2ca02c", 
                                "v2exl_legitlead" = "#d62728"),
                     labels = c("v2exl_legitperf" = "Performance", 
                                "v2exl_legitratio" = "Rational/Legal", 
                                "v2exl_legitideol" = "Ideology", 
                                "v2exl_legitlead" = "Person of the Leader")) +
  scale_linetype_manual(name = '', 
                        values = 'solid',
                        breaks = c("Leadership Change"),
                        labels = c(expression(bold("Leadership Change")))) +  
  theme_minimal(base_family = "Times New Roman", base_size = 12) + 
  theme(plot.title = element_text(hjust = 0.5, size = 13, face = "bold", family = "Times New Roman"),  
        axis.title = element_text(size = 8, face = "italic", margin = margin(t = -0.2, r = 0, b = 0, l = -10), family = "Times New Roman"),  
        legend.position = "bottom",
        legend.box = "vertical",  
        legend.title = element_text(size = 8.5, face = "bold", family = "Times New Roman"),  
        legend.text = element_text(size = 8.5, family = "Times New Roman"),  
        legend.spacing.x = unit(0.2, "cm"),  
        legend.margin = margin(t = 0, r = 0, b = 0, l = 0),  
        panel.grid.major = element_line(color = "grey90"),  
        panel.grid.minor = element_blank(),  
        panel.background = element_rect(fill = "grey98", color = NA),  
        strip.background = element_rect(fill = "grey90", color = "grey50"),  
        strip.text = element_text(size = 10, face = "bold", family = "Times New Roman"),
        plot.background = element_rect(fill = "white")) +  
  guides(color = guide_legend(title.position = "top", title.hjust = 0.55),
         linetype = guide_legend(title.position = "top", title.vjust = 0.5)) +  
  facet_wrap(~ gwf_country, scales = "free_y", ncol = 5)  

combined_plot <- combined_plot + theme(legend.spacing.y = unit(-0.4, "cm"))  

combined_plot
#ggsave("combined_plot.png", plot = combined_plot, width = 10, height = 6, units = "in")



```




```{r appendix}
###### cases

distinct_values <- distinct(data1960, gwf_casename) %>%
  rename("GWF Authoritarian Regimes" = gwf_casename)

#write.csv(distinct_values, file = "distinct_values.csv")

colnames(data1960)

###### descriptives
variables_for_desc <- c("v2exl_legitperf", "v2exl_legitratio", "v2exl_legitideol", "gdppc")

count_NAs <- function(x) {
  sum(is.na(x))
}

desc <- function(x) {
  c(min = round(min(x, na.rm = TRUE), 2),
    max = round(max(x, na.rm = TRUE), 2),
    median = round(median(x, na.rm = TRUE), 2),
    mean = round(mean(x, na.rm = TRUE), 2),
    NAs = count_NAs(x))
}

data1960_desc <- as.data.frame(data1960)


desc_table <- sapply(data1960_desc[variables_for_desc], desc)
desc_table <- as.data.frame(t(desc_table))

#write.csv(desc_table, file = "descriptive_statistics.csv")
desc_table

```

```{r conclusion graph}

data_conc <- data.frame(
  x = c(-0.04, 0.25, -1, -1.2, 1, 1.4, 1.6),
  y = c(0.9, 1.35, 2, 1, 1.9, -0.9, -1.6),
  label = c("Economic Growth", "Election-Year", "↑ Dependence", "↑ Crime?", 
            "↑ Regional Instability?", "↑ Inequality?", "↓ State Capacity?"),
  color = c("Economic Growth", "Election Year", "Increasing Dependence", "Increasing Crime?", 
            "Increased Regional Instability?", "Increasing Inequality?", "Decreasing State Capacity?")
)

color_blind_palette <- c("Economic Growth" = "#0072B2", "Election Year" = "#56B4E9", "Increasing Dependence" = "#009E73", 
                         "Increasing Crime?" = "#984ea3", "Increased Regional Instability?" = "#E69F00", 
                         "Increasing Inequality?" = "#D55E00", "Decreasing State Capacity?" = "#CC79A7")

p <- ggplot(data_conc, aes(x = x, y = y, label = label, color = color)) +
  geom_point(size = 3) +
  geom_text(hjust = 0.5, vjust = -1, size = 3, nudge_y = 0.02, fontface = "bold", family = "Times New Roman") +  
  scale_color_manual(values = color_blind_palette) +
  geom_segment(aes(x = -1.3, xend = 2, y = 0, yend = 0), arrow = arrow(length = unit(0.3, "cm")), color = "black") +
  geom_segment(aes(x = 0, xend = 0, y = -2.2, yend = 3), arrow = arrow(length = unit(0.3, "cm")), color = "black") +
  labs(x = "Δ Ideological\nLegitimation", 
       y = "Δ Pragmatic\nLegitimation") +
  theme_minimal() +
  theme(
    axis.title.x = element_text(size = 8, margin = margin(t = 10), face = "bold", family = "Times New Roman"),  
    axis.title.y = element_text(size = 8, margin = margin(r = 10), face = "bold", family = "Times New Roman"),  
    axis.text = element_text(size = 5, face = "bold", family = "Times New Roman"),  
    legend.position = "none",
    plot.margin = margin(1,1,1,1, "cm"),
    panel.grid = element_blank(),
    axis.text.x = element_blank(),  
    axis.text.y = element_blank()   
  )

p + theme( 
  axis.title.x = element_text(hjust = 0.9, vjust = 81, family = "Times New Roman"),  
  axis.title.y = element_text(hjust = 0.92, vjust = -107, family = "Times New Roman")
)


```