library(tidyverse)
library(plm)
library(stargazer)
library(lmtest)
library(sandwich)
library(magrittr)
library(corrr)
library(knitr)
library(kableExtra)
library(PerformanceAnalytics)
library(hhi)
devtools::install_github("davidsjoberg/ggbump")
library(ggbump)
if(!require(ggflags)) devtools::install_github("rensa/ggflags")
library(ggflags)
pd_df %>%
# filter(pop > 999999) %>%
group_by(year) %>%
mutate(total_annual_budget = sum(budget)) %>%
ungroup() %>%
group_by(un_code, year) %>%
mutate(budget_proportion = budget / total_annual_budget) %>%
ungroup() %>%
group_by(year) %>%
mutate(rank_budget = rank(-budget_proportion, ties.method = "random")) %>%
ungroup() %>%
group_by(country_name) %>%
mutate(mean_rank = mean(rank_budget)) %>%
ungroup() %>%
select(country_name, year, budget_proportion, rank_budget, mean_rank) -> hi
# Only keep top 10 mean recipient countries
bump_df <- hi %>%
filter(country_name %in% (hi %>%
distinct(country_name, mean_rank) %>%
top_n(10, -mean_rank) %>%
pull(country_name)))
bump_df$iso2 <- countrycode::countrycode(bump_df$country_name, "country.name", "iso2c")
bump_df$iso2 <- tolower(bump_df$iso2 )
bump_df$country_name <- with(bump_df, reorder(country_name, mean_rank))
bump_df %>%
filter(year != "2013") %>%
ggplot(aes(year, rank_budget,
group = country_name,
color = country_name,
fill = country_name)) +
geom_bump(aes(smooth = 10),
size = 1.5,
lineend = "round") +
geom_flag(data = bump_df %>%
filter(year == max(year)),
aes(country = iso2),
size = 8,
color = "black") +
scale_y_reverse(breaks = 1:100) +
bbplot::bbc_style() +
theme(legend.text = element_text(size = 10),
panel.grid = element_blank(),
panel.background = element_rect(fill = "#C2C2C2", color = "transparent"),
plot.background = element_rect(fill = "#C2C2C2"),
text = element_text(color = "white"))
We see a very big increase in spending to Russia from low of 21st country to among top 5 countries from 2016 onwards
#Looking at top ranking countries with budget as a percentage of GDP per captita
pd_df %>%
group_by(country_name, year) %>%
mutate(budget_per_GDP = budget / GDP_per_cap_constant_2010 ) %>%
ungroup() %>%
group_by(year) %>%
mutate(rank_budget = rank(-budget_per_GDP, ties.method = "min")) %>%
ungroup() %>%
group_by(country_name) %>%
mutate(mean_rank = mean(-budget_per_GDP)) %>%
ungroup() %>%
select(country_name, year, rank_budget, mean_rank) -> hi2
# Only keep top 10 mean recipient countries
bump_df <- hi2 %>% filter(country_name %in% (hi2 %>%
distinct(country_name, mean_rank) %>%
top_n(10, -mean_rank) %>%
pull(country_name)))
bump_df$iso2 <- countrycode::countrycode(bump_df$country_name, "country.name", "iso2c")
bump_df$iso2 <- tolower(bump_df$iso2)
bump_df$country_name <- with(bump_df, reorder(country_name, rank_budget))
bump_df %>%
filter(year != "2013") %>%
ggplot(aes(year, rank_budget,
group = country_name,
color = country_name,
fill = country_name)) +
geom_bump(aes(smooth = 10),
size = 1.5,
lineend = "round") +
geom_flag(data = bump_df %>%
filter(year == max(year)),
aes(country = iso2),
size = 8,
color = "black") +
scale_y_reverse(breaks = 1:100) +
bbplot::bbc_style() +
theme(legend.text = element_text(size = 10),
panel.grid = element_blank(),
panel.background = element_rect(fill = "#C2C2C2", color = "transparent"),
plot.background = element_rect(fill = "#C2C2C2"),
text = element_text(color = "white"))
pd_df %>%
group_by(country_name, year) %>%
mutate(budget_per_GDP = budget / GDP_per_cap_constant_2010 ) %>%
ungroup() %>%
group_by(year) %>%
mutate(rank_budget = rank(-budget_per_GDP, ties.method = "min")) %>%
ungroup() %>%
group_by(country_name) %>%
mutate(mean_rank = mean(budget_per_GDP)) %>%
ungroup() %>%
select(country_name, year, rank_budget, mean_rank) %>%
# Only keep top 10 mean recipient countries
filter(country_name %in% (hi2 %>%
distinct(country_name, mean_rank) %>%
top_n(10, -mean_rank) %>%
pull(country_name)))
## # A tibble: 70 x 4
## country_name year rank_budget mean_rank
## <chr> <int> <int> <dbl>
## 1 Afghanistan 2017 1 79358.
## 2 Afghanistan 2016 1 79358.
## 3 Afghanistan 2013 1 79358.
## 4 Afghanistan 2014 1 79358.
## 5 Afghanistan 2018 1 79358.
## 6 Afghanistan 2015 1 79358.
## 7 Afghanistan 2019 1 79358.
## 8 Congo, Dem. Rep. 2015 6 7824.
## 9 Congo, Dem. Rep. 2018 3 7824.
## 10 Congo, Dem. Rep. 2016 3 7824.
## # ... with 60 more rows
bump_df$iso2 <- countrycode::countrycode(bump_df$country_name, "country.name", "iso2c")
bump_df$iso2 <- tolower(bump_df$iso2)
bump_df$country_name <- with(bump_df, reorder(country_name, rank_budget))
bump_df %>%
filter(year != "2013") %>%
ggplot(aes(year, rank_budget,
group = country_name,
color = country_name,
fill = country_name)) +
geom_bump(aes(smooth = 10),
size = 1.5,
lineend = "round") +
geom_flag(data = bump_df %>%
filter(year == max(year)),
aes(country = iso2),
size = 8,
color = "black") +
scale_y_reverse(breaks = 1:100) +
bbplot::bbc_style() +
theme(legend.text = element_text(size = 10),
# panel.grid = element_blank(),
panel.background = element_rect(fill = "#C2C2C2", color = "transparent"),
plot.background = element_rect(fill = "#C2C2C2"),
text = element_text(color = "white"))
## Mean spending per country by the US per year
We can see a big drop in the average amount of money that went to each country in 2016! What happened in President Obama's last year? Or did budgets fall after Trump administration?
```r
pd_df %>%
filter(year != 2013) %>%
group_by( year) %>%
summarise(mean_budget = mean(budget, na.rm = TRUE)) %>%
ggplot(aes(x = year, y = mean_budget)) +
geom_line(size = 2, alpha = 0.7) +
# geom_smooth(method = "lm") +
bbplot::bbc_style() +
theme(legend.position = "bottom",
legend.text = element_text(size = 10))
Median spending over the years a bit more volatile!
pd_df %>%
filter(year != 2013) %>%
group_by(year) %>%
summarise(median_budget = median(budget, na.rm = TRUE)) %>%
ggplot(aes(x = year, y = median_budget)) +
geom_line(size = 2, alpha = 0.7) +
# geom_smooth(method = "lm") +
bbplot::bbc_style() +
theme(legend.position = "bottom",
legend.text = element_text(size = 10))
## Extent of recipient concentration
Using the Hirschman - Herfindahl Index of Concentration formula
pd_df %>%
filter(pop > 999999) %>%
# group_by(year) %>%
mutate(total_annual_budget = sum(budget)) %>%
ungroup() %>%
group_by(un_code) %>%
mutate(budget_proportion = budget / total_annual_budget) %>%
ungroup() %>%
# filter(year == 2019) %>%
group_by(country_name) %>%
summarise(country_prop = sum(budget_proportion, na.rm = TRUE)) -> hhi_df
# Square the percentage "market share" for each country
hhi_df_square <- sapply(hhi_df$country_prop, function(x) x^2)
hhi_sum_df <- hhi_df %>%
mutate(sum_square = sum(hhi_df_square, na.rm = TRUE))
# hhi_df_country <- hhi_df %>% pull(country_name)
#
# hhi_df <- data.frame(hhi_df_country, hhi_df_square)
#
# # Sum all squares together
#
# hhi(hhi_df, "country_prop")
#
# View(hhi_sum_square)
pd_df %>%
group_by(country.x) %>%
summarise(country_sum = sum(budget)) %>%
arrange(desc(country_sum)) %>%
head(n = 10) %>% kbl() %>% kable_paper("hover", full_width = F)
| country.x | country_sum |
|---|---|
| Afghanistan | 320227200 |
| Pakistan | 271562912 |
| Iraq | 93217522 |
| Russia | 83698073 |
| India | 79224873 |
| China | 73776214 |
| Japan | 73150831 |
| Kenya | 70679051 |
| Indonesia | 66611463 |
| Brazil | 65248400 |
pd_df %>%
group_by(country.x) %>%
summarise(country_sum = sum(budget)) %>%
arrange(desc(country_sum)) %>%
tail(n = 10) %>% kbl() %>% kable_paper("hover", full_width = F)
| country.x | country_sum |
|---|---|
| Vatican City | 919013 |
| Guyana | 910827 |
| Taiwan | 851513 |
| Cape Verde | 796005 |
| Guinea-Bissau | 627305 |
| Marshall Islands | 429697 |
| Palau | 415811 |
| The Bahamas | 247400 |
| The Gambia | 183550 |
| Curacao | 5880 |
pd_df %>%
filter(pop > 999999) %>%
mutate(budget_pd = budget / pop ) %>%
group_by(country.x) %>%
summarise(country_sum = sum(budget_pd, na.rm = TRUE)) %>%
arrange(desc(country_sum)) %>%
head(n = 10) %>% kbl() %>% kable_paper("hover", full_width = F)
| country.x | country_sum |
|---|---|
| Botswana | 13.436588 |
| Cyprus | 10.089509 |
| Afghanistan | 9.201695 |
| Moldova | 7.840552 |
| Georgia | 7.756889 |
| Swaziland | 6.174381 |
| Palestinian Territories | 6.168170 |
| Bahrain | 5.852899 |
| Bosnia and Herzegovina | 5.382142 |
| Estonia | 4.835109 |
pd_df %>%
filter(pop > 999999) %>%
mutate(budget_pd = budget / pop ) %>%
group_by(country_name) %>%
summarise(country_sum = sum(budget_pd, na.rm = TRUE)) %>%
arrange(desc(country_sum)) %>%
tail(n = 10) %>% kbl() %>% kable_paper("hover", full_width = F)
| country_name | country_sum |
|---|---|
| Indonesia | 0.2582964 |
| Nigeria | 0.2425019 |
| Angola | 0.2399337 |
| Madagascar | 0.2259801 |
| Somalia | 0.2139121 |
| Philippines | 0.2083287 |
| Sudan | 0.1459121 |
| Bangladesh | 0.1242950 |
| India | 0.0603926 |
| China | 0.0537257 |
pd_df %>%
filter(pop > 999999) %>%
mutate(budget_per_GDP = budget / GDP_per_cap_constant_2010 ) %>%
group_by(country_name) %>%
summarise(country_sum = sum(budget_per_GDP, na.rm = TRUE)) %>%
arrange(desc(country_sum)) %>%
head(n = 10) %>% kbl() %>% kable_paper("hover", full_width = F)
| country_name | country_sum |
|---|---|
| Afghanistan | 555506.84 |
| Pakistan | 246565.56 |
| Ethiopia | 106876.77 |
| Mozambique | 73074.41 |
| Kenya | 65855.17 |
| Congo, Dem. Rep. | 54764.73 |
| India | 45287.30 |
| Tanzania | 37525.35 |
| Uganda | 28829.60 |
| Zimbabwe | 27249.68 |
pd_df %>%
filter(pop > 999999) %>%
mutate(budget_per_GDP = budget / GDP_per_cap_constant_2010 ) %>%
group_by(country_name) %>%
summarise(country_sum = sum(budget_per_GDP, na.rm = TRUE)) %>%
arrange(desc(country_sum)) %>%
tail(n = 10) %>% kbl() %>% kable_paper("hover", full_width = F)
| country_name | country_sum |
|---|---|
| Netherlands | 196.28710 |
| Sweden | 160.77004 |
| Equatorial Guinea | 142.52215 |
| Qatar | 135.63607 |
| Denmark | 113.49647 |
| Ireland | 98.87856 |
| Norway | 81.39777 |
| Switzerland | 59.20861 |
| Somalia | 0.00000 |
| Syrian Arab Republic | 0.00000 |
Afghanistan, Pakistan and Iraq account for almost 20% of all PD budget spending
pd_df %>%
filter(pop > 999999) %>%
# group_by(year) %>%
mutate(total_annual_budget = sum(budget)) %>%
ungroup() %>%
group_by(un_code) %>%
mutate(budget_proportion = budget / total_annual_budget) %>%
ungroup() %>%
group_by(country_name) %>%
summarise(country_sum = sum(budget_proportion, na.rm = TRUE)) %>%
arrange(desc(country_sum)) %>%
head(n = 20) %>% kbl() %>% kable_paper("hover", full_width = F)
| country_name | country_sum |
|---|---|
| Afghanistan | 0.0934030 |
| Pakistan | 0.0792088 |
| Iraq | 0.0271894 |
| Russian Federation | 0.0244128 |
| India | 0.0231081 |
| China | 0.0215189 |
| Kenya | 0.0206155 |
| Indonesia | 0.0194291 |
| Brazil | 0.0190315 |
| South Africa | 0.0172208 |
| Ukraine | 0.0169231 |
| Germany | 0.0163838 |
| Japan | 0.0163828 |
| Ethiopia | 0.0143528 |
| Nigeria | 0.0130050 |
| Mexico | 0.0123035 |
| Mozambique | 0.0118875 |
| Israel | 0.0115187 |
| Korea, Rep. | 0.0113068 |
| Egypt, Arab Rep. | 0.0109773 |
pd_df %>%
filter(pop > 999999) %>%
# group_by(year) %>%
mutate(total_annual_budget = sum(budget)) %>%
ungroup() %>%
group_by(un_code) %>%
mutate(budget_proportion = budget / total_annual_budget) %>%
ungroup() %>%
group_by(country_name) %>%
summarise(country_sum = sum(budget_proportion, na.rm = TRUE)) %>%
mutate(rank_budget = rank(country_sum)) %>%
ungroup() %>%
mutate(var = ifelse(rank_budget < 150, "Rest", country_name)) %>%
group_by(var) %>%
summarise(sum_bud = sum(country_sum, na.rm = TRUE)) -> df3
# slice_max(order_by = rank_budget, n = 5)
# df4 <- df3 %>%
# mutate(csum = rev(cumsum(rev(var))),
# pos = sum_bud/2 + lead(csum, 1),
# pos = if_else(is.na(pos), sum_bud/2, pos),
# csum = round(csum, 2),
# pos = round(pos, 2))
df3 %>%
ggplot(aes(x="", y = sum_bud,
fill = fct_inorder(as.factor(var)))) +
geom_bar(stat = "identity") +
coord_polar("y", start = 1) +
scale_fill_brewer(palette = "Set1") +
guides(fill = guide_legend(title = "Country")) +
theme_void() +
theme(legend.position = "bottom") +
ggtitle("Comparing the rest of the countries (purple) with \n
Afghanistan, Iraq, Pakistan proportion of total budget")
# ggrepel::geom_label_repel(data = df4, aes(y = round(pos, 2), label = paste0(pos, "%")), size = 5, nudge_x = 1, show.legend = FALSE) +
Smallest proportion of budget (EXCLUDING countries under 1 million)
pd_df %>%
filter(pop > 999999) %>%
mutate(total_annual_budget = sum(budget)) %>%
ungroup() %>%
group_by(un_code) %>%
mutate(budget_proportion = budget / total_annual_budget) %>%
ungroup() %>%
group_by(country_name) %>%
summarise(country_sum = sum(budget_proportion, na.rm = TRUE)) %>%
arrange(desc(country_sum)) %>%
tail(n = 10) %>% kbl() %>% kable_paper("hover", full_width = F)
| country_name | country_sum |
|---|---|
| Lesotho | 0.0011192 |
| Timor-Leste | 0.0011032 |
| Mauritius | 0.0010210 |
| Somalia | 0.0008679 |
| Congo, Rep. | 0.0008091 |
| Gabon | 0.0007296 |
| Central African Republic | 0.0005856 |
| Equatorial Guinea | 0.0005139 |
| Gambia, The | 0.0005050 |
| Guinea-Bissau | 0.0001830 |
Smallest proportion of budget (including countries under 1 million)
pd_df %>%
# filter(pop > 999999) %>%
group_by(year) %>%
mutate(total_annual_budget = sum(budget)) %>%
ungroup() %>%
group_by(un_code, year) %>%
mutate(budget_proportion = budget / total_annual_budget) %>%
ungroup() %>%
group_by(country_name) %>%
summarise(country_sum = sum(budget_proportion, na.rm = TRUE)) %>%
arrange(desc(country_sum)) %>%
tail(n = 10) %>% kbl() %>% kable_paper("hover", full_width = F)
| country_name | country_sum |
|---|---|
| Gambia, The | 0.0036858 |
| Suriname | 0.0029780 |
| Malta | 0.0028351 |
| Samoa | 0.0024180 |
| Belize | 0.0023101 |
| Micronesia, Fed. Sts. | 0.0020952 |
| Guyana | 0.0018022 |
| Guinea-Bissau | 0.0012978 |
| Palau | 0.0008973 |
| Marshall Islands | 0.0008554 |
library(skimr)
# skim(pd_df)
pd_df %>%
ggplot(aes(x=log(budget))) +
geom_histogram( aes(fill = as.factor(year))) + facet_wrap(~year)
pd_df %>%
ggplot(aes(x=log(budget_pc))) +
geom_histogram( aes(fill = as.factor(year))) + facet_wrap(~year)
Looking at the total budget across different Department of State bureau over the years 2014 to 2019
pd_df %>%
filter(year != 2013) %>%
group_by(bureau, year) %>%
summarise(sum_budget = sum(budget, na.rm = TRUE)) %>%
ggplot(aes(x = year, y = sum_budget)) +
geom_line(aes(color = bureau), size = 2, alpha = 0.7) +
bbplot::bbc_style() +
theme(legend.position = "bottom",
legend.text = element_text(size = 10))
Looking at the budget PER CAPITA across different Department of State bureau over the years 2014 to 2019
pd_df %>%
filter(year != 2013) %>%
group_by(bureau, year) %>%
summarise(sum_budget = sum(budget_pc, na.rm = TRUE)) %>%
ggplot(aes(x = year, y = sum_budget)) +
geom_line(aes(color = bureau), size = 2, alpha = 0.7) +
bbplot::bbc_style() +
theme(legend.position = "bottom",
legend.text = element_text(size = 10))
Looking at the total budget across different economic groups over the years 2014 to 2019
pd_df %>%
filter(year != 2013) %>%
filter(!is.na(economy)) %>%
group_by(economy, year) %>%
summarise(sum_budget = sum(budget, na.rm = TRUE)) %>%
ggplot(aes(x = year, y = sum_budget)) +
geom_line(aes(color = economy), size = 2, alpha = 0.7) +
bbplot::bbc_style() +
theme(legend.position = "bottom",
legend.text = element_text(size = 10))
Looking at the budget PER CAPITA across different economic groups over the years 2014 to 2019
pd_df %>%
filter(year != 2013) %>%
filter(!is.na(economy)) %>%
group_by(economy, year) %>%
summarise(sum_budget = sum(budget_pc, na.rm = TRUE)) %>%
ggplot(aes(x = year, y = sum_budget)) +
geom_line(aes(color = economy), size = 2, alpha = 0.7) +
bbplot::bbc_style() +
theme(legend.position = "bottom",
legend.text = element_text(size = 10))
Looking at sum of budget PER CAPITA across different economic groups over the years 2014 to 2019 with populations OVER 1 MILLION
pd_df %>%
filter(pop > 1000000) %>%
filter(year != 2013) %>%
filter(!is.na(economy)) %>%
group_by(economy, year) %>%
summarise(sum_budget = sum(budget_pc, na.rm = TRUE)) %>%
ggplot(aes(x = year, y = sum_budget)) +
geom_line(aes(color = economy), size = 2, alpha = 0.7) +
bbplot::bbc_style() +
theme(legend.position = "bottom",
legend.text = element_text(size = 10))
Looking at mean total budgets specifically in DEVELOPED regions
plyr::count(pd_df$economy)
## x freq
## 1 1. Developed region: G7 42
## 2 2. Developed region: nonG7 216
## 3 3. Emerging region: BRIC 28
## 4 4. Emerging region: MIKT 28
## 5 5. Emerging region: G20 126
## 6 6. Developing region 484
## 7 7. Least developed region 287
## 8 <NA> 11
pd_df %>%
filter(year != 2013) %>%
filter(!is.na(economy)) %>%
filter(economy == "1. Developed region: G7" | economy == "2. Developed region: nonG7") %>%
group_by(economy, year) %>%
summarise(sum_budget = mean(budget, na.rm = TRUE)) %>%
ggplot(aes(x = year, y = sum_budget)) +
geom_line(aes(color = economy), size = 2, alpha = 0.7) +
bbplot::bbc_style() +
theme(legend.position = "bottom",
legend.text = element_text(size = 10))
Looking at mean total budgets specifically in EMERGING regions
pd_df %>%
filter(year != 2013) %>%
filter(!is.na(economy)) %>%
filter(economy == "3. Emerging region: BRIC" | economy == "4. Emerging region: MIKT" | economy == "5. Emerging region: G20") %>%
group_by(economy, year) %>%
summarise(sum_budget = mean(budget, na.rm = TRUE)) %>%
ggplot(aes(x = year, y = sum_budget)) +
geom_line(aes(color = economy), size = 2, alpha = 0.7) +
bbplot::bbc_style() +
theme(legend.position = "bottom",
legend.text = element_text(size = 10))
Looking at mean total budgets specifically in developed and least developed regions
pd_df %>%
filter(year != 2013) %>%
filter(!is.na(economy)) %>%
filter(economy == "6. Developing region" | economy == "7. Least developed region" ) %>%
group_by(economy, year) %>%
summarise(sum_budget = mean(budget, na.rm = TRUE)) %>%
ggplot(aes(x = year, y = sum_budget)) +
geom_line(aes(color = economy), size = 2, alpha = 0.7) +
bbplot::bbc_style() +
theme(legend.position = "bottom",
legend.text = element_text(size = 10))
Mean total budget across all income groups
pd_df %>%
filter(year != 2013) %>%
filter(!is.na(economy)) %>%
# filter(economy == "6. Developing region" | economy == "7. Least developed region" ) %>%
group_by(economy, year) %>%
summarise(sum_budget = mean(budget, na.rm = TRUE)) %>%
ggplot(aes(x = year, y = sum_budget)) +
geom_line(aes(color = economy), size = 2, alpha = 0.7) +
bbplot::bbc_style() +
theme(legend.position = "bottom",
legend.text = element_text(size = 10))
Mean budget PER CAPITA across all income groups
pd_df %>%
filter(year != 2013) %>%
filter(!is.na(economy)) %>%
# filter(economy == "6. Developing region" | economy == "7. Least developed region" ) %>%
group_by(economy, year) %>%
summarise(sum_budget = mean(budget_pc, na.rm = TRUE)) %>%
ggplot(aes(x = year, y = sum_budget)) +
geom_line(aes(color = economy), size = 2, alpha = 0.7) +
bbplot::bbc_style() +
theme(legend.position = "bottom",
legend.text = element_text(size = 10))
Mean budget PER CAPITA across all income groups with populations ABOVE 1 million
pd_df %>%
filter(pop > 1000000) %>%
filter(year != 2013) %>%
filter(!is.na(economy)) %>%
# filter(economy == "6. Developing region" | economy == "7. Least developed region" ) %>%
group_by(economy, year) %>%
summarise(sum_budget = mean(budget_pc, na.rm = TRUE)) %>%
ggplot(aes(x = year, y = sum_budget)) +
geom_line(aes(color = economy), size = 2, alpha = 0.7) +
bbplot::bbc_style() +
theme(legend.position = "bottom",
legend.text = element_text(size = 10))
Next we can look at total budgets different income groups across the years
Total budgets in each income group
pd_df$income_group_name <- as.factor(pd_df$income_group_name)
pd_df$income_group_name <- factor(pd_df$income_group_name, levels = c("Low Income Country", "Lower Middle Income Country", "Upper Middle Income Country", "High Income Country"))
plyr::count(pd_df$income_group_name)
## x freq
## 1 Low Income Country 195
## 2 Lower Middle Income Country 278
## 3 Upper Middle Income Country 333
## 4 High Income Country 357
## 5 <NA> 59
pd_df %>%
# filter(pop > 1000000) %>%
filter(!is.na(income_group_name)) %>%
filter(income_group_name != "NULL") %>%
filter(year != 2013) %>%
group_by(income_group_name, year) %>%
summarise(sum_budget = sum(budget, na.rm = TRUE)) %>%
ggplot(aes(x = year, y = sum_budget)) +
# geom_point(alpha = 0.8) +
geom_line(aes(color = income_group_name), size = 2, alpha = 0.7) +
bbplot::bbc_style() +
theme(legend.position = "bottom",
legend.text = element_text(size = 10))
Next we can look at the relationship between total budget and GDP with populations above 1 million
pd_df %>%
# filter(year == 2019) %>%
filter(pop > 1000000) %>%
# filter(!is.na(OECD_member)) %>%
filter(!is.na(income_group_name)) %>%
filter(income_group_name != "NULL") %>%
filter(year != 2013) %>%
# group_by(income_group_name, year) %>%
# summarise(sum_budget = sum(budget, na.rm = TRUE)) %>%
ggplot(aes(x = log(GDP_per_cap_constant_2010), y = log(budget), color = as.factor(income_group_name))) +
geom_point( alpha = 0.3) +
geom_smooth(method = "glm", se = FALSE, size = 2) +
facet_wrap(~ year) +
bbplot::bbc_style() +
theme(legend.position = "bottom",
legend.text = element_text(size = 10))
Next we can look at the relationship between budget PER CAPITA and GDP with populations above 1 million
pd_df %>%
filter(pop > 1000000) %>%
# filter(!is.na(OECD_member)) %>%
filter(!is.na(income_group_name)) %>%
filter(income_group_name != "NULL") %>%
filter(year != 2013) %>%
# group_by(income_group_name, year) %>%
# summarise(sum_budget = sum(budget, na.rm = TRUE)) %>%
ggplot(aes(x = log(GDP_per_cap_constant_2010), y = log(budget_pc), color = as.factor(income_group_name))) +
geom_point(aes(color = income_group_name), alpha = 0.4) +
geom_smooth(method = "glm", se = FALSE) +
facet_wrap(~ year) +
bbplot::bbc_style() +
theme(legend.position = "bottom",
legend.text = element_text(size = 10))
library(PerformanceAnalytics)
pd_df %>%
filter(year == 2019) %>%
# filter(pop > 1000000) %>%
select(budget, exports_constant_2010, imports_constant_2010, GDP_constant_2010, pop) %>%
mutate(budget_pc = budget / pop,
exports_pc = exports_constant_2010 / pop,
imports_pc = imports_constant_2010 / pop,
gdp_pc = GDP_constant_2010 / pop) %>%
select(budget_pc, exports_pc, imports_pc, gdp_pc, pop) -> pd_corr
Looking at the correlation between budget per capita and exports per capita, imports per capita and GDP per capita in 2015
pd_corr <- apply(pd_corr, 2, log)
chart.Correlation(pd_corr, histogram=TRUE, pch=19)
pd_df$income_group_name <- as.factor(pd_df$income_group_name)
pd_df$income_group_name <- factor(pd_df$income_group_name, levels = c("Low Income Country", "Lower Middle Income Country", "Upper Middle Income Country", "High Income Country"))
plyr::count(pd_df$income_group_name)
## x freq
## 1 Low Income Country 195
## 2 Lower Middle Income Country 278
## 3 Upper Middle Income Country 333
## 4 High Income Country 357
## 5 <NA> 59
pd_df %>%
# filter(pop > 1000000) %>%
filter(!is.na(income_group_name)) %>%
filter(income_group_name != "NULL") %>%
filter(year != 2013) %>%
group_by(income_group_name, year) %>%
summarise(sum_budget = sum(budget, na.rm = TRUE)) %>%
ggplot(aes(x = year, y = sum_budget)) +
# geom_point(alpha = 0.8) +
geom_line(aes(color = income_group_name), size = 2, alpha = 0.7) +
bbplot::bbc_style() +
theme(legend.position = "bottom",
legend.text = element_text(size = 10))
Next we can look at the relationship between total budget and GDP with populations above 1 million
pd_df %>%
# filter(year == 2019) %>%
filter(pop > 1000000) %>%
# filter(!is.na(OECD_member)) %>%
filter(!is.na(income_group_name)) %>%
filter(income_group_name != "NULL") %>%
filter(year != 2013) %>%
# group_by(income_group_name, year) %>%
# summarise(sum_budget = sum(budget, na.rm = TRUE)) %>%
ggplot(aes(x = log(GDP_per_cap_constant_2010), y = log(budget), color = as.factor(income_group_name))) +
geom_point( alpha = 0.3) +
geom_smooth(method = "glm", se = FALSE, size = 2) +
facet_wrap(~ year) +
bbplot::bbc_style() +
theme(legend.position = "bottom",
legend.text = element_text(size = 10))
Next we can look at the relationship between budget PER CAPITA and GDP with populations above 1 million
pd_df %>%
filter(pop > 1000000) %>%
# filter(!is.na(OECD_member)) %>%
filter(!is.na(income_group_name)) %>%
filter(income_group_name != "NULL") %>%
filter(year != 2013) %>%
# group_by(income_group_name, year) %>%
# summarise(sum_budget = sum(budget, na.rm = TRUE)) %>%
ggplot(aes(x = log(GDP_per_cap_constant_2010), y = log(budget_pc), color = as.factor(income_group_name))) +
geom_point(aes(color = income_group_name), alpha = 0.4) +
geom_smooth(method = "glm", se = FALSE) +
facet_wrap(~ year) +
theme(legend.position = "bottom") + bbplot::bbc_style() +
theme(legend.position = "bottom",
legend.text = element_text(size = 10))
Looking at correlation matrix in 2019
pd_df %>%
filter(year == 2015) %>%
# filter(pop > 1000000) %>%
mutate(budget_pc = budget / pop,
exports_pc = exports_constant_2010 / pop,
imports_pc = imports_constant_2010 / pop,
ln_budget = log(budget),
ln_budget_pc = log(budget_pc),
ln_manufac_exports = log(manufac_exports + 1),
ln_manufac_imports = log(manufac_imports + 1),
ln_energy_exports = log(energy_exports + 1),
ln_energy_imports = log(energy_imports + 1),
ln_agri_exports = log(agri_exports + 1),
ln_agri_imports = log(agri_imports + 1)) %>%
select(ln_budget, ln_budget_pc, ln_manufac_exports, ln_manufac_imports, ln_energy_exports, ln_energy_imports, ln_agri_exports, ln_agri_imports) -> trade_corr
Looking at the correlation between budget per capita and different exports from US and imports to US in 2015
# trade_corr <- apply(trade_corr, 2, log)
chart.Correlation(trade_corr, histogram=TRUE, pch=19)
Looking at the exports from US / imports to US and their correlation to PD budgets in Higher Income Countries in 2015
pd_df %>%
filter(income_group_name == "High Income Country" | income_group_name == "Higher Middle Income Country") %>%
filter(year == 2015) %>%
# filter(pop > 1000000) %>%
mutate(budget_pc = budget / pop,
exports_pc = exports_constant_2010 / pop,
imports_pc = imports_constant_2010 / pop,
ln_budget = log(budget),
ln_budget_pc = log(budget_pc),
ln_manufac_exports = log(manufac_exports + 1),
ln_manufac_imports = log(manufac_imports + 1),
ln_energy_exports = log(energy_exports + 1),
ln_energy_imports = log(energy_imports + 1),
ln_agri_exports = log(agri_exports + 1),
ln_agri_imports = log(agri_imports + 1)) %>%
select(ln_budget, ln_budget_pc, ln_manufac_exports, ln_manufac_imports, ln_energy_exports, ln_energy_imports, ln_agri_exports, ln_agri_imports) -> trade_corr
chart.Correlation(trade_corr, histogram=TRUE, pch=19)
Looking at the exports from US / imports to US and their correlation to PD budgets in Lower Income Countries in 2015
pd_df %>%
filter(income_group_name == "Low Income Country" | income_group_name == "Lower Middle Income Country") %>%
filter(year == 2015) %>%
# filter(pop > 1000000) %>%
mutate(budget_pc = budget / pop,
exports_pc = exports_constant_2010 / pop,
imports_pc = imports_constant_2010 / pop,
ln_budget = log(budget),
ln_budget_pc = log(budget_pc),
ln_manufac_exports = log(manufac_exports + 1),
ln_manufac_imports = log(manufac_imports + 1),
ln_energy_exports = log(energy_exports + 1),
ln_energy_imports = log(energy_imports + 1),
ln_agri_exports = log(agri_exports + 1),
ln_agri_imports = log(agri_imports + 1)) %>%
select(ln_budget, ln_budget_pc, ln_manufac_exports, ln_manufac_imports, ln_energy_exports, ln_energy_imports, ln_agri_exports, ln_agri_imports) -> trade_corr
chart.Correlation(trade_corr, histogram=TRUE, pch=19)
Looking at correlation between total budget and polity score, percentage muslim, percentage christian, distance to DC
pd_df %>%
filter(year == 2015) %>%
# filter(pop > 1000000) %>%
mutate(budget_pc = budget / pop) %>%
mutate(ln_budget_pc = log(budget_pc)) %>%
mutate(ln_budget = log(budget),
ln_oda_pc = log(oda_per_cap_current + 1)) %>%
select(ln_budget, ln_budget_pc, ln_oda_pc, percent_muslim, percent_muslim, percent_christian, polyarchy_demo, cso_repress) -> more_corr
# more_corr <- apply(more_corr, 2, log)
chart.Correlation(more_corr, histogram=TRUE, pch=19)
We can look at 2019 only and look at populations over 1 million
The relationship between budget (per capita and total) and US bases, net FDI inflows, US military military aid
pd_df %>%
filter(year == 2019) %>%
filter(pop > 1000000) %>%
mutate(budget_pc = budget / pop) %>%
mutate(ln_budget_pc = log(budget_pc)) %>%
mutate(ln_budget = log(budget)) %>%
mutate(ln_bases = log(all_bases + 1)) %>%
mutate(ln_fdi_in_net = log(fdi_in_net + 1)) %>%
mutate(ln_sum_mil_aid = log(sum_mil_aid + 1)) %>%
select(ln_budget_pc, ln_budget, ln_bases, ln_fdi_in_net, ln_sum_mil_aid) -> more_corr
# more_corr <- apply(more_corr, 2, log)
chart.Correlation(more_corr, histogram=TRUE, pch=19)
We can look at the relationship between different globalization scores and budgets per capita in Europe (populations over 1 million)
pd_df %>%
# filter(year == 2016) %>%
filter(continent.x == "Europe") %>%
filter(pop > 1000000) %>%
select(budget_pc, kof,
# kof_df, kof_dj,
cultural_kof,
# cultural_kof_df, cultural_kof_dj,
finance_kof,
# finance_kof_df, finance_kof_dj,
info_kof,
# info_kof_df, info_kof_dj,
political_kof,
# political_kof_df, political_kof_dj
) -> kof_corr
kof_corr <- apply(kof_corr, 2, log)
chart.Correlation(kof_corr, histogram=TRUE)
We can look at the relationship between different globalization scores and budgets per capita in Africa (populations over 1 million)
pd_df %>%
filter(continent.x == "Africa") %>%
filter(pop > 1000000) %>%
select(budget_pc, kof,
# kof_df, kof_dj,
cultural_kof,
# cultural_kof_df, cultural_kof_dj,
finance_kof,
# finance_kof_df, finance_kof_dj,
info_kof,
# info_kof_df, info_kof_dj,
political_kof,
# political_kof_df, political_kof_dj
) -> kof_corr
kof_corr <- apply(kof_corr, 2, log)
chart.Correlation(kof_corr, histogram=TRUE)
We can look at the relationship between different globalization scores and budgets per capita in Asia (populations over 1 million)
pd_df %>%
# filter(year == 2016) %>%
filter(continent.x == "Asia") %>%
filter(pop > 1000000) %>%
select(budget_pc, kof,
# kof_df, kof_dj,
cultural_kof,
# cultural_kof_df, cultural_kof_dj,
finance_kof,
# finance_kof_df, finance_kof_dj,
info_kof,
# info_kof_df, info_kof_dj,
political_kof,
# political_kof_df, political_kof_dj
) -> kof_corr
kof_corr <- apply(kof_corr, 2, log)
chart.Correlation(kof_corr, histogram=TRUE)
Null Hyp: that there is NO relationship between US PD spending and GDP, population and democracy scores
Alt Hyp 1: that the US spends more PD budget in countries with higher GDP
Alt Hyp 2: that the US spends more PD budget in countries with lower GDP
pd_df %<>%
distinct(un_code, year, .keep_all = TRUE)
basic_model_fe_ind <- plm(log(budget_pc) ~ log(GDP_per_cap_constant_2010) + log(pop) + polity_score, data = pd_df, index = c("un_code", "year"), model = "within", effect = "individual")
basic_model_fe_time <- plm(log(budget_pc) ~ log(GDP_per_cap_constant_2010) + log(pop) + polity_score, data = pd_df, index = c("un_code", "year"), model = "within", effect = "time")
basic_model_fe_two <- plm(log(budget_pc) ~ log(GDP_per_cap_constant_2010) + log(pop) + polity_score, data = pd_df, index = c("un_code", "year"), model = "within", effect = "twoways")
basic_model_re <- plm(log(budget_pc) ~ log(GDP_per_cap_constant_2010) + log(pop) + polity_score, data = pd_df, index = c("un_code", "year"), model = "random")
stargazer(basic_model_fe_ind, basic_model_fe_time, basic_model_fe_two, basic_model_re,
column.labels =c("FE Individ", "FE Time", "FE Twoway", "RE Swar"),
type = "text")
##
## ====================================================================================================================
## Dependent variable:
## -------------------------------------------------------------------------------------
## log(budget_pc)
## FE Individ FE Time FE Twoway RE Swar
## (1) (2) (3) (4)
## --------------------------------------------------------------------------------------------------------------------
## log(GDP_per_cap_constant_2010) -0.663*** 0.128*** -0.109 0.075*
## (0.211) (0.015) (0.215) (0.040)
##
## log(pop) -7.395*** -0.501*** -5.417*** -0.561***
## (0.373) (0.014) (0.505) (0.038)
##
## polity_score 0.027*** 0.007* 0.025*** 0.013
## (0.011) (0.004) (0.009) (0.008)
##
## Constant 6.444***
## (0.753)
##
## --------------------------------------------------------------------------------------------------------------------
## Observations 1,109 1,109 1,109 1,109
## R2 0.321 0.572 0.125 0.184
## Adjusted R2 0.205 0.568 -0.031 0.181
## F Statistic 149.225*** (df = 3; 946) 488.745*** (df = 3; 1099) 44.837*** (df = 3; 940) 242.825***
## ====================================================================================================================
## Note: *p<0.1; **p<0.05; ***p<0.01
basic_total_fe_ind <- plm(log(budget) ~ log(GDP_per_cap_constant_2010) + log(pop) + polity_score, data = pd_df, index = c("un_code", "year"), model = "within", effect = "individual")
basic_total_fe_time <- plm(log(budget) ~ log(GDP_per_cap_constant_2010) + log(pop) + polity_score, data = pd_df, index = c("un_code", "year"), model = "within", effect = "time")
basic_total_fe_two <- plm(log(budget) ~ log(GDP_per_cap_constant_2010) + log(pop) + polity_score, data = pd_df, index = c("un_code", "year"), model = "within", effect = "twoways")
basic_total_re <- plm(log(budget) ~ log(GDP_per_cap_constant_2010) + log(pop) + polity_score, data = pd_df, index = c("un_code", "year"), model = "random")
stargazer(basic_total_fe_ind, basic_total_fe_time, basic_total_fe_two, basic_total_re,
column.labels =c("FE Individ", "FE Time", "FE Twoway", "RE Swar"),
type = "text")
##
## ====================================================================================================================
## Dependent variable:
## -------------------------------------------------------------------------------------
## log(budget)
## FE Individ FE Time FE Twoway RE Swar
## (1) (2) (3) (4)
## --------------------------------------------------------------------------------------------------------------------
## log(GDP_per_cap_constant_2010) -0.663*** 0.128*** -0.109 0.075*
## (0.211) (0.015) (0.215) (0.040)
##
## log(pop) -6.395*** 0.499*** -4.417*** 0.439***
## (0.373) (0.014) (0.505) (0.038)
##
## polity_score 0.027*** 0.007* 0.025*** 0.013
## (0.011) (0.004) (0.009) (0.008)
##
## Constant 6.444***
## (0.753)
##
## --------------------------------------------------------------------------------------------------------------------
## Observations 1,109 1,109 1,109 1,109
## R2 0.265 0.527 0.089 0.166
## Adjusted R2 0.139 0.523 -0.074 0.164
## F Statistic 113.796*** (df = 3; 946) 407.951*** (df = 3; 1099) 30.461*** (df = 3; 940) 133.374***
## ====================================================================================================================
## Note: *p<0.1; **p<0.05; ***p<0.01
First run the Hausman test.
If the p is smaller than 0.05, use fixed effects
plm::phtest(basic_model_re, basic_model_fe_ind)
##
## Hausman Test
##
## data: log(budget_pc) ~ log(GDP_per_cap_constant_2010) + log(pop) + ...
## chisq = 396.22, df = 3, p-value < 2.2e-16
## alternative hypothesis: one model is inconsistent
The p value is under 0.05 so we can focus on the fixed effects model
Run the Lagrange Multiplier Test for time effects
Breusch-Pagan type for unbalanced data
If the p is smaller than 0.05, we use time-fixed effects (null is NO time fixed effects)
If p is smaller than 0.05, use time-fixed effects (null is NO time fixed effects)
plmtest(basic_model_fe_ind, c("time"), type=("bp"))
##
## Lagrange Multiplier Test - time effects (Breusch-Pagan) for unbalanced
## panels
##
## data: log(budget_pc) ~ log(GDP_per_cap_constant_2010) + log(pop) + ...
## chisq = 2498.1, df = 1, p-value < 2.2e-16
## alternative hypothesis: significant effects
# NB alternatives to Breusch Pagan ("bp")
# honda
# ghm
# kw
The null is that no time-fixed effects are needed
pFtest(basic_model_fe_ind, basic_model_fe_time)
##
## F test for individual effects
##
## data: log(budget_pc) ~ log(GDP_per_cap_constant_2010) + log(pop) + ...
## F = 12.834, df1 = 153, df2 = 946, p-value < 2.2e-16
## alternative hypothesis: significant effects
# If p is smaller than 0.05, use time-fixed effects (null is NO time fixed effects)
basic_total_fe_time <- plm(log(budget) ~ log(GDP_per_cap_constant_2010) + log(pop) + polity_score + I(polity_score^2) + as.factor(subregion), data = pd_df, index = c("un_code", "year"), model = "within", effect = "time")
usaid_total_fe_time <- plm(log(sum_usaid + 1) ~ log(GDP_per_cap_constant_2010) + log(pop) + polity_score + I(polity_score^2) + as.factor(subregion), data = pd_df, index = c("un_code", "year"), model = "within", effect = "time")
plyr::count(pd_df$subregion)
## x freq
## 1 Australia and New Zealand 14
## 2 Caribbean 51
## 3 Central America 56
## 4 Central Asia 35
## 5 Eastern Africa 112
## 6 Eastern Asia 35
## 7 Eastern Europe 70
## 8 Melanesia 14
## 9 Micronesia 21
## 10 Middle Africa 56
## 11 Northern Africa 42
## 12 Northern America 7
## 13 Northern Europe 70
## 14 Polynesia 7
## 15 South-Eastern Asia 77
## 16 South America 84
## 17 Southern Africa 35
## 18 Southern Asia 42
## 19 Southern Europe 91
## 20 Western Africa 112
## 21 Western Asia 125
## 22 Western Europe 48
## 23 <NA> 11
stargazer(basic_total_fe_time,
# coeftest(basic_total_fe_time, vcovHC(basic_total_fe_time, type = "HC0")),
coeftest(basic_total_fe_time, vcovHC(basic_total_fe_time, type = "HC1")),
# coeftest(basic_total_fe_time, vcovHC(basic_total_fe_time, type = "HC2")),
coeftest(basic_total_fe_time, vcovHC(basic_total_fe_time, type = "HC3")),
usaid_total_fe_time,
type = "text")
##
## =============================================================================================================
## Dependent variable:
## ----------------------------------------------------------------------
## log(budget) log(sum_usaid + 1)
## panel coefficient panel
## linear test linear
## (1) (2) (3) (4)
## -------------------------------------------------------------------------------------------------------------
## log(GDP_per_cap_constant_2010) 0.051** 0.051 0.051 -3.087***
## (0.024) (0.048) (0.048) (0.175)
##
## log(pop) 0.490*** 0.490*** 0.490*** 1.130***
## (0.014) (0.034) (0.034) (0.109)
##
## polity_score 0.029*** 0.029*** 0.029*** 0.171***
## (0.005) (0.008) (0.008) (0.034)
##
## I(polity_score2) -0.002*** -0.002 -0.002 -0.035***
## (0.001) (0.002) (0.002) (0.007)
##
## as.factor(subregion)Caribbean -0.177 -0.177 -0.177 8.050***
## (0.201) (0.245) (0.249) (1.745)
##
## as.factor(subregion)Central America -0.269 -0.269* -0.269* 8.904***
## (0.197) (0.147) (0.149) (1.720)
##
## as.factor(subregion)Central Asia 0.569*** 0.569** 0.569** 8.322***
## (0.216) (0.226) (0.230) (1.829)
##
## as.factor(subregion)Eastern Africa -0.470** -0.470* -0.470* 4.030**
## (0.203) (0.256) (0.259) (1.767)
##
## as.factor(subregion)Eastern Asia -0.145 -0.145 -0.145 5.770***
## (0.202) (0.154) (0.157) (1.756)
##
## as.factor(subregion)Eastern Europe -0.097 -0.097 -0.097 5.578***
## (0.188) (0.161) (0.163) (1.692)
##
## as.factor(subregion)Melanesia -0.480* -0.480 -0.480 4.805**
## (0.251) (0.466) (0.492) (2.119)
##
## as.factor(subregion)Middle Africa -0.947*** -0.947*** -0.947*** 4.413**
## (0.208) (0.240) (0.243) (1.807)
##
## as.factor(subregion)Northern Africa -0.510** -0.510* -0.510* 7.172***
## (0.209) (0.287) (0.291) (1.790)
##
## as.factor(subregion)Northern America -0.372 -0.372*** -0.372*** 3.243
## (0.292) (0.054) (0.056) (2.266)
##
## as.factor(subregion)Northern Europe -0.370** -0.370*** -0.370*** 1.016
## (0.185) (0.111) (0.112) (1.696)
##
## as.factor(subregion)South-Eastern Asia -0.352* -0.352* -0.352* 5.632***
## (0.197) (0.181) (0.184) (1.717)
##
## as.factor(subregion)South America -0.332* -0.332* -0.332* 4.940***
## (0.187) (0.177) (0.179) (1.664)
##
## as.factor(subregion)Southern Africa 0.153 0.153 0.153 11.520***
## (0.208) (0.252) (0.257) (1.773)
##
## as.factor(subregion)Southern Asia 0.210 0.210 0.210 4.130**
## (0.213) (0.556) (0.564) (1.816)
##
## as.factor(subregion)Southern Europe -0.016 -0.016 -0.016 4.881***
## (0.184) (0.159) (0.160) (1.645)
##
## as.factor(subregion)Western Africa -0.884*** -0.884*** -0.884*** 3.291*
## (0.199) (0.228) (0.230) (1.743)
##
## as.factor(subregion)Western Asia 0.435** 0.435** 0.435** 8.163***
## (0.188) (0.183) (0.185) (1.672)
##
## as.factor(subregion)Western Europe -0.373* -0.373** -0.373** -0.462
## (0.191) (0.168) (0.170) (1.726)
##
## -------------------------------------------------------------------------------------------------------------
## Observations 1,109 979
## R2 0.650 0.657
## Adjusted R2 0.641 0.646
## F Statistic 87.150*** (df = 23; 1079) 78.973*** (df = 23; 949)
## =============================================================================================================
## Note: *p<0.1; **p<0.05; ***p<0.01
Null Hyp: that there is NO relationship between US PD spending in country j and the cultural affinity between country j and the US
Alt Hyp 1: that the US spends more PD budget in countries with whom it has higher cultural affinity
Alt Hyp 2: that the US spends more PD budget in countries with whom it has lower cultural affinity
culture_model_fe_ind <- plm(log(budget_pc) ~ agree_with_US_score + cultural_kof + log(dist_to_dc) + percent_christian + as.factor(OECD_member), data = pd_df, index = c("un_code", "year"), model = "within", effect = "individual")
culture_model_fe_time <- plm(log(budget_pc) ~ agree_with_US_score + cultural_kof + log(dist_to_dc) + percent_christian + as.factor(OECD_member), data = pd_df, index = c("un_code", "year"), model = "within", effect = "time")
culture_model_fe_two <- plm(log(budget_pc) ~ agree_with_US_score + cultural_kof + log(dist_to_dc) + percent_christian + as.factor(OECD_member), data = pd_df, index = c("un_code", "year"), model = "within", effect = "twoways")
culture_model_re <- plm(log(budget_pc) ~ agree_with_US_score + cultural_kof + log(dist_to_dc) + percent_christian + as.factor(OECD_member), data = pd_df, index = c("un_code", "year"), model = "random")
stargazer(culture_model_fe_ind, culture_model_fe_time, culture_model_fe_two, culture_model_re,
column.labels =c("FE Individ", "FE Time", "FE Twoway", "RE Swar"),
type = "text")
##
## =======================================================================================================
## Dependent variable:
## -------------------------------------------------------------------------------
## log(budget_pc)
## FE Individ FE Time FE Twoway RE Swar
## (1) (2) (3) (4)
## -------------------------------------------------------------------------------------------------------
## agree_with_US_score -1.676*** 1.231*** -1.294*** -1.490***
## (0.222) (0.355) (0.404) (0.218)
##
## cultural_kof -0.004 0.019*** 0.001 0.020***
## (0.010) (0.002) (0.007) (0.004)
##
## log(dist_to_dc) -0.074 -0.089
## (0.073) (0.170)
##
## percent_christian 0.152*** 0.001 0.007 0.003
## (0.054) (0.001) (0.042) (0.002)
##
## as.factor(OECD_member)1 -0.441 -0.985*** 0.084 -0.323
## (0.448) (0.123) (0.337) (0.216)
##
## Constant -1.746
## (1.600)
##
## -------------------------------------------------------------------------------------------------------
## Observations 954 954 954 954
## R2 0.078 0.155 0.013 0.068
## Adjusted R2 -0.111 0.146 -0.197 0.063
## F Statistic 16.668*** (df = 4; 791) 34.699*** (df = 5; 943) 2.617** (df = 4; 786) 69.495***
## =======================================================================================================
## Note: *p<0.1; **p<0.05; ***p<0.01
culture_model_fe_ind <- plm(log(budget) ~ agree_with_US_score + cultural_kof + log(dist_to_dc) + percent_christian + as.factor(OECD_member), data = pd_df, index = c("un_code", "year"), model = "within", effect = "individual")
culture_model_fe_time <- plm(log(budget) ~ agree_with_US_score + cultural_kof + log(dist_to_dc) + percent_christian + as.factor(OECD_member), data = pd_df, index = c("un_code", "year"), model = "within", effect = "time")
culture_model_fe_two <- plm(log(budget) ~ agree_with_US_score + cultural_kof + log(dist_to_dc) + percent_christian + as.factor(OECD_member), data = pd_df, index = c("un_code", "year"), model = "within", effect = "twoways")
culture_model_re <- plm(log(budget) ~ agree_with_US_score + cultural_kof + log(dist_to_dc) + percent_christian + as.factor(OECD_member), data = pd_df, index = c("un_code", "year"), model = "random")
stargazer(culture_model_fe_ind, culture_model_fe_time, culture_model_fe_two, culture_model_re,
column.labels =c("FE Individ", "FE Time", "FE Twoway", "RE Swar"),
type = "text")
##
## =======================================================================================================
## Dependent variable:
## -------------------------------------------------------------------------------
## log(budget)
## FE Individ FE Time FE Twoway RE Swar
## (1) (2) (3) (4)
## -------------------------------------------------------------------------------------------------------
## agree_with_US_score -1.612*** -0.357 -1.216*** -1.529***
## (0.214) (0.362) (0.396) (0.209)
##
## cultural_kof -0.005 0.008*** -0.003 0.009**
## (0.009) (0.002) (0.007) (0.004)
##
## log(dist_to_dc) 0.260*** 0.242
## (0.075) (0.174)
##
## percent_christian 0.152*** -0.004*** 0.016 -0.003
## (0.052) (0.001) (0.041) (0.002)
##
## as.factor(OECD_member)1 -0.467 0.168 0.007 0.308
## (0.432) (0.125) (0.332) (0.216)
##
## Constant 12.250***
## (1.633)
##
## -------------------------------------------------------------------------------------------------------
## Observations 960 960 960 960
## R2 0.078 0.060 0.012 0.060
## Adjusted R2 -0.111 0.050 -0.198 0.056
## F Statistic 16.869*** (df = 4; 796) 12.099*** (df = 5; 949) 2.453** (df = 4; 791) 61.426***
## =======================================================================================================
## Note: *p<0.1; **p<0.05; ***p<0.01
First run the Hausman test.
If the p is smaller than 0.05, use fixed effects
plm::phtest(culture_model_re, culture_model_fe_ind)
##
## Hausman Test
##
## data: log(budget) ~ agree_with_US_score + cultural_kof + log(dist_to_dc) + ...
## chisq = 17.094, df = 4, p-value = 0.001853
## alternative hypothesis: one model is inconsistent
The p value is under 0.05 so we can focus on the fixed effects model
Run the Lagrange Multiplier Test for time effects
Breusch-Pagan type for unbalanced data
If the p is smaller than 0.05, we use time-fixed effects (null is NO time fixed effects)
plmtest(culture_model_fe_ind, c("time"), type=("bp"))
##
## Lagrange Multiplier Test - time effects (Breusch-Pagan) for balanced
## panels
##
## data: log(budget) ~ agree_with_US_score + cultural_kof + log(dist_to_dc) + ...
## chisq = 518.05, df = 1, p-value < 2.2e-16
## alternative hypothesis: significant effects
# NB alternatives to Breusch Pagan ("bp")
# honda
# ghm
# kw
pFtest(culture_model_fe_ind, culture_model_fe_time)
##
## F test for individual effects
##
## data: log(budget) ~ agree_with_US_score + cultural_kof + log(dist_to_dc) + ...
## F = 19.031, df1 = 153, df2 = 796, p-value < 2.2e-16
## alternative hypothesis: significant effects
# If p is smaller than 0.05, use time-fixed effects (null is NO time fixed effects)
phtest(culture_model_re, culture_model_fe_ind)
##
## Hausman Test
##
## data: log(budget) ~ agree_with_US_score + cultural_kof + log(dist_to_dc) + ...
## chisq = 17.094, df = 4, p-value = 0.001853
## alternative hypothesis: one model is inconsistent
Null: that there is relationship between US PD spending in country j and the levels of economic activity between country j and the US
Alt 1: that the US spends more PD budget in countries with whom it has more economic interests
Alt 2: that the US spends more PD budget in countries with whom it has fewer economic interests
pd_df %<>%
mutate(budget_pc = budget / pop) %>%
mutate(energy_imp_pc = energy_imports / pop) %>%
mutate(agri_imp_pc = agri_imports / pop) %>%
mutate(services_imp_pc = services_imports / pop) %>%
mutate(manufac_imp_pc = manufac_imports / pop) %>%
mutate(energy_exp_pc = energy_exports / pop) %>%
mutate(agri_exp_pc = agri_exports / pop) %>%
mutate(services_exp_pc = services_exports / pop) %>%
mutate(manufac_exp_pc = manufac_exports / pop)
Is there a a relationship between PD budget per capita and and energy imports per capita to the US per capita?
pd_df %>%
filter(year < 2017) -> pd_trade
pd_trade %>%
# filter(year == 2013) %>%
filter(!is.na(income_group_name)) %>%
filter(income_group_name != "NULL") %>%
mutate(income_group_name = as.factor(income_group_name)) %>%
ggplot(aes(x = log(budget_pc), y = log(energy_imp_pc), group = as.factor(year))) +
geom_point(aes(color = as.factor(year))) +
geom_smooth(method = "lm", se = FALSE, aes(color = as.factor(year))) +
facet_wrap(~region_name)
Is there a a relationship between PD budget per capita and and manufacturing exports per capita to the US per capita?
pd_df %>%
filter(year < 2017) -> pd_trade
pd_trade %>%
# filter(year == 2013) %>%
filter(!is.na(income_group_name)) %>%
filter(income_group_name != "NULL") %>%
ggplot(aes(x = log(budget_pc), y = log(manufac_exp_pc), group = as.factor(year))) +
geom_point(aes(color = as.factor(year))) +
geom_smooth(method = "lm", se = FALSE, aes(color = as.factor(year))) +
facet_wrap(~income_group_name)
# pd_df %>%
# filter(pop > 1000000) -> pd_small
#
# econ_model_fe_ind <- plm(log(budget_pc) ~ log(trade_openess) + log(GDP_per_cap_constant_2010) + log(pop),
# data = pd_small, index = c("un_code", "year"), model = "within", effect = "individual")
#
# econ_model_fe_time <- plm(log(budget_pc) ~ log(trade_openess) + log(GDP_per_cap_constant_2010) + log(pop),
# data = pd_small, index = c("un_code", "year"), model = "within", effect = "time")
#
# econ_model_fe_two <- plm(log(budget_pc) ~ log(trade_openess) + log(GDP_per_cap_constant_2010) + log(pop),
# data = pd_small, index = c("un_code", "year"), model = "within", effect = "twoways")
#
# econ_model_re <- plm(log(budget_pc) ~ log(trade_openess) + log(GDP_per_cap_constant_2010) + log(pop),
# data = pd_small, index = c("un_code", "year"), model = "random")
#
# stargazer(econ_model_fe_ind, econ_model_fe_time, econ_model_fe_two, econ_model_re,
# column.labels =c("FE Individ", "FE Time", "FE Twoway", "RE Swar"),
# type = "text")
#
#
# plm::phtest(econ_model_re, econ_model_fe_ind)
#
# plmtest(econ_model_fe_ind, c("time"), type=("bp"))
# NB alternatives to Breusch Pagan ("bp")
# honda
# ghm
# kw
# pFtest(econ_model_fe_ind, econ_model_fe_time)
# If p is smaller than 0.05, use time-fixed effects (null is NO time fixed effects)
Economic relations
# transform(INCOME = factor(income_group_name ,levels = c("Low Income Country","Lower Middle Income Country","Upper Middle Income Country", "High Income Country", "NULL")),
# log(energy_exp_pc + 1) + log(agri_exp_pc + 1) + log(manufac_exp_pc + 1) + log(energy_imp_pc + 1) + log(agri_imp_pc + 1) + log(manufac_imp_pc + 1)
econ_model_fe_ind <- plm(log(budget_pc) ~ log(fdi_in_per_cap + 1) + finance_kof + log(energy_imports + 1) + log(manufac_exports + 1),
data = pd_trade, index = c("un_code", "year"), model = "within", effect = "individual")
econ_model_fe_time <- plm(log(budget_pc) ~ log(fdi_in_per_cap + 1) + finance_kof + log(energy_imports + 1) + log(manufac_exports + 1),
data = pd_trade, index = c("un_code", "year"), model = "within", effect = "time")
econ_model_fe_two <- plm(log(budget_pc) ~ log(fdi_in_per_cap + 1) + finance_kof + log(energy_imports + 1) + log(manufac_exports + 1),
data = pd_trade, index = c("un_code", "year"), model = "within", effect = "twoways")
econ_model_re <- plm(log(budget_pc) ~ log(fdi_in_per_cap + 1) + finance_kof + log(energy_imports + 1) + log(manufac_exports + 1),
data = pd_trade, index = c("un_code", "year"), model = "random")
stargazer(econ_model_fe_ind, econ_model_fe_time, econ_model_fe_two, econ_model_re,
column.labels =c("FE Individ", "FE Time", "FE Twoway", "RE Swar"),
type = "text")
##
## ======================================================================================================
## Dependent variable:
## -----------------------------------------------------------------------------
## log(budget_pc)
## FE Individ FE Time FE Twoway RE Swar
## (1) (2) (3) (4)
## ------------------------------------------------------------------------------------------------------
## log(fdi_in_per_cap + 1) 0.090* 0.007 0.034 0.070
## (0.052) (0.042) (0.035) (0.047)
##
## finance_kof -0.056*** 0.034*** -0.009 0.017***
## (0.012) (0.003) (0.008) (0.005)
##
## log(energy_imports + 1) 0.069** -0.013 0.027 0.003
## (0.034) (0.016) (0.023) (0.024)
##
## log(manufac_exports + 1) 0.270*** -0.206*** 0.029 -0.117***
## (0.075) (0.021) (0.054) (0.035)
##
## Constant -2.102***
## (0.277)
##
## ------------------------------------------------------------------------------------------------------
## Observations 624 624 624 624
## R2 0.093 0.298 0.009 0.035
## Adjusted R2 -0.239 0.290 -0.363 0.028
## F Statistic 11.738*** (df = 4; 456) 65.288*** (df = 4; 616) 0.996 (df = 4; 453) 26.636***
## ======================================================================================================
## Note: *p<0.1; **p<0.05; ***p<0.01
First run the Hausman test.
If the p is smaller than 0.05, use fixed effects
plm::phtest(econ_model_re, econ_model_fe_ind)
##
## Hausman Test
##
## data: log(budget_pc) ~ log(fdi_in_per_cap + 1) + finance_kof + log(energy_imports + ...
## chisq = 101.03, df = 4, p-value < 2.2e-16
## alternative hypothesis: one model is inconsistent
The p value is under 0.05 so we can focus on the fixed effects model
Run the Lagrange Multiplier Test for time effects
Breusch-Pagan type for unbalanced data
If the p is smaller than 0.05, we use time-fixed effects (null is NO time fixed effects)
plmtest(econ_model_fe_ind, c("time"), type=("bp"))
##
## Lagrange Multiplier Test - time effects (Breusch-Pagan) for unbalanced
## panels
##
## data: log(budget_pc) ~ log(fdi_in_per_cap + 1) + finance_kof + log(energy_imports + ...
## chisq = 1596, df = 1, p-value < 2.2e-16
## alternative hypothesis: significant effects
# NB alternatives to Breusch Pagan ("bp")
# honda for balanced data
# kw
pFtest(culture_model_fe_ind, culture_model_fe_time)
##
## F test for individual effects
##
## data: log(budget) ~ agree_with_US_score + cultural_kof + log(dist_to_dc) + ...
## F = 19.031, df1 = 153, df2 = 796, p-value < 2.2e-16
## alternative hypothesis: significant effects
Dependent variable : ODA per capita received
Independent variables : * FDI inflows per capita * Financial globalization
oda_model_fe_ind <- plm(log(oda_per_cap_current + 1) ~ log(fdi_in_per_cap + 1) + finance_kof, data = pd_trade, index = c("un_code", "year"), model = "within", effect = "twoways")
oda_model_fe_time <- plm(log(oda_per_cap_current + 1) ~ log(fdi_in_per_cap + 1) + finance_kof, data = pd_trade, index = c("un_code", "year"), model = "within", effect = "time")
oda_model_fe_two <- plm(log(oda_per_cap_current + 1) ~ log(fdi_in_per_cap + 1) + finance_kof, data = pd_trade, index = c("un_code", "year"), model = "within", effect = "twoways")
oda_model_re <- plm(log(oda_per_cap_current + 1) ~ log(fdi_in_per_cap + 1) + finance_kof, data = pd_trade, index = c("un_code", "year"), model = "random")
stargazer(econ_model_fe_ind, econ_model_fe_time, econ_model_fe_two, econ_model_re,
oda_model_fe_ind, oda_model_fe_time, oda_model_fe_two, oda_model_re, type = "text")
##
## ==============================================================================================================================================================================
## Dependent variable:
## -----------------------------------------------------------------------------------------------------------------------------------------------------
## log(budget_pc) log(oda_per_cap_current + 1)
## (1) (2) (3) (4) (5) (6) (7) (8)
## ------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
## log(fdi_in_per_cap + 1) 0.090* 0.007 0.034 0.070 0.001 0.279*** 0.001 0.034
## (0.052) (0.042) (0.035) (0.047) (0.033) (0.079) (0.033) (0.033)
##
## finance_kof -0.056*** 0.034*** -0.009 0.017*** 0.003 0.005 0.003 0.004
## (0.012) (0.003) (0.008) (0.005) (0.007) (0.005) (0.007) (0.006)
##
## log(energy_imports + 1) 0.069** -0.013 0.027 0.003
## (0.034) (0.016) (0.023) (0.024)
##
## log(manufac_exports + 1) 0.270*** -0.206*** 0.029 -0.117***
## (0.075) (0.021) (0.054) (0.035)
##
## Constant -2.102*** 3.392***
## (0.277) (0.311)
##
## ------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
## Observations 624 624 624 624 448 448 448 448
## R2 0.093 0.298 0.009 0.035 0.001 0.037 0.001 0.009
## Adjusted R2 -0.239 0.290 -0.363 0.028 -0.370 0.026 -0.370 0.005
## F Statistic 11.738*** (df = 4; 456) 65.288*** (df = 4; 616) 0.996 (df = 4; 453) 26.636*** 0.088 (df = 2; 326) 8.441*** (df = 2; 442) 0.088 (df = 2; 326) 1.904
## ==============================================================================================================================================================================
## Note: *p<0.1; **p<0.05; ***p<0.01
pd_df %>%
# filter(year == 2015) %>%
ggplot(aes(x = cso_repress, y = log(budget + 1))) +
geom_smooth(method = "lm") +
geom_text(aes(label = country.x, color = income_group_name), size = 2) +
facet_wrap(~income_group_name, scales = "free")
pd_df %<>%
mutate(is_Africa = ifelse(continent.x == "Africa", 1, 0))
china_influence <- plm(log(budget_pc + 1) ~ common_relig_china +
as.factor(rta_china) +
log(tradeflow_with_china)*as.factor(is_Africa) +
# log(import_to_china_cow) +
# log(export_from_china_cow) +
int_org_with_china,
# log(dist_to_beijing),
# as.factor(contiguity_china),
data = pd_df,
index = c("un_code", "year"),
model = "random")
stargazer(china_influence, type = "text")
##
## ===========================================================================
## Dependent variable:
## ---------------------------
## log(budget_pc + 1)
## ---------------------------------------------------------------------------
## common_relig_china 1.688
## (3.594)
##
## as.factor(rta_china)1 -0.128*
## (0.077)
##
## log(tradeflow_with_china) -0.060***
## (0.014)
##
## as.factor(is_Africa)1 -0.708***
## (0.254)
##
## int_org_with_china 0.005
## (0.004)
##
## log(tradeflow_with_china):as.factor(is_Africa)1 0.042**
## (0.020)
##
## Constant 0.956***
## (0.170)
##
## ---------------------------------------------------------------------------
## Observations 194
## R2 0.140
## Adjusted R2 0.112
## F Statistic 30.436***
## ===========================================================================
## Note: *p<0.1; **p<0.05; ***p<0.01
pd_df %<>%
mutate(is_Africa = ifelse(continent.x == "Africa", 1, 0))
russia_influence <- plm(log(budget_pc + 1) ~
common_relig_russia +
log(dist_to_moscow) +
# as.factor(rta_russia) +
as.factor(common_ethn_lang_russia) +
as.factor(colonial_relationship_russia) +
log(tradeflow_with_russia) +
log(int_org_with_russia),
# log(import_to_china_cow) +
# log(export_from_china_cow) +
# log(dist_to_beijing),
# as.factor(contiguity_china),
data = pd_df,
index = c("un_code", "year"),
model = "within", effect = "time")
stargazer(russia_influence, type = "text")
##
## ====================================================================
## Dependent variable:
## ---------------------------
## log(budget_pc + 1)
## --------------------------------------------------------------------
## common_relig_russia -1.189**
## (0.572)
##
## log(dist_to_moscow) -0.127***
## (0.030)
##
## as.factor(common_ethn_lang_russia)1 0.124
## (0.095)
##
## as.factor(colonial_relationship_russia)1 0.016
## (0.185)
##
## log(tradeflow_with_russia) -0.016**
## (0.007)
##
## log(int_org_with_russia) -0.215
## (0.137)
##
## --------------------------------------------------------------------
## Observations 191
## R2 0.159
## Adjusted R2 0.127
## F Statistic 5.755*** (df = 6; 183)
## ====================================================================
## Note: *p<0.1; **p<0.05; ***p<0.01
Scatterplot to look at the relationship between total budget and number of US military bases in each region across the years
pd_df %>%
filter(year != 2013) %>%
filter(country.x != 'Syria') %>%
# filter(all_bases != 0) %>%
ggplot(aes(x = log(budget), y = log(all_bases))) +
geom_point(aes(color = bureau), alpha = 0.4) +
geom_smooth(method = "lm", se = FALSE, aes(color = bureau)) +
# geom_text(aes(label = country.x, color = bureau)) +
facet_wrap(~year) + bbplot::bbc_style()
Look at number of US military bases in the country and budget PER CAPITA
pd_df %>%
filter(year != 2013) %>%
filter(country.x != 'Syria') %>%
# filter(all_bases != 0) %>%
ggplot(aes(x = log(budget_pc), y = log(all_bases))) +
geom_point(aes(color = bureau), alpha = 0.4) +
geom_smooth(method = "lm", se = FALSE, aes(color = bureau)) +
# geom_text(aes(label = country.x, color = bureau)) +
facet_wrap(~year) + bbplot::bbc_style()
They look very similar, so I look at 2014 more closely
We can look only at countries with bases (remove all countries from the scatterplot with no military bases) in 2014
pd_df %>%
filter(country.x != 'Syria') %>%
filter(year == 2014) %>%
filter(all_bases != 0) %>%
ggplot(aes(x = log(budget), y = log(all_bases))) +
# geom_point(aes(color = bureau), alpha = 0.4) +
geom_smooth(method = "lm", se = FALSE, aes(color = bureau), alpha = 0.2) +
geom_text(aes(label = country.x, color = bureau)) +
facet_wrap(~year)
pd_df %<>%
distinct(un_code, year, .keep_all = TRUE)
military_model_time <- plm(log(budget) ~ log(military_exp + 1) + all_bases*as.factor(al_qaeda_main) + log(troops + 1) ,
data = pd_df,
index = c("un_code", "year"),
model = "within", effect = "time")
military_model<- plm(log(budget) ~ log(military_exp + 1) + all_bases*as.factor(al_qaeda_main) + log(troops + 1)
,
data = pd_df,
index = c("un_code", "year"),
model = "within")
military_model_two <- plm(log(budget) ~log(military_exp + 1) + all_bases*as.factor(al_qaeda_main) + log(troops + 1)
,
data = pd_df,
index = c("un_code", "year"),
model = "within", effect = "twoway")
military_model_re <- plm(log(budget) ~ log(military_exp + 1) + all_bases*as.factor(al_qaeda_main) + log(troops + 1)
,
data = pd_df,
index = c("un_code", "year"),
model = "random")
stargazer(military_model, military_model_time, military_model_two, military_model_re,
type = "text")
##
## ================================================================================================================
## Dependent variable:
## ----------------------------------------------------------------------------
## log(budget)
## (1) (2) (3) (4)
## ----------------------------------------------------------------------------------------------------------------
## log(military_exp + 1) -0.117* 0.078*** -0.084 0.095***
## (0.067) (0.009) (0.069) (0.018)
##
## all_bases -0.003 0.032**
## (0.007) (0.013)
##
## as.factor(al_qaeda_main)1 3.535*** 2.033
## (1.267) (2.868)
##
## log(troops + 1) -0.041** 0.156*** -0.029 0.007
## (0.021) (0.017) (0.023) (0.019)
##
## all_bases:as.factor(al_qaeda_main)1 -0.102 0.143
## (0.167) (0.375)
##
## Constant 12.079***
## (0.353)
##
## ----------------------------------------------------------------------------------------------------------------
## Observations 888 888 888 888
## R2 0.010 0.377 0.004 0.333
## Adjusted R2 -0.200 0.370 -0.215 0.329
## F Statistic 3.673** (df = 2; 732) 106.313*** (df = 5; 877) 1.447 (df = 2; 727) 80.097***
## ================================================================================================================
## Note: *p<0.1; **p<0.05; ***p<0.01
plm::phtest(military_model_re, military_model)
##
## Hausman Test
##
## data: log(budget) ~ log(military_exp + 1) + all_bases * as.factor(al_qaeda_main) + ...
## chisq = 39.146, df = 2, p-value = 3.159e-09
## alternative hypothesis: one model is inconsistent
pFtest( military_model, military_model_time)
##
## F test for individual effects
##
## data: log(budget) ~ log(military_exp + 1) + all_bases * as.factor(al_qaeda_main) + ...
## F = 23.957, df1 = 145, df2 = 732, p-value < 2.2e-16
## alternative hypothesis: significant effects
plmtest(military_model_two, c("ind"), type=("bp"))
##
## Lagrange Multiplier Test - (Breusch-Pagan) for unbalanced panels
##
## data: log(budget) ~ log(military_exp + 1) + all_bases * as.factor(al_qaeda_main) + ...
## chisq = 1226.2, df = 1, p-value < 2.2e-16
## alternative hypothesis: significant effects
Looking at the role of information globalization and level of government censorship of the media.
Do we see a negative relatiosni
censor_model_time <- plm(log(budget_pc) ~ gov_censor_media + info_kof,
data = pd_df,
index = c("un_code", "year"),
model = "within", effect = "time")
censor_model_ind <- plm(log(budget_pc) ~ gov_censor_media + info_kof,
data = pd_df,
index = c("un_code", "year"),
model = "within", effect = "individual")
censor_model_two <- plm(log(budget_pc) ~ gov_censor_media + info_kof,
data = pd_df,
index = c("un_code", "year"),
model = "within", effect = "twoways")
censor_model_swar <- plm(log(budget_pc) ~ gov_censor_media + info_kof,
data = pd_df,
index = c("un_code", "year"),
model = "random")
censor_model_wal <- plm(log(budget_pc) ~ gov_censor_media + info_kof,
data = pd_df,
index = c("un_code", "year"),
model = "random",
random.method = "walhus")
censor_model_ner <- plm(log(budget_pc) ~ gov_censor_media + info_kof,
data = pd_df,
index = c("un_code", "year"),
model = "random",
random.method = "nerlove")
censor_model_amemiya <- plm(log(budget_pc) ~ gov_censor_media + info_kof,
data = pd_df,
index = c("un_code", "year"),
model = "random",
random.method = "amemiya")
stargazer(censor_model_time, censor_model_ind,censor_model_two, censor_model_swar,
column.labels = c("Time FE", "Individual FE", "Twoway FE","Swamy RE"), type = "text")
##
## ==============================================================================================
## Dependent variable:
## -----------------------------------------------------------------------------
## log(budget_pc)
## Time FE Individual FE Twoway FE Swamy RE
## (1) (2) (3) (4)
## ----------------------------------------------------------------------------------------------
## gov_censor_media 0.002 0.213*** 0.073** 0.159***
## (0.025) (0.049) (0.037) (0.039)
##
## info_kof 0.024*** -0.056*** 0.001 0.002
## (0.003) (0.010) (0.008) (0.005)
##
## Constant -2.168***
## (0.368)
##
## ----------------------------------------------------------------------------------------------
## Observations 959 959 959 959
## R2 0.103 0.056 0.005 0.020
## Adjusted R2 0.097 -0.135 -0.204 0.018
## F Statistic 54.717*** (df = 2; 951) 23.488*** (df = 2; 797) 1.988 (df = 2; 792) 19.356***
## ==============================================================================================
## Note: *p<0.1; **p<0.05; ***p<0.01
phtest(censor_model_swar, censor_model_ind)
##
## Hausman Test
##
## data: log(budget_pc) ~ gov_censor_media + info_kof
## chisq = 48.81, df = 2, p-value = 2.518e-11
## alternative hypothesis: one model is inconsistent
pFtest(censor_model_ind, censor_model_time)
##
## F test for individual effects
##
## data: log(budget_pc) ~ gov_censor_media + info_kof
## F = 17.346, df1 = 154, df2 = 797, p-value < 2.2e-16
## alternative hypothesis: significant effects
plm::pooltest(log(budget_pc) ~ gov_censor_media + info_kof,
data = pd_df,
index = c("un_code", "year"),
model = "within")
##
## F statistic
##
## data: log(budget_pc) ~ gov_censor_media + info_kof
## F = 1.4723, df1 = 318, df2 = 479, p-value = 6.571e-05
## alternative hypothesis: unstability
pFtest(log(budget_pc) ~ gov_censor_media + info_kof,
data = pd_df,
index = c("un_code", "year"),
model = "within", effect = "time", type = "ghm")
##
## F test for time effects
##
## data: log(budget_pc) ~ gov_censor_media + info_kof
## F = 23.937, df1 = 5, df2 = 951, p-value < 2.2e-16
## alternative hypothesis: significant effects
## PUBLIC DIPLOMACY
# pd_df %>%
# filter(year == 2014) %>%
# # filter(!is.na(region_name)) %>%
# filter(income_group_name != "NULL") %>%
# mutate(mil_aid_pc = sum_mil_aid / pop) %>%
# mutate(mil_exp_pc = military_exp / pop) %>%
# ggplot(aes(x = log(budget_pc), y = trade_kof_dj)) + geom_point(aes(color = income_group_name)) + geom_smooth(method = "lm") + facet_wrap(~region_name)
pd_df %>%
filter(year != 2019) %>%
filter(income_group_name != "NULL") %>%
filter(!is.na(income_group_name)) %>%
ggplot(aes(y = log(budget_pc), x = trade_kof_dj)) + geom_point(aes(color = income_group_name)) + geom_smooth(method = "lm") + facet_wrap(~year)
## ODA
pd_df %>%
filter(year == 2016) %>%
ggplot(aes(x = log(oda_per_cap_current), y = trade_kof_dj)) + geom_point(aes(color = region_name)) + geom_smooth(method = "lm") + facet_wrap(~bureau)
pd_df %>%
# filter(year == 2014) %>%
# filter(region_name == "Sub-Saharan Africa") %>%
filter(!is.na(region_name)) %>%
# filter(continent.x != "Pacific") %>%
ggplot(aes(x = log(budget_pc), y = log(dist_to_dc))) +
geom_point(aes(color = region_name)) +
geom_smooth(method = "lm") +
facet_wrap(~year) +
bbplot::bbc_style() +
theme(legend.position = "bottom",
legend.text = element_text(size = 10))
pd_df %>%
# filter(year == 2014) %>%
# filter(region_name == "Sub-Saharan Africa") %>%
filter(!is.na(region_name)) %>%
# filter(continent.x != "Pacific") %>%
ggplot(aes(x = log(budget_pc), y = log(dist_to_beijing))) +
geom_point(aes(color = region_name)) +
geom_smooth(method = "lm") +
facet_wrap(~year) +
bbplot::bbc_style() +
theme(legend.position = "bottom",
legend.text = element_text(size = 10))
pd_df %>%
filter(year == 2014) %>%
# filter(region_name == "Sub-Saharan Africa") %>%
filter(!is.na(region_name)) %>%
# filter(continent.x != "Pacific") %>%
ggplot(aes(x = log(budget_pc), y = log(dist_to_dc))) + geom_point(aes(color = region_name)) + geom_smooth(method = "lm") + facet_wrap(~bureau) +
bbplot::bbc_style() +
theme(legend.position = "bottom",
legend.text = element_text(size = 10))
censor_model_time <- plm(log(budget_pc) ~ rule_of_law,
data = pd_df,
index = c("un_code", "year"),
model = "within", effect = "time")
censor_model_ind <- plm(log(budget_pc) ~rule_of_law,
data = pd_df,
index = c("un_code", "year"),
model = "within", effect = "individual")
censor_model_two <- plm(log(budget_pc) ~ rule_of_law,
data = pd_df,
index = c("un_code", "year"),
model = "within", effect = "twoways")
censor_model_swar <- plm(log(budget_pc) ~ rule_of_law,
data = pd_df,
index = c("un_code", "year"),
model = "random")
censor_model_wal <- plm(log(budget_pc) ~ rule_of_law,
data = pd_df,
index = c("un_code", "year"),
model = "random",
random.method = "walhus")
censor_model_ner <- plm(log(budget_pc) ~ rule_of_law,
data = pd_df,
index = c("un_code", "year"),
model = "random",
random.method = "nerlove")
censor_model_amemiya <- plm(log(budget_pc) ~ rule_of_law,
data = pd_df,
index = c("un_code", "year"),
model = "random",
random.method = "amemiya")
stargazer(censor_model_time, censor_model_ind,censor_model_two, censor_model_swar, censor_model_wal, censor_model_ner, censor_model_amemiya,
column.labels = c( "Time FE", "Individual FE", "Twoway FE","Swamy RE", "Walhus RE", "Nerlove RE","Amemiya RE"),
type = "text") %>% kbl() %>% kable_paper()
##
## =======================================================================================================================
## Dependent variable:
## ----------------------------------------------------------------------------------------------------------
## log(budget_pc)
## Time FE Individual FE Twoway FE Swamy RE Walhus RE Nerlove RE Amemiya RE
## (1) (2) (3) (4) (5) (6) (7)
## -----------------------------------------------------------------------------------------------------------------------
## rule_of_law 0.967*** 0.148 0.390 0.658*** 0.659*** 0.609*** 0.649***
## (0.104) (0.329) (0.243) (0.205) (0.204) (0.216) (0.207)
##
## Constant -2.342*** -2.342*** -2.314*** -2.337***
## (0.138) (0.138) (0.148) (0.140)
##
## -----------------------------------------------------------------------------------------------------------------------
## Observations 1,118 1,118 1,118 1,118 1,118 1,118 1,118
## R2 0.073 0.0002 0.003 0.010 0.010 0.008 0.009
## Adjusted R2 0.067 -0.167 -0.171 0.009 0.009 0.007 0.008
## F Statistic 86.895*** (df = 1; 1110) 0.204 (df = 1; 957) 2.580 (df = 1; 951) 10.345*** 10.396*** 7.922*** 9.845***
## =======================================================================================================================
## Note: *p<0.1; **p<0.05; ***p<0.01
| x |
|---|
| ======================================================================================================================= |
| Dependent variable: |
| ———————————————————————————————————- |
| log(budget_pc) |
| Time FE Individual FE Twoway FE Swamy RE Walhus RE Nerlove RE Amemiya RE |
|
| ———————————————————————————————————————– |
| rule_of_law 0.967*** 0.148 0.390 0.658*** 0.659*** 0.609*** 0.649*** |
| (0.104) (0.329) (0.243) (0.205) (0.204) (0.216) (0.207) |
| Constant -2.342*** -2.342*** -2.314*** -2.337*** |
| (0.138) (0.138) (0.148) (0.140) |
| ———————————————————————————————————————– |
| Observations 1,118 1,118 1,118 1,118 1,118 1,118 1,118 |
| R2 0.073 0.0002 0.003 0.010 0.010 0.008 0.009 |
| Adjusted R2 0.067 -0.167 -0.171 0.009 0.009 0.007 0.008 |
| F Statistic 86.895*** (df = 1; 1110) 0.204 (df = 1; 957) 2.580 (df = 1; 951) 10.345*** 10.396*** 7.922*** 9.845*** |
| ======================================================================================================================= |
| Note: p<0.1; p<0.05; p<0.01 |
Is there a correlation between PD budgets per capita and level of financial globalization?
Is there a correlation between PD budgets per capita and manufacturing exports from the US per capita?
We can break it down by region
pd_df %>%
mutate(budget_pc = budget / pop) %>%
mutate(energy_imp_pc = energy_imports / pop) %>%
mutate(agri_imp_pc = agri_imports / pop) %>%
mutate(services_imp_pc = services_imports / pop) %>%
mutate(manufac_imp_pc = manufac_imports / pop) %>%
mutate(energy_exp_pc = energy_exports / pop) %>%
mutate(agri_exp_pc = agri_exports / pop) %>%
mutate(services_exp_pc = services_exports / pop) %>%
mutate(manufac_exp_pc = manufac_exports / pop) %>%
filter(year == 2016) %>%
ggplot(aes(x = log(budget_pc),
y = log(cultural_kof_dj ))) +
geom_point(aes(color = bureau)) +
facet_wrap(~bureau)