COVID19 Models

Load Packages

library(tidyverse)
library(codebookr)
library(summarytools)
library(broom) 
library(performance)
library(gt)
library(gtsummary)
library(janitor)
library(forcats)
library(margins)
library(ggplot2)
library(expss)
library(glmtoolbox)
library(DescTools)

Import Data

load("~/Desktop/R-Code/SDOH_Vax/covid_21_later_1.8.22.rda")
View(covid_21_later_1.8.22)

Include only patients from Michigan

covid_2021_MI <- covid_21_later_1.8.22[ which(covid_21_later_1.8.22$PATIENT_STATE_CD=='MI'), ]

Make fully vaccinated and fully vaccinated + boosted variables

covid_2021_MI %>% 
mutate(full_vax = case_when(total_cov_vax>= 2 ~ '1',TRUE ~ "0")) %>% 
  mutate(vax_boost = case_when(total_cov_vax>= 3 ~ '1',TRUE ~ "0")) -> cov_clean1

cov_clean1$vax_boost = as.numeric(cov_clean1$vax_boost)

Baseline characteristics

cov_clean1 %>% 
  dplyr::select(ibd_3, age_yrs, gender, race_5, ethnic_3, lang_3, relig_affil, mstat_5, act_tob, max_ch, insurance, IC, insurance, pop_dens,r_pct, RPL_THEMES, RPL_4, RPL_THEME1, RPL_THEME2, RPL_THEME3, RPL_THEME4, full_vax, vax_boost) -> cov_baseline
cov_baseline %>% tbl_summary(
        statistic = list(all_continuous() ~ "{mean} ({sd})"),
        missing_text = "(Missing)")

Characteristic	N = 8,8561
IBD Diagnosis
CD	4,386 (50%)
UC	4,421 (50%)
IC	2 (<0.1%)
Unknown	47 (0.5%)
Age	48 (19)
Gender
Male	4,017 (45%)
Female	4,839 (55%)
Race
White	7,621 (86%)
Black	596 (6.7%)
Asian or Pacific Islander	237 (2.7%)
American Indian or Alaska Native	35 (0.4%)
Other	367 (4.1%)
Ethnicity
Hispanic	188 (2.2%)
Non-Hispanic	8,406 (98%)
(Missing)	262
Preferred Language
English	8,753 (99%)
Other	103 (1.2%)
Any Religious Affiliation	4,911 (58%)
(Missing)	383
Marital Status
Married	3,548 (53%)
Single	2,708 (41%)
Divorced/Separated	268 (4.0%)
Widowed	149 (2.2%)
(Missing)	2,183
Active Tobacco Use	993 (11%)
(Missing)	130
Charlson Comorbidity Index	3.5 (5.1)
Insurance Type
Private Insurance	5,765 (65%)
Medicaid	1,338 (15%)
Medicare	1,626 (18%)
Other Governmental	69 (0.8%)
Other	43 (0.5%)
(Missing)	15
Immunocompromised	4,935 (63%)
(Missing)	1,026
Population Density of Census Tract	2,194 (2,314)
(Missing)	171
Percentage Republican Voters in Census Tract	44 (18)
(Missing)	1,214
Social Vulnerability Index	0.35 (0.25)
(Missing)	28
SVI by Quartile
First	3,729 (42%)
Second	2,642 (30%)
Third	1,668 (19%)
Fourth	789 (8.9%)
(Missing)	28
Socioeconomic Status	0.33 (0.25)
(Missing)	69
Household Composition	0.37 (0.26)
(Missing)	28
Minority and Language Status	0.49 (0.29)
(Missing)	24
Housing and Transportation	0.42 (0.28)
(Missing)	45
full_vax
0	2,440 (28%)
1	6,416 (72%)
vax_boost	4,536 (51%)
1 n (%); Mean (SD)

Bivariate Analysis

Fully vaccinated

cov_clean1$full_vax = as.numeric(cov_clean1$full_vax)
tbl_full_biv <-
  tbl_uvregression(
    cov_clean1[c("full_vax", "ibd_3", "age_yrs", "gender", "race_5", "ethnic_3", "lang_3", "mstat_5", "relig_affil", "act_tob", "max_ch", "IC", "pop_dens", "r_pct", "RPL_THEMES", "RPL_THEME1", "RPL_THEME2", "RPL_THEME3", "RPL_THEME4")],
    method = glm,
    y = full_vax,
    method.args = list(family = binomial),
    exponentiate = TRUE)
print(tbl_full_biv, method = render)

`...` must be empty.
✖ Problematic argument:
• method = render

Characteristic	N	OR1	95% CI1	p-value
IBD Diagnosis	8,856
CD		—	—
UC		1.20	1.09, 1.32	<0.001
IC		0.42	0.02, 10.5	0.5
Unknown		0.98	0.53, 1.90	>0.9
Age	8,856	1.02	1.02, 1.02	<0.001
Gender	8,856
Male		—	—
Female		1.02	0.93, 1.12	0.6
Race	8,856
White		—	—
Black		0.74	0.62, 0.89	<0.001
Asian or Pacific Islander		2.91	1.98, 4.46	<0.001
American Indian or Alaska Native		0.72	0.36, 1.49	0.4
Other		0.72	0.57, 0.90	0.003
Ethnicity	8,594
Hispanic		—	—
Non-Hispanic		1.30	0.95, 1.76	0.10
Preferred Language	8,856
English		—	—
Other		0.74	0.49, 1.12	0.14
Marital Status	6,673
Married		—	—
Single		0.60	0.53, 0.67	<0.001
Divorced/Separated		0.69	0.52, 0.91	0.008
Widowed		1.44	0.94, 2.30	0.11
Any Religious Affiliation	8,473
Yes		—	—
No		0.93	0.84, 1.02	0.12
Active Tobacco Use	8,726
No		—	—
Yes		0.44	0.39, 0.51	<0.001
Charlson Comorbidity Index	8,856	1.03	1.02, 1.04	<0.001
Immunocompromised	7,830	0.87	0.78, 0.97	0.010
Population Density of Census Tract	8,685	1.00	1.00, 1.00	<0.001
Percentage Republican Voters in Census Tract	7,642	0.98	0.98, 0.98	<0.001
Social Vulnerability Index	8,828	0.32	0.27, 0.39	<0.001
Socioeconomic Status	8,787	0.25	0.21, 0.30	<0.001
Household Composition	8,828	0.29	0.24, 0.34	<0.001
Minority and Language Status	8,832	1.68	1.43, 1.97	<0.001
Housing and Transportation	8,811	0.56	0.48, 0.66	<0.001
1 OR = Odds Ratio, CI = Confidence Interval

NULL

Vaccinated + Boosted

cov_clean1$vax_boost = as.numeric(cov_clean1$vax_boost)
tbl_boost_biv <-
  tbl_uvregression(
    cov_clean1[c("vax_boost", "ibd_3", "age_yrs", "gender", "race_5", "ethnic_3", "lang_3", "mstat_5", "relig_affil", "act_tob", "max_ch", "IC", "pop_dens", "r_pct", "RPL_THEMES", "RPL_THEME1", "RPL_THEME2", "RPL_THEME3", "RPL_THEME4")],
    method = glm,
    y = vax_boost,
    method.args = list(family = binomial),
    exponentiate = TRUE)
print(tbl_boost_biv, method = render)

`...` must be empty.
✖ Problematic argument:
• method = render

Characteristic	N	OR1	95% CI1	p-value
IBD Diagnosis	8,856
CD		—	—
UC		1.19	1.10, 1.29	<0.001
IC		1.04	0.04, 26.3	>0.9
Unknown		0.84	0.47, 1.49	0.6
Age	8,856	1.02	1.02, 1.03	<0.001
Gender	8,856
Male		—	—
Female		1.08	0.99, 1.17	0.084
Race	8,856
White		—	—
Black		0.64	0.54, 0.75	<0.001
Asian or Pacific Islander		2.16	1.64, 2.88	<0.001
American Indian or Alaska Native		0.62	0.31, 1.20	0.2
Other		0.66	0.53, 0.82	<0.001
Ethnicity	8,594
Hispanic		—	—
Non-Hispanic		1.49	1.11, 2.00	0.008
Preferred Language	8,856
English		—	—
Other		0.63	0.42, 0.93	0.021
Marital Status	6,673
Married		—	—
Single		0.55	0.50, 0.61	<0.001
Divorced/Separated		0.67	0.53, 0.86	0.002
Widowed		1.09	0.78, 1.53	0.6
Any Religious Affiliation	8,473
Yes		—	—
No		0.94	0.86, 1.02	0.2
Active Tobacco Use	8,726
No		—	—
Yes		0.48	0.42, 0.55	<0.001
Charlson Comorbidity Index	8,856	1.05	1.04, 1.05	<0.001
Immunocompromised	7,830	0.85	0.77, 0.93	<0.001
Population Density of Census Tract	8,685	1.00	1.00, 1.00	<0.001
Percentage Republican Voters in Census Tract	7,642	0.98	0.98, 0.98	<0.001
Social Vulnerability Index	8,828	0.32	0.27, 0.38	<0.001
Socioeconomic Status	8,787	0.25	0.21, 0.30	<0.001
Household Composition	8,828	0.29	0.25, 0.34	<0.001
Minority and Language Status	8,832	1.63	1.41, 1.89	<0.001
Housing and Transportation	8,811	0.60	0.51, 0.69	<0.001
1 OR = Odds Ratio, CI = Confidence Interval

NULL

Multivariable Models: Full Vaccination

SVI Continuous

covid_SVI <- glm(full_vax ~ ibd_3 + age_yrs + gender + race_5 + ethnic_3 + mstat_5 + relig_affil
                      + lang_3 + act_tob + max_ch + IC + insurance
                      + pop_dens + RPL_THEMES,
              family = "binomial", 
              data = cov_clean1)
summary(covid_SVI )

Call:
glm(formula = full_vax ~ ibd_3 + age_yrs + gender + race_5 + 
    ethnic_3 + mstat_5 + relig_affil + lang_3 + act_tob + max_ch + 
    IC + insurance + pop_dens + RPL_THEMES, family = "binomial", 
    data = cov_clean1)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-2.5928   0.1930   0.5884   0.7478   1.8613  

Coefficients:
                                         Estimate Std. Error z value Pr(>|z|)    
(Intercept)                             0.2702567  0.3055591   0.884   0.3764    
ibd_3UC                                -0.0167532  0.0711872  -0.235   0.8139    
ibd_3Unknown                            0.6353846  0.7855093   0.809   0.4186    
age_yrs                                 0.0278880  0.0027693  10.071  < 2e-16 ***
genderFemale                            0.0785449  0.0676024   1.162   0.2453    
race_5Black                             0.0058227  0.1328929   0.044   0.9651    
race_5Asian or Pacific Islander         1.3896421  0.3546876   3.918 8.93e-05 ***
race_5American Indian or Alaska Native  0.3361174  0.5910656   0.569   0.5696    
race_5Other                            -0.1109427  0.2023429  -0.548   0.5835    
ethnic_3Non-Hispanic                   -0.3536417  0.2520787  -1.403   0.1606    
mstat_5Single                           0.1619629  0.0870996   1.860   0.0630 .  
mstat_5Divorced/Separated               0.0289382  0.1734601   0.167   0.8675    
mstat_5Widowed                          0.2102497  0.2742062   0.767   0.4432    
relig_affilNo                           0.1532223  0.0693411   2.210   0.0271 *  
lang_3Other                             0.2017215  0.4186213   0.482   0.6299    
act_tobYes                             -0.5236502  0.0981832  -5.333 9.64e-08 ***
max_ch                                  0.0030914  0.0076423   0.405   0.6858    
IC                                      0.1664620  0.0773338   2.153   0.0314 *  
insuranceMedicaid                      -1.0101336  0.0935137 -10.802  < 2e-16 ***
insuranceMedicare                      -0.5033601  0.1094899  -4.597 4.28e-06 ***
insuranceOther Governmental            -1.7199228  0.4084509  -4.211 2.54e-05 ***
insuranceOther                         -0.1871227  0.4499810  -0.416   0.6775    
pop_dens                                0.0001202  0.0000167   7.193 6.33e-13 ***
RPL_THEMES                             -0.9791784  0.1393509  -7.027 2.11e-12 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 6060.3  on 5388  degrees of freedom
Residual deviance: 5515.1  on 5365  degrees of freedom
  (3467 observations deleted due to missingness)
AIC: 5563.1

Number of Fisher Scoring iterations: 5

broom::glance(covid_SVI )

broom::tidy(covid_SVI , exponentiate = TRUE)

tbl_regression(covid_SVI, exponentiate = TRUE)

Characteristic	OR1	95% CI1	p-value
IBD Diagnosis
CD	—	—
UC	0.98	0.86, 1.13	0.8
Unknown	1.89	0.49, 12.5	0.4
Age	1.03	1.02, 1.03	<0.001
Gender
Male	—	—
Female	1.08	0.95, 1.23	0.2
Race
White	—	—
Black	1.01	0.78, 1.31	>0.9
Asian or Pacific Islander	4.01	2.12, 8.64	<0.001
American Indian or Alaska Native	1.40	0.48, 5.15	0.6
Other	0.89	0.61, 1.34	0.6
Ethnicity
Hispanic	—	—
Non-Hispanic	0.70	0.42, 1.14	0.2
Marital Status
Married	—	—
Single	1.18	0.99, 1.40	0.063
Divorced/Separated	1.03	0.74, 1.46	0.9
Widowed	1.23	0.74, 2.17	0.4
Any Religious Affiliation
Yes	—	—
No	1.17	1.02, 1.34	0.027
Preferred Language
English	—	—
Other	1.22	0.56, 2.90	0.6
Active Tobacco Use
No	—	—
Yes	0.59	0.49, 0.72	<0.001
Charlson Comorbidity Index	1.00	0.99, 1.02	0.7
Immunocompromised	1.18	1.01, 1.37	0.031
Insurance Type
Private Insurance	—	—
Medicaid	0.36	0.30, 0.44	<0.001
Medicare	0.60	0.49, 0.75	<0.001
Other Governmental	0.18	0.08, 0.40	<0.001
Other	0.83	0.36, 2.15	0.7
Population Density of Census Tract	1.00	1.00, 1.00	<0.001
Social Vulnerability Index	0.38	0.29, 0.49	<0.001
1 OR = Odds Ratio, CI = Confidence Interval

# Hosmer-Lemeshow Goodness-of-Fit Test
hltest(covid_SVI)

   The Hosmer-Lemeshow goodness-of-fit test

         Statistic =  4.53295 
degrees of freedom =  8 
           p-value =  0.80613 

# C-Statistic/AUROC 
Cstat(covid_SVI)

[1] 0.6983068

# Model performance 
model_performance(covid_SVI)

# Indices of model performance

AIC      |     AICc |      BIC | Tjur's R2 |  RMSE | Sigma | Log_loss | Score_log |   PCP
-----------------------------------------------------------------------------------------
5563.082 | 5563.306 | 5721.293 |     0.105 | 0.409 | 1.014 |    0.512 |      -Inf | 0.665

performance::check_model(covid_SVI)

Variable `Component` is not in your data frame :/

# Margins 

cplot(covid_SVI, "RPL_THEMES", what = "prediction", main = "Predicted Likelihood of COVID Vaccine Given SVI")

SVI Quartiles

covid_SVI4 <- glm(full_vax ~ ibd_3 + age_yrs + gender + race_5 + ethnic_3 + mstat_5 + relig_affil
                      + lang_3 + act_tob + max_ch + IC + insurance
                      + pop_dens + RPL_4,
              family = "binomial", 
              data = cov_clean1)
summary(covid_SVI4 )

Call:
glm(formula = full_vax ~ ibd_3 + age_yrs + gender + race_5 + 
    ethnic_3 + mstat_5 + relig_affil + lang_3 + act_tob + max_ch + 
    IC + insurance + pop_dens + RPL_4, family = "binomial", data = cov_clean1)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-2.5780   0.1844   0.5909   0.7445   1.8788  

Coefficients:
                                         Estimate Std. Error z value Pr(>|z|)    
(Intercept)                             1.134e-01  3.045e-01   0.372   0.7096    
ibd_3UC                                -2.412e-02  7.126e-02  -0.339   0.7350    
ibd_3Unknown                            6.419e-01  7.883e-01   0.814   0.4155    
age_yrs                                 2.771e-02  2.767e-03  10.014  < 2e-16 ***
genderFemale                            8.031e-02  6.762e-02   1.188   0.2350    
race_5Black                            -6.728e-04  1.336e-01  -0.005   0.9960    
race_5Asian or Pacific Islander         1.386e+00  3.545e-01   3.909 9.28e-05 ***
race_5American Indian or Alaska Native  2.843e-01  5.933e-01   0.479   0.6318    
race_5Other                            -1.143e-01  2.028e-01  -0.564   0.5730    
ethnic_3Non-Hispanic                   -3.555e-01  2.525e-01  -1.408   0.1592    
mstat_5Single                           1.531e-01  8.722e-02   1.755   0.0792 .  
mstat_5Divorced/Separated               1.158e-02  1.739e-01   0.067   0.9469    
mstat_5Widowed                          1.959e-01  2.745e-01   0.714   0.4755    
relig_affilNo                           1.578e-01  6.938e-02   2.275   0.0229 *  
lang_3Other                             2.294e-01  4.202e-01   0.546   0.5851    
act_tobYes                             -5.296e-01  9.830e-02  -5.387 7.15e-08 ***
max_ch                                  3.369e-03  7.651e-03   0.440   0.6597    
IC                                      1.587e-01  7.726e-02   2.054   0.0400 *  
insuranceMedicaid                      -1.021e+00  9.351e-02 -10.923  < 2e-16 ***
insuranceMedicare                      -5.113e-01  1.095e-01  -4.669 3.03e-06 ***
insuranceOther Governmental            -1.710e+00  4.090e-01  -4.180 2.91e-05 ***
insuranceOther                         -2.067e-01  4.498e-01  -0.460   0.6458    
pop_dens                                1.202e-04  1.675e-05   7.177 7.12e-13 ***
RPL_4Second                            -5.144e-02  8.141e-02  -0.632   0.5275    
RPL_4Third                             -4.209e-01  9.309e-02  -4.522 6.14e-06 ***
RPL_4Fourth                            -7.624e-01  1.241e-01  -6.146 7.96e-10 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 6060.3  on 5388  degrees of freedom
Residual deviance: 5514.4  on 5363  degrees of freedom
  (3467 observations deleted due to missingness)
AIC: 5566.4

Number of Fisher Scoring iterations: 5

broom::glance(covid_SVI4 )

broom::tidy(covid_SVI4 , exponentiate = TRUE)

tbl_regression(covid_SVI4, exponentiate = TRUE)

Characteristic	OR1	95% CI1	p-value
IBD Diagnosis
CD	—	—
UC	0.98	0.85, 1.12	0.7
Unknown	1.90	0.49, 12.6	0.4
Age	1.03	1.02, 1.03	<0.001
Gender
Male	—	—
Female	1.08	0.95, 1.24	0.2
Race
White	—	—
Black	1.00	0.77, 1.30	>0.9
Asian or Pacific Islander	4.00	2.11, 8.60	<0.001
American Indian or Alaska Native	1.33	0.45, 4.91	0.6
Other	0.89	0.60, 1.34	0.6
Ethnicity
Hispanic	—	—
Non-Hispanic	0.70	0.42, 1.14	0.2
Marital Status
Married	—	—
Single	1.17	0.98, 1.38	0.079
Divorced/Separated	1.01	0.72, 1.43	>0.9
Widowed	1.22	0.73, 2.15	0.5
Any Religious Affiliation
Yes	—	—
No	1.17	1.02, 1.34	0.023
Preferred Language
English	—	—
Other	1.26	0.57, 2.99	0.6
Active Tobacco Use
No	—	—
Yes	0.59	0.49, 0.71	<0.001
Charlson Comorbidity Index	1.00	0.99, 1.02	0.7
Immunocompromised	1.17	1.01, 1.36	0.040
Insurance Type
Private Insurance	—	—
Medicaid	0.36	0.30, 0.43	<0.001
Medicare	0.60	0.48, 0.74	<0.001
Other Governmental	0.18	0.08, 0.40	<0.001
Other	0.81	0.35, 2.11	0.6
Population Density of Census Tract	1.00	1.00, 1.00	<0.001
SVI by Quartile
First	—	—
Second	0.95	0.81, 1.11	0.5
Third	0.66	0.55, 0.79	<0.001
Fourth	0.47	0.37, 0.60	<0.001
1 OR = Odds Ratio, CI = Confidence Interval

# Hosmer-Lemeshow Goodness-of-Fit Test
hltest(covid_SVI4)

   The Hosmer-Lemeshow goodness-of-fit test

         Statistic =  9.10998 
degrees of freedom =  8 
           p-value =  0.3331 

# C-Statistic/AUROC 
Cstat(covid_SVI4)

[1] 0.6973355

# Model performance 
model_performance(covid_SVI4)

# Indices of model performance

AIC      |     AICc |      BIC | Tjur's R2 |  RMSE | Sigma | Log_loss | Score_log |   PCP
-----------------------------------------------------------------------------------------
5566.356 | 5566.617 | 5737.751 |     0.106 | 0.409 | 1.014 |    0.512 |      -Inf | 0.665

performance::check_model(covid_SVI4)

Variable `Component` is not in your data frame :/

# Margins 
cplot(covid_SVI4, "RPL_4", what = "prediction", main = "Predicted Likelihood of COVID Vaccine Given SVI Quartile")

All Themes

covid_all_themes <- glm(full_vax ~ ibd_3 + age_yrs + gender + race_5 + ethnic_3 + 
                      lang_3 + mstat_5 + relig_affil + act_tob + max_ch + IC
                      + pop_dens + RPL_THEME1 + RPL_THEME2 + RPL_THEME3
                      + RPL_THEME4,
              family = "binomial",
              data = cov_clean1)
summary(covid_all_themes )

Call:
glm(formula = full_vax ~ ibd_3 + age_yrs + gender + race_5 + 
    ethnic_3 + lang_3 + mstat_5 + relig_affil + act_tob + max_ch + 
    IC + pop_dens + RPL_THEME1 + RPL_THEME2 + RPL_THEME3 + RPL_THEME4, 
    family = "binomial", data = cov_clean1)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-2.4875  -0.8094   0.5922   0.7832   1.6428  

Coefficients:
                                         Estimate Std. Error z value Pr(>|z|)    
(Intercept)                             3.239e-01  3.081e-01   1.052 0.292992    
ibd_3UC                                -2.166e-03  7.057e-02  -0.031 0.975510    
ibd_3Unknown                            6.241e-01  7.831e-01   0.797 0.425484    
age_yrs                                 2.555e-02  2.512e-03  10.170  < 2e-16 ***
genderFemale                            3.179e-02  6.708e-02   0.474 0.635626    
race_5Black                            -1.012e-01  1.322e-01  -0.766 0.443945    
race_5Asian or Pacific Islander         1.198e+00  3.562e-01   3.363 0.000771 ***
race_5American Indian or Alaska Native  2.899e-01  5.962e-01   0.486 0.626836    
race_5Other                            -3.473e-01  2.001e-01  -1.736 0.082586 .  
ethnic_3Non-Hispanic                   -3.755e-01  2.509e-01  -1.496 0.134550    
lang_3Other                            -1.760e-01  4.079e-01  -0.432 0.666104    
mstat_5Single                          -4.669e-02  8.373e-02  -0.558 0.577045    
mstat_5Divorced/Separated              -2.302e-01  1.692e-01  -1.361 0.173653    
mstat_5Widowed                          5.660e-02  2.720e-01   0.208 0.835143    
relig_affilNo                           1.251e-01  6.869e-02   1.821 0.068644 .  
act_tobYes                             -6.156e-01  9.631e-02  -6.392 1.64e-10 ***
max_ch                                 -2.452e-03  7.532e-03  -0.326 0.744737    
IC                                      2.270e-01  7.692e-02   2.952 0.003161 ** 
pop_dens                                9.617e-05  1.783e-05   5.392 6.95e-08 ***
RPL_THEME1                             -1.279e+00  2.022e-01  -6.324 2.54e-10 ***
RPL_THEME2                             -5.463e-01  1.746e-01  -3.129 0.001757 ** 
RPL_THEME3                              4.282e-01  1.294e-01   3.310 0.000933 ***
RPL_THEME4                              2.214e-01  1.415e-01   1.565 0.117618    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 6024  on 5344  degrees of freedom
Residual deviance: 5552  on 5322  degrees of freedom
  (3511 observations deleted due to missingness)
AIC: 5598

Number of Fisher Scoring iterations: 5

broom::glance(covid_all_themes )

broom::tidy(covid_all_themes , exponentiate = TRUE)

tbl_regression(covid_all_themes, exponentiate = TRUE)

Characteristic	OR1	95% CI1	p-value
IBD Diagnosis
CD	—	—
UC	1.00	0.87, 1.15	>0.9
Unknown	1.87	0.49, 12.3	0.4
Age	1.03	1.02, 1.03	<0.001
Gender
Male	—	—
Female	1.03	0.90, 1.18	0.6
Race
White	—	—
Black	0.90	0.70, 1.17	0.4
Asian or Pacific Islander	3.31	1.74, 7.15	<0.001
American Indian or Alaska Native	1.34	0.45, 4.94	0.6
Other	0.71	0.48, 1.05	0.083
Ethnicity
Hispanic	—	—
Non-Hispanic	0.69	0.41, 1.11	0.13
Preferred Language
English	—	—
Other	0.84	0.39, 1.95	0.7
Marital Status
Married	—	—
Single	0.95	0.81, 1.12	0.6
Divorced/Separated	0.79	0.57, 1.11	0.2
Widowed	1.06	0.64, 1.86	0.8
Any Religious Affiliation
Yes	—	—
No	1.13	0.99, 1.30	0.069
Active Tobacco Use
No	—	—
Yes	0.54	0.45, 0.65	<0.001
Charlson Comorbidity Index	1.00	0.98, 1.01	0.7
Immunocompromised	1.25	1.08, 1.46	0.003
Population Density of Census Tract	1.00	1.00, 1.00	<0.001
Socioeconomic Status	0.28	0.19, 0.41	<0.001
Household Composition	0.58	0.41, 0.82	0.002
Minority and Language Status	1.53	1.19, 1.98	<0.001
Housing and Transportation	1.25	0.95, 1.65	0.12
1 OR = Odds Ratio, CI = Confidence Interval

# Hosmer-Lemeshow Goodness-of-Fit Test
hltest(covid_all_themes)

   The Hosmer-Lemeshow goodness-of-fit test

         Statistic =  12.25271 
degrees of freedom =  8 
           p-value =  0.14028 

# C-Statistic/AUROC 
Cstat(covid_all_themes)

[1] 0.6892033

# Model performance 
model_performance(covid_all_themes)

# Indices of model performance

AIC      |     AICc |      BIC | Tjur's R2 |  RMSE | Sigma | Log_loss | Score_log |   PCP
-----------------------------------------------------------------------------------------
5597.959 | 5598.166 | 5749.389 |     0.089 | 0.414 | 1.021 |    0.519 |      -Inf | 0.657

performance::check_model(covid_all_themes, panel = TRUE)

Variable `Component` is not in your data frame :/

# Margins 
cplot(covid_all_themes, "RPL_THEME1", what = "prediction", main = "Predicted Likelihood of COVID Vaccine Given THEME1")

cplot(covid_all_themes, "RPL_THEME2", what = "prediction", main = "Predicted Likelihood of COVID Vaccine Given THEME2")

cplot(covid_all_themes, "RPL_THEME3", what = "prediction", main = "Predicted Likelihood of COVID Vaccine Given THEME3")

cplot(covid_all_themes, "RPL_THEME4", what = "prediction", main = "Predicted Likelihood of COVID Vaccine Given THEME4")

Sensitivity Analysis with r_pct

sens_covid_SVI4 <- glm(full_vax ~ ibd_3 + age_yrs + gender + race_5 + ethnic_3 + mstat_5 + relig_affil
                      + lang_3 + act_tob + max_ch + IC + insurance 
                      + pop_dens + r_pct + RPL_4,
              family = "binomial", 
              data = cov_clean1)
summary(sens_covid_SVI4 )

Call:
glm(formula = full_vax ~ ibd_3 + age_yrs + gender + race_5 + 
    ethnic_3 + mstat_5 + relig_affil + lang_3 + act_tob + max_ch + 
    IC + insurance + pop_dens + r_pct + RPL_4, family = "binomial", 
    data = cov_clean1)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-2.8439   0.2206   0.5586   0.7490   1.8961  

Coefficients:
                                         Estimate Std. Error z value Pr(>|z|)    
(Intercept)                             1.389e+00  3.625e-01   3.831 0.000128 ***
ibd_3UC                                -2.152e-02  7.709e-02  -0.279 0.780170    
ibd_3Unknown                            1.445e+00  1.067e+00   1.354 0.175738    
age_yrs                                 2.804e-02  2.988e-03   9.386  < 2e-16 ***
genderFemale                            1.057e-01  7.284e-02   1.452 0.146570    
race_5Black                            -3.387e-01  1.481e-01  -2.287 0.022219 *  
race_5Asian or Pacific Islander         1.142e+00  3.787e-01   3.016 0.002565 ** 
race_5American Indian or Alaska Native  3.053e-01  6.815e-01   0.448 0.654205    
race_5Other                            -1.127e-01  2.163e-01  -0.521 0.602441    
ethnic_3Non-Hispanic                   -4.003e-01  2.755e-01  -1.453 0.146286    
mstat_5Single                           1.206e-01  9.442e-02   1.277 0.201633    
mstat_5Divorced/Separated              -1.391e-02  1.863e-01  -0.075 0.940479    
mstat_5Widowed                          2.486e-01  3.114e-01   0.798 0.424677    
relig_affilNo                           1.264e-01  7.479e-02   1.689 0.091161 .  
lang_3Other                             3.290e-01  4.386e-01   0.750 0.453203    
act_tobYes                             -5.146e-01  1.067e-01  -4.824 1.41e-06 ***
max_ch                                  8.447e-03  8.332e-03   1.014 0.310712    
IC                                      2.094e-01  8.380e-02   2.499 0.012463 *  
insuranceMedicaid                      -9.717e-01  1.009e-01  -9.631  < 2e-16 ***
insuranceMedicare                      -4.748e-01  1.189e-01  -3.992 6.55e-05 ***
insuranceOther Governmental            -1.872e+00  4.295e-01  -4.359 1.31e-05 ***
insuranceOther                         -3.217e-01  4.630e-01  -0.695 0.487113    
pop_dens                                2.215e-05  2.074e-05   1.068 0.285671    
r_pct                                  -2.329e-02  2.812e-03  -8.283  < 2e-16 ***
RPL_4Second                            -3.233e-02  8.878e-02  -0.364 0.715744    
RPL_4Third                             -4.128e-01  1.008e-01  -4.095 4.23e-05 ***
RPL_4Fourth                            -9.167e-01  1.314e-01  -6.979 2.97e-12 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 5327.9  on 4753  degrees of freedom
Residual deviance: 4746.6  on 4727  degrees of freedom
  (4102 observations deleted due to missingness)
AIC: 4800.6

Number of Fisher Scoring iterations: 5

broom::glance(sens_covid_SVI4 )

broom::tidy(sens_covid_SVI4 , exponentiate = TRUE)

tbl_regression(sens_covid_SVI4, exponentiate = TRUE)

Characteristic	OR1	95% CI1	p-value
IBD Diagnosis
CD	—	—
UC	0.98	0.84, 1.14	0.8
Unknown	4.24	0.78, 79.4	0.2
Age	1.03	1.02, 1.03	<0.001
Gender
Male	—	—
Female	1.11	0.96, 1.28	0.15
Race
White	—	—
Black	0.71	0.53, 0.96	0.022
Asian or Pacific Islander	3.13	1.59, 7.13	0.003
American Indian or Alaska Native	1.36	0.40, 6.30	0.7
Other	0.89	0.59, 1.38	0.6
Ethnicity
Hispanic	—	—
Non-Hispanic	0.67	0.38, 1.13	0.15
Marital Status
Married	—	—
Single	1.13	0.94, 1.36	0.2
Divorced/Separated	0.99	0.69, 1.43	>0.9
Widowed	1.28	0.72, 2.45	0.4
Any Religious Affiliation
Yes	—	—
No	1.13	0.98, 1.31	0.091
Preferred Language
English	—	—
Other	1.39	0.61, 3.44	0.5
Active Tobacco Use
No	—	—
Yes	0.60	0.49, 0.74	<0.001
Charlson Comorbidity Index	1.01	0.99, 1.03	0.3
Immunocompromised	1.23	1.05, 1.45	0.012
Insurance Type
Private Insurance	—	—
Medicaid	0.38	0.31, 0.46	<0.001
Medicare	0.62	0.49, 0.79	<0.001
Other Governmental	0.15	0.06, 0.36	<0.001
Other	0.72	0.31, 1.93	0.5
Population Density of Census Tract	1.00	1.00, 1.00	0.3
Percentage Republican Voters in Census Tract	0.98	0.97, 0.98	<0.001
SVI by Quartile
First	—	—
Second	0.97	0.81, 1.15	0.7
Third	0.66	0.54, 0.81	<0.001
Fourth	0.40	0.31, 0.52	<0.001
1 OR = Odds Ratio, CI = Confidence Interval

# Hosmer-Lemeshow Goodness-of-Fit Test
hltest(sens_covid_SVI4)

   The Hosmer-Lemeshow goodness-of-fit test

         Statistic =  18.40193 
degrees of freedom =  9 
           p-value =  0.030787 

# C-Statistic/AUROC 
Cstat(sens_covid_SVI4)

[1] 0.7205842

# Model performance 
model_performance(sens_covid_SVI4)

# Indices of model performance

AIC      |     AICc |      BIC | Tjur's R2 |  RMSE | Sigma | Log_loss | Score_log |   PCP
-----------------------------------------------------------------------------------------
4800.633 | 4800.953 | 4975.235 |     0.125 | 0.404 | 1.002 |    0.499 |      -Inf | 0.673

performance::check_model(sens_covid_SVI4)

Variable `Component` is not in your data frame :/

# Margins 
cplot(sens_covid_SVI4, "RPL_4", what = "prediction", main = "Predicted Likelihood of COVID Vaccine Given SVI Quartile")

Multivariable Models: Vaccination + Boosted

SVI Continuous

covid_SVI_boost <- glm(vax_boost ~ ibd_3 + age_yrs + gender + race_5 + ethnic_3 + mstat_5 + relig_affil
                      + lang_3 + act_tob + max_ch + IC + insurance
                      + pop_dens + RPL_THEMES,
              family = "binomial", 
              data = cov_clean1)
summary(covid_SVI_boost )

Call:
glm(formula = vax_boost ~ ibd_3 + age_yrs + gender + race_5 + 
    ethnic_3 + mstat_5 + relig_affil + lang_3 + act_tob + max_ch + 
    IC + insurance + pop_dens + RPL_THEMES, family = "binomial", 
    data = cov_clean1)

Deviance Residuals: 
   Min      1Q  Median      3Q     Max  
-2.235  -1.124   0.696   1.006   2.232  

Coefficients:
                                         Estimate Std. Error z value Pr(>|z|)    
(Intercept)                            -8.412e-01  2.691e-01  -3.126  0.00177 ** 
ibd_3UC                                 3.144e-04  6.220e-02   0.005  0.99597    
ibd_3Unknown                            5.371e-02  5.569e-01   0.096  0.92318    
age_yrs                                 2.865e-02  2.406e-03  11.904  < 2e-16 ***
genderFemale                            1.070e-01  5.937e-02   1.802  0.07153 .  
race_5Black                            -9.653e-02  1.231e-01  -0.784  0.43300    
race_5Asian or Pacific Islander         1.178e+00  2.269e-01   5.193 2.07e-07 ***
race_5American Indian or Alaska Native -2.154e-01  4.845e-01  -0.445  0.65658    
race_5Other                            -2.788e-01  1.848e-01  -1.509  0.13129    
ethnic_3Non-Hispanic                   -2.635e-01  2.227e-01  -1.183  0.23662    
mstat_5Single                           8.165e-02  7.519e-02   1.086  0.27748    
mstat_5Divorced/Separated              -1.172e-01  1.548e-01  -0.757  0.44888    
mstat_5Widowed                         -2.395e-01  2.073e-01  -1.155  0.24805    
relig_affilNo                           1.459e-01  6.095e-02   2.394  0.01667 *  
lang_3Other                            -6.628e-01  3.603e-01  -1.840  0.06584 .  
act_tobYes                             -4.195e-01  9.659e-02  -4.343 1.41e-05 ***
max_ch                                  1.038e-02  6.487e-03   1.600  0.10952    
IC                                      1.899e-01  6.719e-02   2.827  0.00470 ** 
insuranceMedicaid                      -9.635e-01  9.635e-02 -10.000  < 2e-16 ***
insuranceMedicare                      -2.966e-01  9.289e-02  -3.193  0.00141 ** 
insuranceOther Governmental            -1.072e+00  4.251e-01  -2.521  0.01170 *  
insuranceOther                         -2.286e-01  3.884e-01  -0.589  0.55617    
pop_dens                                8.752e-05  1.378e-05   6.351 2.13e-10 ***
RPL_THEMES                             -1.065e+00  1.260e-01  -8.451  < 2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 7428.8  on 5388  degrees of freedom
Residual deviance: 6757.6  on 5365  degrees of freedom
  (3467 observations deleted due to missingness)
AIC: 6805.6

Number of Fisher Scoring iterations: 4

broom::glance(covid_SVI_boost )

broom::tidy(covid_SVI_boost , exponentiate = TRUE)

tbl_regression(covid_SVI_boost, exponentiate = TRUE)

Characteristic	OR1	95% CI1	p-value
IBD Diagnosis
CD	—	—
UC	1.00	0.89, 1.13	>0.9
Unknown	1.06	0.36, 3.33	>0.9
Age	1.03	1.02, 1.03	<0.001
Gender
Male	—	—
Female	1.11	0.99, 1.25	0.072
Race
White	—	—
Black	0.91	0.71, 1.16	0.4
Asian or Pacific Islander	3.25	2.11, 5.16	<0.001
American Indian or Alaska Native	0.81	0.31, 2.11	0.7
Other	0.76	0.53, 1.09	0.13
Ethnicity
Hispanic	—	—
Non-Hispanic	0.77	0.50, 1.19	0.2
Marital Status
Married	—	—
Single	1.09	0.94, 1.26	0.3
Divorced/Separated	0.89	0.66, 1.21	0.4
Widowed	0.79	0.53, 1.19	0.2
Any Religious Affiliation
Yes	—	—
No	1.16	1.03, 1.30	0.017
Preferred Language
English	—	—
Other	0.52	0.25, 1.04	0.066
Active Tobacco Use
No	—	—
Yes	0.66	0.54, 0.79	<0.001
Charlson Comorbidity Index	1.01	1.00, 1.02	0.11
Immunocompromised	1.21	1.06, 1.38	0.005
Insurance Type
Private Insurance	—	—
Medicaid	0.38	0.32, 0.46	<0.001
Medicare	0.74	0.62, 0.89	0.001
Other Governmental	0.34	0.14, 0.77	0.012
Other	0.80	0.37, 1.72	0.6
Population Density of Census Tract	1.00	1.00, 1.00	<0.001
Social Vulnerability Index	0.34	0.27, 0.44	<0.001
1 OR = Odds Ratio, CI = Confidence Interval

# Hosmer-Lemeshow Goodness-of-Fit Test
hltest(covid_SVI_boost)

   The Hosmer-Lemeshow goodness-of-fit test

         Statistic =  8.13924 
degrees of freedom =  8 
           p-value =  0.41999 

# C-Statistic/AUROC 
Cstat(covid_SVI_boost)

[1] 0.6936377

# Model performance 
model_performance(covid_SVI_boost)

# Indices of model performance

AIC      |     AICc |      BIC | Tjur's R2 |  RMSE | Sigma | Log_loss | Score_log | Score_spherical |   PCP
-----------------------------------------------------------------------------------------------------------
6805.559 | 6805.783 | 6963.770 |     0.119 | 0.468 | 1.122 |    0.627 |      -Inf |       1.856e-04 | 0.563

performance::check_model(covid_SVI_boost)

Variable `Component` is not in your data frame :/

# Margins 
cplot(covid_SVI_boost, "RPL_THEMES", what = "prediction", main = "Predicted Likelihood of COVID Vaccine Given SVI")

SVI Quartiles

covid_SVI4_boost <- glm(vax_boost ~ ibd_3 + age_yrs + gender + race_5 + ethnic_3 + mstat_5 + relig_affil
                      + lang_3 + act_tob + max_ch + IC + insurance
                      + pop_dens + RPL_4,
              family = "binomial", 
              data = cov_clean1)
summary(covid_SVI4_boost )

Call:
glm(formula = vax_boost ~ ibd_3 + age_yrs + gender + race_5 + 
    ethnic_3 + mstat_5 + relig_affil + lang_3 + act_tob + max_ch + 
    IC + insurance + pop_dens + RPL_4, family = "binomial", data = cov_clean1)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-2.1959  -1.1273   0.7119   1.0031   2.2797  

Coefficients:
                                         Estimate Std. Error z value Pr(>|z|)    
(Intercept)                            -9.962e-01  2.687e-01  -3.708 0.000209 ***
ibd_3UC                                -5.285e-03  6.222e-02  -0.085 0.932311    
ibd_3Unknown                            4.466e-02  5.571e-01   0.080 0.936098    
age_yrs                                 2.851e-02  2.405e-03  11.853  < 2e-16 ***
genderFemale                            1.091e-01  5.937e-02   1.837 0.066229 .  
race_5Black                            -1.006e-01  1.239e-01  -0.812 0.416558    
race_5Asian or Pacific Islander         1.176e+00  2.268e-01   5.185 2.16e-07 ***
race_5American Indian or Alaska Native -2.747e-01  4.861e-01  -0.565 0.571983    
race_5Other                            -2.880e-01  1.851e-01  -1.556 0.119696    
ethnic_3Non-Hispanic                   -2.632e-01  2.232e-01  -1.180 0.238124    
mstat_5Single                           7.485e-02  7.522e-02   0.995 0.319714    
mstat_5Divorced/Separated              -1.369e-01  1.550e-01  -0.883 0.377259    
mstat_5Widowed                         -2.555e-01  2.076e-01  -1.231 0.218449    
relig_affilNo                           1.507e-01  6.096e-02   2.471 0.013464 *  
lang_3Other                            -6.467e-01  3.616e-01  -1.789 0.073686 .  
act_tobYes                             -4.264e-01  9.664e-02  -4.413 1.02e-05 ***
max_ch                                  1.071e-02  6.491e-03   1.649 0.099069 .  
IC                                      1.823e-01  6.712e-02   2.716 0.006611 ** 
insuranceMedicaid                      -9.716e-01  9.631e-02 -10.088  < 2e-16 ***
insuranceMedicare                      -3.053e-01  9.289e-02  -3.287 0.001014 ** 
insuranceOther Governmental            -1.056e+00  4.256e-01  -2.482 0.013073 *  
insuranceOther                         -2.406e-01  3.879e-01  -0.620 0.535147    
pop_dens                                8.648e-05  1.380e-05   6.266 3.69e-10 ***
RPL_4Second                            -1.065e-01  6.916e-02  -1.540 0.123636    
RPL_4Third                             -4.837e-01  8.355e-02  -5.789 7.07e-09 ***
RPL_4Fourth                            -8.558e-01  1.216e-01  -7.041 1.91e-12 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 7428.8  on 5388  degrees of freedom
Residual deviance: 6758.8  on 5363  degrees of freedom
  (3467 observations deleted due to missingness)
AIC: 6810.8

Number of Fisher Scoring iterations: 4

broom::glance(covid_SVI4_boost )

broom::tidy(covid_SVI4_boost , exponentiate = TRUE)

tbl_regression(covid_SVI4_boost, exponentiate = TRUE)

Characteristic	OR1	95% CI1	p-value
IBD Diagnosis
CD	—	—
UC	1.0	0.88, 1.12	>0.9
Unknown	1.05	0.36, 3.31	>0.9
Age	1.03	1.02, 1.03	<0.001
Gender
Male	—	—
Female	1.12	0.99, 1.25	0.066
Race
White	—	—
Black	0.90	0.71, 1.15	0.4
Asian or Pacific Islander	3.24	2.11, 5.14	<0.001
American Indian or Alaska Native	0.76	0.29, 1.99	0.6
Other	0.75	0.52, 1.08	0.12
Ethnicity
Hispanic	—	—
Non-Hispanic	0.77	0.50, 1.19	0.2
Marital Status
Married	—	—
Single	1.08	0.93, 1.25	0.3
Divorced/Separated	0.87	0.64, 1.18	0.4
Widowed	0.77	0.52, 1.17	0.2
Any Religious Affiliation
Yes	—	—
No	1.16	1.03, 1.31	0.013
Preferred Language
English	—	—
Other	0.52	0.25, 1.06	0.074
Active Tobacco Use
No	—	—
Yes	0.65	0.54, 0.79	<0.001
Charlson Comorbidity Index	1.01	1.00, 1.02	0.10
Immunocompromised	1.20	1.05, 1.37	0.007
Insurance Type
Private Insurance	—	—
Medicaid	0.38	0.31, 0.46	<0.001
Medicare	0.74	0.61, 0.88	0.001
Other Governmental	0.35	0.14, 0.78	0.013
Other	0.79	0.37, 1.70	0.5
Population Density of Census Tract	1.00	1.00, 1.00	<0.001
SVI by Quartile
First	—	—
Second	0.90	0.79, 1.03	0.12
Third	0.62	0.52, 0.73	<0.001
Fourth	0.42	0.33, 0.54	<0.001
1 OR = Odds Ratio, CI = Confidence Interval

# Hosmer-Lemeshow Goodness-of-Fit Test
hltest(covid_SVI4_boost)

   The Hosmer-Lemeshow goodness-of-fit test

         Statistic =  7.17048 
degrees of freedom =  8 
           p-value =  0.51836 

# C-Statistic/AUROC 
Cstat(covid_SVI4_boost)

[1] 0.6934182

# Model performance 
model_performance(covid_SVI4_boost)

# Indices of model performance

AIC      |     AICc |      BIC | Tjur's R2 |  RMSE | Sigma | Log_loss | Score_log | Score_spherical |   PCP
-----------------------------------------------------------------------------------------------------------
6810.841 | 6811.103 | 6982.236 |     0.118 | 0.468 | 1.123 |    0.627 |      -Inf |       1.856e-04 | 0.563

performance::check_model(covid_SVI4_boost)

Variable `Component` is not in your data frame :/

# Margins 
cplot(covid_SVI4_boost, "RPL_4", what = "prediction", main = "Predicted Likelihood of COVID Vaccine Given SVI")

All Themes

covid_all_themes_boost <- glm(vax_boost ~ ibd_3 + age_yrs + gender + race_5 + ethnic_3 + 
                      lang_3 + mstat_5 + relig_affil + act_tob + max_ch + IC
                      + pop_dens + RPL_THEME1 + RPL_THEME2 + RPL_THEME3
                      + RPL_THEME4,
              family = "binomial",
              data = cov_clean1)
summary(covid_all_themes_boost )

Call:
glm(formula = vax_boost ~ ibd_3 + age_yrs + gender + race_5 + 
    ethnic_3 + lang_3 + mstat_5 + relig_affil + act_tob + max_ch + 
    IC + pop_dens + RPL_THEME1 + RPL_THEME2 + RPL_THEME3 + RPL_THEME4, 
    family = "binomial", data = cov_clean1)

Deviance Residuals: 
   Min      1Q  Median      3Q     Max  
-2.273  -1.095   0.663   1.028   2.227  

Coefficients:
                                         Estimate Std. Error z value Pr(>|z|)    
(Intercept)                            -8.931e-01  2.730e-01  -3.272 0.001067 ** 
ibd_3UC                                -1.352e-03  6.220e-02  -0.022 0.982653    
ibd_3Unknown                            3.721e-02  5.590e-01   0.067 0.946931    
age_yrs                                 2.787e-02  2.199e-03  12.675  < 2e-16 ***
genderFemale                            8.148e-02  5.942e-02   1.371 0.170308    
race_5Black                            -1.700e-01  1.233e-01  -1.378 0.168069    
race_5Asian or Pacific Islander         9.548e-01  2.296e-01   4.159 3.20e-05 ***
race_5American Indian or Alaska Native -3.139e-01  4.951e-01  -0.634 0.526123    
race_5Other                            -5.096e-01  1.826e-01  -2.791 0.005258 ** 
ethnic_3Non-Hispanic                   -2.598e-01  2.218e-01  -1.171 0.241614    
lang_3Other                            -1.034e+00  3.566e-01  -2.899 0.003748 ** 
mstat_5Single                          -8.760e-02  7.312e-02  -1.198 0.230879    
mstat_5Divorced/Separated              -3.271e-01  1.522e-01  -2.149 0.031611 *  
mstat_5Widowed                         -3.408e-01  2.074e-01  -1.644 0.100254    
relig_affilNo                           1.182e-01  6.093e-02   1.939 0.052475 .  
act_tobYes                             -4.926e-01  9.546e-02  -5.160 2.47e-07 ***
max_ch                                  7.238e-03  6.462e-03   1.120 0.262709    
IC                                      2.405e-01  6.740e-02   3.568 0.000359 ***
pop_dens                                6.641e-05  1.501e-05   4.425 9.63e-06 ***
RPL_THEME1                             -1.330e+00  1.805e-01  -7.366 1.76e-13 ***
RPL_THEME2                             -5.759e-01  1.558e-01  -3.695 0.000220 ***
RPL_THEME3                              4.490e-01  1.126e-01   3.988 6.67e-05 ***
RPL_THEME4                              2.642e-01  1.252e-01   2.110 0.034827 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 7372.0  on 5344  degrees of freedom
Residual deviance: 6727.1  on 5322  degrees of freedom
  (3511 observations deleted due to missingness)
AIC: 6773.1

Number of Fisher Scoring iterations: 4

broom::glance(covid_all_themes_boost )

broom::tidy(covid_all_themes_boost , exponentiate = TRUE)

tbl_regression(covid_all_themes_boost, exponentiate = TRUE)

Characteristic	OR1	95% CI1	p-value
IBD Diagnosis
CD	—	—
UC	1.00	0.88, 1.13	>0.9
Unknown	1.04	0.35, 3.28	>0.9
Age	1.03	1.02, 1.03	<0.001
Gender
Male	—	—
Female	1.08	0.97, 1.22	0.2
Race
White	—	—
Black	0.84	0.66, 1.07	0.2
Asian or Pacific Islander	2.60	1.68, 4.14	<0.001
American Indian or Alaska Native	0.73	0.27, 1.95	0.5
Other	0.60	0.42, 0.86	0.005
Ethnicity
Hispanic	—	—
Non-Hispanic	0.77	0.50, 1.19	0.2
Preferred Language
English	—	—
Other	0.36	0.17, 0.71	0.004
Marital Status
Married	—	—
Single	0.92	0.79, 1.06	0.2
Divorced/Separated	0.72	0.53, 0.97	0.032
Widowed	0.71	0.48, 1.07	0.10
Any Religious Affiliation
Yes	—	—
No	1.13	1.00, 1.27	0.052
Active Tobacco Use
No	—	—
Yes	0.61	0.51, 0.74	<0.001
Charlson Comorbidity Index	1.01	1.0, 1.02	0.3
Immunocompromised	1.27	1.11, 1.45	<0.001
Population Density of Census Tract	1.00	1.00, 1.00	<0.001
Socioeconomic Status	0.26	0.19, 0.38	<0.001
Household Composition	0.56	0.41, 0.76	<0.001
Minority and Language Status	1.57	1.26, 1.95	<0.001
Housing and Transportation	1.30	1.02, 1.67	0.035
1 OR = Odds Ratio, CI = Confidence Interval

# Hosmer-Lemeshow Goodness-of-Fit Test
hltest(covid_all_themes_boost)

   The Hosmer-Lemeshow goodness-of-fit test

         Statistic =  16.87449 
degrees of freedom =  8 
           p-value =  0.031443 

# C-Statistic/AUROC 
Cstat(covid_all_themes_boost)

[1] 0.6933198

# Model performance 
model_performance(covid_all_themes_boost)

# Indices of model performance

AIC      |     AICc |      BIC | Tjur's R2 |  RMSE | Sigma | Log_loss | Score_log | Score_spherical |   PCP
-----------------------------------------------------------------------------------------------------------
6773.148 | 6773.355 | 6924.578 |     0.115 | 0.469 | 1.124 |    0.629 |      -Inf |       1.871e-04 | 0.561

performance::check_model(covid_all_themes_boost, panel = TRUE)

Variable `Component` is not in your data frame :/

# Margins 
cplot(covid_all_themes_boost, "RPL_THEME1", what = "prediction", main = "Predicted Likelihood of COVID Vaccine Given THEME1")

cplot(covid_all_themes_boost, "RPL_THEME2", what = "prediction", main = "Predicted Likelihood of COVID Vaccine Given THEME2")

cplot(covid_all_themes_boost, "RPL_THEME3", what = "prediction", main = "Predicted Likelihood of COVID Vaccine Given THEME3")

cplot(covid_all_themes_boost, "RPL_THEME4", what = "prediction", main = "Predicted Likelihood of COVID Vaccine Given THEME4")

Sensitivity Analysis with r_pct

sens_covid_boost <- glm(vax_boost ~ ibd_3 + age_yrs + gender + race_5 + ethnic_3 + mstat_5 + relig_affil
                      + lang_3 + act_tob + max_ch + IC + insurance 
                      + pop_dens + r_pct + RPL_4,
              family = "binomial", 
              data = cov_clean1)
summary(sens_covid_boost )

Call:
glm(formula = vax_boost ~ ibd_3 + age_yrs + gender + race_5 + 
    ethnic_3 + mstat_5 + relig_affil + lang_3 + act_tob + max_ch + 
    IC + insurance + pop_dens + r_pct + RPL_4, family = "binomial", 
    data = cov_clean1)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-2.3149  -1.0816   0.6091   0.9816   2.3756  

Coefficients:
                                         Estimate Std. Error z value Pr(>|z|)    
(Intercept)                             4.098e-01  3.130e-01   1.309 0.190521    
ibd_3UC                                -1.073e-02  6.763e-02  -0.159 0.873955    
ibd_3Unknown                            3.123e-01  5.959e-01   0.524 0.600177    
age_yrs                                 2.890e-02  2.611e-03  11.067  < 2e-16 ***
genderFemale                            1.218e-01  6.427e-02   1.896 0.057987 .  
race_5Black                            -5.338e-01  1.377e-01  -3.878 0.000105 ***
race_5Asian or Pacific Islander         8.774e-01  2.437e-01   3.600 0.000318 ***
race_5American Indian or Alaska Native -6.365e-02  5.515e-01  -0.115 0.908130    
race_5Other                            -3.972e-01  1.970e-01  -2.016 0.043764 *  
ethnic_3Non-Hispanic                   -2.297e-01  2.394e-01  -0.960 0.337190    
mstat_5Single                           3.467e-02  8.176e-02   0.424 0.671533    
mstat_5Divorced/Separated              -1.366e-01  1.681e-01  -0.813 0.416401    
mstat_5Widowed                         -2.994e-01  2.292e-01  -1.306 0.191449    
relig_affilNo                           1.436e-01  6.614e-02   2.171 0.029952 *  
lang_3Other                            -6.683e-01  3.798e-01  -1.759 0.078496 .  
act_tobYes                             -4.139e-01  1.052e-01  -3.934 8.35e-05 ***
max_ch                                  1.565e-02  7.056e-03   2.217 0.026595 *  
IC                                      2.692e-01  7.335e-02   3.670 0.000242 ***
insuranceMedicaid                      -9.260e-01  1.042e-01  -8.890  < 2e-16 ***
insuranceMedicare                      -2.681e-01  1.012e-01  -2.650 0.008041 ** 
insuranceOther Governmental            -1.203e+00  4.495e-01  -2.676 0.007452 ** 
insuranceOther                         -2.111e-01  4.077e-01  -0.518 0.604559    
pop_dens                               -2.479e-05  1.720e-05  -1.441 0.149480    
r_pct                                  -2.772e-02  2.379e-03 -11.650  < 2e-16 ***
RPL_4Second                            -1.442e-01  7.568e-02  -1.906 0.056713 .  
RPL_4Third                             -4.590e-01  9.094e-02  -5.047 4.48e-07 ***
RPL_4Fourth                            -9.812e-01  1.303e-01  -7.531 5.04e-14 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 6549.3  on 4753  degrees of freedom
Residual deviance: 5792.1  on 4727  degrees of freedom
  (4102 observations deleted due to missingness)
AIC: 5846.1

Number of Fisher Scoring iterations: 4

broom::glance(sens_covid_boost )

broom::tidy(sens_covid_boost , exponentiate = TRUE)

tbl_regression(sens_covid_boost, exponentiate = TRUE)

Characteristic	OR1	95% CI1	p-value
IBD Diagnosis
CD	—	—
UC	0.99	0.87, 1.13	0.9
Unknown	1.37	0.44, 4.77	0.6
Age	1.03	1.02, 1.03	<0.001
Gender
Male	—	—
Female	1.13	1.00, 1.28	0.058
Race
White	—	—
Black	0.59	0.45, 0.77	<0.001
Asian or Pacific Islander	2.40	1.51, 3.95	<0.001
American Indian or Alaska Native	0.94	0.32, 2.88	>0.9
Other	0.67	0.46, 0.99	0.044
Ethnicity
Hispanic	—	—
Non-Hispanic	0.79	0.50, 1.27	0.3
Marital Status
Married	—	—
Single	1.04	0.88, 1.22	0.7
Divorced/Separated	0.87	0.63, 1.21	0.4
Widowed	0.74	0.48, 1.17	0.2
Any Religious Affiliation
Yes	—	—
No	1.15	1.01, 1.31	0.030
Preferred Language
English	—	—
Other	0.51	0.24, 1.07	0.078
Active Tobacco Use
No	—	—
Yes	0.66	0.54, 0.81	<0.001
Charlson Comorbidity Index	1.02	1.00, 1.03	0.027
Immunocompromised	1.31	1.13, 1.51	<0.001
Insurance Type
Private Insurance	—	—
Medicaid	0.40	0.32, 0.49	<0.001
Medicare	0.76	0.63, 0.93	0.008
Other Governmental	0.30	0.12, 0.71	0.007
Other	0.81	0.36, 1.82	0.6
Population Density of Census Tract	1.00	1.00, 1.00	0.15
Percentage Republican Voters in Census Tract	0.97	0.97, 0.98	<0.001
SVI by Quartile
First	—	—
Second	0.87	0.75, 1.00	0.057
Third	0.63	0.53, 0.76	<0.001
Fourth	0.37	0.29, 0.48	<0.001
1 OR = Odds Ratio, CI = Confidence Interval

# Hosmer-Lemeshow Goodness-of-Fit Test
hltest(sens_covid_boost)

   The Hosmer-Lemeshow goodness-of-fit test

         Statistic =  10.56831 
degrees of freedom =  9 
           p-value =  0.30646 

# C-Statistic/AUROC 
Cstat(sens_covid_boost)

[1] 0.7210684

# Model performance 
model_performance(sens_covid_boost)

# Indices of model performance

AIC      |     AICc |      BIC | Tjur's R2 |  RMSE | Sigma | Log_loss | Score_log | Score_spherical |   PCP
-----------------------------------------------------------------------------------------------------------
5846.120 | 5846.439 | 6020.722 |     0.149 | 0.459 | 1.107 |    0.609 |      -Inf |       2.103e-04 | 0.578

performance::check_model(sens_covid_boost)

Variable `Component` is not in your data frame :/

# Margins 
cplot(sens_covid_boost, "RPL_4", what = "prediction", main = "Predicted Likelihood of COVID Vaccine Given SVI Quartile")