Question 1


All Deaths Columns 1-4 

#column 1
rdd_cd<-mutate(rdd_cd, da = ifelse(agecell >= 21,1,0))
rc.1 = summary(lm(all ~ agecell + da, data=rdd_cd))
rc.1
## 
## Call:
## lm(formula = all ~ agecell + da, data = rdd_cd)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -5.0559 -1.8483  0.1149  1.4909  5.8043 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 112.3097    12.6681   8.866 1.96e-11 ***
## agecell      -0.9747     0.6325  -1.541     0.13    
## da            7.6627     1.4403   5.320 3.15e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.493 on 45 degrees of freedom
##   (2 observations deleted due to missingness)
## Multiple R-squared:  0.5946, Adjusted R-squared:  0.5765 
## F-statistic: 32.99 on 2 and 45 DF,  p-value: 1.508e-09
#column 2
rdd_cd <- mutate(rdd_cd, age = agecell - 21)
rc.2 = summary(lm(all ~ age + da + I(age^2) + age*da + I(age^2)*da, data=rdd_cd))
rc.2 
## 
## Call:
## lm(formula = all ~ age + da + I(age^2) + age * da + I(age^2) * 
##     da, data = rdd_cd)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -4.3343 -1.3946  0.1849  1.2848  5.0817 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  93.0729     1.4038  66.301  < 2e-16 ***
## age          -0.8306     3.2901  -0.252    0.802    
## da            9.5478     1.9853   4.809 1.97e-05 ***
## I(age^2)     -0.8403     1.6153  -0.520    0.606    
## age:da       -6.0170     4.6529  -1.293    0.203    
## da:I(age^2)   2.9042     2.2843   1.271    0.211    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.285 on 42 degrees of freedom
##   (2 observations deleted due to missingness)
## Multiple R-squared:  0.6821, Adjusted R-squared:  0.6442 
## F-statistic: 18.02 on 5 and 42 DF,  p-value: 1.624e-09
#column 3
rdd_cd <- filter(rdd_cd, (agecell >= 20 & agecell < 22))
rc.3 = summary(lm(all ~ agecell + da, data=rdd_cd))
rc.3
## 
## Call:
## lm(formula = all ~ agecell + da, data = rdd_cd)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -4.2577 -1.1736  0.0033  1.2509  3.9358 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  160.094     34.860   4.592 0.000158 ***
## agecell       -3.256      1.700  -1.916 0.069094 .  
## da             9.753      1.934   5.043 5.41e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.362 on 21 degrees of freedom
##   (2 observations deleted due to missingness)
## Multiple R-squared:  0.7029, Adjusted R-squared:  0.6746 
## F-statistic: 24.84 on 2 and 21 DF,  p-value: 2.924e-06
#column 4
rdd_cd <- filter(rdd_cd, (agecell >= 20 & agecell < 22))
rdd_cd <- mutate(rdd_cd, age = agecell - 21)
rc.4 = summary(lm(all ~ age + da + I(age^2) + age*da + I(age^2)*da, data=rdd_cd))
rc.4 
## 
## Call:
## lm(formula = all ~ age + da + I(age^2) + age * da + I(age^2) * 
##     da, data = rdd_cd)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -4.3270 -1.1898  0.1008  1.1612  3.7022 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  94.3403     2.0464  46.101   <2e-16 ***
## age           9.3988     9.6193   0.977   0.3415    
## da            9.6111     2.8940   3.321   0.0038 ** 
## I(age^2)     11.1633     9.4514   1.181   0.2529    
## age:da      -24.4478    13.6038  -1.797   0.0891 .  
## da:I(age^2)  -0.8742    13.3663  -0.065   0.9486    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.333 on 18 degrees of freedom
##   (2 observations deleted due to missingness)
## Multiple R-squared:  0.7517, Adjusted R-squared:  0.6827 
## F-statistic:  10.9 on 5 and 18 DF,  p-value: 6.056e-05


Motor Vehicle Deaths Columns 1-4 

#column 1
rdd_cd<-mutate(rdd_cd, da = ifelse(agecell >= 21,1,0))
rc.5 = summary(lm(mva ~ agecell + da, data=rdd_cd))
rc.5
## 
## Call:
## lm(formula = mva ~ agecell + da, data = rdd_cd)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -2.5318 -0.8494 -0.1800  0.7577  3.3094 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  95.4814     6.7549  14.135  < 2e-16 ***
## agecell      -3.1488     0.3372  -9.337 4.26e-12 ***
## da            4.5340     0.7680   5.904 4.34e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.329 on 45 degrees of freedom
##   (2 observations deleted due to missingness)
## Multiple R-squared:  0.7025, Adjusted R-squared:  0.6893 
## F-statistic: 53.14 on 2 and 45 DF,  p-value: 1.419e-12
#column 2
rdd_cd <- mutate(rdd_cd, age = agecell - 21)
rc.6 = summary(lm(mva ~ age + da + I(age^2) + age*da + I(age^2)*da, data=rdd_cd))
rc.6
## 
## Call:
## lm(formula = mva ~ age + da + I(age^2) + age * da + I(age^2) * 
##     da, data = rdd_cd)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -2.4646 -0.7761 -0.2646  0.8613  3.2418 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  29.8090     0.8166  36.505  < 2e-16 ***
## age          -2.9330     1.9138  -1.533 0.132877    
## da            4.6629     1.1548   4.038 0.000224 ***
## I(age^2)     -0.1852     0.9396  -0.197 0.844661    
## age:da       -0.8231     2.7065  -0.304 0.762527    
## da:I(age^2)   0.1985     1.3288   0.149 0.881980    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.329 on 42 degrees of freedom
##   (2 observations deleted due to missingness)
## Multiple R-squared:  0.7224, Adjusted R-squared:  0.6894 
## F-statistic: 21.86 on 5 and 42 DF,  p-value: 1.017e-10
#column 3
rdd_cd <- filter(rdd_cd, (agecell >= 20 & agecell < 22))
rc.7 = summary(lm(mva ~ agecell + da, data=rdd_cd))
rc.7
## 
## Call:
## lm(formula = mva ~ agecell + da, data = rdd_cd)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -2.1088 -0.8810 -0.3918  0.7157  3.1297 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 102.4206    19.7941   5.174 3.98e-05 ***
## agecell      -3.4683     0.9651  -3.594 0.001708 ** 
## da            4.7593     1.0981   4.334 0.000292 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.341 on 21 degrees of freedom
##   (2 observations deleted due to missingness)
## Multiple R-squared:  0.4736, Adjusted R-squared:  0.4234 
## F-statistic: 9.445 on 2 and 21 DF,  p-value: 0.001186
#column 4
rdd_cd <- filter(rdd_cd, (agecell >= 20 & agecell < 22))
rdd_cd <- mutate(rdd_cd, age = agecell - 21)
rc.8 = summary(lm(mva ~ age + da + I(age^2) + age*da + I(age^2)*da, data=rdd_cd))
rc.8
## 
## Call:
## lm(formula = mva ~ age + da + I(age^2) + age * da + I(age^2) * 
##     da, data = rdd_cd)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.89264 -0.64184 -0.03073  0.64433  2.39390 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  30.1883     1.1040  27.345 4.11e-16 ***
## age           0.6801     5.1894   0.131  0.89718    
## da            5.8925     1.5613   3.774  0.00139 ** 
## I(age^2)      4.4599     5.0989   0.875  0.39327    
## age:da      -15.1667     7.3390  -2.067  0.05347 .  
## da:I(age^2)   6.9652     7.2109   0.966  0.34688    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.258 on 18 degrees of freedom
##   (2 observations deleted due to missingness)
## Multiple R-squared:  0.6029, Adjusted R-squared:  0.4926 
## F-statistic: 5.465 on 5 and 18 DF,  p-value: 0.003125


Suicide Deaths Columns 1-4 

rdd_cd <- read_csv("rdd_cd.csv")
#column 1
rdd_cd<-mutate(rdd_cd, da = ifelse(agecell >= 21,1,0))
rc.9 = summary(lm(suicide ~ agecell + da, data=rdd_cd))
rc.9
## 
## Call:
## lm(formula = suicide ~ agecell + da, data = rdd_cd)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -1.4256 -0.5982 -0.0939  0.4860  1.6818 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  15.2644     3.9999   3.816 0.000411 ***
## agecell      -0.1814     0.1997  -0.908 0.368501    
## da            1.7943     0.4548   3.946 0.000276 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.7872 on 45 degrees of freedom
##   (2 observations deleted due to missingness)
## Multiple R-squared:  0.4754, Adjusted R-squared:  0.4521 
## F-statistic: 20.39 on 2 and 45 DF,  p-value: 4.961e-07
#column 2
rdd_cd <- mutate(rdd_cd, age = agecell - 21)
rc.10 = summary(lm(suicide ~ age + da + I(age^2) + age*da + I(age^2)*da, data=rdd_cd))
rc.10
## 
## Call:
## lm(formula = suicide ~ age + da + I(age^2) + age * da + I(age^2) * 
##     da, data = rdd_cd)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -1.3870 -0.5570 -0.0451  0.5336  1.8212 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 11.69811    0.49417  23.672   <2e-16 ***
## age          0.13823    1.15817   0.119   0.9056    
## da           1.81433    0.69886   2.596   0.0129 *  
## I(age^2)     0.05552    0.56861   0.098   0.9227    
## age:da      -0.70018    1.63791  -0.427   0.6712    
## da:I(age^2)  0.03088    0.80414   0.038   0.9696    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.8044 on 42 degrees of freedom
##   (2 observations deleted due to missingness)
## Multiple R-squared:  0.4887, Adjusted R-squared:  0.4279 
## F-statistic:  8.03 on 5 and 42 DF,  p-value: 2.203e-05
#column 3
rdd_cd <- filter(rdd_cd, (agecell >= 20 & agecell < 22))
rc.11 = summary(lm(suicide ~ agecell + da, data=rdd_cd))
rc.11
## 
## Call:
## lm(formula = suicide ~ agecell + da, data = rdd_cd)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.54946 -0.57179  0.09345  0.56743  1.46789 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)  
## (Intercept) 11.754110  12.491886   0.941   0.3574  
## agecell     -0.005422   0.609040  -0.009   0.9930  
## da           1.724426   0.693023   2.488   0.0213 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.8466 on 21 degrees of freedom
##   (2 observations deleted due to missingness)
## Multiple R-squared:  0.5409, Adjusted R-squared:  0.4972 
## F-statistic: 12.37 on 2 and 21 DF,  p-value: 0.0002819
#column 4
rdd_cd <- filter(rdd_cd, (agecell >= 20 & agecell < 22))
rdd_cd <- mutate(rdd_cd, age = agecell - 21)
rc.12 = summary(lm(suicide ~ age + da + I(age^2) + age*da + I(age^2)*da, data=rdd_cd))
rc.12
## 
## Call:
## lm(formula = suicide ~ age + da + I(age^2) + age * da + I(age^2) * 
##     da, data = rdd_cd)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.61676 -0.61178 -0.01709  0.50744  1.49905 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   11.919      0.796  14.974 1.33e-11 ***
## age            1.484      3.742   0.397    0.696    
## da             1.297      1.126   1.152    0.264    
## I(age^2)       1.407      3.676   0.383    0.706    
## age:da        -0.385      5.292  -0.073    0.943    
## da:I(age^2)   -2.630      5.199  -0.506    0.619    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.9073 on 18 degrees of freedom
##   (2 observations deleted due to missingness)
## Multiple R-squared:  0.548,  Adjusted R-squared:  0.4224 
## F-statistic: 4.364 on 5 and 18 DF,  p-value: 0.008896


QUESTION 2 

coef1<-rc.1$coefficients[3,1]
coef2<-rc.2$coefficients[3,1]
coef3<-rc.3$coefficients[3,1]
coef4<-rc.4$coefficients[3,1]
coef5<-rc.5$coefficients[3,1]
coef6<-rc.6$coefficients[3,1]
coef7<-rc.7$coefficients[3,1]
coef8<-rc.8$coefficients[3,1]
coef9<-rc.9$coefficients[3,1]
coef10<-rc.10$coefficients[3,1]
coef11<-rc.11$coefficients[3,1]
coef12<-rc.12$coefficients[3,1]

sdcoef1<-rc.1$coefficients[3,2]
sdcoef2<-rc.2$coefficients[3,2]
sdcoef3<-rc.3$coefficients[3,2]
sdcoef4<-rc.4$coefficients[3,2]
sdcoef5<-rc.5$coefficients[3,2]
sdcoef6<-rc.6$coefficients[3,2]
sdcoef7<-rc.7$coefficients[3,2]
sdcoef8<-rc.8$coefficients[3,2]
sdcoef9<-rc.9$coefficients[3,2]
sdcoef10<-rc.10$coefficients[3,2]
sdcoef11<-rc.11$coefficients[3,2]
sdcoef12<-rc.12$coefficients[3,2]


all_deaths_coef <- c(coef1,coef2,coef3,coef4)
all_deaths_sd <- c(sdcoef1,sdcoef2,sdcoef3, sdcoef4)

mva_deaths_coef <- c(coef5,coef6,coef7,coef8)
mva_deaths_sd <- c(sdcoef5,sdcoef6,sdcoef7, sdcoef8)

suicide_deaths_coef <- c(coef9,coef10,coef11,coef12)
suicide_deaths_sd <- c(sdcoef9,sdcoef10,sdcoef11, sdcoef12)

age19_22_1 <- c(coef1, sdcoef1, coef5,sdcoef5, coef9, sdcoef9)
age19_22_2 <- c(coef2, sdcoef2, coef6,sdcoef6, coef10, sdcoef10)

age20_21_1 <- c(coef3, sdcoef3, coef7,sdcoef7, coef11, sdcoef11)
age20_21_2 <- c(coef4, sdcoef4, coef8,sdcoef8, coef12, sdcoef12)

mortality21.data <- matrix(c(coef1,coef2, coef3, coef4, sdcoef1, sdcoef2, sdcoef3, sdcoef4, coef5, coef6, coef7, coef8, sdcoef5, sdcoef6, 
                  sdcoef7, sdcoef8, coef9, coef10, coef11, coef12, sdcoef9, sdcoef10, sdcoef11, sdcoef12),ncol=4,byrow=TRUE)

colnames(mortality21.data) <- c("age19_22_1","age19_22_2","age20_21_1", "age20_21_2")
rownames(mortality21.data) <- c("all_deaths_coef", "all_deaths_sd", "mva_deaths_coef", "mva_deaths_sd", "suicide_deaths_coef", "suicide_deaths_sd")
mortality21.data<- as.table(mortality21.data)

mortality21.data %>%
  kbl(caption = "Sharp RD estimates of MLDA on Mortality") %>%
    kable_paper("hover", full_width = F)  %>%
  kable_material_dark()
Sharp RD estimates of MLDA on Mortality
age19_22_1 age19_22_2 age20_21_1 age20_21_2
all_deaths_coef 7.6627089 9.5477885 9.7533106 9.611077
all_deaths_sd 1.4402859 1.9852773 1.9339773 2.894032
mva_deaths_coef 4.5340329 4.6628585 4.7592835 5.892489
mva_deaths_sd 0.7679953 1.1547995 1.0981324 1.561275
suicide_deaths_coef 1.7942890 1.8143320 1.7244260 1.296599
suicide_deaths_sd 0.4547684 0.6988606 0.6930231 1.125703


QUESTION 3 

all_dealths_coef_plot <- RDestimate(all~ + agecell, data = rdd_cd, cutpoint = 21)
plot(all_dealths_coef_plot)
title(main="Do All Death Rates Grow When people Turn 21?", xlab = "Age", ylab = "Death Rates per (100,000)")

mva_dealths_coef_plot <- RDestimate(mva ~ + agecell, data = rdd_cd, cutpoint = 21)
plot(all_dealths_coef_plot)
title(main="Do Motor Vehicle Deaths Grow When people Turn 21?", xlab = "Age", ylab = "Motor Vehicle Death Rates per (100,000)")

suicide_dealths_coef_plot <- RDestimate(suicide ~ + agecell, data = rdd_cd, cutpoint = 21)
plot(all_dealths_coef_plot)
title(main="Do Sucide Deaths Grow When people Turn 21?", xlab = "Age", ylab = "Suicide Death Rates per (100,000)")


QUESTION 4 

    rdd_cd %>% 
    select(agecell, mva) %>% 
    mutate(D = as.factor(ifelse(agecell >= 21, 1, 0))) %>% 
    ggplot(aes(x = agecell, y = mva, color = D)) +
  xlab("Age") +
  ylab("Motor Vehicle Death Rates Per 100,000") +
  ggtitle("Do Motor Vehicle Deaths Grow When people Turn 21?") +
          geom_vline(xintercept=21, linetype="longdash") +
    geom_point(show.legend = FALSE) + 
  scale_x_continuous(breaks = c(19, 20, 21, 22)) +
    geom_smooth(method = "lm",show.legend = FALSE)

In context of RDD model, the figure displays a sharp cut off at age 21 - also the MLDA. Through the graph we can see a clear jump in MVA death rates at and after the cut off. Adding an interaction term between running variable and cut off gives us different coefficients. The assignment of treatment and control is based on some clear-cut threshold of observed variables. The causual inference is made by comparing the coefficient of interest on both sides of the cutoff point, the more we can infer this causual effect. Everything around the cut off point should be very similar on average in observed and unobserved chracteristics. Which means that distribution of observed and unobservedcharacterstics should be continous around threshold, thereby creating a relitively small window for any difference in outcome, and as such implying causual effect. The jump at the cut off, implies a causual relationship between MLDA and MVA death rates.