Synth-report.Rmd

---
title: "Synth reporting"
author: "Andrew Baxter"
date: "10/09/2019"
output: 
  html_document:
    toc: true
    toc_depth: 2
    number_sections: true
    toc_float: 
      collapsed: false
  # word_document:
  #   fig_height: 6
  #   fig_width: 8
  #   reference_docx: report_template.docx
always_allow_html: yes
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = FALSE, 
                      message = FALSE, 
                      warning = FALSE,
                      error = FALSE,
                      fig.width = 6,
                      fig.height = 5,
                      dpi = 150,
                      out.width = "100%")
library(readxl)
library(tidyverse)
library(stringr)
library(broom)
library(magrittr)
library(Synth)
library(Rcpp)
library(purrr)
library(svglite)
library(SPHSUgraphs)
library(plotly)
library(knitr)
`-.gg` <- function(e1, e2) e2(e1)
```
```{r load_functions, include=FALSE}
load('Data/synth_data.rdata') # outputted from 'Synth_data.R'
load('Data/special_preds.rdata')  # outputted from 'Synth_create_sps.R'
load("Data/placebos_country.rdata")  # outputted from 'Synth_create_sps.R'
load("Data/time_placebos.rdata")  # outputted from 'Synth_time_pbs.R'
source('R/Synth_functions.R')
```
```{r results='asis'}
cat("
<style>
caption {
      color: black;
      font-weight: bold;
      font-size: 1.0em;
    }
</style>
")

```

# Setting up data

Data from Human Mortality Database and Human Fertility Database was used to calculate birth rates. This was tidied to give Under-18 rates:

```{r tidy_u18}
kable(head(synthData_u18 %>% mutate_if(is.numeric, ~ round(.,2))))
```

Abortion estimates were added to give Under-20 rates:

```{r tidy_u_19}
kable(head(synthData_u20 %>% mutate_if(is.numeric, ~ round(.,2))))
```

# Iterating through year combinations

For each comparison, I iterated through all combinations of years as special predictors to minimise MSPE (whilst prioritising fewest groupings). For example, for the under-18 basic model with years as special predictors:

```{r example_iterations, eval = FALSE, echo = TRUE}
it_u18_rateSp <- testSynthIterations(
  yrs = 1985:1998,
  pred = "rate",
  data = synthData_u18[,1:4],
  ccodes = u_18_ccodes,
  n = 4,
  predictors = NULL,
  time.optimise = 1985:1998
) %>%
  arrange(groups, mspe)
```

# Generating Synthetic Control models for all Under-18 comparisons

## Model 1: Rate only as predictor

### England vs Synthetic Control
```{r md_u18_rateSp, include = FALSE}
opt_y1 <- 1985
synthPrep(synthData_u18[,1:4],
          grp = "u18_rateSp",
          dependent = "rate",
          special.predictors = sp_u18_rateSp,
          time.optimise.ssr = opt_y1:1998)
```

```{r md_u18_rateSp_graph}
gr_u18_rateSp <- gg_synth(md = md_u18_rateSp, post = TRUE) 


gr_u18_rateSp + labs(title = sPredText(dp_u18_rateSp))

#  ggsave("graphs/Synth-u18-sp.svg", width = 200, height = 150, units = "mm")

```

### Weights and balance

```{r st_u18_rateSp_w}
st_u18_rateSp$tab.w %>% arrange(desc(w.weights)) %>% select(Country = unit.names, Weight = w.weights) %>% 
  kable(caption = "Country weights")

st_u18_rateSp$tab.v %>% kable(caption = "Predictor weights")

st_u18_rateSp$tab.pred %>% kable(caption = "Predictor balance between synthetic and treated units")
```


### Placebo testing by country and time

```{r plot_u18_pb_nat}
pb_plot <- 
  gg_gaps(md_u18_rateSp, pl_u18_rateSp) + 
  # theme(text = element_text(size = 12),
  #       title = element_text(size = 12)) +
  labs(title = "England vs Synthetic Control, no added predictors - country placebo")

ggplotly(pb_plot)  # %>% 
#  api_create(filename = "Plot 43")
```
```{r pp_u18_sp}

pp_u18_sp <- gg_pre_postMSPE(md_u18_rateSp, pl_u18_rateSp) + labs(title = "Post/pre-TPS mean squared prediction error")
pp_u18_sp

```

```{r save_u18_pb_nat, include=FALSE}
#  ggsave(plot = pp_u18_sp, "graphs/Synth-u18-sp-pp-ratio.svg", width = 200, height = 150, units = "mm")
#  ggsave(plot = pb_plot, "graphs/Synth-u18-sp-gaps.svg", width = 200, height = 150, units = "mm")
```

```{r u18_pb_time, eval=FALSE, include=FALSE}
super_md_u18_sp <- tibble()

for (y in 1995:2005){
    sp <-   purrr::map(sp_u18_rateSp, 
    function(x) list(var = 
                       ifelse(min(x$yrs)<y,
                       x[[1]],
                       NA),
                     yrs = 
                       x$yrs[x$yrs<y],
                     op = 
                       ifelse(min(x$yrs)<y,
                       x[[3]],
                       NA))) %>% 
  purrr::modify_if(anyNA, function(x) NULL)

dp <- dataprep(
  foo = synthData_u18[,1:4],
  special.predictors = sp,
  predictors.op = "mean",
  time.predictors.prior = opt_y1:y,
  dependent = "rate",
  unit.variable = "Code",
  unit.names.variable = "Country",
  time.variable = "Year",
  treatment.identifier = u_18_ccodes$Code[u_18_ccodes$Country =="England and Wales"],
  controls.identifier = u_18_ccodes$Code[u_18_ccodes$Country !="England and Wales"],
  time.optimize.ssr = 1985:y,
  time.plot = 1985:2013
)
  md <- predvalues_synth(dp, synth_outputs = FALSE, yr = y)
  md <- md %>% 
  mutate(IntYr = y,
         mspe = pre_MSPE(md))
  super_md_u18_sp <- bind_rows(super_md_u18_sp, md)
}
```

```{r u18_pb_time_graph}
tb_u18_sp <- super_md_u18_sp %>% 
  ggplot(aes(Year, Rate, col = Group, linetype = Group, frame = IntYr)) +
  geom_segment(aes(x = IntYr-0.5, xend = IntYr-0.5, y = 0, yend = 20), linetype = "dotted", size = 1, col = "grey", alpha = 0.7) +
  geom_line(size = 1.5) +
  theme_sphsu_light() +
  ylab(paste0("Under-18 birth rate (per 1,000 women)")) +
  theme(legend.title = element_blank(),
        panel.grid = element_blank(),
        line = element_blank()) +
  scale_linetype_manual(name = "Data", values = c("Synthetic" = "dashed", "Treated" = "solid")) +
  scale_colour_manual(name = "Data", values = c("Synthetic" = sphsu_cols("Turquoise", names = FALSE), "Treated" = sphsu_cols("Thistle", names = FALSE))) +
  labs(title = "England vs Synthetic Control, no added predictors - time placebo") +
  scale_y_continuous(limits = c(0, 20))

# tb_u18_sp + facet_wrap(super_md_u18_sp$IntYr, ncol = 3)

tp_u18_sp_plotly <- ggplotly(tb_u18_sp) %>% 
  animation_slider(currentvalue = list(prefix = "Year ", font = list(color="black"))) %>% 
  animation_opts(easing = 'sin-in-out')

tp_u18_sp_plotly
#  api_create(tp_u18_sp_plotly, filename = "England under-18 - special predictors - placebo time")

```

## Model 2: GDP as predictor

### England vs Synthetic Control
```{r u18_gdp, include = FALSE}

# dataprep and synth
# opt_y1 <- 1994
# synthPrep(synthData_u18[,1:5] %>% filter(!Country %in% exclude_u18_gdp),
#           "u18_gdp",
#           dependent = "rate",
#           predictors = "GDPperCap",
#           special.predictors = sp_u18_gdp,
#           time.optimise.ssr = opt_y1:1998)

# dataprep and synth - 1990 onwards method
opt_y1 <- 1990
synthPrep(synthData_u18[,1:5] %>% filter(!Country %in% exclude_u18_gdp_1990),
          "u18_gdp",
          dependent = "rate",
          predictors = "GDPperCap",
          special.predictors = sp_u18_sp_1990,
          # special.predictors = sp_u18_gdp_1990,
          time.optimise.ssr = opt_y1:1998,
          time.predictors.prior = 1990:1998,
          time.plot = 1990:2013)
```
```{r md_u18_gdp_graph}
# graph
gr_u18_gdp <- gg_synth(md = md_u18_gdp, post = TRUE) 


gr_u18_gdp + labs(title = sPredText(dp_u18_gdp))


#  ggsave("graphs/Synth-u18-gdp-1990.svg", width = 200, height = 150, units = "mm")
```

### Weights and balance
```{r st_u18_gdp_w}
# weights kable
st_u18_gdp$tab.w %>% arrange(desc(w.weights)) %>% select(Country = unit.names, Weight = w.weights) %>% 
  kable(caption = "Country weights")

st_u18_gdp$tab.v %>% kable(caption = "Predictor weights")

st_u18_gdp$tab.pred %>% kable(caption = "Predictor balance between synthetic and treated units")
```

### Placebo testing by country and time
```{r plot_u18_gdp_pb_nat}
# gaps graph
pb_plot <- gg_gaps(md_u18_gdp, pl_u18_gdp %>% filter(Year>1989)) +
  labs(title = "England vs Synthetic Control, GDP as predictor - country placebo")

ggplotly(pb_plot)  # %>% 
#  api_create(filename = "Plot 45")
```
```{r pp_u18_gdp}

pp_u18_gdp <- gg_pre_postMSPE(md_u18_gdp, pl_u18_gdp) + labs(title = "Post/pre-TPS mean squared prediction error")
pp_u18_gdp

```
```{r save_u18_gdp_pb_nat, include=FALSE}
#  ggsave(plot = pp_u18_gdp, "graphs/Synth-u18-gdp-pp-ratio.svg", width = 200, height = 150, units = "mm")
#  ggsave(plot = pb_plot, "graphs/Synth-u18-gdp-gaps.svg", width = 200, height = 150, units = "mm")
```

```{r u18_gdp_pb_time, eval=FALSE, include = FALSE}
# placebo construct
super_md_u18_gdp <- tibble()

for (y in 1995:2005){
    sp <-   purrr::map(sp_u18_sp_1990, 
    function(x) list(var = 
                       ifelse(min(x$yrs)<y,
                       x[[1]],
                       NA),
                     yrs = 
                       x$yrs[x$yrs<y],
                     op = 
                       ifelse(min(x$yrs)<y,
                       x[[3]],
                       NA))) %>% 
  purrr::modify_if(anyNA, function(x) NULL)

dp <- dataprep(
  foo = synthData_u18[,1:5] %>% filter(!Country %in% exclude_u18_gdp_1990),
  special.predictors = sp,
  predictors = "GDPperCap",
  predictors.op = "mean",
  time.predictors.prior = opt_y1:y,
  dependent = "rate",
  unit.variable = "Code",
  unit.names.variable = "Country",
  time.variable = "Year",
  treatment.identifier = u_18_ccodes$Code[u_18_ccodes$Country =="England and Wales"],
  controls.identifier = u_18_ccodes %>% filter(!Country %in% exclude_u18_gdp_1990, Country != "England and Wales") %>% pull(Code),
  time.optimize.ssr = opt_y1:y,
  time.plot = opt_y1:2013
)
  md <- predvalues_synth(dp, synth_outputs = FALSE, yr = y)
  md <- md %>% 
  mutate(IntYr = y,
         mspe = pre_MSPE(md))
  super_md_u18_gdp <- bind_rows(super_md_u18_gdp, md)
}
```

```{r u18_gdp_pb_time_graph}
# placebo graph
tp_u18_gdp_plotly <- super_md_u18_gdp %>% 
  ggplot(aes(Year, Rate, col = Group, linetype = Group, frame = IntYr)) +
  geom_segment(aes(x = IntYr-0.5, xend = IntYr-0.5, y = 0, yend = 20), linetype = "dotted", size = 1, col = "grey", alpha = 0.7) +
  geom_line(size = 1.5) +
  theme_sphsu_light() +
  ylab(paste0("Under-18 birth rate (per 1,000 women)")) +
  theme(legend.title = element_blank(),
        panel.grid = element_blank(),
        line = element_blank()) +
  scale_linetype_manual(name = "Data", values = c("Synthetic" = "dashed", "Treated" = "solid")) +
  scale_colour_manual(name = "Data", values = c("Synthetic" = sphsu_cols("Turquoise", names = FALSE), "Treated" = sphsu_cols("Thistle", names = FALSE))) +
    labs(title = "England vs Synthetic Control, GDP per capita as predictor - time placebo") +
  scale_y_continuous(limits = c(0, 20)) - ggplotly

tp_u18_gdp_plotly <- tp_u18_gdp_plotly %>% 
  animation_slider(currentvalue = list(prefix = "Year ", font = list(color="black"))) %>% 
  animation_opts(easing = 'sin-in-out')

tp_u18_gdp_plotly
#  api_create(tp_u18_gdp_plotly, filename = "England under-18 - GDP per capita as predictor - placebo time")

```

## Model 3: All predictors

### England vs Synthetic Control
```{r u18_all, include = FALSE}


# dataprep and synth
# opt_y1 <- 1994
# synthPrep(synthData_u18 %>% filter(!Country %in% exclude_u18_all),
#           "u18_all",
#           dependent = "rate",
#           predictors = c("GDPperCap", "MobilePhones", "UrbanPop", "MF_ratio"),
#           special.predictors = sp_u18_all,
#           time.optimise.ssr = opt_y1:1998)


# dataprep and synth - excluding pre-1990 data
opt_y1 <- 1990
synthPrep(synthData_u18 %>% filter(!Country %in% exclude_u18_all_1990),
          "u18_all",
          dependent = "rate",
          predictors = c("GDPperCap", "MobilePhones", "UrbanPop", "MF_ratio"),
          special.predictors = sp_u18_sp_1990,
          time.optimise.ssr = opt_y1:1998,
          time.predictors.prior = opt_y1:1998,
          time.plot = opt_y1:2013)
```
```{r md_u18_all_graph}
# graph
gr_u18_all <- gg_synth(md = md_u18_all, post = TRUE) 


gr_u18_all + labs(title = sPredText(dp_u18_all))


#  ggsave("graphs/Synth-u18-all-1990.svg", width = 200, height = 150, units = "mm")
```

### Weights and balance
```{r st_u18_all_w}
# weights kable
st_u18_all$tab.w %>% arrange(desc(w.weights)) %>% select(Country = unit.names, Weight = w.weights) %>% 
  kable(caption = "Country weights")

st_u18_all$tab.v %>% kable(caption = "Predictor weights")

st_u18_all$tab.pred %>% kable(caption = "Predictor balance between synthetic and treated units")
```

### Placebo testing by country and time
```{r plot_u18_all_pb_nat, include = FALSE}
# gaps graph
pb_plot <- 
gg_gaps(md_u18_all, pl_u18_all %>% filter(Year>1989)) + 
  labs(title = "England vs Synthetic Control, all predictors - country placebo")

ggplotly(pb_plot)  # %>% 
#  api_create(filename = "Plot 47")
```
```{r pp_u18_all, include = FALSE}

pp_u18_all <- gg_pre_postMSPE(md_u18_all, pl_u18_all) + labs(title = "Post/pre-TPS mean squared prediction error")
pp_u18_all

```
```{r save_u18_all_pb_nat, include=FALSE}
#  ggsave(plot = pp_u18_all, "graphs/Synth-u18-all-pp-ratio.svg", width = 200, height = 150, units = "mm")
#  ggsave(plot = pb_plot,"graphs/Synth-u18-all-gaps.svg", width = 200, height = 150, units = "mm")
```

```{r u18_all_pb_time, eval=FALSE, include = FALSE}
# placebo construct
super_md_u18_all <- tibble()

for (y in 1995:2005){
  
  sp <-   purrr::map(sp_u18_sp_1990, 
    function(x) list(var = 
                       ifelse(min(x$yrs)<y,
                       x[[1]],
                       NA),
                     yrs = 
                       x$yrs[x$yrs<y],
                     op = 
                       ifelse(min(x$yrs)<y,
                       x[[3]],
                       NA))) %>% 
  purrr::modify_if(anyNA, function(x) NULL)

dp <- dataprep(
  foo = synthData_u18 %>% filter(Year>1989, !Country %in% exclude_u18_all_1990),
  special.predictors = sp,
  predictors = c("GDPperCap", "MobilePhones", "UrbanPop", "MF_ratio"),
  predictors.op = "mean",
  time.predictors.prior = opt_y1:y,
  dependent = "rate",
  unit.variable = "Code",
  unit.names.variable = "Country",
  time.variable = "Year",
  treatment.identifier = u_18_ccodes$Code[u_18_ccodes$Country =="England and Wales"],
  controls.identifier = u_18_ccodes %>% filter(!Country %in% exclude_u18_all_1990, Country != "England and Wales") %>% pull(Code),
  time.optimize.ssr = opt_y1:y,
  time.plot = opt_y1:2013
)
  md <- predvalues_synth(dp, synth_outputs = FALSE, yr = y)
  md <- md %>% 
  mutate(IntYr = y,
         mspe = pre_MSPE(md))
  super_md_u18_all <- bind_rows(super_md_u18_all, md)
}
```

```{r u18_all_pb_time_graph, include = FALSE}
# placebo graph
tp_u18_all_plotly <- super_md_u18_all %>% 
  ggplot(aes(Year, Rate, col = Group, linetype = Group, frame = IntYr)) +
  geom_segment(aes(x = IntYr-0.5, xend = IntYr-0.5, y = 0, yend = 20), linetype = "dotted", size = 1, col = "grey", alpha = 0.7) +
  geom_line(size = 1.5) +
  theme_sphsu_light() +
  ylab(paste0("Under-18 birth rate (per 1,000 women)")) +
  theme(legend.title = element_blank(),
        panel.grid = element_blank(),
        line = element_blank()) +
  scale_linetype_manual(name = "Data", values = c("Synthetic" = "dashed", "Treated" = "solid")) +
  scale_colour_manual(name = "Data", values = c("Synthetic" = sphsu_cols("Turquoise", names = FALSE), "Treated" = sphsu_cols("Thistle", names = FALSE))) +
    labs(title = "England vs Synthetic Control, all predictors - time placebo") +
  scale_y_continuous(limits = c(0, 20)) - ggplotly

tp_u18_all_plotly <- tp_u18_all_plotly %>% 
  animation_slider(currentvalue = list(prefix = "Year ", font = list(color="black"))) %>% 
  animation_opts(easing = 'sin-in-out')
  
  tp_u18_all_plotly
  #  api_create(tp_u18_all_plotly, filename = "England under-18 - all predictors - placebo time")
```

# Generating Synthetic Control models for all Under-20 comparisons

## Model 4: England vs Synthetic Control

### England vs Synthetic Control
```{r u20_sp, include = FALSE}
opt_y1 <- 1990
# dataprep and synth
synthPrep(synthData_u20[,1:4],
          "u20_sp",
          dependent = "pRate",
          special.predictors = sp_u20_sp,
          time.optimise.ssr = opt_y1:1998)
```
```{r md_u20_sp_graph}
# graph
gr_u20_sp <- gg_synth(md = md_u20_sp, post = TRUE) 


gr_u20_sp + labs(title = sPredText(dp_u20_sp))
#  ggsave("graphs/Synth-u20-sp.svg", width = 200, height = 150, units = "mm")
```

### Weights and balance
```{r st_u20_sp_w}
# weights kable
st_u20_sp$tab.w %>% arrange(desc(w.weights)) %>% select(Country = unit.names, Weight = w.weights) %>% 
  kable(caption = "Country weights")

st_u20_sp$tab.v %>% kable(caption = "Predictor weights")

st_u20_sp$tab.pred %>% kable(caption = "Predictor balance between synthetic and treated units")
```

### Placebo testing by country and time
```{r plot_u20_sp_pb_nat, include = FALSE}
# gaps graph
pb_plot <- 
gg_gaps(md_u20_sp, pl_u20_sp) + 
  labs(title = "England vs Synthetic Control, no added predictors - country placebo")

ggplotly(pb_plot)  #%>% 
#  api_create(filename = "Plot 49")
```
```{r pp_u20_sp, include = FALSE}

pp_u20_sp <- gg_pre_postMSPE(md_u20_sp, pl_u20_sp) + labs(title = "Post/pre-TPS mean squared prediction error")
pp_u20_sp

```
```{r save_u20_sp_pb_nat, include=FALSE}
#  ggsave(plot = pp_u20_sp, "graphs/Synth-u20-sp-pp-ratio.svg", width = 200, height = 150, units = "mm")
#  ggsave(plot = pb_plot, "graphs/Synth-u20-sp-gaps.svg", width = 200, height = 150, units = "mm")
```

```{r u20_sp_pb_time, eval=FALSE, include = FALSE}
# placebo construct
super_md_u20_sp <- tibble()

for (y in 1995:2005){
    sp <-   purrr::map(sp_u20_sp, 
    function(x) list(var = 
                       ifelse(min(x$yrs)<y,
                       x[[1]],
                       NA),
                     yrs = 
                       x$yrs[x$yrs<y],
                     op = 
                       ifelse(min(x$yrs)<y,
                       x[[3]],
                       NA))) %>% 
  purrr::modify_if(anyNA, function(x) NULL)

dp <- dataprep(
  foo = synthData_u20[,1:4],
  special.predictors = sp,
  predictors.op = "mean",
  time.predictors.prior = 1990:y,
  dependent = "pRate",
  unit.variable = "Code",
  unit.names.variable = "Country",
  time.variable = "Year",
  treatment.identifier = u_20_ccodes$Code[u_20_ccodes$Country =="England and Wales"],
  controls.identifier = u_20_ccodes %>% filter(Country != "England and Wales") %>% pull(Code),
  time.optimize.ssr = opt_y1:y,
  time.plot = 1990:2013
)
  md <- predvalues_synth(dp, synth_outputs = FALSE, yr = y)
  md <- md %>% 
  mutate(IntYr = y,
         mspe = pre_MSPE(md))
  super_md_u20_sp <- bind_rows(super_md_u20_sp, md)
}
```

```{r u20_sp_pb_time_graph, include = FALSE}
# placebo graph
tp_u20_sp_plotly <- super_md_u20_sp %>% 
  ggplot(aes(Year, Rate, col = Group, linetype = Group, frame = IntYr)) +
  geom_segment(aes(x = IntYr-0.5, xend = IntYr-0.5, y = 0, yend = 80), linetype = "dotted", size = 1, col = "grey", alpha = 0.7) +
  geom_line(size = 1.5) +
  theme_sphsu_light() +
  ylab(paste0("Under-20 pregnancy rate (per 1,000 women)")) +
  theme(legend.title = element_blank(),
        panel.grid = element_blank(),
        line = element_blank()) +
  scale_linetype_manual(name = "Data", values = c("Synthetic" = "dashed", "Treated" = "solid")) +
  scale_colour_manual(name = "Data", values = c("Synthetic" = sphsu_cols("Turquoise", names = FALSE), "Treated" = sphsu_cols("Thistle", names = FALSE))) +
    labs(title = "England vs Synthetic Control, no added predictors - time placebo") +
  scale_y_continuous(limits = c(0, 80)) - ggplotly

tp_u20_sp_plotly <- tp_u20_sp_plotly %>% 
  animation_slider(currentvalue = list(prefix = "Year ", font = list(color="black"))) %>% 
  animation_opts(easing = 'sin-in-out')
  
  tp_u20_sp_plotly
  
  #  api_create(tp_u20_sp_plotly, filename = "England under-20 - special predictors - placebo time")
```

## Model 5: GDP as predictor

### England vs Synthetic Control
```{r u20_gdp, include = FALSE}
opt_y1 <- 1990
# dataprep and synth
synthPrep(synthData_u20[,c(1:4, 6)] %>% filter(!Country %in% exclude_u20_gdp),
          "u20_gdp",
          dependent = "pRate",
          predictors = "GDPperCap",
          special.predictors = sp_u20_gdp, 
          time.optimise.ssr = opt_y1:1998)
```
```{r md_u20_gdp_graph}
# graph
gr_u20_gdp <- gg_synth(md = md_u20_gdp, post = TRUE) 


gr_u20_gdp + labs(title = sPredText(dp_u20_gdp))
#  ggsave("graphs/Synth-u20-gdp.svg", width = 200, height = 150, units = "mm")
```

### Weights and balance
```{r st_u20_gdp_w}
# weights kable
st_u20_gdp$tab.w %>% arrange(desc(w.weights)) %>% select(Country = unit.names, Weight = w.weights) %>% 
  kable(caption = "Country weights")

st_u20_gdp$tab.v %>% kable(caption = "Predictor weights")

st_u20_gdp$tab.pred %>% kable(caption = "Predictor balance between synthetic and treated units")
```

### Placebo testing by country and time
```{r plot_u20_gdp_pb_nat, include = FALSE}
# gaps graph
pb_plot <- 
gg_gaps(md_u20_gdp, pl_u20_gdp) +
  labs(title = "England vs Synthetic Control, GDP as predictor - country placebo")

ggplotly(pb_plot)  # %>% 
#  api_create(filename = "Plot 51")
```
```{r pp_u20_gdp, include = FALSE}

pp_u20_gdp <- gg_pre_postMSPE(md_u20_gdp, pl_u20_gdp) + labs(title = "Post/pre-TPS mean squared prediction error")
pp_u20_gdp

```
```{r save_u20_gdp_pb_nat, include=FALSE}
#  ggsave(plot = pp_u20_gdp, "graphs/Synth-u20-gdp-pp-ratio.svg", width = 200, height = 150, units = "mm")
#  ggsave(plot = pb_plot, "graphs/Synth-u20-gdp-gaps.svg", width = 200, height = 150, units = "mm")
```

```{r u20_gdp_pb_time, eval=FALSE, include = FALSE}
# placebo construct
super_md_u20_gdp <- tibble()

for (y in 1995:2005){
  
sp <-   purrr::map(sp_u20_gdp, 
    function(x) list(var = 
                       ifelse(min(x$yrs)<y,
                       x[[1]],
                       NA),
                     yrs = 
                       x$yrs[x$yrs<y],
                     op = 
                       ifelse(min(x$yrs)<y,
                       x[[3]],
                       NA))) %>% 
  purrr::modify_if(anyNA, function(x) NULL)

dp <- dataprep(
  foo = synthData_u20[,c(1:4, 6)] %>% filter(!Country %in% exclude_u20_gdp),
  special.predictors = sp,
  predictors = "GDPperCap",
  predictors.op = "mean",
  time.predictors.prior = 1990:y,
  dependent = "pRate",
  unit.variable = "Code",
  unit.names.variable = "Country",
  time.variable = "Year",
  treatment.identifier = u_20_ccodes$Code[u_20_ccodes$Country =="England and Wales"],
  controls.identifier = u_20_ccodes %>% filter(!Country %in% exclude_u20_gdp, Country != "England and Wales") %>% pull(Code),
  time.optimize.ssr = opt_y1:y,
  time.plot = 1990:2013
)
  md <- predvalues_synth(dp, synth_outputs = FALSE, yr = y)
  md <- md %>% 
  mutate(IntYr = y,
         mspe = pre_MSPE(md))
  super_md_u20_gdp <- bind_rows(super_md_u20_gdp, md)
}
```

```{r u20_gdp_pb_time_graph, include = FALSE}
# placebo graph
tp_u20_gdp_plotly <- super_md_u20_gdp %>% 
  ggplot(aes(Year, Rate, col = Group, linetype = Group, frame = IntYr)) +
  geom_segment(aes(x = IntYr-0.5, xend = IntYr-0.5, y = 0, yend = 80), linetype = "dotted", size = 1, col = "grey", alpha = 0.7) +
  geom_line(size = 1.5) +
  theme_sphsu_light() +
  ylab(paste0("Under-20 pregnancy rate (per 1,000 women)")) +
  theme(legend.title = element_blank(),
        panel.grid = element_blank(),
        line = element_blank()) +
  scale_linetype_manual(name = "Data", values = c("Synthetic" = "dashed", "Treated" = "solid")) +
  scale_colour_manual(name = "Data", values = c("Synthetic" = sphsu_cols("Turquoise", names = FALSE), "Treated" = sphsu_cols("Thistle", names = FALSE))) +
    labs(title = "England vs Synthetic Control, GDP per capita predictor - time placebo") +
  scale_y_continuous(limits = c(0, 80)) - ggplotly

tp_u20_gdp_plotly <- tp_u20_gdp_plotly %>% 
  animation_slider(currentvalue = list(prefix = "Year ", font = list(color="black"))) %>% 
  animation_opts(easing = 'sin-in-out')
  
  tp_u20_gdp_plotly
  #  api_create(tp_u20_gdp_plotly, filename = "England under-20 - GDP per capita as predictor - placebo time")
```

## Model 6: All predictors

### England vs Synthetic Control
```{r u20_all, include = FALSE}
opt_y1 <- 1990
# dataprep and synth
synthPrep(synthData_u20[,-5] %>% filter(!Country %in% exclude_u20_all),
          "u20_all",
          dependent = "pRate",
          predictors = c("GDPperCap", "MobilePhones", "UrbanPop", "MF_ratio"),
          special.predictors = sp_u20_all, 
          time.optimise.ssr = opt_y1:1998)
```
```{r md_u20_all_graph}
# graph
gr_u20_all <- gg_synth(md = md_u20_all, post = TRUE) 


gr_u20_all + labs(title = sPredText(dp_u20_all))

#  ggsave("graphs/Synth-u20-all.svg", width = 200, height = 150, units = "mm")
```

### Weights and balance
```{r st_u20_all_w}
# weights kable
st_u20_all$tab.w %>% arrange(desc(w.weights)) %>% select(Country = unit.names, Weight = w.weights) %>% 
  kable(caption = "Country weights")

st_u20_all$tab.v %>% kable(caption = "Predictor weights")

st_u20_all$tab.pred %>% kable(caption = "Predictor balance between synthetic and treated units")
```

### Placebo testing by country and time
```{r plot_u20_all_pb_nat, include = FALSE}
# gaps graph
pb_plot <- 
gg_gaps(md_u20_all, pl_u20_all) + 
  labs(title = "England vs Synthetic Control, all predictors - country placebo")

ggplotly(pb_plot)  # %>% 
#  api_create(filename = "Plot 53")
```
```{r pp_u20_all, include = FALSE}

pp_u20_all <- gg_pre_postMSPE(md_u20_all, pl_u20_all) + labs(title = "Post/pre-TPS mean squared prediction error")
pp_u20_all

```
```{r save_u20_all_pb_nat, include=FALSE}
#  ggsave(plot = pp_u20_all, "graphs/Synth-u20-all-pp-ratio.svg", width = 200, height = 150, units = "mm")
#  ggsave(plot = pb_plot, "graphs/Synth-u20-all-gaps.svg", width = 200, height = 150, units = "mm")
```

```{r u20_all_pb_time, eval=FALSE, include = FALSE}
# placebo construct
super_md_u20_all <- tibble()

for (y in 1995:2005){
      sp <-   purrr::map(sp_u20_all, 
    function(x) list(var = 
                       ifelse(min(x$yrs)<y,
                       x[[1]],
                       NA),
                     yrs = 
                       x$yrs[x$yrs<y],
                     op = 
                       ifelse(min(x$yrs)<y,
                       x[[3]],
                       NA))) %>% 
  purrr::modify_if(anyNA, function(x) NULL)

dp <- dataprep(
  foo = synthData_u20[,-5] %>% filter(!Country %in% exclude_u20_all),
  special.predictors = sp,
  predictors = c("GDPperCap", "MobilePhones", "UrbanPop", "MF_ratio"),
  predictors.op = "mean",
  time.predictors.prior = 1990:y,
  dependent = "pRate",
  unit.variable = "Code",
  unit.names.variable = "Country",
  time.variable = "Year",
  treatment.identifier = u_20_ccodes$Code[u_20_ccodes$Country =="England and Wales"],
  controls.identifier = u_20_ccodes %>% filter(!Country %in% exclude_u20_all, Country != "England and Wales") %>% pull(Code),
  time.optimize.ssr = opt_y1:1998,
  time.plot = 1990:2013
)
  md <- predvalues_synth(dp, synth_outputs = FALSE, yr = y)
  md <- md %>% 
  mutate(IntYr = y,
         mspe = pre_MSPE(md))
  super_md_u20_all <- bind_rows(super_md_u20_all, md)
}
```

```{r u20_all_pb_time_graph, include = FALSE}
# placebo graph
tp_u20_all_plotly <- super_md_u20_all %>% 
  ggplot(aes(Year, Rate, col = Group, linetype = Group, frame = IntYr)) +
  geom_segment(aes(x = IntYr-0.5, xend = IntYr-0.5, y = 0, yend = 80), linetype = "dotted", size = 1, col = "grey", alpha = 0.7) +
  geom_line(size = 1.5) +
  theme_sphsu_light() +
  ylab(paste0("Under-20 pregnancy rate (per 1,000 women)")) +
  theme(legend.title = element_blank(),
        panel.grid = element_blank(),
        line = element_blank()) +
  scale_linetype_manual(name = "Data", values = c("Synthetic" = "dashed", "Treated" = "solid")) +
  scale_colour_manual(name = "Data", values = c("Synthetic" = sphsu_cols("Turquoise", names = FALSE), "Treated" = sphsu_cols("Thistle", names = FALSE))) +
    labs(title = "England vs Synthetic Control, all predictors - time placebo") +
  scale_y_continuous(limits = c(0, 80)) - ggplotly

tp_u20_all_plotly <- tp_u20_all_plotly %>% 
  animation_slider(currentvalue = list(prefix = "Year ", font = list(color="black"))) %>% 
  animation_opts(easing = 'sin-in-out')
  
  tp_u20_all_plotly
  #  api_create(tp_u20_all_plotly, filename = "England under-20 - all predictors - placebo time")
```