Conceptions to under-18s were extracted from the dataset and coded into three datasets:
I first constructed a model with a single intervention at 1999 to test England data:
EngMod <- lm(Value ~ Time + Cat1 + Trend1, data = EngData)
summary(EngMod)$coefficients
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 40.800000 2.4243796 16.829047 1.147315e-13
## Time 0.750000 0.5421078 1.383489 1.810468e-01
## Cat1 4.016667 2.4104749 1.666338 1.104935e-01
## Trend1 -2.295614 0.5575524 -4.117306 4.908456e-04
A visual inspection of the model and a test for autocorrelation raise concerns about this model:
## [1] "Durbin Watson Test"
## lag Autocorrelation D-W Statistic p-value
## 1 0.77051355 0.3725459 0.000
## 2 0.50992122 0.8417204 0.002
## 3 0.17667011 1.4588131 0.236
## 4 -0.02768092 1.8009124 0.970
## 5 -0.25166863 2.2352735 0.208
## 6 -0.42708533 2.5828226 0.018
## 7 -0.56415048 2.8136162 0.000
## 8 -0.56958102 2.6501264 0.000
## 9 -0.43563602 2.1664036 0.018
## 10 -0.25030550 1.5869374 0.256
## 11 -0.09536613 1.2216464 0.760
## 12 0.01388605 0.9759968 0.812
## 13 0.07924590 0.8377393 0.768
## 14 0.12616731 0.7167512 0.682
## 15 0.14731591 0.6190418 0.612
## Alternative hypothesis: rho[lag] != 0
A second intervention at 2007 makes a better-fit model. This corresponds with a reevaluation and adaptation of the model after 2005 or potentially a phase-in period:
## [1] "Durbin Watson Test"
## lag Autocorrelation D-W Statistic p-value
## 1 0.18335719 1.4090132 0.004
## 2 -0.21037978 2.1848548 0.954
## 3 -0.40649345 2.4644492 0.244
## 4 -0.18427433 1.6585720 0.766
## 5 -0.18777350 1.5415165 0.874
## 6 0.05529295 1.0427883 0.250
## 7 0.13015466 0.8708738 0.162
## 8 0.07780085 0.9297698 0.312
## 9 0.07435073 0.9254829 0.326
## 10 0.05638279 0.9025666 0.380
## 11 0.14285767 0.7206965 0.196
## 12 -0.02674649 1.0574924 0.954
## 13 -0.11503662 1.2260738 0.396
## 14 -0.21481723 1.4232226 0.046
## 15 -0.14384408 1.2723560 0.086
## Alternative hypothesis: rho[lag] != 0
EngMod2_q3 <- gls(
Value ~ Time + Cat1 + Trend1 + Cat2 + Trend2,
data = EngData,
correlation = corARMA(q = 3, form = ~ Time),
method = "ML"
)
printCoefficients(EngMod2_q3)
## # A tibble: 6 x 4
## Coefficient Value Std.Error `p-value`
## <chr> <chr> <chr> <chr>
## 1 (Intercept) 40.931 0.749 0.000
## 2 Time 0.693 0.178 0.001
## 3 Cat1 -1.143 0.944 0.241
## 4 Trend1 -1.140 0.191 0.000
## 5 Cat2 1.230 0.869 0.173
## 6 Trend2 -2.246 0.136 0.000
modScot99_07_null <- gls(
Value ~ Time +
England +
Time_Eng +
Cat1 +
Trend1 +
Cat1_Eng +
Trend1_Eng +
Cat2 +
Trend2 +
Cat2_Eng +
Trend2_Eng,
data = EngScotContro,
correlation = NULL,
method = "ML"
)
printCoefficients(modScot99_07_null)
## # A tibble: 12 x 4
## Coefficient Value Std.Error `p-value`
## <chr> <chr> <chr> <chr>
## 1 (Intercept) 38.330 1.786 0.000
## 2 Time 0.950 0.344 0.009
## 3 England 2.470 2.009 0.227
## 4 Time_Eng -0.200 0.400 0.621
## 5 Cat1 -4.372 1.154 0.001
## 6 Trend1 -0.885 0.371 0.023
## 7 Cat1_Eng 2.799 1.583 0.085
## 8 Trend1_Eng -0.298 0.447 0.509
## 9 Cat2 0.849 1.035 0.417
## 10 Trend2 -2.805 0.198 0.000
## 11 Cat2_Eng 0.503 1.463 0.733
## 12 Trend2_Eng 0.544 0.281 0.061
modScotPI_null <-
gls(
Value ~ Time + England + Time_Eng + Cat2 + Trend2 + Cat2_Eng + Trend2_Eng,
data = EngScotContro %>% filter(Year < 1999 | Year > 2007),
correlation = NULL,
method = "ML"
)
printCoefficients(modScotPI_null)
## # A tibble: 8 x 4
## Coefficient Value Std.Error `p-value`
## <chr> <chr> <chr> <chr>
## 1 (Intercept) 38.330 1.912 0.000
## 2 Time 0.950 0.368 0.017
## 3 England 2.470 2.150 0.263
## 4 Time_Eng -0.200 0.429 0.645
## 5 Cat2 -11.487 4.168 0.012
## 6 Trend2 -3.690 0.398 0.000
## 7 Cat2_Eng 0.618 5.025 0.903
## 8 Trend2_Eng 0.245 0.479 0.613
modScotPI_99_07_null <- gls(
Value ~ Time +
England +
Time_Eng +
Cat1 +
Trend1 +
Cat1_Eng +
Trend1_Eng +
Cat2 +
Trend2 +
Cat2_Eng +
Trend2_Eng,
data = {EngScotContro %>%
filter(Year < 1999 | Year > 2000) %>%
mutate(Trend1=ifelse(Cat1==0,0,Trend1-2),
Trend1_Eng=Trend1*England)},
correlation = NULL,
method = "ML"
)
printCoefficients(modScotPI_99_07_null)
## # A tibble: 12 x 4
## Coefficient Value Std.Error `p-value`
## <chr> <chr> <chr> <chr>
## 1 (Intercept) 38.330 1.623 0.000
## 2 Time 0.950 0.312 0.005
## 3 England 2.470 1.825 0.186
## 4 Time_Eng -0.200 0.364 0.586
## 5 Cat1 -7.966 1.566 0.000
## 6 Trend1 -0.504 0.364 0.176
## 7 Cat1_Eng 3.387 2.040 0.107
## 8 Trend1_Eng -0.546 0.450 0.233
## 9 Cat2 0.003 0.984 0.997
## 10 Trend2 -3.187 0.226 0.000
## 11 Cat2_Eng 1.056 1.392 0.454
## 12 Trend2_Eng 0.792 0.320 0.019
modWal99_07_null <- gls(
Value ~ Time +
England +
Time_Eng +
Cat1 +
Trend1 +
Cat1_Eng +
Trend1_Eng +
Cat2 +
Trend2 +
Cat2_Eng +
Trend2_Eng,
data = EngWalContro,
correlation = NULL,
method = "ML"
)
printCoefficients(modWal99_07_null)
## # A tibble: 12 x 4
## Coefficient Value Std.Error `p-value`
## <chr> <chr> <chr> <chr>
## 1 (Intercept) 44.729 1.031 0.000
## 2 Time 1.339 0.231 0.000
## 3 England -3.929 1.458 0.010
## 4 Time_Eng -0.589 0.326 0.079
## 5 Cat1 -4.692 1.215 0.000
## 6 Trend1 -1.973 0.279 0.000
## 7 Cat1_Eng 3.120 1.719 0.077
## 8 Trend1_Eng 0.789 0.395 0.053
## 9 Cat2 1.553 1.161 0.189
## 10 Trend2 -2.135 0.223 0.000
## 11 Cat2_Eng -0.200 1.642 0.904
## 12 Trend2_Eng -0.127 0.315 0.690
modWalPI_null <-
gls(
Value ~ Time + England + Time_Eng + Cat2 + Trend2 + Cat2_Eng + Trend2_Eng,
data = EngWalContro %>% filter(Year < 1999 | Year > 2007),
correlation = NULL,
method = "ML"
)
printCoefficients(modWalPI_null)
## # A tibble: 8 x 4
## Coefficient Value Std.Error `p-value`
## <chr> <chr> <chr> <chr>
## 1 (Intercept) 44.729 1.120 0.000
## 2 Time 1.339 0.250 0.000
## 3 England -3.929 1.583 0.020
## 4 Time_Eng -0.589 0.354 0.109
## 5 Cat2 -20.893 3.194 0.000
## 6 Trend2 -4.108 0.303 0.000
## 7 Cat2_Eng 10.024 4.517 0.036
## 8 Trend2_Eng 0.663 0.429 0.135
modWalPI_99_07_null <- gls(
Value ~ Time +
England +
Time_Eng +
Cat1 +
Trend1 +
Cat1_Eng +
Trend1_Eng +
Cat2 +
Trend2 +
Cat2_Eng +
Trend2_Eng,
data = {EngWalContro %>%
filter(Year < 1999 | Year > 2000) %>%
mutate(Trend1=ifelse(Cat1==0,0,Trend1-2),
Trend1_Eng=Trend1*England)},
correlation = NULL,
method = "ML"
)
printCoefficients(modWalPI_99_07_null)
## # A tibble: 12 x 4
## Coefficient Value Std.Error `p-value`
## <chr> <chr> <chr> <chr>
## 1 (Intercept) 44.729 0.972 0.000
## 2 Time 1.339 0.217 0.000
## 3 England -3.929 1.375 0.007
## 4 Time_Eng -0.589 0.307 0.064
## 5 Cat1 -10.525 1.522 0.000
## 6 Trend1 -1.579 0.307 0.000
## 7 Cat1_Eng 5.946 2.152 0.009
## 8 Trend1_Eng 0.529 0.435 0.232
## 9 Cat2 0.682 1.146 0.556
## 10 Trend2 -2.529 0.263 0.000
## 11 Cat2_Eng 0.377 1.620 0.817
## 12 Trend2_Eng 0.134 0.372 0.721