From fe6a535bf4fa5f535348bed024a0c7401c3d0cec Mon Sep 17 00:00:00 2001 From: juanitorduz Date: Sun, 7 Nov 2021 13:36:38 +0100 Subject: [PATCH 1/3] rename folder --- R/synthetic_control/rep_for_btsa.R | 474 +++++++++++++++++++++ R/synthetic_control/rep_original.r | 640 +++++++++++++++++++++++++++++ 2 files changed, 1114 insertions(+) create mode 100644 R/synthetic_control/rep_for_btsa.R create mode 100755 R/synthetic_control/rep_original.r diff --git a/R/synthetic_control/rep_for_btsa.R b/R/synthetic_control/rep_for_btsa.R new file mode 100644 index 0000000..745fc44 --- /dev/null +++ b/R/synthetic_control/rep_for_btsa.R @@ -0,0 +1,474 @@ +## Replication Code for +# A. Abadie, A. Diamond, and J. Hainmueller. 2014. +# Comparative Politics and the Synthetic Control Method +# American Journal of Political Science. +# Original source for the code and the data: +# https://dataverse.harvard.edu/dataset.xhtml?persistentId=doi:10.7910/DVN/24714 +# This version for the Berlin Time Series Analysis Meet-up + +library(foreign) +library("Synth") +library(xtable) +library(ggplot2) +library(scales) +library(tidyverse) +# Load Data +d <- read.dta("repgermany.dta") + + +# Data set-up + countrytosynthetise <- 7 + dataprep.out <- + dataprep( + foo = d, + predictors = c("gdp","trade","infrate"), + dependent = "gdp", + unit.variable = 1, + time.variable = 3, + special.predictors = list( + list("industry", 1971:1980, c("mean")), + list("schooling",c(1970,1975), c("mean")), + list("invest70" ,1980, c("mean")) + ), + treatment.identifier = countrytosynthetise, + controls.identifier = unique(d$index)[-which(unique(d$index) %in% countrytosynthetise)], + time.predictors.prior = 1971:1989, + time.optimize.ssr = 1971:1990, + unit.names.variable = 2, + time.plot = 1960:2003 + ) + +# fit training model +synth.out <- + synth( + data.prep.obj=dataprep.out, + Margin.ipop=.005,Sigf.ipop=7,Bound.ipop=6 + ) + +# data prep for main model +dataprep.out <- + dataprep( + foo = d, + predictors = c("gdp","trade","infrate"), + dependent = "gdp", + unit.variable = 1, + time.variable = 3, + special.predictors = list( + list("industry" ,1981:1990, c("mean")), + list("schooling",c(1980,1985), c("mean")), + list("invest80" ,1980, c("mean")) + ), + treatment.identifier = countrytosynthetise, + controls.identifier = unique(d$index)[-which(unique(d$index) %in% countrytosynthetise)], + time.predictors.prior = 1960:1989, + time.optimize.ssr = 1960:1990, + unit.names.variable = 2, + time.plot = 1960:2003 +) + +# fit main model with v from training model +synth.out <- synth( + data.prep.obj=dataprep.out, + custom.v=as.numeric(synth.out$solution.v) + ) + +#### Table 2 +synth.tables <- synth.tab( + dataprep.res = dataprep.out, + synth.res = synth.out + ); synth.tables + +# Replace means for OECD sample (computed externally using proper pop weighting) +synth.tables$tab.pred[,3] <- c(8021.1,31.9,7.4,34.2,44.1,25.9) +colnames(synth.tables$tab.pred)[3] <- "Rest of OECD Sample" +rownames(synth.tables$tab.pred) <- c("GDP per-capita","Trade openness", + "Inflation rate","Industry share", + "Schooling","Investment rate") + +xtable(round(synth.tables$tab.pred,1),digits=1) + +#### Table 1 +# synth weights +tab1 <- data.frame(synth.tables$tab.w) +tab1[,1] <- round(tab1[,1],2) +# regression weights +X0 <- cbind(1,t(dataprep.out$X0)) +X1 <- as.matrix(c(1,dataprep.out$X1)) +W <- X0%*%solve(t(X0)%*%X0)%*%X1 +Wdat <- data.frame(unit.numbers=as.numeric(rownames(X0)), + regression.w=round(W,2)) +tab1 <- merge(tab1,Wdat,by="unit.numbers") +tab1 <- tab1[order(tab1[,3]),] + + +# Table 1 +xtable(cbind(tab1[1:9,c(3,2,4)], + tab1[10:18,c(3,2,4)] + ) + ) + +#### Figures 1, 2 and 3 +data.to.plot <- data.frame(year = rownames(dataprep.out$Y1plot), + west.germany = dataprep.out$Y1plot[, 1]) + +data.to.plot <- as.data.frame(d[,c("gdp","country", "year" )]) +data.to.plot$west.germany <- ifelse(data.to.plot$country == "West Germany", "West Germany", "Not West Germany") +data.to.plot$west.germany.size <- ifelse(data.to.plot$west.germany == "West Germany", "1","2") +graphcolors <- c("#999999", "#ef8a62") + + +table(data.to.plot$country, data.to.plot$west.germany.size) + + +## PLOT WITH ONLY WEST GERMANY +ggplot(data = data.to.plot[data.to.plot$west.germany == "West Germany",], aes(x = year, y = gdp, + group = country, + color = west.germany, + linetype = west.germany.size)) + + geom_line(size = 1.5) + + labs(color = "West Germany?")+ + guides(linetype = FALSE, color = FALSE) + + theme_bw(base_size = 15) + + theme(legend.position = "bottom") + + geom_vline(xintercept = 1990, + color = "#88419d") + + annotate("text", x = 1985, + label="\nReunification", y = 1000, + colour = "#88419d", alpha = 0.8, angle = 0) + + scale_y_continuous(name="per-capita GDP (PPP, 2002 USD)", + labels = comma_format(big.mark = ",", + decimal.mark = ".")) + + scale_color_manual(values = c("#ef8a62")) + + coord_cartesian(ylim = c(0, 30000)) #+ + #ggsave(paste0("figures/germany.pdf"), + # width = 15, height = 10, units = 'cm') + + + +## PLOT WITH WEST GERMANY AND ALL THE OTHER COUNTRIES +ggplot(data = data.to.plot, aes(x = year, y = gdp, + group = country, + color = west.germany, + linetype = west.germany.size)) + + geom_line(size = 1.5) + + labs(color = "West Germany?")+ + guides(linetype = FALSE) + + theme_bw(base_size = 10) + + theme(legend.position = "bottom") + + geom_vline(xintercept = 1990, + color = "#88419d") + + annotate("text", x = 1985, + label="Reunification", y = 1000, + colour = "#88419d", alpha = 0.8, angle = 0) + + scale_y_continuous(name="per-capita GDP (PPP, 2002 USD)", + labels = comma_format(big.mark = ",", + decimal.mark = ".")) + + scale_color_manual(values = graphcolors) + + coord_cartesian(ylim = c(0, 30000)) #+ + #ggsave(paste0("figures/germanynotgermany.pdf"), + # width = 15, height = 10, units = 'cm') + + +## PLOT WITH WEST GERMANY AND THE AVERAGE FROM OTHER COUNTRIES +data.to.plot.average <- data.to.plot %>% + group_by(west.germany, year, west.germany.size) %>% + summarise(GDP = mean(gdp)) + +ggplot(data = data.to.plot.average, + aes(x = year, y = GDP, + group = west.germany, + color = west.germany)) + + geom_line(size = 1.5) + + labs(color = "West Germany?")+ + theme_bw(base_size = 10) + + theme(legend.position = "bottom") + + geom_vline(xintercept = 1990, + color = "#88419d") + + annotate("text", x = 1985, + label="Reunification", y = 1000, + colour = "#88419d", alpha = 0.8, angle = 0) + + scale_y_continuous(name="per-capita GDP (PPP, 2002 USD)", + labels = comma_format(big.mark = ",", + decimal.mark = ".")) + + scale_color_manual(values = graphcolors) + + coord_cartesian(ylim = c(0, 30000)) #+ + #ggsave(paste0("figures/germanyaverage.pdf"), + # width = 15, height = 10, units = 'cm') + + + +data.to.plot.average$west.germany[data.to.plot.average$west.germany == "Not West Germany"] <- "Average of Other OECD Countries" + +ggplot(data = data.to.plot.average, + aes(x = year, y = GDP, + group = west.germany, + color = west.germany)) + + geom_line(size = 1.5) + + labs(color = "West Germany?")+ + theme_bw(base_size = 10) + + theme(legend.position = "bottom") + + geom_vline(xintercept = 1990, + color = "#88419d") + + annotate("text", x = 1985, + label="Reunification", y = 1000, + colour = "#88419d", alpha = 0.8, angle = 0) + + scale_y_continuous(name="per-capita GDP (PPP, 2002 USD)", + labels = comma_format(big.mark = ",", + decimal.mark = ".")) + + scale_color_manual(values = graphcolors) + + coord_cartesian(ylim = c(0, 30000)) #+ + #ggsave(paste0("figures/germanyaverage2.pdf"), + # width = 15, height = 10, units = 'cm') + + + + +# SYNTHETIC UNIT - Weights +# Figure 4 +data.to.plot.synthetic <- data.to.plot.average + +synthY0 <- (dataprep.out$Y0%*%synth.out$solution.w) +data.to.plot.synthetic$GDP[data.to.plot.synthetic$west.germany == "Average of Other OECD Countries"] <- (dataprep.out$Y0%*%synth.out$solution.w) +data.to.plot.synthetic$west.germany[data.to.plot.synthetic$west.germany == "Average of Other OECD Countries"] <- "Synthetic West Germany" + + +ggplot(data = data.to.plot.synthetic, + aes(x = year, y = GDP, + group = west.germany, + color = west.germany)) + + geom_line(size = 1.5) + + labs(color = "West Germany?")+ + theme_bw(base_size = 10) + + theme(legend.position = "bottom") + + geom_vline(xintercept = 1990, + color = "#88419d") + + annotate("text", x = 1985, + label="Reunification", y = 1000, + colour = "#88419d", alpha = 0.8, angle = 0) + + scale_y_continuous(name="per-capita GDP (PPP, 2002 USD)", + labels = comma_format(big.mark = ",", + decimal.mark = ".")) + + scale_color_manual(values = graphcolors) + + coord_cartesian(ylim = c(0, 30000)) #+ + #ggsave(paste0("figures/germanysynth.pdf"), + # width = 15, height = 10, units = 'cm') + + + +# SYNTHETIC UNIT - Regression +# Figure 8 +data.to.plot.synth.reg <- data.to.plot.average + +synthY0.reg <- (dataprep.out$Y0%*%W) +data.to.plot.synth.reg$GDP[data.to.plot.synth.reg$west.germany == "Average of Other OECD Countries"] <- (dataprep.out$Y0%*%W) +data.to.plot.synth.reg$west.germany[data.to.plot.synth.reg$west.germany == "Average of Other OECD Countries"] <- "Synthetic West Germany" + + +ggplot(data = data.to.plot.synth.reg, + aes(x = year, y = GDP, + group = west.germany, + color = west.germany)) + + geom_line(size = 1.5) + + labs(color = "West Germany?")+ + theme_bw(base_size = 10) + + theme(legend.position = "bottom") + + geom_vline(xintercept = 1990, + color = "#88419d") + + annotate("text", x = 1985, + label="Reunification", y = 1000, + colour = "#88419d", alpha = 0.8, angle = 0) + + scale_y_continuous(name="per-capita GDP (PPP, 2002 USD)", + labels = comma_format(big.mark = ",", + decimal.mark = ".")) + + scale_color_manual(values = graphcolors) + + coord_cartesian(ylim = c(0, 30000)) #+ + #ggsave(paste0("figures/germanysynthreg.pdf"), + # width = 15, height = 10, units = 'cm') + + + +# Gap +# Figure 5 +gap <- dataprep.out$Y1-(dataprep.out$Y0%*%synth.out$solution.w) + +data.to.plot.gap <- data.to.plot.synthetic +data.to.plot.gap$GDP[data.to.plot.gap$west.germany == "West Germany"] <- dataprep.out$Y1-(dataprep.out$Y0%*%synth.out$solution.w) + +data.to.plot.gap$GDP[data.to.plot.gap$west.germany == "Synthetic West Germany"] <- 0 + + +ggplot(data = data.to.plot.gap, + aes(x = year, y = GDP, + group = west.germany, + color = west.germany)) + + geom_line(size = 1.5) + + labs(color = "West Germany?")+ + theme_bw(base_size = 10) + + theme(legend.position = "bottom") + + geom_vline(xintercept = 1990, + color = "#88419d") + + annotate("text", x = 1985, + label="Reunification", y = 1000, + colour = "#88419d", alpha = 0.8, angle = 0) + + scale_y_continuous(name="Gap in per-capita GDP (PPP, 2002 USD)", + labels = comma_format(big.mark = ",", + decimal.mark = ".")) + + scale_color_manual(values = graphcolors) + + coord_cartesian(ylim = c(-5000, 5000)) #+ + #ggsave(paste0("figures/germanygap.pdf"), + # width = 15, height = 10, units = 'cm') + + + + +# loop across control units +storegaps <- + matrix(NA, + length(1960:2003), + length(unique(d$index))-1 + ) +rownames(storegaps) <- 1960:2003 +i <- 1 +co <- unique(d$index) + +for(k in unique(d$index)[-7]){ + + # data prep for training model + dataprep.out <- + dataprep( + foo = d, + predictors = c("gdp","trade","infrate"), + dependent = "gdp", + unit.variable = 1, + time.variable = 3, + special.predictors = list( + list("industry",1971:1980, c("mean")), + list("schooling" ,c(1970,1975), c("mean")), + list("invest70" ,1980, c("mean")) + ), + treatment.identifier = k, + controls.identifier = co[-which(co==k)], + time.predictors.prior = 1971:1989, + time.optimize.ssr = 1971:1990, + unit.names.variable = 2, + time.plot = 1960:2003 + ) + + # fit training model + synth.out <- + synth( + data.prep.obj=dataprep.out, + Margin.ipop=.005,Sigf.ipop=7,Bound.ipop=6 + ) + + # data prep for main model +dataprep.out <- + dataprep( + foo = d, + predictors = c("gdp","trade","infrate"), + dependent = "gdp", + unit.variable = 1, + time.variable = 3, + special.predictors = list( + list("industry" ,1981:1990, c("mean")), + list("schooling",c(1980,1985), c("mean")), + list("invest80" ,1980, c("mean")) + ), + treatment.identifier = k, + controls.identifier = co[-which(co==k)], + time.predictors.prior = 1960:1990, + time.optimize.ssr = 1960:1989, + unit.names.variable = 2, + time.plot = 1960:2003 + ) + +# fit main model +synth.out <- synth( + data.prep.obj=dataprep.out, + custom.v=as.numeric(synth.out$solution.v) + ) + + storegaps[,i] <- + dataprep.out$Y1- + (dataprep.out$Y0%*%synth.out$solution.w) + i <- i + 1 +} # close loop over control units +d <- d[order(d$index,d$year),] +colnames(storegaps) <- unique(d$country)[-7] +storegaps <- cbind(gap,storegaps) +colnames(storegaps)[1] <- c("West Germany") + + + +# Placebo gaps +# Figure 6 +sgaps <- as.data.frame(storegaps) +sgaps$year <- rownames(sgaps) +sgaps <- reshape2::melt(sgaps,id.vars = "year") + + +sgaps$west.germany <- ifelse(sgaps$variable == "West Germany", "West Germany", "Not West Germany") +sgaps$west.germany.size <- ifelse(sgaps$west.germany == "West Germany", "1","2") + +sgaps$year <- as.numeric(sgaps$year) + +data.to.plot$west.germany <- ifelse(data.to.plot$country == "West Germany", "West Germany", "Not West Germany") +ggplot(data = sgaps, aes(x = year, y = value, + group = variable, + color = west.germany, + linetype = west.germany.size)) + + geom_line() + + labs(color = "West Germany?")+ + guides(linetype = FALSE) + + theme_bw(base_size = 10) + + theme(legend.position = "bottom") + + geom_vline(xintercept = 1990, + color = "#88419d") + + annotate("text", x = 1985, + label="Reunification", y = 4000, + colour = "#88419d", alpha = 0.8, angle = 0) + + scale_y_continuous(name="Gap in per-capita GDP (PPP, 2002 USD)", + labels = comma_format(big.mark = ",", + decimal.mark = ".")) + + scale_color_manual(values = graphcolors) + + coord_cartesian(ylim = c(-5000, 5000)) #+ + #ggsave(paste0("figures/allrmspe.pdf"), + # width = 15, height = 10, units = 'cm') + + + + + + + +# Figure 7 +# compute ratio of post-reunification RMSPE +# to pre-reunification RMSPE +rmse <- function(x){sqrt(mean(x^2))} +preloss <- apply(storegaps[1:30,],2,rmse) +postloss <- apply(storegaps[31:44,],2,rmse) +prepost <- sort(postloss/preloss) +prepost <- t(t(prepost)) + +prepost <- data.frame(Country = rownames(prepost), RMSPE = prepost) +rownames(prepost) <- NULL +prepost$Country <- factor(prepost$Country, + levels = prepost$Country[order(prepost$RMSPE)]) + + +prepost$west.germany <- ifelse(prepost$Country== "West Germany", "West Germany", "Not West Germany") + +ggplot(prepost, aes(x = Country, y = RMSPE)) + + geom_point(stat = 'identity', size = 6, aes(color = west.germany)) + + theme_bw(base_size = 10) + + coord_flip() + + ylab("Post-Period RMSE / Pre-Period RMSE") + + scale_color_manual(values = graphcolors) + + guides(color = FALSE) #+ + #ggsave(paste0("figures/RMSPERatio.pdf"), + # width = 15, height = 10, units = 'cm') + + + + + diff --git a/R/synthetic_control/rep_original.r b/R/synthetic_control/rep_original.r new file mode 100755 index 0000000..f1b77a1 --- /dev/null +++ b/R/synthetic_control/rep_original.r @@ -0,0 +1,640 @@ +## Replication Code for +# A. Abadie, A. Diamond, and J. Hainmueller. 2014. +# Comparative Politics and the Synthetic Control Method +# American Journal of Political Science. + +rm(list=ls()) +library(foreign) +library("Synth") +library(xtable) + +# Load Data +d <- read.dta("repgermany.dta") + +## Table 1 & 2, Figure 1, 2, & 3 + +## pick v by cross-validation +# data setup for training model +dataprep.out <- + dataprep( + foo = d, + predictors = c("gdp","trade","infrate"), + dependent = "gdp", + unit.variable = 1, + time.variable = 3, + special.predictors = list( + list("industry", 1971:1980, c("mean")), + list("schooling",c(1970,1975), c("mean")), + list("invest70" ,1980, c("mean")) + ), + treatment.identifier = 7, + controls.identifier = unique(d$index)[-7], + time.predictors.prior = 1971:1980, + time.optimize.ssr = 1981:1990, + unit.names.variable = 2, + time.plot = 1960:2003 + ) + +# fit training model +synth.out <- + synth( + data.prep.obj=dataprep.out, + Margin.ipop=.005,Sigf.ipop=7,Bound.ipop=6 + ) + +# data prep for main model +dataprep.out <- + dataprep( + foo = d, + predictors = c("gdp","trade","infrate"), + dependent = "gdp", + unit.variable = 1, + time.variable = 3, + special.predictors = list( + list("industry" ,1981:1990, c("mean")), + list("schooling",c(1980,1985), c("mean")), + list("invest80" ,1980, c("mean")) + ), + treatment.identifier = 7, + controls.identifier = unique(d$index)[-7], + time.predictors.prior = 1981:1990, + time.optimize.ssr = 1960:1989, + unit.names.variable = 2, + time.plot = 1960:2003 +) + +# fit main model with v from training model +synth.out <- synth( + data.prep.obj=dataprep.out, + custom.v=as.numeric(synth.out$solution.v) + ) + +#### Table 2 +synth.tables <- synth.tab( + dataprep.res = dataprep.out, + synth.res = synth.out + ); synth.tables + +# Replace means for OECD sample (computed externally using proper pop weighting) +synth.tables$tab.pred[,3] <- c(8021.1,31.9,7.4,34.2,44.1,25.9) +colnames(synth.tables$tab.pred)[3] <- "Rest of OECD Sample" +rownames(synth.tables$tab.pred) <- c("GDP per-capita","Trade openness", + "Inflation rate","Industry share", + "Schooling","Investment rate") + +xtable(round(synth.tables$tab.pred,1),digits=1) + +#### Table 1 +# synth weights +tab1 <- data.frame(synth.tables$tab.w) +tab1[,1] <- round(tab1[,1],2) +# regression weights +X0 <- cbind(1,t(dataprep.out$X0)) +X1 <- as.matrix(c(1,dataprep.out$X1)) +W <- X0%*%solve(t(X0)%*%X0)%*%X1 +Wdat <- data.frame(unit.numbers=as.numeric(rownames(X0)), + regression.w=round(W,2)) +tab1 <- merge(tab1,Wdat,by="unit.numbers") +tab1 <- tab1[order(tab1[,3]),] + +xtable(cbind(tab1[1:9,c(3,2,4)], + tab1[10:18,c(3,2,4)] + ) + ) + +#### Figure 1: Trends in Per-Capita GDP: West Germany vs. Rest of the OECD Sample +Text.height <- 23000 +Cex.set <- .8 +#pdf(file = "ger_vs_oecd.pdf", width = 5.5, height = 5.5, family = "Times",pointsize = 12) +plot(1960:2003,dataprep.out$Y1plot, + type="l",ylim=c(0,33000),col="black",lty="solid", + ylab ="per-capita GDP (PPP, 2002 USD)", + xlab ="year", + xaxs = "i", yaxs = "i", + lwd=2 + ) +lines(1960:2003,aggregate(d[,c("gdp")],by=list(d$year),mean,na.rm=T)[,2] + ,col="black",lty="dashed",lwd=2) # mean 2 +abline(v=1990,lty="dotted") +legend(x="bottomright", + legend=c("West Germany","rest of the OECD sample") + ,lty=c("solid","dashed"),col=c("black","black") + ,cex=.8,bg="white",lwd=c(2,2)) +arrows(1987,Text.height,1989,Text.height,col="black",length=.1) +text(1982.5,Text.height,"reunification",cex=Cex.set) +#dev.off() + +#### Figure 2: Trends in Per-Capita GDP: West Germany vs. Synthetic West Germany +#pdf(file = "ger_vs_synthger2.pdf", width = 5.5, height = 5.5, family = "Times",pointsize = 12) +synthY0 <- (dataprep.out$Y0%*%synth.out$solution.w) +plot(1960:2003,dataprep.out$Y1plot, + type="l",ylim=c(0,33000),col="black",lty="solid", + ylab ="per-capita GDP (PPP, 2002 USD)", + xlab ="year", + xaxs = "i", yaxs = "i", + lwd=2 + ) +lines(1960:2003,synthY0,col="black",lty="dashed",lwd=2) +abline(v=1990,lty="dotted") +legend(x="bottomright", + legend=c("West Germany","synthetic West Germany") + ,lty=c("solid","dashed"),col=c("black","black") + ,cex=.8,bg="white",lwd=c(2,2)) +arrows(1987,Text.height,1989,Text.height,col="black",length=.1) +text(1982.5,Text.height,"reunification",cex=Cex.set) +#dev.off() + +### Figure 3: Per-Capita GDP Gap Between West Germany and Synthetic West Germany +#pdf(file = "ger_vs_synthger_gaps2.pdf", width = 5.5, height = 5.5, family = "Times",pointsize = 12) +gap <- dataprep.out$Y1-(dataprep.out$Y0%*%synth.out$solution.w) +plot(1960:2003,gap, + type="l",ylim=c(-4500,4500),col="black",lty="solid", + ylab =c("gap in per-capita GDP (PPP, 2002 USD)"), + xlab ="year", + xaxs = "i", yaxs = "i", + lwd=2 + ) +abline(v=1990,lty="dotted") +abline(h=0,lty="dotted") +arrows(1987,1000,1989,1000,col="black",length=.1) +text(1982.5,1000,"reunification",cex=Cex.set) +#dev.off() + +### Figure 4: Placebo Reunification 1975 - Trends in Per-Capita GDP: West Germany vs. Synthetic West Germany + +# data prep for training model +dataprep.out <- + dataprep( + foo = d, + predictors = c("gdp","trade","infrate"), + dependent = "gdp", + unit.variable = 1, + time.variable = 3, + special.predictors = list( + list("industry",1971, c("mean")), + list("schooling",c(1960,1965), c("mean")), + list("invest60" ,1980, c("mean")) + ), + treatment.identifier = 7, + controls.identifier = unique(d$index)[-7], + time.predictors.prior = 1960:1964, + time.optimize.ssr = 1965:1975, + unit.names.variable = 2, + time.plot = 1960:1990 + ) + +# fit training model +synth.out <- synth( + data.prep.obj=dataprep.out, + Margin.ipop=.005,Sigf.ipop=7,Bound.ipop=6 +) + + +# data prep for main model +dataprep.out <- + dataprep( + foo = d, + predictors = c("gdp","trade","infrate"), + dependent = "gdp", + unit.variable = 1, + time.variable = 3, + special.predictors = list( + list("industry" ,1971:1975, c("mean")), + list("schooling",c(1970,1975), c("mean")), + list("invest70" ,1980, c("mean")) + ), + treatment.identifier = 7, + controls.identifier = unique(d$index)[-7], + time.predictors.prior = 1965:1975, + time.optimize.ssr = 1960:1975, + unit.names.variable = 2, + time.plot = 1960:1990 + ) + +# fit main model +synth.out <- synth( + data.prep.obj=dataprep.out, + custom.v=as.numeric(synth.out$solution.v) +) + +Cex.set <- 1 +#pdf(file = "2intimeplacebo1975.pdf", width = 5.5, height = 5.5, family = "Times",pointsize = 12) +plot(1960:1990,dataprep.out$Y1plot, + type="l",ylim=c(0,33000),col="black",lty="solid", + ylab ="per-capita GDP (PPP, 2002 USD)", + xlab ="year", + xaxs = "i", yaxs = "i", + lwd=2 + ) +lines(1960:1990,(dataprep.out$Y0%*%synth.out$solution.w),col="black",lty="dashed",lwd=2) +abline(v=1975,lty="dotted") +legend(x="bottomright", + legend=c("West Germany","synthetic West Germany") + ,lty=c("solid","dashed"),col=c("black","black") + ,cex=.8,bg="white",lwd=c(2,2)) +arrows(1973,20000,1974.5,20000,col="black",length=.1) +text(1967.5,20000,"placebo reunification",cex=Cex.set) +#dev.off() + + +### Figure 5: Ratio of post-reunification RMSPE to pre-reunification RMSPE: West Germany and control countries. + +# loop across control units +storegaps <- + matrix(NA, + length(1960:2003), + length(unique(d$index))-1 + ) +rownames(storegaps) <- 1960:2003 +i <- 1 +co <- unique(d$index) + +for(k in unique(d$index)[-7]){ + + # data prep for training model + dataprep.out <- + dataprep( + foo = d, + predictors = c("gdp","trade","infrate"), + dependent = "gdp", + unit.variable = 1, + time.variable = 3, + special.predictors = list( + list("industry",1971:1980, c("mean")), + list("schooling" ,c(1970,1975), c("mean")), + list("invest70" ,1980, c("mean")) + ), + treatment.identifier = k, + controls.identifier = co[-which(co==k)], + time.predictors.prior = 1971:1980, + time.optimize.ssr = 1981:1990, + unit.names.variable = 2, + time.plot = 1960:2003 + ) + + # fit training model + synth.out <- + synth( + data.prep.obj=dataprep.out, + Margin.ipop=.005,Sigf.ipop=7,Bound.ipop=6 + ) + + # data prep for main model +dataprep.out <- + dataprep( + foo = d, + predictors = c("gdp","trade","infrate"), + dependent = "gdp", + unit.variable = 1, + time.variable = 3, + special.predictors = list( + list("industry" ,1981:1990, c("mean")), + list("schooling",c(1980,1985), c("mean")), + list("invest80" ,1980, c("mean")) + ), + treatment.identifier = k, + controls.identifier = co[-which(co==k)], + time.predictors.prior = 1981:1990, + time.optimize.ssr = 1960:1989, + unit.names.variable = 2, + time.plot = 1960:2003 + ) + +# fit main model +synth.out <- synth( + data.prep.obj=dataprep.out, + custom.v=as.numeric(synth.out$solution.v) + ) + + storegaps[,i] <- + dataprep.out$Y1- + (dataprep.out$Y0%*%synth.out$solution.w) + i <- i + 1 +} # close loop over control units +d <- d[order(d$index,d$year),] +colnames(storegaps) <- unique(d$country)[-7] +storegaps <- cbind(gap,storegaps) +colnames(storegaps)[1] <- c("West Germany") + +# compute ratio of post-reunification RMSPE +# to pre-reunification RMSPE +rmse <- function(x){sqrt(mean(x^2))} +preloss <- apply(storegaps[1:30,],2,rmse) +postloss <- apply(storegaps[31:44,],2,rmse) + +#pdf("2ratio_post_to_preperiod_rmse2a.pdf") +dotchart(sort(postloss/preloss), + xlab="Post-Period RMSE / Pre-Period RMSE", + pch=19) +#dev.off() + +### Figure 6: Leave-one-out distribution of the synthetic control for West Germany + +# loop over leave one outs +storegaps <- + matrix(NA, + length(1960:2003), + 5) +colnames(storegaps) <- c(1,3,9,12,14) +co <- unique(d$index)[-7] + +for(k in 1:5){ + +# data prep for training model +omit <- c(1,3,9,12,14)[k] + dataprep.out <- + dataprep( + foo = d, + predictors = c("gdp","trade","infrate"), + dependent = "gdp", + unit.variable = 1, + time.variable = 3, + special.predictors = list( + list("industry",1971:1980, c("mean")), + list("schooling" ,c(1970,1975), c("mean")), + list("invest70" ,1980, c("mean")) + ), + treatment.identifier = 7, + controls.identifier = co[-which(co==omit)], + time.predictors.prior = 1971:1980, + time.optimize.ssr = 1981:1990, + unit.names.variable = 2, + time.plot = 1960:2003 + ) + + # fit training model + synth.out <- synth( + data.prep.obj=dataprep.out, + Margin.ipop=.005,Sigf.ipop=7,Bound.ipop=6 + ) + +# data prep for main model +dataprep.out <- + dataprep( + foo = d, + predictors = c("gdp","trade","infrate"), + dependent = "gdp", + unit.variable = 1, + time.variable = 3, + special.predictors = list( + list("industry" ,1981:1990, c("mean")), + list("schooling",c(1980,1985), c("mean")), + list("invest80" ,1980, c("mean")) + ), + treatment.identifier = 7, + controls.identifier = co[-which(co==omit)], + time.predictors.prior = 1981:1990, + time.optimize.ssr = 1960:1989, + unit.names.variable = 2, + time.plot = 1960:2003 + ) + + # fit main model + synth.out <- synth( + data.prep.obj=dataprep.out, + custom.v=as.numeric(synth.out$solution.v) + ) + storegaps[,k] <- (dataprep.out$Y0%*%synth.out$solution.w) +} # close loop over leave one outs + +Text.height <- 23000 +Cex.set <- .8 +#pdf(file = "1jackknife2.pdf", width = 5.5, height = 5.5, family = "Times",pointsize = 12) +plot(1960:2003,dataprep.out$Y1plot, + type="l",ylim=c(0,33000),col="black",lty="solid", + ylab ="per-capita GDP (PPP, 2002 USD)", + xlab ="year", + xaxs = "i", yaxs = "i",lwd=2 + ) + +abline(v=1990,lty="dotted") +arrows(1987,23000,1989,23000,col="black",length=.1) + for(i in 1:5){ + lines(1960:2003,storegaps[,i],col="darkgrey",lty="solid") + } +lines(1960:2003,synthY0,col="black",lty="dashed",lwd=2) +lines(1960:2003,dataprep.out$Y1plot,col="black",lty="solid",lwd=2) +text(1982.5,23000,"reunification",cex=.8) +legend(x="bottomright", + legend=c("West Germany", + "synthetic West Germany", + "synthetic West Germany (leave-one-out)") + ,lty=c("solid","dashed","solid"), + col=c("black","black","darkgrey") + ,cex=.8,bg="white",lwdc(2,2,1)) +#dev.off() + + +### Table 3: Synthetic Weights from Combinations of Control Countries +rm(list=ls()) +library(gtools) +library(kernlab) + +# data prep for training model +d <- read.dta("repgermany.dta") +dataprep.out <- + dataprep( + foo = d, + predictors = c("gdp","trade","infrate"), + dependent = "gdp", + unit.variable = 1, + time.variable = 3, + special.predictors = list( + list("industry", 1971:1980, c("mean")), + list("schooling",c(1970,1975), c("mean")), + list("invest70" ,1980, c("mean")) + ), + treatment.identifier = 7, + controls.identifier = unique(d$index)[-7], + time.predictors.prior = 1971:1980, + time.optimize.ssr = 1981:1990, + unit.names.variable = 2, + time.plot = 1960:2003 + ) + +# fit training model +synth.out <- + synth( + data.prep.obj=dataprep.out, + Margin.ipop=.005,Sigf.ipop=7,Bound.ipop=6 + ) + +# data prep for main model +dataprep.out <- + dataprep( + foo = d, + predictors = c("gdp","trade","infrate"), + dependent = "gdp", + unit.variable = 1, + time.variable = 3, + special.predictors = list( + list("industry" ,1981:1990, c("mean")), + list("schooling",c(1980,1985), c("mean")), + list("invest80" ,1980, c("mean")) + ), + treatment.identifier = 7, + controls.identifier = unique(d$index)[-7], + time.predictors.prior = 1981:1990, + time.optimize.ssr = 1960:1989, + unit.names.variable = 2, + time.plot = 1960:2003 + ) + +# fit main model with v from training model +synth.out <- synth( + data.prep.obj=dataprep.out, + custom.v=as.numeric(synth.out$solution.v) +) + +synth.tables <- synth.tab( + dataprep.res = dataprep.out, + synth.res = synth.out +) + +table3 <- list() +synth.tables$tab.w[,1] <- round(synth.tables$tab.w[,1],2) +table3[[5]] <-synth.tables$tab.w[order(-1*synth.tables$tab.w[,1]),2:1][1:5,] + +# compute loss for all combinations +# of 4, 3, 2, 1 sized donor pools + +# get W and v +solution.w <- round(synth.out$solution.w,3) +V <- diag(as.numeric(synth.out$solution.v)) + +# compute scaled Xs +nvarsV <- dim(dataprep.out$X0)[1] +big.dataframe <- cbind(dataprep.out$X0, dataprep.out$X1) +divisor <- sqrt(apply(big.dataframe, 1, var)) +scaled.matrix <- + t(t(big.dataframe) %*% ( 1/(divisor) * + diag(rep(dim(big.dataframe)[1], 1)) )) +X0.scaled <- scaled.matrix[,c(1:(dim(dataprep.out$X0)[2]))] +X1.scaled <- as.matrix(scaled.matrix[,dim(scaled.matrix)[2]]) + +dn <- d[d$year==1970,c("country","index")] +dn <- dn[order(dn$index),] +dn <- dn[-7,] + +table2store <- matrix(NA,nrow(dataprep.out$X1),4) +fig7store <- matrix(NA,length(1960:2003),4) + +# loop through number of controls +for(pp in 4:1){ + store <- combinations(length(unique(d$index)[-7]), + r=pp, v=unique(d$index)[-7]) + store.loss <- matrix(NA,nrow=nrow(store),1) + store.w <- matrix(NA,nrow=nrow(store),pp) + store.c <- store.w + +# loop through combinations + for(k in 1:nrow(store)){ + # index positions of control units + posvector <- c() + for(i in 1:pp){ + posvector <- c(posvector,which(dn$index==store[k,i])) + } + + # run quad optimization + X0temp <- X0.scaled[ , posvector ] + H <- t(X0temp) %*% V %*% (X0temp) + c <- -1*c(t(X1.scaled) %*% V %*% (X0temp) ) + + if(pp == 1){ + solution.w <- matrix(1) + } else { + res <- ipop(c = c, H = H, A = t(rep(1, length(c))), + b = 1, l = rep(0, length(c)), + u = rep(1, length(c)), r = 0, + margin = 0.005,sigf = 7, bound = 6) + solution.w <- as.matrix(primal(res)) + } + loss.w <- t(X1.scaled - X0temp %*% solution.w) %*% V %*% (X1.scaled - X0temp %*% solution.w) + + store.loss[k] <- loss.w + store.w[k,] <- t(solution.w) + store.c[k,] <- dn$country[posvector] + } # close loop over combinations + +# get best fitting combination + dat <- data.frame(store.loss, + store, + store.c, + store.w + ) + colnames(dat) <- c("loss", + paste("CNo",1:pp,sep=""), + paste("CNa",1:pp,sep=""), + paste("W",1:pp,sep="") + ) + dat <- dat[order(dat$loss),] + Countries <- dat[1,paste("CNo",1:pp,sep="")] + Cweights <- as.numeric(dat[1,paste("W",1:pp,sep="")]) + + outdat <- data.frame(unit.names=as.vector( + (t(as.vector(dat[1,paste("CNa",1:pp,sep="")])))), + w.weights=round(Cweights,2)) + +table3[[pp]]<- outdat[order(-1*outdat$w.weights),] + + # get posvector for fitting + posvector <- c() + if(pp == 1 ){ + posvector <- c(posvector,which(dn$index==Countries)) + } else { + for(i in 1:pp){ + posvector <- c(posvector,which(dn$index==Countries[1,i])) + } + } + + X0t <- as.matrix(dataprep.out$X0[,posvector])%*% as.matrix(Cweights) + table2store[,(4:1)[pp]] <- X0t + + fig7store[,(4:1)[pp]] <- + dataprep.out$Y0[,posvector]%*%as.matrix(Cweights) + +} # close loop over number of countries + +# Table 3 +table3 + +# Table 4 +synth.tables$tab.pred[,3] <- c(8021.1,31.9,7.4,34.2,44.1,25.9) +table4 <- round( + cbind(synth.tables$tab.pred[,1:2], + table2store, + synth.tables$tab.pred[,3]),1) +rownames(table4) <- c("GDP per-capita","Trade openness", + "Inflation rate","Industry share", + "Schooling","Investment rate") +colnames(table4)[2:7] <- c(5:1,"OECD Sample") +table4 + +## Figure 7: Per-Capita GDP Gaps Between West Germany and Sparse Synthetic Controls +Text.height <- 23000 +Cex.set <- .8 + +par(mfrow=c(2,2)) +for(pp in 4:1){ +#pdf(file = paste("2ger_vs_synth","CValt",pp,".pdf",sep=""), width = 5.5, height = 5.5, family = "Times",pointsize = 12) + plot(1960:2003,dataprep.out$Y1, + type="l",ylim=c(0,33000),col="black",lty="solid", + ylab ="per-capita GDP (PPP, 2002 USD)", + xlab ="year", + xaxs = "i", yaxs = "i", + lwd=2, + main=paste("No. of control countries: ",pp,sep="") + ) + lines(1960:2003,fig7store[,c(4:1)[pp]],col="black",lty="dashed",lwd=2) + abline(v=1990,lty="dotted") + legend(x="bottomright", + legend=c("West Germany","synthetic West Germany") + ,lty=c("solid","dashed"),col=c("black","black") + ,cex=.8,bg="white",lwd=c(2,2)) + arrows(1987,Text.height,1989,Text.height,col="black",length=.1) + text(1982.5,Text.height,"reunification",cex=Cex.set) + #dev.off() +} + + + From 12549fe21346d6efd855be16a241c0a03b5fa044 Mon Sep 17 00:00:00 2001 From: juanitorduz Date: Sun, 7 Nov 2021 13:43:46 +0100 Subject: [PATCH 2/3] update renv --- renv.lock | 528 +++++++++++++++++++++++++++++++++--------------- renv/.gitignore | 2 + renv/activate.R | 427 ++++++++++++++++++++++++++++++++++----- 3 files changed, 748 insertions(+), 209 deletions(-) diff --git a/renv.lock b/renv.lock index 52ece59..4608524 100644 --- a/renv.lock +++ b/renv.lock @@ -1,6 +1,6 @@ { "R": { - "Version": "3.5.1", + "Version": "4.1.2", "Repositories": [ { "Name": "CRAN", @@ -11,31 +11,31 @@ "Packages": { "DBI": { "Package": "DBI", - "Version": "1.1.0", + "Version": "1.1.1", "Source": "Repository", "Repository": "CRAN", - "Hash": "4744be45519d675af66c28478720fce5" + "Hash": "030aaec5bc6553f35347cbb1e70b1a17" }, "MASS": { "Package": "MASS", - "Version": "7.3-50", + "Version": "7.3-54", "Source": "Repository", "Repository": "CRAN", - "Hash": "b24dc52a732ee76079252f793c210b37" + "Hash": "0e59129db205112e3963904db67fd0dc" }, "Matrix": { "Package": "Matrix", - "Version": "1.2-14", + "Version": "1.3-4", "Source": "Repository", "Repository": "CRAN", - "Hash": "e84dcbc31943fe880bbb80fd59baca6d" + "Hash": "4ed05e9c9726267e4a5872e09c04587c" }, "R6": { "Package": "R6", - "Version": "2.4.1", + "Version": "2.5.1", "Source": "Repository", "Repository": "CRAN", - "Hash": "292b54f8f4b94669b08f94e5acce6be2" + "Hash": "470851b6d5d0ac559e9d01bb352b4021" }, "RColorBrewer": { "Package": "RColorBrewer", @@ -46,10 +46,17 @@ }, "Rcpp": { "Package": "Rcpp", - "Version": "1.0.5", + "Version": "1.0.7", "Source": "Repository", "Repository": "CRAN", - "Hash": "125dc7a0ed375eb68c0ce533b48d291f" + "Hash": "dab19adae4440ae55aa8a9d238b246bb" + }, + "Synth": { + "Package": "Synth", + "Version": "1.1-5", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "bd6eb9b3b5b37e14f7cec093a48312bd" }, "askpass": { "Package": "askpass", @@ -67,10 +74,10 @@ }, "backports": { "Package": "backports", - "Version": "1.1.8", + "Version": "1.3.0", "Source": "Repository", "Repository": "CRAN", - "Hash": "3ef0eac19317fd03c0c854aed581d473" + "Hash": "194ad71f8ed59393272a0c4be2eac215" }, "base64enc": { "Package": "base64enc", @@ -79,19 +86,40 @@ "Repository": "CRAN", "Hash": "543776ae6848fde2f48ff3816d0628bc" }, + "bit": { + "Package": "bit", + "Version": "4.0.4", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "f36715f14d94678eea9933af927bc15d" + }, + "bit64": { + "Package": "bit64", + "Version": "4.0.5", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "9fe98599ca456d6552421db0d6772d8f" + }, + "blob": { + "Package": "blob", + "Version": "1.2.2", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "dc5f7a6598bb025d20d66bb758f12879" + }, "broom": { "Package": "broom", - "Version": "0.5.5", + "Version": "0.7.10", "Source": "Repository", "Repository": "CRAN", - "Hash": "70ec3aef8a0df13661bf6cc2ff877d61" + "Hash": "ddf8bc55ea050f984835dd2d23cd6828" }, "callr": { "Package": "callr", - "Version": "3.4.2", + "Version": "3.7.0", "Source": "Repository", "Repository": "CRAN", - "Hash": "d16f6e202179e69c1400cc1f845bfc09" + "Hash": "461aa75a11ce2400245190ef5d3995df" }, "cellranger": { "Package": "cellranger", @@ -102,66 +130,94 @@ }, "cli": { "Package": "cli", - "Version": "2.0.2", + "Version": "3.1.0", "Source": "Repository", "Repository": "CRAN", - "Hash": "ff0becff7bfdfe3f75d29aff8f3172dd" + "Hash": "66a3834e54593c89d8beefb312347e58" }, "clipr": { "Package": "clipr", - "Version": "0.7.0", + "Version": "0.7.1", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "ebaa97ac99cc2daf04e77eecc7b781d7" + }, + "codetools": { + "Package": "codetools", + "Version": "0.2-18", "Source": "Repository", "Repository": "CRAN", - "Hash": "08cf4045c149a0f0eaf405324c7495bd" + "Hash": "019388fc48e48b3da0d3a76ff94608a8" }, "colorspace": { "Package": "colorspace", - "Version": "1.4-1", + "Version": "2.0-2", "Source": "Repository", "Repository": "CRAN", - "Hash": "6b436e95723d1f0e861224dd9b094dfb" + "Hash": "6baccb763ee83c0bd313460fdb8b8a84" + }, + "cpp11": { + "Package": "cpp11", + "Version": "0.4.1", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "70976176dfd7f179f212783aab2547b1" }, "crayon": { "Package": "crayon", - "Version": "1.3.4", + "Version": "1.4.2", "Source": "Repository", "Repository": "CRAN", - "Hash": "0d57bc8e27b7ba9e45dba825ebc0de6b" + "Hash": "0a6a65d92bd45b47b94b84244b528d17" }, "curl": { "Package": "curl", - "Version": "4.3", + "Version": "4.3.2", "Source": "Repository", "Repository": "CRAN", - "Hash": "2b7d10581cc730804e9ed178c8374bd6" + "Hash": "022c42d49c28e95d69ca60446dbabf88" + }, + "data.table": { + "Package": "data.table", + "Version": "1.14.2", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "36b67b5adf57b292923f5659f5f0c853" }, "dbplyr": { "Package": "dbplyr", - "Version": "1.4.2", + "Version": "2.1.1", "Source": "Repository", "Repository": "CRAN", - "Hash": "022c8630abecb00f22740d021ab89595" + "Hash": "1f37fa4ab2f5f7eded42f78b9a887182" }, "digest": { "Package": "digest", - "Version": "0.6.25", + "Version": "0.6.28", "Source": "Repository", "Repository": "CRAN", - "Hash": "f697db7d92b7028c4b3436e9603fb636" + "Hash": "49b5c6e230bfec487b8917d5a0c77cca" }, "dplyr": { "Package": "dplyr", - "Version": "0.8.5", + "Version": "1.0.7", "Source": "Repository", "Repository": "CRAN", - "Hash": "57a42ddf80f429764ff7987128c3fd0a" + "Hash": "36f1ae62f026c8ba9f9b5c9a08c03297" + }, + "dtplyr": { + "Package": "dtplyr", + "Version": "1.1.0", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "1e14e4c5b2814de5225312394bc316da" }, "ellipsis": { "Package": "ellipsis", - "Version": "0.3.0", + "Version": "0.3.2", "Source": "Repository", "Repository": "CRAN", - "Hash": "7067d90c1c780bfe80c0d497e3d7b49d" + "Hash": "bb0eec2fe32e88d9e2836c2f73ea2077" }, "evaluate": { "Package": "evaluate", @@ -172,59 +228,108 @@ }, "fansi": { "Package": "fansi", - "Version": "0.4.1", + "Version": "0.5.0", "Source": "Repository", "Repository": "CRAN", - "Hash": "7fce217eaaf8016e72065e85c73027b5" + "Hash": "d447b40982c576a72b779f0a3b3da227" }, "farver": { "Package": "farver", - "Version": "2.0.3", + "Version": "2.1.0", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "c98eb5133d9cb9e1622b8691487f11bb" + }, + "fastmap": { + "Package": "fastmap", + "Version": "1.1.0", "Source": "Repository", "Repository": "CRAN", - "Hash": "dad6793a5a1f73c8e91f1a1e3e834b05" + "Hash": "77bd60a6157420d4ffa93b27cf6a58b8" }, "forcats": { "Package": "forcats", - "Version": "0.5.0", + "Version": "0.5.1", "Source": "Repository", "Repository": "CRAN", - "Hash": "1cb4279e697650f0bd78cd3601ee7576" + "Hash": "81c3244cab67468aac4c60550832655d" }, "foreign": { "Package": "foreign", - "Version": "0.8-70", + "Version": "0.8-81", "Source": "Repository", "Repository": "CRAN", - "Hash": "3d11f0b67c9c71cd3fde607063f70635" + "Hash": "74628ea7a3be5ee8a7b5bb0a8e84882e" }, "fs": { "Package": "fs", - "Version": "1.3.2", + "Version": "1.5.0", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "44594a07a42e5f91fac9f93fda6d0109" + }, + "furrr": { + "Package": "furrr", + "Version": "0.2.3", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "2aba4ab06e8707ac054c6683cb6fed56" + }, + "future": { + "Package": "future", + "Version": "1.23.0", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "7bf6fbed7f00cae876901fd70c04f3a4" + }, + "gargle": { + "Package": "gargle", + "Version": "1.2.0", "Source": "Repository", "Repository": "CRAN", - "Hash": "77cefb9b0bb109ec5ce526e6948c91aa" + "Hash": "9d234e6a87a6f8181792de6dc4a00e39" }, "generics": { "Package": "generics", - "Version": "0.0.2", + "Version": "0.1.1", "Source": "Repository", "Repository": "CRAN", - "Hash": "b8cff1d1391fd1ad8b65877f4c7f2e53" + "Hash": "3f6bcfb0ee5d671d9fd1893d2faa79cb" }, "ggplot2": { "Package": "ggplot2", - "Version": "3.3.0", + "Version": "3.3.5", "Source": "Repository", "Repository": "CRAN", - "Hash": "911561e07da928345f1ae2d69f97f3ea" + "Hash": "d7566c471c7b17e095dd023b9ef155ad" + }, + "globals": { + "Package": "globals", + "Version": "0.14.0", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "eca8023ed5ca6372479ebb9b3207f5ae" }, "glue": { "Package": "glue", - "Version": "1.4.1", + "Version": "1.4.2", "Source": "Repository", "Repository": "CRAN", - "Hash": "f43e0d5e85ccb0a4045670c0607ee504" + "Hash": "6efd734b14c6471cfe443345f3e35e29" + }, + "googledrive": { + "Package": "googledrive", + "Version": "2.0.0", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "c3a25adbbfbb03f12e6f88c5fb1f3024" + }, + "googlesheets4": { + "Package": "googlesheets4", + "Version": "1.0.0", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "9a6564184dc4a81daea4f1d7ce357c6a" }, "gtable": { "Package": "gtable", @@ -233,131 +338,152 @@ "Repository": "CRAN", "Hash": "ac5c6baf7822ce8732b343f14c072c4d" }, - "haven": { - "Package": "haven", - "Version": "2.2.0", + "gtools": { + "Package": "gtools", + "Version": "3.9.2", "Source": "Repository", "Repository": "CRAN", - "Hash": "e3f662e125e9fdffd1ee4e94baea3451" + "Hash": "2ace6c4a06297d0b364e0444384a2b82" }, - "here": { - "Package": "here", - "Version": "0.1", + "haven": { + "Package": "haven", + "Version": "2.4.3", "Source": "Repository", "Repository": "CRAN", - "Hash": "2c0406b8e0a4c3516ab37be62da74e3c" + "Hash": "10bec8a8264f3eb59531e8c4c0303f96" }, "highr": { "Package": "highr", - "Version": "0.8", + "Version": "0.9", "Source": "Repository", "Repository": "CRAN", - "Hash": "4dc5bb88961e347a0f4d8aad597cbfac" + "Hash": "8eb36c8125038e648e5d111c0d7b2ed4" }, "hms": { "Package": "hms", - "Version": "0.5.3", + "Version": "1.1.1", "Source": "Repository", "Repository": "CRAN", - "Hash": "726671f634529d470545f9fd1a9d1869" + "Hash": "5b8a2dd0fdbe2ab4f6081e6c7be6dfca" }, "htmltools": { "Package": "htmltools", - "Version": "0.5.0", + "Version": "0.5.2", "Source": "Repository", "Repository": "CRAN", - "Hash": "7d651b7131794fe007b1ad6f21aaa401" + "Hash": "526c484233f42522278ab06fb185cb26" }, "httr": { "Package": "httr", - "Version": "1.4.1", + "Version": "1.4.2", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "a525aba14184fec243f9eaec62fbed43" + }, + "ids": { + "Package": "ids", + "Version": "1.0.1", "Source": "Repository", "Repository": "CRAN", - "Hash": "7146fea4685b4252ebf478978c75f597" + "Hash": "99df65cfef20e525ed38c3d2577f7190" }, "isoband": { "Package": "isoband", - "Version": "0.2.0", + "Version": "0.2.5", "Source": "Repository", "Repository": "CRAN", - "Hash": "15f6d57a664cd953a31ae4ea61e5e60e" + "Hash": "7ab57a6de7f48a8dc84910d1eca42883" + }, + "jquerylib": { + "Package": "jquerylib", + "Version": "0.1.4", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "5aab57a3bd297eee1c1d862735972182" }, "jsonlite": { "Package": "jsonlite", - "Version": "1.7.0", + "Version": "1.7.2", "Source": "Repository", "Repository": "CRAN", - "Hash": "2657f20b9a74c996c602e74ebe540b06" + "Hash": "98138e0994d41508c7a6b84a0600cfcb" + }, + "kernlab": { + "Package": "kernlab", + "Version": "0.9-29", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "46b2f61dadae64579e4cbe8af3c088d6" }, "knitr": { "Package": "knitr", - "Version": "1.29", + "Version": "1.36", "Source": "Repository", "Repository": "CRAN", - "Hash": "e5f4c41c17df8cdf7b0df12117c0d99a" + "Hash": "46344b93f8854714cdf476433a59ed10" }, "labeling": { "Package": "labeling", - "Version": "0.3", + "Version": "0.4.2", "Source": "Repository", "Repository": "CRAN", - "Hash": "73832978c1de350df58108c745ed0e3e" + "Hash": "3d5108641f47470611a32d0bdf357a72" }, "lattice": { "Package": "lattice", - "Version": "0.20-35", + "Version": "0.20-45", "Source": "Repository", "Repository": "CRAN", - "Hash": "10a22a9a66fbe7944e9ef98985d7c927" + "Hash": "b64cdbb2b340437c4ee047a1f4c4377b" }, "lifecycle": { "Package": "lifecycle", - "Version": "0.2.0", + "Version": "1.0.1", "Source": "Repository", "Repository": "CRAN", - "Hash": "361811f31f71f8a617a9a68bf63f1f42" + "Hash": "a6b6d352e3ed897373ab19d8395c98d0" + }, + "listenv": { + "Package": "listenv", + "Version": "0.8.0", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "0bde42ee282efb18c7c4e63822f5b4f7" }, "lubridate": { "Package": "lubridate", - "Version": "1.7.4", + "Version": "1.8.0", "Source": "Repository", "Repository": "CRAN", - "Hash": "796afeea047cda6bdb308d374a33eeb6" + "Hash": "2ff5eedb6ee38fb1b81205c73be1be5a" }, "magrittr": { "Package": "magrittr", - "Version": "1.5", + "Version": "2.0.1", "Source": "Repository", "Repository": "CRAN", - "Hash": "1bb58822a20301cee84a41678e25d9b7" - }, - "markdown": { - "Package": "markdown", - "Version": "1.1", - "Source": "Repository", - "Repository": "CRAN", - "Hash": "61e4a10781dd00d7d81dd06ca9b94e95" + "Hash": "41287f1ac7d28a92f0a286ed507928d3" }, "mgcv": { "Package": "mgcv", - "Version": "1.8-24", + "Version": "1.8-38", "Source": "Repository", "Repository": "CRAN", - "Hash": "b71c00b21c099b3f4beca9d784ad7ec4" + "Hash": "be3c61ffbb1e3d3b3df214d192ac5444" }, "mime": { "Package": "mime", - "Version": "0.9", + "Version": "0.12", "Source": "Repository", "Repository": "CRAN", - "Hash": "e87a35ec73b157552814869f45a63aa3" + "Hash": "18e9c28c1d3ca1560ce30658b22ce104" }, "modelr": { "Package": "modelr", - "Version": "0.1.6", + "Version": "0.1.8", "Source": "Repository", "Repository": "CRAN", - "Hash": "30a2db10e829133fa0a1e6eaed25bc73" + "Hash": "9fd59716311ee82cba83dc2826fc5577" }, "munsell": { "Package": "munsell", @@ -368,24 +494,45 @@ }, "nlme": { "Package": "nlme", - "Version": "3.1-137", + "Version": "3.1-153", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "2d632e0d963a653a0329756ce701ecdd" + }, + "numDeriv": { + "Package": "numDeriv", + "Version": "2016.8-1.1", "Source": "Repository", "Repository": "CRAN", - "Hash": "4320ab595f7bbff5458bc6a044a57fc0" + "Hash": "df58958f293b166e4ab885ebcad90e02" }, "openssl": { "Package": "openssl", - "Version": "1.4.1", + "Version": "1.4.5", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "5406fd37ef0bf9b88c8a4f264d6ec220" + }, + "optimx": { + "Package": "optimx", + "Version": "2021-10.12", "Source": "Repository", "Repository": "CRAN", - "Hash": "49f7258fd86ebeaea1df24d9ded00478" + "Hash": "e384e8590262bb9f6c2c8c67569c4b1a" + }, + "parallelly": { + "Package": "parallelly", + "Version": "1.28.1", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "5300c9fc71841550bdca64d39e82af0e" }, "pillar": { "Package": "pillar", - "Version": "1.4.3", + "Version": "1.6.4", "Source": "Repository", "Repository": "CRAN", - "Hash": "fa3ed60396b6998d0427c57dab90fba4" + "Hash": "60200b6aa32314ac457d3efbb5ccbd98" }, "pkgconfig": { "Package": "pkgconfig", @@ -401,33 +548,54 @@ "Repository": "CRAN", "Hash": "ec0e5ab4e5f851f6ef32cd1d1984957f" }, + "prettyunits": { + "Package": "prettyunits", + "Version": "1.1.1", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "95ef9167b75dde9d2ccc3c7528393e7e" + }, "processx": { "Package": "processx", - "Version": "3.4.2", + "Version": "3.5.2", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "0cbca2bc4d16525d009c4dbba156b37c" + }, + "progress": { + "Package": "progress", + "Version": "1.2.2", "Source": "Repository", "Repository": "CRAN", - "Hash": "20a082f2bde0ffcd8755779fd476a274" + "Hash": "14dc9f7a3c91ebb14ec5bb9208a07061" }, "ps": { "Package": "ps", - "Version": "1.3.2", + "Version": "1.6.0", "Source": "Repository", "Repository": "CRAN", - "Hash": "98777535b61c57d1749344345e2a4ccd" + "Hash": "32620e2001c1dce1af49c49dccbb9420" }, "purrr": { "Package": "purrr", + "Version": "0.3.4", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "97def703420c8ab10d8f0e6c72101e02" + }, + "rappdirs": { + "Package": "rappdirs", "Version": "0.3.3", "Source": "Repository", "Repository": "CRAN", - "Hash": "22aca7d1181718e927d403a8c2d69d62" + "Hash": "5e3c5dc0b071b21fa128676560dbe94d" }, "readr": { "Package": "readr", - "Version": "1.3.1", + "Version": "2.0.2", "Source": "Repository", "Repository": "CRAN", - "Hash": "af8ab99cd936773a148963905736907b" + "Hash": "7cb2c3ecfbc2c6786221d2c0c1f6ed68" }, "readxl": { "Package": "readxl", @@ -443,68 +611,75 @@ "Repository": "CRAN", "Hash": "c66b930d20bb6d858cd18e1cebcfae5c" }, + "rematch2": { + "Package": "rematch2", + "Version": "2.1.2", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "76c9e04c712a05848ae7a23d2f170a40" + }, "renv": { "Package": "renv", - "Version": "0.10.0", + "Version": "0.14.0", "Source": "Repository", "Repository": "CRAN", - "Hash": "0d3f1c2b63e70b3d918246b4e2ca59b7" + "Hash": "30e5eba91b67f7f4d75d31de14bbfbdc" }, "reprex": { "Package": "reprex", - "Version": "0.3.0", + "Version": "2.0.1", "Source": "Repository", "Repository": "CRAN", - "Hash": "b06bfb3504cc8a4579fd5567646f745b" + "Hash": "911d101becedc0fde495bd910984bdc8" }, "reshape2": { "Package": "reshape2", - "Version": "1.4.3", + "Version": "1.4.4", "Source": "Repository", "Repository": "CRAN", - "Hash": "15a23ad30f51789188e439599559815c" + "Hash": "bb5996d0bd962d214a11140d77589917" }, "rlang": { "Package": "rlang", - "Version": "0.4.7", + "Version": "0.4.12", "Source": "Repository", "Repository": "CRAN", - "Hash": "c06d2a6887f4b414f8e927afd9ee976a" + "Hash": "0879f5388fe6e4d56d7ef0b7ccb031e5" }, "rmarkdown": { "Package": "rmarkdown", - "Version": "2.3", + "Version": "2.11", "Source": "Repository", "Repository": "CRAN", - "Hash": "202260e1b2c410edc086d5b8f1ed946e" + "Hash": "320017b52d05a943981272b295750388" }, - "rprojroot": { - "Package": "rprojroot", - "Version": "1.3-2", + "rsample": { + "Package": "rsample", + "Version": "0.1.0", "Source": "Repository", "Repository": "CRAN", - "Hash": "f6a407ae5dd21f6f80a6708bbb6eb3ae" + "Hash": "9e4f4f3b91998715bcc740f88ea328a3" }, "rstudioapi": { "Package": "rstudioapi", - "Version": "0.11", + "Version": "0.13", "Source": "Repository", "Repository": "CRAN", - "Hash": "33a5b27a03da82ac4b1d43268f80088a" + "Hash": "06c85365a03fdaf699966cc1d3cf53ea" }, "rvest": { "Package": "rvest", - "Version": "0.3.5", + "Version": "1.0.2", "Source": "Repository", "Repository": "CRAN", - "Hash": "6a20c2cdf133ebc7ac45888c9ccc052b" + "Hash": "bb099886deffecd6f9b298b7d4492943" }, "scales": { "Package": "scales", - "Version": "1.1.0", + "Version": "1.1.1", "Source": "Repository", "Repository": "CRAN", - "Hash": "a1c68369c629ea3188d0676e37069c65" + "Hash": "6f76f71042411426ec8df6c54f34e6dd" }, "selectr": { "Package": "selectr", @@ -513,12 +688,19 @@ "Repository": "CRAN", "Hash": "3838071b66e0c566d55cc26bd6e27bf4" }, + "slider": { + "Package": "slider", + "Version": "0.2.2", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "5237bd176dc0c4dd7eb8dcdafe514de3" + }, "stringi": { "Package": "stringi", - "Version": "1.4.6", + "Version": "1.7.5", "Source": "Repository", "Repository": "CRAN", - "Hash": "e99d8d656980d2dd416a962ae55aec90" + "Hash": "cd50dc9b449de3d3b47cdc9976886999" }, "stringr": { "Package": "stringr", @@ -529,94 +711,122 @@ }, "sys": { "Package": "sys", - "Version": "3.3", + "Version": "3.4", "Source": "Repository", "Repository": "CRAN", - "Hash": "507f3116a38d37ad330a038b3be07b66" + "Hash": "b227d13e29222b4574486cfcbde077fa" }, "tibble": { "Package": "tibble", - "Version": "2.1.3", + "Version": "3.1.5", "Source": "Repository", "Repository": "CRAN", - "Hash": "8248ee35d1e15d1e506f05f5a5d46a75" + "Hash": "36eb05ad4cfdfeaa56f5a9b2a1311efd" }, "tidyr": { "Package": "tidyr", - "Version": "1.0.2", + "Version": "1.1.4", "Source": "Repository", "Repository": "CRAN", - "Hash": "fb73a010ace00d6c584c2b53a21b969c" + "Hash": "c8fbdbd9fcac223d6c6fe8e406f368e1" }, "tidyselect": { "Package": "tidyselect", - "Version": "1.0.0", + "Version": "1.1.1", "Source": "Repository", "Repository": "CRAN", - "Hash": "7d4b0f1ab542d8cb7a40c593a4de2f36" + "Hash": "7243004a708d06d4716717fa1ff5b2fe" }, "tidyverse": { "Package": "tidyverse", - "Version": "1.3.0", + "Version": "1.3.1", "Source": "Repository", "Repository": "CRAN", - "Hash": "bd51be662f359fa99021f3d51e911490" + "Hash": "fc4c72b6ae9bb283416bd59a3303bbab" }, "tinytex": { "Package": "tinytex", - "Version": "0.25", + "Version": "0.35", "Source": "Repository", "Repository": "CRAN", - "Hash": "a4b9662282097d1033c60420dcb83350" + "Hash": "1d7220fe46159fb9f5c99a44354a2bff" + }, + "tzdb": { + "Package": "tzdb", + "Version": "0.2.0", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "5e069fb033daf2317bd628d3100b75c5" }, "utf8": { "Package": "utf8", - "Version": "1.1.4", + "Version": "1.2.2", "Source": "Repository", "Repository": "CRAN", - "Hash": "4a5081acfb7b81a572e4384a7aaf2af1" + "Hash": "c9c462b759a5cc844ae25b5942654d13" + }, + "uuid": { + "Package": "uuid", + "Version": "1.0-3", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "2097822ba5e4440b81a0c7525d0315ce" }, "vctrs": { "Package": "vctrs", - "Version": "0.2.4", + "Version": "0.3.8", "Source": "Repository", "Repository": "CRAN", - "Hash": "6c839a149a30cb4ffc70443efa74c197" + "Hash": "ecf749a1b39ea72bd9b51b76292261f1" }, "viridisLite": { "Package": "viridisLite", - "Version": "0.3.0", + "Version": "0.4.0", "Source": "Repository", "Repository": "CRAN", - "Hash": "ce4f6271baa94776db692f1cb2055bee" + "Hash": "55e157e2aa88161bdb0754218470d204" }, - "whisker": { - "Package": "whisker", - "Version": "0.4", + "vroom": { + "Package": "vroom", + "Version": "1.5.5", "Source": "Repository", "Repository": "CRAN", - "Hash": "ca970b96d894e90397ed20637a0c1bbe" + "Hash": "9c3b3a3f947c7936cea7485349247e5b" + }, + "warp": { + "Package": "warp", + "Version": "0.2.0", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "2982481615756e24e79fee95bdc95daa" }, "withr": { "Package": "withr", - "Version": "2.1.2", + "Version": "2.4.2", "Source": "Repository", "Repository": "CRAN", - "Hash": "aa57ed55ff2df4bea697a07df528993d" + "Hash": "ad03909b44677f930fa156d47d7a3aeb" }, "xfun": { "Package": "xfun", - "Version": "0.16", + "Version": "0.28", "Source": "Repository", "Repository": "CRAN", - "Hash": "b4106139b90981a8bfea9c10bab0baf1" + "Hash": "f7f3a61ab62cd046d307577a8ae12999" }, "xml2": { "Package": "xml2", - "Version": "1.2.5", + "Version": "1.3.2", + "Source": "Repository", + "Repository": "CRAN", + "Hash": "d4d71a75dd3ea9eb5fa28cc21f9585e2" + }, + "xtable": { + "Package": "xtable", + "Version": "1.8-4", "Source": "Repository", "Repository": "CRAN", - "Hash": "13efdfe4c591af40c1c29dc6408281ed" + "Hash": "b8acdf8af494d9ec19ccb2481a9b11c2" }, "yaml": { "Package": "yaml", diff --git a/renv/.gitignore b/renv/.gitignore index 82740ba..993aebf 100644 --- a/renv/.gitignore +++ b/renv/.gitignore @@ -1,3 +1,5 @@ +local/ +lock/ library/ python/ staging/ diff --git a/renv/activate.R b/renv/activate.R index c2fe5e9..304fd90 100644 --- a/renv/activate.R +++ b/renv/activate.R @@ -2,13 +2,27 @@ local({ # the requested version of renv - version <- "0.10.0" + version <- "0.14.0" # the project directory project <- getwd() + # allow environment variable to control activation + activate <- Sys.getenv("RENV_ACTIVATE_PROJECT") + if (!nzchar(activate)) { + + # don't auto-activate when R CMD INSTALL is running + if (nzchar(Sys.getenv("R_INSTALL_PKG"))) + return(FALSE) + + } + + # bail if activation was explicitly disabled + if (tolower(activate) %in% c("false", "f", "0")) + return(FALSE) + # avoid recursion - if (!is.na(Sys.getenv("RENV_R_INITIALIZING", unset = NA))) + if (nzchar(Sys.getenv("RENV_R_INITIALIZING"))) return(invisible(TRUE)) # signal that we're loading renv during R startup @@ -39,12 +53,6 @@ local({ # load bootstrap tools bootstrap <- function(version, library) { - # fix up repos - repos <- getOption("repos") - on.exit(options(repos = repos), add = TRUE) - repos[repos == "@CRAN@"] <- "https://cloud.r-project.org" - options(repos = repos) - # attempt to download renv tarball <- tryCatch(renv_bootstrap_download(version), error = identity) if (inherits(tarball, "error")) @@ -57,34 +65,58 @@ local({ } - renv_bootstrap_download_impl <- function(url, destfile) { + renv_bootstrap_tests_running <- function() { + getOption("renv.tests.running", default = FALSE) + } - mode <- "wb" + renv_bootstrap_repos <- function() { - # https://bugs.r-project.org/bugzilla/show_bug.cgi?id=17715 - fixup <- - Sys.info()[["sysname"]] == "Windows" && - substring(url, 1L, 5L) == "file:" + # check for repos override + repos <- Sys.getenv("RENV_CONFIG_REPOS_OVERRIDE", unset = NA) + if (!is.na(repos)) + return(repos) - if (fixup) - mode <- "w+b" + # if we're testing, re-use the test repositories + if (renv_bootstrap_tests_running()) + return(getOption("renv.tests.repos")) - download.file( - url = url, - destfile = destfile, - mode = mode, - quiet = TRUE + # retrieve current repos + repos <- getOption("repos") + + # ensure @CRAN@ entries are resolved + repos[repos == "@CRAN@"] <- getOption( + "renv.repos.cran", + "https://cloud.r-project.org" ) + # add in renv.bootstrap.repos if set + default <- c(FALLBACK = "https://cloud.r-project.org") + extra <- getOption("renv.bootstrap.repos", default = default) + repos <- c(repos, extra) + + # remove duplicates that might've snuck in + dupes <- duplicated(repos) | duplicated(names(repos)) + repos[!dupes] + } renv_bootstrap_download <- function(version) { - methods <- list( - renv_bootstrap_download_cran_latest, - renv_bootstrap_download_cran_archive, - renv_bootstrap_download_github - ) + # if the renv version number has 4 components, assume it must + # be retrieved via github + nv <- numeric_version(version) + components <- unclass(nv)[[1]] + + methods <- if (length(components) == 4L) { + list( + renv_bootstrap_download_github + ) + } else { + list( + renv_bootstrap_download_cran_latest, + renv_bootstrap_download_cran_archive + ) + } for (method in methods) { path <- tryCatch(method(version), error = identity) @@ -96,21 +128,44 @@ local({ } + renv_bootstrap_download_impl <- function(url, destfile) { + + mode <- "wb" + + # https://bugs.r-project.org/bugzilla/show_bug.cgi?id=17715 + fixup <- + Sys.info()[["sysname"]] == "Windows" && + substring(url, 1L, 5L) == "file:" + + if (fixup) + mode <- "w+b" + + utils::download.file( + url = url, + destfile = destfile, + mode = mode, + quiet = TRUE + ) + + } + renv_bootstrap_download_cran_latest <- function(version) { - # check for renv on CRAN matching this version - db <- as.data.frame(available.packages(), stringsAsFactors = FALSE) - if (!"renv" %in% rownames(db)) - stop("renv is not available on your declared package repositories") + spec <- renv_bootstrap_download_cran_latest_find(version) - entry <- db["renv", ] - if (!identical(entry$Version, version)) - stop("renv is not available on your declared package repositories") + message("* Downloading renv ", version, " ... ", appendLF = FALSE) - message("* Downloading renv ", version, " from CRAN ... ", appendLF = FALSE) + type <- spec$type + repos <- spec$repos info <- tryCatch( - download.packages("renv", destdir = tempdir()), + utils::download.packages( + pkgs = "renv", + destdir = tempdir(), + repos = repos, + type = type, + quiet = TRUE + ), condition = identity ) @@ -119,19 +174,65 @@ local({ return(FALSE) } - message("OK") + # report success and return + message("OK (downloaded ", type, ")") info[1, 2] } + renv_bootstrap_download_cran_latest_find <- function(version) { + + # check whether binaries are supported on this system + binary <- + getOption("renv.bootstrap.binary", default = TRUE) && + !identical(.Platform$pkgType, "source") && + !identical(getOption("pkgType"), "source") && + Sys.info()[["sysname"]] %in% c("Darwin", "Windows") + + types <- c(if (binary) "binary", "source") + + # iterate over types + repositories + for (type in types) { + for (repos in renv_bootstrap_repos()) { + + # retrieve package database + db <- tryCatch( + as.data.frame( + utils::available.packages(type = type, repos = repos), + stringsAsFactors = FALSE + ), + error = identity + ) + + if (inherits(db, "error")) + next + + # check for compatible entry + entry <- db[db$Package %in% "renv" & db$Version %in% version, ] + if (nrow(entry) == 0) + next + + # found it; return spec to caller + spec <- list(entry = entry, type = type, repos = repos) + return(spec) + + } + } + + # if we got here, we failed to find renv + fmt <- "renv %s is not available from your declared package repositories" + stop(sprintf(fmt, version)) + + } + renv_bootstrap_download_cran_archive <- function(version) { name <- sprintf("renv_%s.tar.gz", version) - repos <- getOption("repos") + repos <- renv_bootstrap_repos() urls <- file.path(repos, "src/contrib/Archive/renv", name) destfile <- file.path(tempdir(), name) - message("* Downloading renv ", version, " from CRAN archive ... ", appendLF = FALSE) + message("* Downloading renv ", version, " ... ", appendLF = FALSE) for (url in urls) { @@ -190,7 +291,7 @@ local({ return(FALSE) } - message("Done!") + message("OK") return(destfile) } @@ -222,7 +323,7 @@ local({ } - renv_bootstrap_prefix <- function() { + renv_bootstrap_platform_prefix <- function() { # construct version prefix version <- paste(R.version$major, R.version$minor, sep = ".") @@ -241,8 +342,8 @@ local({ components <- c(prefix, R.version$platform) # include prefix if provided by user - prefix <- Sys.getenv("RENV_PATHS_PREFIX") - if (nzchar(prefix)) + prefix <- renv_bootstrap_platform_prefix_impl() + if (!is.na(prefix) && nzchar(prefix)) components <- c(prefix, components) # build prefix @@ -250,6 +351,152 @@ local({ } + renv_bootstrap_platform_prefix_impl <- function() { + + # if an explicit prefix has been supplied, use it + prefix <- Sys.getenv("RENV_PATHS_PREFIX", unset = NA) + if (!is.na(prefix)) + return(prefix) + + # if the user has requested an automatic prefix, generate it + auto <- Sys.getenv("RENV_PATHS_PREFIX_AUTO", unset = NA) + if (auto %in% c("TRUE", "True", "true", "1")) + return(renv_bootstrap_platform_prefix_auto()) + + # empty string on failure + "" + + } + + renv_bootstrap_platform_prefix_auto <- function() { + + prefix <- tryCatch(renv_bootstrap_platform_os(), error = identity) + if (inherits(prefix, "error") || prefix %in% "unknown") { + + msg <- paste( + "failed to infer current operating system", + "please file a bug report at https://github.com/rstudio/renv/issues", + sep = "; " + ) + + warning(msg) + + } + + prefix + + } + + renv_bootstrap_platform_os <- function() { + + sysinfo <- Sys.info() + sysname <- sysinfo[["sysname"]] + + # handle Windows + macOS up front + if (sysname == "Windows") + return("windows") + else if (sysname == "Darwin") + return("macos") + + # check for os-release files + for (file in c("/etc/os-release", "/usr/lib/os-release")) + if (file.exists(file)) + return(renv_bootstrap_platform_os_via_os_release(file, sysinfo)) + + # check for redhat-release files + if (file.exists("/etc/redhat-release")) + return(renv_bootstrap_platform_os_via_redhat_release()) + + "unknown" + + } + + renv_bootstrap_platform_os_via_os_release <- function(file, sysinfo) { + + # read /etc/os-release + release <- utils::read.table( + file = file, + sep = "=", + quote = c("\"", "'"), + col.names = c("Key", "Value"), + comment.char = "#", + stringsAsFactors = FALSE + ) + + vars <- as.list(release$Value) + names(vars) <- release$Key + + # get os name + os <- tolower(sysinfo[["sysname"]]) + + # read id + id <- "unknown" + for (field in c("ID", "ID_LIKE")) { + if (field %in% names(vars) && nzchar(vars[[field]])) { + id <- vars[[field]] + break + } + } + + # read version + version <- "unknown" + for (field in c("UBUNTU_CODENAME", "VERSION_CODENAME", "VERSION_ID", "BUILD_ID")) { + if (field %in% names(vars) && nzchar(vars[[field]])) { + version <- vars[[field]] + break + } + } + + # join together + paste(c(os, id, version), collapse = "-") + + } + + renv_bootstrap_platform_os_via_redhat_release <- function() { + + # read /etc/redhat-release + contents <- readLines("/etc/redhat-release", warn = FALSE) + + # infer id + id <- if (grepl("centos", contents, ignore.case = TRUE)) + "centos" + else if (grepl("redhat", contents, ignore.case = TRUE)) + "redhat" + else + "unknown" + + # try to find a version component (very hacky) + version <- "unknown" + + parts <- strsplit(contents, "[[:space:]]")[[1L]] + for (part in parts) { + + nv <- tryCatch(numeric_version(part), error = identity) + if (inherits(nv, "error")) + next + + version <- nv[1, 1] + break + + } + + paste(c("linux", id, version), collapse = "-") + + } + + renv_bootstrap_library_root_name <- function(project) { + + # use project name as-is if requested + asis <- Sys.getenv("RENV_PATHS_LIBRARY_ROOT_ASIS", unset = "FALSE") + if (asis) + return(basename(project)) + + # otherwise, disambiguate based on project's path + id <- substring(renv_bootstrap_hash_text(project), 1L, 8L) + paste(basename(project), id, sep = "-") + + } + renv_bootstrap_library_root <- function(project) { path <- Sys.getenv("RENV_PATHS_LIBRARY", unset = NA) @@ -257,10 +504,13 @@ local({ return(path) path <- Sys.getenv("RENV_PATHS_LIBRARY_ROOT", unset = NA) - if (!is.na(path)) - return(file.path(path, basename(project))) + if (!is.na(path)) { + name <- renv_bootstrap_library_root_name(project) + return(file.path(path, name)) + } - file.path(project, "renv/library") + prefix <- renv_bootstrap_profile_prefix() + paste(c(project, prefix, "renv/library"), collapse = "/") } @@ -279,8 +529,8 @@ local({ paste("renv", loadedversion, sep = "@") fmt <- paste( - "renv %1$s was loaded from project library, but renv %2$s is recorded in lockfile.", - "Use `renv::record(\"%3$s\")` to record this version in the lockfile.", + "renv %1$s was loaded from project library, but this project is configured to use renv %2$s.", + "Use `renv::record(\"%3$s\")` to record renv %1$s in the lockfile.", "Use `renv::restore(packages = \"renv\")` to install renv %2$s into the project library.", sep = "\n" ) @@ -292,6 +542,16 @@ local({ } + renv_bootstrap_hash_text <- function(text) { + + hashfile <- tempfile("renv-hash-") + on.exit(unlink(hashfile), add = TRUE) + + writeLines(text, con = hashfile) + tools::md5sum(hashfile) + + } + renv_bootstrap_load <- function(project, libpath, version) { # try to load renv from the project library @@ -307,12 +567,69 @@ local({ TRUE } + + renv_bootstrap_profile_load <- function(project) { + + # if RENV_PROFILE is already set, just use that + profile <- Sys.getenv("RENV_PROFILE", unset = NA) + if (!is.na(profile) && nzchar(profile)) + return(profile) + + # check for a profile file (nothing to do if it doesn't exist) + path <- file.path(project, "renv/local/profile") + if (!file.exists(path)) + return(NULL) + + # read the profile, and set it if it exists + contents <- readLines(path, warn = FALSE) + if (length(contents) == 0L) + return(NULL) + + # set RENV_PROFILE + profile <- contents[[1L]] + if (nzchar(profile)) + Sys.setenv(RENV_PROFILE = profile) + + profile + + } + + renv_bootstrap_profile_prefix <- function() { + profile <- renv_bootstrap_profile_get() + if (!is.null(profile)) + return(file.path("renv/profiles", profile)) + } + + renv_bootstrap_profile_get <- function() { + profile <- Sys.getenv("RENV_PROFILE", unset = "") + renv_bootstrap_profile_normalize(profile) + } + + renv_bootstrap_profile_set <- function(profile) { + profile <- renv_bootstrap_profile_normalize(profile) + if (is.null(profile)) + Sys.unsetenv("RENV_PROFILE") + else + Sys.setenv(RENV_PROFILE = profile) + } + + renv_bootstrap_profile_normalize <- function(profile) { + + if (is.null(profile) || profile %in% c("", "default")) + return(NULL) + + profile + + } + + # load the renv profile, if any + renv_bootstrap_profile_load(project) # construct path to library root root <- renv_bootstrap_library_root(project) # construct library prefix for platform - prefix <- renv_bootstrap_prefix() + prefix <- renv_bootstrap_platform_prefix() # construct full libpath libpath <- file.path(root, prefix) @@ -321,12 +638,22 @@ local({ if (renv_bootstrap_load(project, libpath, version)) return(TRUE) - # load failed; attempt to bootstrap + # load failed; inform user we're about to bootstrap + prefix <- paste("# Bootstrapping renv", version) + postfix <- paste(rep.int("-", 77L - nchar(prefix)), collapse = "") + header <- paste(prefix, postfix) + message(header) + + # perform bootstrap bootstrap(version, libpath) + # exit early if we're just testing bootstrap + if (!is.na(Sys.getenv("RENV_BOOTSTRAP_INSTALL_ONLY", unset = NA))) + return(TRUE) + # try again to load if (requireNamespace("renv", lib.loc = libpath, quietly = TRUE)) { - message("Successfully installed and loaded renv ", version, ".") + message("* Successfully installed and loaded renv ", version, ".") return(renv::load()) } From 93e81c9568859571a103f4da37fa33f471354dde Mon Sep 17 00:00:00 2001 From: juanitorduz Date: Sun, 7 Nov 2021 13:46:22 +0100 Subject: [PATCH 3/3] rename folder --- R/2020-07/rep_for_btsa.R | 474 ----------------------------- R/2020-07/rep_original.r | 640 --------------------------------------- 2 files changed, 1114 deletions(-) delete mode 100644 R/2020-07/rep_for_btsa.R delete mode 100755 R/2020-07/rep_original.r diff --git a/R/2020-07/rep_for_btsa.R b/R/2020-07/rep_for_btsa.R deleted file mode 100644 index 745fc44..0000000 --- a/R/2020-07/rep_for_btsa.R +++ /dev/null @@ -1,474 +0,0 @@ -## Replication Code for -# A. Abadie, A. Diamond, and J. Hainmueller. 2014. -# Comparative Politics and the Synthetic Control Method -# American Journal of Political Science. -# Original source for the code and the data: -# https://dataverse.harvard.edu/dataset.xhtml?persistentId=doi:10.7910/DVN/24714 -# This version for the Berlin Time Series Analysis Meet-up - -library(foreign) -library("Synth") -library(xtable) -library(ggplot2) -library(scales) -library(tidyverse) -# Load Data -d <- read.dta("repgermany.dta") - - -# Data set-up - countrytosynthetise <- 7 - dataprep.out <- - dataprep( - foo = d, - predictors = c("gdp","trade","infrate"), - dependent = "gdp", - unit.variable = 1, - time.variable = 3, - special.predictors = list( - list("industry", 1971:1980, c("mean")), - list("schooling",c(1970,1975), c("mean")), - list("invest70" ,1980, c("mean")) - ), - treatment.identifier = countrytosynthetise, - controls.identifier = unique(d$index)[-which(unique(d$index) %in% countrytosynthetise)], - time.predictors.prior = 1971:1989, - time.optimize.ssr = 1971:1990, - unit.names.variable = 2, - time.plot = 1960:2003 - ) - -# fit training model -synth.out <- - synth( - data.prep.obj=dataprep.out, - Margin.ipop=.005,Sigf.ipop=7,Bound.ipop=6 - ) - -# data prep for main model -dataprep.out <- - dataprep( - foo = d, - predictors = c("gdp","trade","infrate"), - dependent = "gdp", - unit.variable = 1, - time.variable = 3, - special.predictors = list( - list("industry" ,1981:1990, c("mean")), - list("schooling",c(1980,1985), c("mean")), - list("invest80" ,1980, c("mean")) - ), - treatment.identifier = countrytosynthetise, - controls.identifier = unique(d$index)[-which(unique(d$index) %in% countrytosynthetise)], - time.predictors.prior = 1960:1989, - time.optimize.ssr = 1960:1990, - unit.names.variable = 2, - time.plot = 1960:2003 -) - -# fit main model with v from training model -synth.out <- synth( - data.prep.obj=dataprep.out, - custom.v=as.numeric(synth.out$solution.v) - ) - -#### Table 2 -synth.tables <- synth.tab( - dataprep.res = dataprep.out, - synth.res = synth.out - ); synth.tables - -# Replace means for OECD sample (computed externally using proper pop weighting) -synth.tables$tab.pred[,3] <- c(8021.1,31.9,7.4,34.2,44.1,25.9) -colnames(synth.tables$tab.pred)[3] <- "Rest of OECD Sample" -rownames(synth.tables$tab.pred) <- c("GDP per-capita","Trade openness", - "Inflation rate","Industry share", - "Schooling","Investment rate") - -xtable(round(synth.tables$tab.pred,1),digits=1) - -#### Table 1 -# synth weights -tab1 <- data.frame(synth.tables$tab.w) -tab1[,1] <- round(tab1[,1],2) -# regression weights -X0 <- cbind(1,t(dataprep.out$X0)) -X1 <- as.matrix(c(1,dataprep.out$X1)) -W <- X0%*%solve(t(X0)%*%X0)%*%X1 -Wdat <- data.frame(unit.numbers=as.numeric(rownames(X0)), - regression.w=round(W,2)) -tab1 <- merge(tab1,Wdat,by="unit.numbers") -tab1 <- tab1[order(tab1[,3]),] - - -# Table 1 -xtable(cbind(tab1[1:9,c(3,2,4)], - tab1[10:18,c(3,2,4)] - ) - ) - -#### Figures 1, 2 and 3 -data.to.plot <- data.frame(year = rownames(dataprep.out$Y1plot), - west.germany = dataprep.out$Y1plot[, 1]) - -data.to.plot <- as.data.frame(d[,c("gdp","country", "year" )]) -data.to.plot$west.germany <- ifelse(data.to.plot$country == "West Germany", "West Germany", "Not West Germany") -data.to.plot$west.germany.size <- ifelse(data.to.plot$west.germany == "West Germany", "1","2") -graphcolors <- c("#999999", "#ef8a62") - - -table(data.to.plot$country, data.to.plot$west.germany.size) - - -## PLOT WITH ONLY WEST GERMANY -ggplot(data = data.to.plot[data.to.plot$west.germany == "West Germany",], aes(x = year, y = gdp, - group = country, - color = west.germany, - linetype = west.germany.size)) + - geom_line(size = 1.5) + - labs(color = "West Germany?")+ - guides(linetype = FALSE, color = FALSE) + - theme_bw(base_size = 15) + - theme(legend.position = "bottom") + - geom_vline(xintercept = 1990, - color = "#88419d") + - annotate("text", x = 1985, - label="\nReunification", y = 1000, - colour = "#88419d", alpha = 0.8, angle = 0) + - scale_y_continuous(name="per-capita GDP (PPP, 2002 USD)", - labels = comma_format(big.mark = ",", - decimal.mark = ".")) + - scale_color_manual(values = c("#ef8a62")) + - coord_cartesian(ylim = c(0, 30000)) #+ - #ggsave(paste0("figures/germany.pdf"), - # width = 15, height = 10, units = 'cm') - - - -## PLOT WITH WEST GERMANY AND ALL THE OTHER COUNTRIES -ggplot(data = data.to.plot, aes(x = year, y = gdp, - group = country, - color = west.germany, - linetype = west.germany.size)) + - geom_line(size = 1.5) + - labs(color = "West Germany?")+ - guides(linetype = FALSE) + - theme_bw(base_size = 10) + - theme(legend.position = "bottom") + - geom_vline(xintercept = 1990, - color = "#88419d") + - annotate("text", x = 1985, - label="Reunification", y = 1000, - colour = "#88419d", alpha = 0.8, angle = 0) + - scale_y_continuous(name="per-capita GDP (PPP, 2002 USD)", - labels = comma_format(big.mark = ",", - decimal.mark = ".")) + - scale_color_manual(values = graphcolors) + - coord_cartesian(ylim = c(0, 30000)) #+ - #ggsave(paste0("figures/germanynotgermany.pdf"), - # width = 15, height = 10, units = 'cm') - - -## PLOT WITH WEST GERMANY AND THE AVERAGE FROM OTHER COUNTRIES -data.to.plot.average <- data.to.plot %>% - group_by(west.germany, year, west.germany.size) %>% - summarise(GDP = mean(gdp)) - -ggplot(data = data.to.plot.average, - aes(x = year, y = GDP, - group = west.germany, - color = west.germany)) + - geom_line(size = 1.5) + - labs(color = "West Germany?")+ - theme_bw(base_size = 10) + - theme(legend.position = "bottom") + - geom_vline(xintercept = 1990, - color = "#88419d") + - annotate("text", x = 1985, - label="Reunification", y = 1000, - colour = "#88419d", alpha = 0.8, angle = 0) + - scale_y_continuous(name="per-capita GDP (PPP, 2002 USD)", - labels = comma_format(big.mark = ",", - decimal.mark = ".")) + - scale_color_manual(values = graphcolors) + - coord_cartesian(ylim = c(0, 30000)) #+ - #ggsave(paste0("figures/germanyaverage.pdf"), - # width = 15, height = 10, units = 'cm') - - - -data.to.plot.average$west.germany[data.to.plot.average$west.germany == "Not West Germany"] <- "Average of Other OECD Countries" - -ggplot(data = data.to.plot.average, - aes(x = year, y = GDP, - group = west.germany, - color = west.germany)) + - geom_line(size = 1.5) + - labs(color = "West Germany?")+ - theme_bw(base_size = 10) + - theme(legend.position = "bottom") + - geom_vline(xintercept = 1990, - color = "#88419d") + - annotate("text", x = 1985, - label="Reunification", y = 1000, - colour = "#88419d", alpha = 0.8, angle = 0) + - scale_y_continuous(name="per-capita GDP (PPP, 2002 USD)", - labels = comma_format(big.mark = ",", - decimal.mark = ".")) + - scale_color_manual(values = graphcolors) + - coord_cartesian(ylim = c(0, 30000)) #+ - #ggsave(paste0("figures/germanyaverage2.pdf"), - # width = 15, height = 10, units = 'cm') - - - - -# SYNTHETIC UNIT - Weights -# Figure 4 -data.to.plot.synthetic <- data.to.plot.average - -synthY0 <- (dataprep.out$Y0%*%synth.out$solution.w) -data.to.plot.synthetic$GDP[data.to.plot.synthetic$west.germany == "Average of Other OECD Countries"] <- (dataprep.out$Y0%*%synth.out$solution.w) -data.to.plot.synthetic$west.germany[data.to.plot.synthetic$west.germany == "Average of Other OECD Countries"] <- "Synthetic West Germany" - - -ggplot(data = data.to.plot.synthetic, - aes(x = year, y = GDP, - group = west.germany, - color = west.germany)) + - geom_line(size = 1.5) + - labs(color = "West Germany?")+ - theme_bw(base_size = 10) + - theme(legend.position = "bottom") + - geom_vline(xintercept = 1990, - color = "#88419d") + - annotate("text", x = 1985, - label="Reunification", y = 1000, - colour = "#88419d", alpha = 0.8, angle = 0) + - scale_y_continuous(name="per-capita GDP (PPP, 2002 USD)", - labels = comma_format(big.mark = ",", - decimal.mark = ".")) + - scale_color_manual(values = graphcolors) + - coord_cartesian(ylim = c(0, 30000)) #+ - #ggsave(paste0("figures/germanysynth.pdf"), - # width = 15, height = 10, units = 'cm') - - - -# SYNTHETIC UNIT - Regression -# Figure 8 -data.to.plot.synth.reg <- data.to.plot.average - -synthY0.reg <- (dataprep.out$Y0%*%W) -data.to.plot.synth.reg$GDP[data.to.plot.synth.reg$west.germany == "Average of Other OECD Countries"] <- (dataprep.out$Y0%*%W) -data.to.plot.synth.reg$west.germany[data.to.plot.synth.reg$west.germany == "Average of Other OECD Countries"] <- "Synthetic West Germany" - - -ggplot(data = data.to.plot.synth.reg, - aes(x = year, y = GDP, - group = west.germany, - color = west.germany)) + - geom_line(size = 1.5) + - labs(color = "West Germany?")+ - theme_bw(base_size = 10) + - theme(legend.position = "bottom") + - geom_vline(xintercept = 1990, - color = "#88419d") + - annotate("text", x = 1985, - label="Reunification", y = 1000, - colour = "#88419d", alpha = 0.8, angle = 0) + - scale_y_continuous(name="per-capita GDP (PPP, 2002 USD)", - labels = comma_format(big.mark = ",", - decimal.mark = ".")) + - scale_color_manual(values = graphcolors) + - coord_cartesian(ylim = c(0, 30000)) #+ - #ggsave(paste0("figures/germanysynthreg.pdf"), - # width = 15, height = 10, units = 'cm') - - - -# Gap -# Figure 5 -gap <- dataprep.out$Y1-(dataprep.out$Y0%*%synth.out$solution.w) - -data.to.plot.gap <- data.to.plot.synthetic -data.to.plot.gap$GDP[data.to.plot.gap$west.germany == "West Germany"] <- dataprep.out$Y1-(dataprep.out$Y0%*%synth.out$solution.w) - -data.to.plot.gap$GDP[data.to.plot.gap$west.germany == "Synthetic West Germany"] <- 0 - - -ggplot(data = data.to.plot.gap, - aes(x = year, y = GDP, - group = west.germany, - color = west.germany)) + - geom_line(size = 1.5) + - labs(color = "West Germany?")+ - theme_bw(base_size = 10) + - theme(legend.position = "bottom") + - geom_vline(xintercept = 1990, - color = "#88419d") + - annotate("text", x = 1985, - label="Reunification", y = 1000, - colour = "#88419d", alpha = 0.8, angle = 0) + - scale_y_continuous(name="Gap in per-capita GDP (PPP, 2002 USD)", - labels = comma_format(big.mark = ",", - decimal.mark = ".")) + - scale_color_manual(values = graphcolors) + - coord_cartesian(ylim = c(-5000, 5000)) #+ - #ggsave(paste0("figures/germanygap.pdf"), - # width = 15, height = 10, units = 'cm') - - - - -# loop across control units -storegaps <- - matrix(NA, - length(1960:2003), - length(unique(d$index))-1 - ) -rownames(storegaps) <- 1960:2003 -i <- 1 -co <- unique(d$index) - -for(k in unique(d$index)[-7]){ - - # data prep for training model - dataprep.out <- - dataprep( - foo = d, - predictors = c("gdp","trade","infrate"), - dependent = "gdp", - unit.variable = 1, - time.variable = 3, - special.predictors = list( - list("industry",1971:1980, c("mean")), - list("schooling" ,c(1970,1975), c("mean")), - list("invest70" ,1980, c("mean")) - ), - treatment.identifier = k, - controls.identifier = co[-which(co==k)], - time.predictors.prior = 1971:1989, - time.optimize.ssr = 1971:1990, - unit.names.variable = 2, - time.plot = 1960:2003 - ) - - # fit training model - synth.out <- - synth( - data.prep.obj=dataprep.out, - Margin.ipop=.005,Sigf.ipop=7,Bound.ipop=6 - ) - - # data prep for main model -dataprep.out <- - dataprep( - foo = d, - predictors = c("gdp","trade","infrate"), - dependent = "gdp", - unit.variable = 1, - time.variable = 3, - special.predictors = list( - list("industry" ,1981:1990, c("mean")), - list("schooling",c(1980,1985), c("mean")), - list("invest80" ,1980, c("mean")) - ), - treatment.identifier = k, - controls.identifier = co[-which(co==k)], - time.predictors.prior = 1960:1990, - time.optimize.ssr = 1960:1989, - unit.names.variable = 2, - time.plot = 1960:2003 - ) - -# fit main model -synth.out <- synth( - data.prep.obj=dataprep.out, - custom.v=as.numeric(synth.out$solution.v) - ) - - storegaps[,i] <- - dataprep.out$Y1- - (dataprep.out$Y0%*%synth.out$solution.w) - i <- i + 1 -} # close loop over control units -d <- d[order(d$index,d$year),] -colnames(storegaps) <- unique(d$country)[-7] -storegaps <- cbind(gap,storegaps) -colnames(storegaps)[1] <- c("West Germany") - - - -# Placebo gaps -# Figure 6 -sgaps <- as.data.frame(storegaps) -sgaps$year <- rownames(sgaps) -sgaps <- reshape2::melt(sgaps,id.vars = "year") - - -sgaps$west.germany <- ifelse(sgaps$variable == "West Germany", "West Germany", "Not West Germany") -sgaps$west.germany.size <- ifelse(sgaps$west.germany == "West Germany", "1","2") - -sgaps$year <- as.numeric(sgaps$year) - -data.to.plot$west.germany <- ifelse(data.to.plot$country == "West Germany", "West Germany", "Not West Germany") -ggplot(data = sgaps, aes(x = year, y = value, - group = variable, - color = west.germany, - linetype = west.germany.size)) + - geom_line() + - labs(color = "West Germany?")+ - guides(linetype = FALSE) + - theme_bw(base_size = 10) + - theme(legend.position = "bottom") + - geom_vline(xintercept = 1990, - color = "#88419d") + - annotate("text", x = 1985, - label="Reunification", y = 4000, - colour = "#88419d", alpha = 0.8, angle = 0) + - scale_y_continuous(name="Gap in per-capita GDP (PPP, 2002 USD)", - labels = comma_format(big.mark = ",", - decimal.mark = ".")) + - scale_color_manual(values = graphcolors) + - coord_cartesian(ylim = c(-5000, 5000)) #+ - #ggsave(paste0("figures/allrmspe.pdf"), - # width = 15, height = 10, units = 'cm') - - - - - - - -# Figure 7 -# compute ratio of post-reunification RMSPE -# to pre-reunification RMSPE -rmse <- function(x){sqrt(mean(x^2))} -preloss <- apply(storegaps[1:30,],2,rmse) -postloss <- apply(storegaps[31:44,],2,rmse) -prepost <- sort(postloss/preloss) -prepost <- t(t(prepost)) - -prepost <- data.frame(Country = rownames(prepost), RMSPE = prepost) -rownames(prepost) <- NULL -prepost$Country <- factor(prepost$Country, - levels = prepost$Country[order(prepost$RMSPE)]) - - -prepost$west.germany <- ifelse(prepost$Country== "West Germany", "West Germany", "Not West Germany") - -ggplot(prepost, aes(x = Country, y = RMSPE)) + - geom_point(stat = 'identity', size = 6, aes(color = west.germany)) + - theme_bw(base_size = 10) + - coord_flip() + - ylab("Post-Period RMSE / Pre-Period RMSE") + - scale_color_manual(values = graphcolors) + - guides(color = FALSE) #+ - #ggsave(paste0("figures/RMSPERatio.pdf"), - # width = 15, height = 10, units = 'cm') - - - - - diff --git a/R/2020-07/rep_original.r b/R/2020-07/rep_original.r deleted file mode 100755 index f1b77a1..0000000 --- a/R/2020-07/rep_original.r +++ /dev/null @@ -1,640 +0,0 @@ -## Replication Code for -# A. Abadie, A. Diamond, and J. Hainmueller. 2014. -# Comparative Politics and the Synthetic Control Method -# American Journal of Political Science. - -rm(list=ls()) -library(foreign) -library("Synth") -library(xtable) - -# Load Data -d <- read.dta("repgermany.dta") - -## Table 1 & 2, Figure 1, 2, & 3 - -## pick v by cross-validation -# data setup for training model -dataprep.out <- - dataprep( - foo = d, - predictors = c("gdp","trade","infrate"), - dependent = "gdp", - unit.variable = 1, - time.variable = 3, - special.predictors = list( - list("industry", 1971:1980, c("mean")), - list("schooling",c(1970,1975), c("mean")), - list("invest70" ,1980, c("mean")) - ), - treatment.identifier = 7, - controls.identifier = unique(d$index)[-7], - time.predictors.prior = 1971:1980, - time.optimize.ssr = 1981:1990, - unit.names.variable = 2, - time.plot = 1960:2003 - ) - -# fit training model -synth.out <- - synth( - data.prep.obj=dataprep.out, - Margin.ipop=.005,Sigf.ipop=7,Bound.ipop=6 - ) - -# data prep for main model -dataprep.out <- - dataprep( - foo = d, - predictors = c("gdp","trade","infrate"), - dependent = "gdp", - unit.variable = 1, - time.variable = 3, - special.predictors = list( - list("industry" ,1981:1990, c("mean")), - list("schooling",c(1980,1985), c("mean")), - list("invest80" ,1980, c("mean")) - ), - treatment.identifier = 7, - controls.identifier = unique(d$index)[-7], - time.predictors.prior = 1981:1990, - time.optimize.ssr = 1960:1989, - unit.names.variable = 2, - time.plot = 1960:2003 -) - -# fit main model with v from training model -synth.out <- synth( - data.prep.obj=dataprep.out, - custom.v=as.numeric(synth.out$solution.v) - ) - -#### Table 2 -synth.tables <- synth.tab( - dataprep.res = dataprep.out, - synth.res = synth.out - ); synth.tables - -# Replace means for OECD sample (computed externally using proper pop weighting) -synth.tables$tab.pred[,3] <- c(8021.1,31.9,7.4,34.2,44.1,25.9) -colnames(synth.tables$tab.pred)[3] <- "Rest of OECD Sample" -rownames(synth.tables$tab.pred) <- c("GDP per-capita","Trade openness", - "Inflation rate","Industry share", - "Schooling","Investment rate") - -xtable(round(synth.tables$tab.pred,1),digits=1) - -#### Table 1 -# synth weights -tab1 <- data.frame(synth.tables$tab.w) -tab1[,1] <- round(tab1[,1],2) -# regression weights -X0 <- cbind(1,t(dataprep.out$X0)) -X1 <- as.matrix(c(1,dataprep.out$X1)) -W <- X0%*%solve(t(X0)%*%X0)%*%X1 -Wdat <- data.frame(unit.numbers=as.numeric(rownames(X0)), - regression.w=round(W,2)) -tab1 <- merge(tab1,Wdat,by="unit.numbers") -tab1 <- tab1[order(tab1[,3]),] - -xtable(cbind(tab1[1:9,c(3,2,4)], - tab1[10:18,c(3,2,4)] - ) - ) - -#### Figure 1: Trends in Per-Capita GDP: West Germany vs. Rest of the OECD Sample -Text.height <- 23000 -Cex.set <- .8 -#pdf(file = "ger_vs_oecd.pdf", width = 5.5, height = 5.5, family = "Times",pointsize = 12) -plot(1960:2003,dataprep.out$Y1plot, - type="l",ylim=c(0,33000),col="black",lty="solid", - ylab ="per-capita GDP (PPP, 2002 USD)", - xlab ="year", - xaxs = "i", yaxs = "i", - lwd=2 - ) -lines(1960:2003,aggregate(d[,c("gdp")],by=list(d$year),mean,na.rm=T)[,2] - ,col="black",lty="dashed",lwd=2) # mean 2 -abline(v=1990,lty="dotted") -legend(x="bottomright", - legend=c("West Germany","rest of the OECD sample") - ,lty=c("solid","dashed"),col=c("black","black") - ,cex=.8,bg="white",lwd=c(2,2)) -arrows(1987,Text.height,1989,Text.height,col="black",length=.1) -text(1982.5,Text.height,"reunification",cex=Cex.set) -#dev.off() - -#### Figure 2: Trends in Per-Capita GDP: West Germany vs. Synthetic West Germany -#pdf(file = "ger_vs_synthger2.pdf", width = 5.5, height = 5.5, family = "Times",pointsize = 12) -synthY0 <- (dataprep.out$Y0%*%synth.out$solution.w) -plot(1960:2003,dataprep.out$Y1plot, - type="l",ylim=c(0,33000),col="black",lty="solid", - ylab ="per-capita GDP (PPP, 2002 USD)", - xlab ="year", - xaxs = "i", yaxs = "i", - lwd=2 - ) -lines(1960:2003,synthY0,col="black",lty="dashed",lwd=2) -abline(v=1990,lty="dotted") -legend(x="bottomright", - legend=c("West Germany","synthetic West Germany") - ,lty=c("solid","dashed"),col=c("black","black") - ,cex=.8,bg="white",lwd=c(2,2)) -arrows(1987,Text.height,1989,Text.height,col="black",length=.1) -text(1982.5,Text.height,"reunification",cex=Cex.set) -#dev.off() - -### Figure 3: Per-Capita GDP Gap Between West Germany and Synthetic West Germany -#pdf(file = "ger_vs_synthger_gaps2.pdf", width = 5.5, height = 5.5, family = "Times",pointsize = 12) -gap <- dataprep.out$Y1-(dataprep.out$Y0%*%synth.out$solution.w) -plot(1960:2003,gap, - type="l",ylim=c(-4500,4500),col="black",lty="solid", - ylab =c("gap in per-capita GDP (PPP, 2002 USD)"), - xlab ="year", - xaxs = "i", yaxs = "i", - lwd=2 - ) -abline(v=1990,lty="dotted") -abline(h=0,lty="dotted") -arrows(1987,1000,1989,1000,col="black",length=.1) -text(1982.5,1000,"reunification",cex=Cex.set) -#dev.off() - -### Figure 4: Placebo Reunification 1975 - Trends in Per-Capita GDP: West Germany vs. Synthetic West Germany - -# data prep for training model -dataprep.out <- - dataprep( - foo = d, - predictors = c("gdp","trade","infrate"), - dependent = "gdp", - unit.variable = 1, - time.variable = 3, - special.predictors = list( - list("industry",1971, c("mean")), - list("schooling",c(1960,1965), c("mean")), - list("invest60" ,1980, c("mean")) - ), - treatment.identifier = 7, - controls.identifier = unique(d$index)[-7], - time.predictors.prior = 1960:1964, - time.optimize.ssr = 1965:1975, - unit.names.variable = 2, - time.plot = 1960:1990 - ) - -# fit training model -synth.out <- synth( - data.prep.obj=dataprep.out, - Margin.ipop=.005,Sigf.ipop=7,Bound.ipop=6 -) - - -# data prep for main model -dataprep.out <- - dataprep( - foo = d, - predictors = c("gdp","trade","infrate"), - dependent = "gdp", - unit.variable = 1, - time.variable = 3, - special.predictors = list( - list("industry" ,1971:1975, c("mean")), - list("schooling",c(1970,1975), c("mean")), - list("invest70" ,1980, c("mean")) - ), - treatment.identifier = 7, - controls.identifier = unique(d$index)[-7], - time.predictors.prior = 1965:1975, - time.optimize.ssr = 1960:1975, - unit.names.variable = 2, - time.plot = 1960:1990 - ) - -# fit main model -synth.out <- synth( - data.prep.obj=dataprep.out, - custom.v=as.numeric(synth.out$solution.v) -) - -Cex.set <- 1 -#pdf(file = "2intimeplacebo1975.pdf", width = 5.5, height = 5.5, family = "Times",pointsize = 12) -plot(1960:1990,dataprep.out$Y1plot, - type="l",ylim=c(0,33000),col="black",lty="solid", - ylab ="per-capita GDP (PPP, 2002 USD)", - xlab ="year", - xaxs = "i", yaxs = "i", - lwd=2 - ) -lines(1960:1990,(dataprep.out$Y0%*%synth.out$solution.w),col="black",lty="dashed",lwd=2) -abline(v=1975,lty="dotted") -legend(x="bottomright", - legend=c("West Germany","synthetic West Germany") - ,lty=c("solid","dashed"),col=c("black","black") - ,cex=.8,bg="white",lwd=c(2,2)) -arrows(1973,20000,1974.5,20000,col="black",length=.1) -text(1967.5,20000,"placebo reunification",cex=Cex.set) -#dev.off() - - -### Figure 5: Ratio of post-reunification RMSPE to pre-reunification RMSPE: West Germany and control countries. - -# loop across control units -storegaps <- - matrix(NA, - length(1960:2003), - length(unique(d$index))-1 - ) -rownames(storegaps) <- 1960:2003 -i <- 1 -co <- unique(d$index) - -for(k in unique(d$index)[-7]){ - - # data prep for training model - dataprep.out <- - dataprep( - foo = d, - predictors = c("gdp","trade","infrate"), - dependent = "gdp", - unit.variable = 1, - time.variable = 3, - special.predictors = list( - list("industry",1971:1980, c("mean")), - list("schooling" ,c(1970,1975), c("mean")), - list("invest70" ,1980, c("mean")) - ), - treatment.identifier = k, - controls.identifier = co[-which(co==k)], - time.predictors.prior = 1971:1980, - time.optimize.ssr = 1981:1990, - unit.names.variable = 2, - time.plot = 1960:2003 - ) - - # fit training model - synth.out <- - synth( - data.prep.obj=dataprep.out, - Margin.ipop=.005,Sigf.ipop=7,Bound.ipop=6 - ) - - # data prep for main model -dataprep.out <- - dataprep( - foo = d, - predictors = c("gdp","trade","infrate"), - dependent = "gdp", - unit.variable = 1, - time.variable = 3, - special.predictors = list( - list("industry" ,1981:1990, c("mean")), - list("schooling",c(1980,1985), c("mean")), - list("invest80" ,1980, c("mean")) - ), - treatment.identifier = k, - controls.identifier = co[-which(co==k)], - time.predictors.prior = 1981:1990, - time.optimize.ssr = 1960:1989, - unit.names.variable = 2, - time.plot = 1960:2003 - ) - -# fit main model -synth.out <- synth( - data.prep.obj=dataprep.out, - custom.v=as.numeric(synth.out$solution.v) - ) - - storegaps[,i] <- - dataprep.out$Y1- - (dataprep.out$Y0%*%synth.out$solution.w) - i <- i + 1 -} # close loop over control units -d <- d[order(d$index,d$year),] -colnames(storegaps) <- unique(d$country)[-7] -storegaps <- cbind(gap,storegaps) -colnames(storegaps)[1] <- c("West Germany") - -# compute ratio of post-reunification RMSPE -# to pre-reunification RMSPE -rmse <- function(x){sqrt(mean(x^2))} -preloss <- apply(storegaps[1:30,],2,rmse) -postloss <- apply(storegaps[31:44,],2,rmse) - -#pdf("2ratio_post_to_preperiod_rmse2a.pdf") -dotchart(sort(postloss/preloss), - xlab="Post-Period RMSE / Pre-Period RMSE", - pch=19) -#dev.off() - -### Figure 6: Leave-one-out distribution of the synthetic control for West Germany - -# loop over leave one outs -storegaps <- - matrix(NA, - length(1960:2003), - 5) -colnames(storegaps) <- c(1,3,9,12,14) -co <- unique(d$index)[-7] - -for(k in 1:5){ - -# data prep for training model -omit <- c(1,3,9,12,14)[k] - dataprep.out <- - dataprep( - foo = d, - predictors = c("gdp","trade","infrate"), - dependent = "gdp", - unit.variable = 1, - time.variable = 3, - special.predictors = list( - list("industry",1971:1980, c("mean")), - list("schooling" ,c(1970,1975), c("mean")), - list("invest70" ,1980, c("mean")) - ), - treatment.identifier = 7, - controls.identifier = co[-which(co==omit)], - time.predictors.prior = 1971:1980, - time.optimize.ssr = 1981:1990, - unit.names.variable = 2, - time.plot = 1960:2003 - ) - - # fit training model - synth.out <- synth( - data.prep.obj=dataprep.out, - Margin.ipop=.005,Sigf.ipop=7,Bound.ipop=6 - ) - -# data prep for main model -dataprep.out <- - dataprep( - foo = d, - predictors = c("gdp","trade","infrate"), - dependent = "gdp", - unit.variable = 1, - time.variable = 3, - special.predictors = list( - list("industry" ,1981:1990, c("mean")), - list("schooling",c(1980,1985), c("mean")), - list("invest80" ,1980, c("mean")) - ), - treatment.identifier = 7, - controls.identifier = co[-which(co==omit)], - time.predictors.prior = 1981:1990, - time.optimize.ssr = 1960:1989, - unit.names.variable = 2, - time.plot = 1960:2003 - ) - - # fit main model - synth.out <- synth( - data.prep.obj=dataprep.out, - custom.v=as.numeric(synth.out$solution.v) - ) - storegaps[,k] <- (dataprep.out$Y0%*%synth.out$solution.w) -} # close loop over leave one outs - -Text.height <- 23000 -Cex.set <- .8 -#pdf(file = "1jackknife2.pdf", width = 5.5, height = 5.5, family = "Times",pointsize = 12) -plot(1960:2003,dataprep.out$Y1plot, - type="l",ylim=c(0,33000),col="black",lty="solid", - ylab ="per-capita GDP (PPP, 2002 USD)", - xlab ="year", - xaxs = "i", yaxs = "i",lwd=2 - ) - -abline(v=1990,lty="dotted") -arrows(1987,23000,1989,23000,col="black",length=.1) - for(i in 1:5){ - lines(1960:2003,storegaps[,i],col="darkgrey",lty="solid") - } -lines(1960:2003,synthY0,col="black",lty="dashed",lwd=2) -lines(1960:2003,dataprep.out$Y1plot,col="black",lty="solid",lwd=2) -text(1982.5,23000,"reunification",cex=.8) -legend(x="bottomright", - legend=c("West Germany", - "synthetic West Germany", - "synthetic West Germany (leave-one-out)") - ,lty=c("solid","dashed","solid"), - col=c("black","black","darkgrey") - ,cex=.8,bg="white",lwdc(2,2,1)) -#dev.off() - - -### Table 3: Synthetic Weights from Combinations of Control Countries -rm(list=ls()) -library(gtools) -library(kernlab) - -# data prep for training model -d <- read.dta("repgermany.dta") -dataprep.out <- - dataprep( - foo = d, - predictors = c("gdp","trade","infrate"), - dependent = "gdp", - unit.variable = 1, - time.variable = 3, - special.predictors = list( - list("industry", 1971:1980, c("mean")), - list("schooling",c(1970,1975), c("mean")), - list("invest70" ,1980, c("mean")) - ), - treatment.identifier = 7, - controls.identifier = unique(d$index)[-7], - time.predictors.prior = 1971:1980, - time.optimize.ssr = 1981:1990, - unit.names.variable = 2, - time.plot = 1960:2003 - ) - -# fit training model -synth.out <- - synth( - data.prep.obj=dataprep.out, - Margin.ipop=.005,Sigf.ipop=7,Bound.ipop=6 - ) - -# data prep for main model -dataprep.out <- - dataprep( - foo = d, - predictors = c("gdp","trade","infrate"), - dependent = "gdp", - unit.variable = 1, - time.variable = 3, - special.predictors = list( - list("industry" ,1981:1990, c("mean")), - list("schooling",c(1980,1985), c("mean")), - list("invest80" ,1980, c("mean")) - ), - treatment.identifier = 7, - controls.identifier = unique(d$index)[-7], - time.predictors.prior = 1981:1990, - time.optimize.ssr = 1960:1989, - unit.names.variable = 2, - time.plot = 1960:2003 - ) - -# fit main model with v from training model -synth.out <- synth( - data.prep.obj=dataprep.out, - custom.v=as.numeric(synth.out$solution.v) -) - -synth.tables <- synth.tab( - dataprep.res = dataprep.out, - synth.res = synth.out -) - -table3 <- list() -synth.tables$tab.w[,1] <- round(synth.tables$tab.w[,1],2) -table3[[5]] <-synth.tables$tab.w[order(-1*synth.tables$tab.w[,1]),2:1][1:5,] - -# compute loss for all combinations -# of 4, 3, 2, 1 sized donor pools - -# get W and v -solution.w <- round(synth.out$solution.w,3) -V <- diag(as.numeric(synth.out$solution.v)) - -# compute scaled Xs -nvarsV <- dim(dataprep.out$X0)[1] -big.dataframe <- cbind(dataprep.out$X0, dataprep.out$X1) -divisor <- sqrt(apply(big.dataframe, 1, var)) -scaled.matrix <- - t(t(big.dataframe) %*% ( 1/(divisor) * - diag(rep(dim(big.dataframe)[1], 1)) )) -X0.scaled <- scaled.matrix[,c(1:(dim(dataprep.out$X0)[2]))] -X1.scaled <- as.matrix(scaled.matrix[,dim(scaled.matrix)[2]]) - -dn <- d[d$year==1970,c("country","index")] -dn <- dn[order(dn$index),] -dn <- dn[-7,] - -table2store <- matrix(NA,nrow(dataprep.out$X1),4) -fig7store <- matrix(NA,length(1960:2003),4) - -# loop through number of controls -for(pp in 4:1){ - store <- combinations(length(unique(d$index)[-7]), - r=pp, v=unique(d$index)[-7]) - store.loss <- matrix(NA,nrow=nrow(store),1) - store.w <- matrix(NA,nrow=nrow(store),pp) - store.c <- store.w - -# loop through combinations - for(k in 1:nrow(store)){ - # index positions of control units - posvector <- c() - for(i in 1:pp){ - posvector <- c(posvector,which(dn$index==store[k,i])) - } - - # run quad optimization - X0temp <- X0.scaled[ , posvector ] - H <- t(X0temp) %*% V %*% (X0temp) - c <- -1*c(t(X1.scaled) %*% V %*% (X0temp) ) - - if(pp == 1){ - solution.w <- matrix(1) - } else { - res <- ipop(c = c, H = H, A = t(rep(1, length(c))), - b = 1, l = rep(0, length(c)), - u = rep(1, length(c)), r = 0, - margin = 0.005,sigf = 7, bound = 6) - solution.w <- as.matrix(primal(res)) - } - loss.w <- t(X1.scaled - X0temp %*% solution.w) %*% V %*% (X1.scaled - X0temp %*% solution.w) - - store.loss[k] <- loss.w - store.w[k,] <- t(solution.w) - store.c[k,] <- dn$country[posvector] - } # close loop over combinations - -# get best fitting combination - dat <- data.frame(store.loss, - store, - store.c, - store.w - ) - colnames(dat) <- c("loss", - paste("CNo",1:pp,sep=""), - paste("CNa",1:pp,sep=""), - paste("W",1:pp,sep="") - ) - dat <- dat[order(dat$loss),] - Countries <- dat[1,paste("CNo",1:pp,sep="")] - Cweights <- as.numeric(dat[1,paste("W",1:pp,sep="")]) - - outdat <- data.frame(unit.names=as.vector( - (t(as.vector(dat[1,paste("CNa",1:pp,sep="")])))), - w.weights=round(Cweights,2)) - -table3[[pp]]<- outdat[order(-1*outdat$w.weights),] - - # get posvector for fitting - posvector <- c() - if(pp == 1 ){ - posvector <- c(posvector,which(dn$index==Countries)) - } else { - for(i in 1:pp){ - posvector <- c(posvector,which(dn$index==Countries[1,i])) - } - } - - X0t <- as.matrix(dataprep.out$X0[,posvector])%*% as.matrix(Cweights) - table2store[,(4:1)[pp]] <- X0t - - fig7store[,(4:1)[pp]] <- - dataprep.out$Y0[,posvector]%*%as.matrix(Cweights) - -} # close loop over number of countries - -# Table 3 -table3 - -# Table 4 -synth.tables$tab.pred[,3] <- c(8021.1,31.9,7.4,34.2,44.1,25.9) -table4 <- round( - cbind(synth.tables$tab.pred[,1:2], - table2store, - synth.tables$tab.pred[,3]),1) -rownames(table4) <- c("GDP per-capita","Trade openness", - "Inflation rate","Industry share", - "Schooling","Investment rate") -colnames(table4)[2:7] <- c(5:1,"OECD Sample") -table4 - -## Figure 7: Per-Capita GDP Gaps Between West Germany and Sparse Synthetic Controls -Text.height <- 23000 -Cex.set <- .8 - -par(mfrow=c(2,2)) -for(pp in 4:1){ -#pdf(file = paste("2ger_vs_synth","CValt",pp,".pdf",sep=""), width = 5.5, height = 5.5, family = "Times",pointsize = 12) - plot(1960:2003,dataprep.out$Y1, - type="l",ylim=c(0,33000),col="black",lty="solid", - ylab ="per-capita GDP (PPP, 2002 USD)", - xlab ="year", - xaxs = "i", yaxs = "i", - lwd=2, - main=paste("No. of control countries: ",pp,sep="") - ) - lines(1960:2003,fig7store[,c(4:1)[pp]],col="black",lty="dashed",lwd=2) - abline(v=1990,lty="dotted") - legend(x="bottomright", - legend=c("West Germany","synthetic West Germany") - ,lty=c("solid","dashed"),col=c("black","black") - ,cex=.8,bg="white",lwd=c(2,2)) - arrows(1987,Text.height,1989,Text.height,col="black",length=.1) - text(1982.5,Text.height,"reunification",cex=Cex.set) - #dev.off() -} - - -