Load data
clean_data <- read_feather(path=here::here("assets/clean_data/clean_data.feather"))
clean_data %<>% sample_frac(size=0.5)
expglm <- speedglm(dead ~ age_grp + sex +
employ_status +
generation + education + loinc_decile + rur_urb +
province + mar_stat,
family=binomial(), data=clean_data)
summary(expglm)
## Generalized Linear Model of class 'speedglm':
##
## Call: speedglm(formula = dead ~ age_grp + sex + employ_status + generation + education + loinc_decile + rur_urb + province + mar_stat, data = clean_data, family = binomial())
##
## Coefficients:
## ------------------------------------------------------------------
## Estimate Std. Error z value
## (Intercept) -5.99317 0.02211 -271.04
## age_grp25-29 0.11430 0.02738 4.17
## age_grp30-34 0.24090 0.02664 9.04
## age_grp35-39 0.68427 0.02439 28.05
## age_grp40-44 1.14069 0.02250 50.70
## age_grp45-49 1.70439 0.02141 79.62
## age_grp50-54 1.83550 0.02141 85.73
## age_grp55-59 2.52395 0.02077 121.50
## age_grp60-64 2.96274 0.02070 143.12
## age_grp65-69 2.98354 0.02086 143.02
## age_grp70-74 3.09258 0.02094 147.67
## age_grp75-79 3.29382 0.02104 156.59
## age_grp80-84 3.56057 0.02132 167.03
## age_grp85-89 3.91021 0.02232 175.20
## age_grp90 plus 4.19239 0.02505 167.34
## sexMale 0.68077 0.00436 156.21
## employ_statusPart-time 0.22512 0.00893 25.20
## employ_statusDidn't work 1.16630 0.00589 197.98
## generation1st Generation -0.26225 0.00536 -48.95
## generation2nd Generation 0.09857 0.00572 17.23
## educationNone 0.19136 0.00542 35.31
## educationTrades Cert/Dip 0.04130 0.00842 4.91
## educationApprenticeship Cert/Dip 0.00592 0.01003 0.59
## educationNon-Uni Cert/Dip < 1 year -0.07689 0.01498 -5.13
## educationNon-Uni Cert/Dip 1-2 years -0.16197 0.00957 -16.93
## educationNon-Uni Cert/Dip > 2years -0.16186 0.00993 -16.30
## educationUni Cert/Dip below Bachelor's -0.13662 0.01079 -12.66
## educationBachelor's degree -0.36841 0.00946 -38.94
## educationUni Cert/Dip above Bachelor's -0.24855 0.01767 -14.07
## educationMedicine etc. -0.10077 0.03143 -3.21
## educationMaster's degree -0.44269 0.01499 -29.53
## educationDoctorate -0.16592 0.02503 -6.63
## loinc_decile1 - lowest 0.39083 0.01253 31.19
## loinc_decile2 0.39347 0.00986 39.92
## loinc_decile3 0.42100 0.00947 44.44
## loinc_decile4 0.33117 0.00957 34.59
## loinc_decile5 0.23907 0.00973 24.56
## loinc_decile6 0.16404 0.00991 16.55
## loinc_decile7 0.12514 0.01006 12.44
## loinc_decile8 0.08228 0.01019 8.08
## loinc_decile9 0.04414 0.01033 4.27
## rur_urbRural -0.02882 0.00493 -5.85
## provinceNewfoundland and Labrador -0.11480 0.01542 -7.45
## provincePrince Edward Island 0.15375 0.02880 5.34
## provinceNova Scotia 0.01595 0.01152 1.38
## provinceNew Brunswick -0.08670 0.01306 -6.64
## provinceQuebec -0.14852 0.00556 -26.72
## provinceManitoba 0.10111 0.01012 9.99
## provinceSaskatchewan 0.10982 0.01086 10.11
## provinceAlberta 0.00923 0.00769 1.20
## provinceBritish Columbia -0.03613 0.00646 -5.59
## provinceYukon 0.32780 0.04600 7.13
## provinceNorthwest Territories 0.22313 0.03746 5.96
## provinceNunavut 0.42089 0.04205 10.01
## mar_statDivorced 0.18579 0.00705 26.35
## mar_statSeparated 0.23323 0.01182 19.74
## mar_statNever married 0.28332 0.00678 41.82
## mar_statWidowed 0.74883 0.00571 131.10
## Pr(>|z|)
## (Intercept) 0.00e+00 ***
## age_grp25-29 2.99e-05 ***
## age_grp30-34 1.55e-19 ***
## age_grp35-39 3.58e-173 ***
## age_grp40-44 0.00e+00 ***
## age_grp45-49 0.00e+00 ***
## age_grp50-54 0.00e+00 ***
## age_grp55-59 0.00e+00 ***
## age_grp60-64 0.00e+00 ***
## age_grp65-69 0.00e+00 ***
## age_grp70-74 0.00e+00 ***
## age_grp75-79 0.00e+00 ***
## age_grp80-84 0.00e+00 ***
## age_grp85-89 0.00e+00 ***
## age_grp90 plus 0.00e+00 ***
## sexMale 0.00e+00 ***
## employ_statusPart-time 3.61e-140 ***
## employ_statusDidn't work 0.00e+00 ***
## generation1st Generation 0.00e+00 ***
## generation2nd Generation 1.63e-66 ***
## educationNone 4.60e-273 ***
## educationTrades Cert/Dip 9.30e-07 ***
## educationApprenticeship Cert/Dip 5.55e-01
## educationNon-Uni Cert/Dip < 1 year 2.85e-07 ***
## educationNon-Uni Cert/Dip 1-2 years 2.63e-64 ***
## educationNon-Uni Cert/Dip > 2years 1.07e-59 ***
## educationUni Cert/Dip below Bachelor's 1.00e-36 ***
## educationBachelor's degree 0.00e+00 ***
## educationUni Cert/Dip above Bachelor's 6.17e-45 ***
## educationMedicine etc. 1.35e-03 **
## educationMaster's degree 1.10e-191 ***
## educationDoctorate 3.40e-11 ***
## loinc_decile1 - lowest 1.49e-213 ***
## loinc_decile2 0.00e+00 ***
## loinc_decile3 0.00e+00 ***
## loinc_decile4 3.38e-262 ***
## loinc_decile5 3.36e-133 ***
## loinc_decile6 1.54e-61 ***
## loinc_decile7 1.51e-35 ***
## loinc_decile8 6.68e-16 ***
## loinc_decile9 1.92e-05 ***
## rur_urbRural 4.93e-09 ***
## provinceNewfoundland and Labrador 9.61e-14 ***
## provincePrince Edward Island 9.38e-08 ***
## provinceNova Scotia 1.66e-01
## provinceNew Brunswick 3.19e-11 ***
## provinceQuebec 2.52e-157 ***
## provinceManitoba 1.70e-23 ***
## provinceSaskatchewan 4.85e-24 ***
## provinceAlberta 2.30e-01
## provinceBritish Columbia 2.27e-08 ***
## provinceYukon 1.03e-12 ***
## provinceNorthwest Territories 2.58e-09 ***
## provinceNunavut 1.40e-23 ***
## mar_statDivorced 4.65e-153 ***
## mar_statSeparated 9.72e-87 ***
## mar_statNever married 0.00e+00 ***
## mar_statWidowed 0.00e+00 ***
##
## -------------------------------------------------------------------
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---
## null df: 4346648; null deviance: 2381703;
## residuals df: 4346591; residuals deviance: 1738094;
## # obs.: 4346649; # non-zero weighted obs.: 4346649;
## AIC: 1738210; log Likelihood: -869047;
## RSS: 4695749; dispersion: 1; iterations: 7;
## rank: 58; max tolerance: 2.38e-10; convergence: TRUE.
interactionglm <- glm(dead ~ age_grp + sex +
employ_status +
generation + education + loinc_decile + rur_urb +
province + mar_stat + sex*mar_stat,
family=binomial(), data=clean_data)
summary(interactionglm)
anova(expglm, interactionglm, test = "Chi")
expglm_df <- tidy(expglm) %>%
bind_cols(tibble(or=exp(coef(expglm)))) %>%
bind_cols(as.data.frame(exp(confint.default(expglm)))) %>%
rename(ci95_lower = "2.5 %", ci95_upper = "97.5 %")
kable(expglm_df)
| (Intercept) |
-5.993 |
0.022 |
-271.04 |
0.00e+00 |
0.002 |
0.002 |
0.003 |
| age_grp25-29 |
0.114 |
0.027 |
4.17 |
2.99e-05 |
1.121 |
1.063 |
1.183 |
| age_grp30-34 |
0.241 |
0.027 |
9.04 |
1.55e-19 |
1.272 |
1.208 |
1.341 |
| age_grp35-39 |
0.684 |
0.024 |
28.05 |
3.58e-173 |
1.982 |
1.890 |
2.079 |
| age_grp40-44 |
1.141 |
0.022 |
50.70 |
0.00e+00 |
3.129 |
2.994 |
3.270 |
| age_grp45-49 |
1.704 |
0.021 |
79.62 |
0.00e+00 |
5.498 |
5.272 |
5.734 |
| age_grp50-54 |
1.835 |
0.021 |
85.73 |
0.00e+00 |
6.268 |
6.011 |
6.537 |
| age_grp55-59 |
2.524 |
0.021 |
121.50 |
0.00e+00 |
12.478 |
11.980 |
12.996 |
| age_grp60-64 |
2.963 |
0.021 |
143.12 |
0.00e+00 |
19.351 |
18.582 |
20.152 |
| age_grp65-69 |
2.984 |
0.021 |
143.03 |
0.00e+00 |
19.758 |
18.966 |
20.582 |
| age_grp70-74 |
3.093 |
0.021 |
147.67 |
0.00e+00 |
22.034 |
21.148 |
22.957 |
| age_grp75-79 |
3.294 |
0.021 |
156.59 |
0.00e+00 |
26.946 |
25.857 |
28.080 |
| age_grp80-84 |
3.561 |
0.021 |
167.03 |
0.00e+00 |
35.183 |
33.744 |
36.684 |
| age_grp85-89 |
3.910 |
0.022 |
175.20 |
0.00e+00 |
49.909 |
47.773 |
52.141 |
| age_grp90 plus |
4.192 |
0.025 |
167.34 |
0.00e+00 |
66.181 |
63.010 |
69.512 |
| sexMale |
0.681 |
0.004 |
156.21 |
0.00e+00 |
1.975 |
1.959 |
1.992 |
| employ_statusPart-time |
0.225 |
0.009 |
25.20 |
3.61e-140 |
1.252 |
1.231 |
1.275 |
| employ_statusDidn’t work |
1.166 |
0.006 |
197.98 |
0.00e+00 |
3.210 |
3.173 |
3.247 |
| generation1st Generation |
-0.262 |
0.005 |
-48.95 |
0.00e+00 |
0.769 |
0.761 |
0.777 |
| generation2nd Generation |
0.099 |
0.006 |
17.23 |
1.63e-66 |
1.104 |
1.091 |
1.116 |
| educationNone |
0.191 |
0.005 |
35.31 |
4.60e-273 |
1.211 |
1.198 |
1.224 |
| educationTrades Cert/Dip |
0.041 |
0.008 |
4.91 |
9.30e-07 |
1.042 |
1.025 |
1.060 |
| educationApprenticeship Cert/Dip |
0.006 |
0.010 |
0.59 |
5.55e-01 |
1.006 |
0.986 |
1.026 |
| educationNon-Uni Cert/Dip < 1 year |
-0.077 |
0.015 |
-5.13 |
2.85e-07 |
0.926 |
0.899 |
0.954 |
| educationNon-Uni Cert/Dip 1-2 years |
-0.162 |
0.010 |
-16.93 |
2.63e-64 |
0.850 |
0.835 |
0.867 |
| educationNon-Uni Cert/Dip > 2years |
-0.162 |
0.010 |
-16.30 |
1.07e-59 |
0.851 |
0.834 |
0.867 |
| educationUni Cert/Dip below Bachelor’s |
-0.137 |
0.011 |
-12.66 |
1.00e-36 |
0.872 |
0.854 |
0.891 |
| educationBachelor’s degree |
-0.368 |
0.009 |
-38.94 |
0.00e+00 |
0.692 |
0.679 |
0.705 |
| educationUni Cert/Dip above Bachelor’s |
-0.249 |
0.018 |
-14.07 |
6.17e-45 |
0.780 |
0.753 |
0.807 |
| educationMedicine etc. |
-0.101 |
0.031 |
-3.21 |
1.35e-03 |
0.904 |
0.850 |
0.962 |
| educationMaster’s degree |
-0.443 |
0.015 |
-29.53 |
1.10e-191 |
0.642 |
0.624 |
0.661 |
| educationDoctorate |
-0.166 |
0.025 |
-6.63 |
3.40e-11 |
0.847 |
0.807 |
0.890 |
| loinc_decile1 - lowest |
0.391 |
0.013 |
31.19 |
1.49e-213 |
1.478 |
1.442 |
1.515 |
| loinc_decile2 |
0.393 |
0.010 |
39.92 |
0.00e+00 |
1.482 |
1.454 |
1.511 |
| loinc_decile3 |
0.421 |
0.009 |
44.44 |
0.00e+00 |
1.523 |
1.495 |
1.552 |
| loinc_decile4 |
0.331 |
0.010 |
34.59 |
3.38e-262 |
1.393 |
1.367 |
1.419 |
| loinc_decile5 |
0.239 |
0.010 |
24.56 |
3.36e-133 |
1.270 |
1.246 |
1.295 |
| loinc_decile6 |
0.164 |
0.010 |
16.55 |
1.54e-61 |
1.178 |
1.156 |
1.201 |
| loinc_decile7 |
0.125 |
0.010 |
12.44 |
1.51e-35 |
1.133 |
1.111 |
1.156 |
| loinc_decile8 |
0.082 |
0.010 |
8.08 |
6.68e-16 |
1.086 |
1.064 |
1.108 |
| loinc_decile9 |
0.044 |
0.010 |
4.27 |
1.92e-05 |
1.045 |
1.024 |
1.067 |
| rur_urbRural |
-0.029 |
0.005 |
-5.85 |
4.93e-09 |
0.972 |
0.962 |
0.981 |
| provinceNewfoundland and Labrador |
-0.115 |
0.015 |
-7.45 |
9.61e-14 |
0.892 |
0.865 |
0.919 |
| provincePrince Edward Island |
0.154 |
0.029 |
5.34 |
9.38e-08 |
1.166 |
1.102 |
1.234 |
| provinceNova Scotia |
0.016 |
0.012 |
1.38 |
1.66e-01 |
1.016 |
0.993 |
1.039 |
| provinceNew Brunswick |
-0.087 |
0.013 |
-6.64 |
3.19e-11 |
0.917 |
0.894 |
0.941 |
| provinceQuebec |
-0.149 |
0.006 |
-26.72 |
2.52e-157 |
0.862 |
0.853 |
0.871 |
| provinceManitoba |
0.101 |
0.010 |
9.99 |
1.70e-23 |
1.106 |
1.085 |
1.129 |
| provinceSaskatchewan |
0.110 |
0.011 |
10.11 |
4.85e-24 |
1.116 |
1.093 |
1.140 |
| provinceAlberta |
0.009 |
0.008 |
1.20 |
2.30e-01 |
1.009 |
0.994 |
1.025 |
| provinceBritish Columbia |
-0.036 |
0.006 |
-5.59 |
2.27e-08 |
0.965 |
0.952 |
0.977 |
| provinceYukon |
0.328 |
0.046 |
7.13 |
1.03e-12 |
1.388 |
1.268 |
1.519 |
| provinceNorthwest Territories |
0.223 |
0.037 |
5.96 |
2.58e-09 |
1.250 |
1.161 |
1.345 |
| provinceNunavut |
0.421 |
0.042 |
10.01 |
1.40e-23 |
1.523 |
1.403 |
1.654 |
| mar_statDivorced |
0.186 |
0.007 |
26.35 |
4.65e-153 |
1.204 |
1.188 |
1.221 |
| mar_statSeparated |
0.233 |
0.012 |
19.74 |
9.72e-87 |
1.263 |
1.234 |
1.292 |
| mar_statNever married |
0.283 |
0.007 |
41.82 |
0.00e+00 |
1.328 |
1.310 |
1.345 |
| mar_statWidowed |
0.749 |
0.006 |
131.10 |
0.00e+00 |
2.115 |
2.091 |
2.138 |
intglm_df <- tidy(interactionglm) %>%
bind_cols(tibble(or=exp(coef(interactionglm)))) %>%
bind_cols(as.data.frame(exp(confint.default(interactionglm)))) %>%
rename(ci95_lower = "2.5 %", ci95_upper = "97.5 %")
kable(intglm_df)
# save(expglm, file = here::here("assets/clean_data/expglm.rds"))
feather::write_feather(expglm_df, path=here::here("assets/clean_data/expglm_df.feather"))
# save(interactionglm, file = here::here("assets/clean_data/intglm.rds"))
#feather::write_feather(intglm_df, path=here::here("assets/clean_data/intglm_df.feather"))
Effect plots
sm_effs <- Effect(c("sex","mar_stat"), interactionglm)
plot(sm_effs)
plot_model(interactionglm, type = "pred", terms=c("sex", "mar_stat"), axis.title=c("Sex", "Predicted probability of death"), title="Sex and Marital Status interaction", legend.title="Marital Status" )
ge_effs <- Effect(c("education","generation"), interactionglm)
plot(ge_effs, axes=list(x=list(labels=letters[1:13]), y=list(lab="Probability of death")))
plot_model(interactionglm, type = "pred", terms=c("education", "generation"), axis.title=c("Generation", "Predicted probability of death"), title="Education and Generation interaction", legend.title="Marital Status" )