Skip to content

Instantly share code, notes, and snippets.

@chrishanretty
Created September 10, 2022 12:01
Show Gist options
  • Save chrishanretty/ffb900c6ebc09ea7d9592062f6fcc1bc to your computer and use it in GitHub Desktop.
Save chrishanretty/ffb900c6ebc09ea7d9592062f6fcc1bc to your computer and use it in GitHub Desktop.
Update results of Kim (2021) to include binary variable for hereditary head of state
### Replicate
###
### Kim, Nam Kyu. "Previous Military Rule and Democratic Survival." Journal of Conflict Resolution 65, no. 2-3 (2021): 534-562
###
### including an additional covariate for a hereditary head of state
###
library(rio)
library(tidyverse)
library(vdemdata)
library(countrycode)
library(brms)
library(modelsummary)
set.seed(1247)
### Cut-down the vdem data
data("vdem")
vdem <- vdem %>%
dplyr::select(country_name,
iso3c = country_text_id,
year,
v2expathhs)
### Read in the Kim data
dat <- rio::import("main_data.dta")
dat <- dat %>%
mutate(iso3c = countrycode(cowcode, "cown", "iso3c"))
### Merge
dat <- merge(dat, vdem,
by.x = c("iso3c", "year"),
by.y = c("iso3c", "year"),
all.x = TRUE,
all.y = FALSE)
table(dat$v2expathhs, useNA = "always")
dat <- dat %>%
mutate(is_monarchy = as.numeric(v2expathhs == 4))
### Replicate the original analysis
### specifically model one
f1 <- gwf_fail ~ leadermil_prior +
reg_Lv2x_polyarchy + log_demyears + Le_migdppcln + Le_migdpgro +
postcold + poly(gwf_duration, 3)
f2 <- update(f1, . ~ . + is_monarchy)
mod <- glm(f1, data = dat %>% filter(gwf_dem == 1),
family = binomial(link = "probit"))
confint(mod)
### Restrict it just countries with non-missing values on monarchy
mod2 <- update(mod,
data = dat %>%
filter(gwf_dem == 1) %>%
filter(!is.na(is_monarchy)))
### update it now to include is_monarchy
mod3 <- update(mod2,
formula = f2)
### Does it look any different in a Bayesian model which might deal
### with perfect separation?
mod4 <- brm(f2,
data = dat %>%
filter(gwf_dem == 1) %>%
filter(!is.na(is_monarchy)),
family = bernoulli(link = "probit"),
prior = set_prior("normal(0, 1)", coef = "leadermil_prior") +
set_prior("normal(0, 1)", coef = "Le_migdpgro") +
set_prior("normal(0, 1)", coef = "Le_migdppcln") +
set_prior("normal(0, 1)", coef = "log_demyears") +
set_prior("normal(0, 1)", coef = "postcold") +
set_prior("normal(0, 1)", coef = "reg_Lv2x_polyarchy"),
chains = 4,
cores = 4)
cm <- c('Intercept' = 'Intercept',
'leadermil_prior' = 'Prior military rule',
'Le_migdpgro' = 'GDP growth',
'Le_migdppcln' = 'log(GDP per capita)',
'log_demyears' = 'log(Previous democracy years)',
'postcold' = 'Post-Cold War',
'reg_Lv2x_polyarchy' = 'Regional democracy',
'b_leadermil_prior' = 'Prior military rule',
'b_Le_migdpgro' = 'GDP growth',
'b_Le_migdppcln' = 'log(GDP per capita)',
'b_log_demyears' = 'log(Previous democracy years)',
'b_postcold' = 'Post-Cold War',
'b_reg_Lv2x_polyarchy' = 'Regional democracy',
'is_monarchy' = 'Hereditary HoS',
'b_is_monarchy' = 'Hereditary HoS'
)
modelsummary(list("Original Kim model" = mod,
"Model w/ NA removed" = mod2,
"Model w/ monarchy" = mod3,
"Bayesian model" = mod4),
coef_map = cm,
statistic = "conf.int")
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment