Healthcare Utilization

you need to hand-load vax_dataset, ibd_first_last_encounter,admit_count, ED_visits_11_22_22, IBD_surgery_11_17_22 prior to running this code

Load Packages

library(tidyverse)

── Attaching packages ───────────────────────────────────── tidyverse 1.3.2 ──✔ ggplot2 3.3.6      ✔ purrr   0.3.5 
✔ tibble  3.1.8      ✔ dplyr   1.0.10
✔ tidyr   1.2.1      ✔ stringr 1.4.1 
✔ readr   2.1.3      ✔ forcats 0.5.2 ── Conflicts ──────────────────────────────────────── tidyverse_conflicts() ──
✖ lubridate::as.difftime() masks base::as.difftime()
✖ lubridate::date()        masks base::date()
✖ dplyr::filter()          masks stats::filter()
✖ lubridate::intersect()   masks base::intersect()
✖ dplyr::lag()             masks stats::lag()
✖ lubridate::setdiff()     masks base::setdiff()
✖ lubridate::union()       masks base::union()

library(codebookr)
library(broom)
library(performance)
library(gt)
library(gtsummary)

#StandWithUkraine

library(janitor)

Attaching package: ‘janitor’

The following objects are masked from ‘package:stats’:

    chisq.test, fisher.test

library(forcats)
library(here)

here() starts at /nfs/isilon/path-ibdx/jsheehan

view(vax_dataset)
view(ibd_first_last_encounter)

Data Cleaning

Make person years variable, exclude those with less than 1 year of data

ibd_first_last_encounter$pdays <- ibd_first_last_encounter$last_enc - ibd_first_last_encounter$first_enc

ibd_first_last_encounter$pdays = as.numeric(ibd_first_last_encounter$pdays)

ibd_first_last_encounter %>% 
  mutate(pyears = pdays/365) -> ibd_pyears

pyears_ibd <- subset(ibd_pyears, pyears >= 1)

pyears_ibd %>% 
  dplyr::select(sbj_id, pyears) -> sbj_pyears

Merge with vax_data_set

sbj_pyears$sbj_id = as.numeric(sbj_pyears$sbj_id)
use_1 <- left_join(x= vax_dataset, y = sbj_pyears, by = "sbj_id")

  mutate(admit_2 = case_when(IBD_hosp>= 1 ~ '1',TRUE ~ "0")) %>% 
  dplyr::select(sbj_id, admit_2) -> admit_long
ddply(admit_long,.(sbj_id, admit_2),nrow) -> admit_count

Make new hospitalization variable

IBD_hospitalizations_11_30_22 %>% 
  dplyr::select(sbj_id, IBD_hosp) -> sbj_admits
sbj_admits[sbj_admits == 0] <- NA
sbj_admits.omit <- na.omit(sbj_admits)

ddply(sbj_admits.omit,.(sbj_id, IBD_hosp),nrow) -> admit_count
  

Merge admit_count data

admit_count$sbj_id = as.numeric(admit_count$sbj_id)
use_1$sbj_id = as.numeric(use_1$sbj_id)
left_join(x = use_1, y = admit_count, by = "sbj_id") -> utilize_2

Rename V1 to admit_count

colnames(utilize_2)[colnames(utilize_2) == "V1"] <- "num_admits"
utilize_2$num_admits = as.numeric(utilize_2$num_admits)

utilize_2 %>% 
  mutate(num_admits = coalesce(num_admits, 0)) -> utilize_3

Make ED Count Variable

ED_visits_11_22_22$sbj_id = as.factor(ED_visits_11_22_22$sbj_id)
ED_visits_11_22_22 %>% 
  mutate(ed_2 = case_when(ed_days>= 0 ~ '1',TRUE ~ "0")) %>% 
  dplyr::select(sbj_id, ed_2) -> ED_visits
ddply(ED_visits,.(sbj_id, ed_2),nrow) -> ED_count

Merge ED count data

utilize_3$sbj_id = as.numeric(utilize_3$sbj_id)
ED_count$sbj_id = as.numeric(ED_count$sbj_id)
left_join(x = utilize_3, y = ED_count, by = "sbj_id") -> use_4

Rename ED visit variable

colnames(use_4)[colnames(use_4) == "V1"] <- "ed_count"
use_4$ed_count = as.numeric(use_4$ed_count)
use_4 %>% 
  mutate(ed_count = coalesce(ed_count, 0)) -> use_5

Make surgery count variable

IBD_surgery_11_17_22 %>% 
  mutate(surg_2 = case_when(sbj_id>= 0 ~ '1',TRUE ~ "0")) %>% 
  dplyr::select(sbj_id, surg_2) -> surgeries
ddply(surgeries,.(sbj_id, surg_2),nrow) -> surg_count

Merge surgery data

use_5$sbj_id = as.numeric(use_5$sbj_id)
surg_count$sbj_id = as.numeric(surg_count$sbj_id)
left_join(x = use_5, y = surg_count, by = "sbj_id") -> use_6

Rename surg count variable and get change NA to zeros

colnames(use_6)[colnames(use_6) == "V1"] <- "surg_count"
use_6$surg_count = as.numeric(use_6$surg_count)
use_6 %>% 
  mutate(surg_count = coalesce(surg_count, 0)) -> use_7
      

Make readmission variable

IBD_hospitalizations_11_30_22 %>% 
  dplyr::select(sbj_id, readmit) -> readmits
readmits$readmit<-ifelse(readmits$readmit=="Y",1,0)
readmits[readmits == 0] <- NA
readmits.omit <- na.omit(readmits)
ddply(readmits.omit,.(sbj_id, readmit),nrow) -> readmit_count

Merge readmit data

readmit_count$sbj_id = as.numeric(readmit_count$sbj_id)
use_7$sbj_id = as.numeric(use_7$sbj_id)
left_join(x = use_7, y = readmit_count, by = "sbj_id") -> use_8

Rename readmission variable

colnames(use_8)[colnames(use_8) == "V1"] <- "readmit_count"
use_8$readmit_count = as.numeric(use_8$readmit_count)
use_8 %>% 
  mutate(readmit_count = coalesce(readmit_count, 0)) -> use_9

Make LOS long and average

IBD_hospitalizations_11_30_22 %>% 
  dplyr::select(sbj_id, inpatient_days) -> sbj_id_los
sbj_id_los %>% 
    # build grouping by combination of variables
    dplyr::group_by(sbj_id) %>%
    # add row number which works per group due to prior grouping
    dplyr::mutate(numbers = dplyr::row_number()) %>%
    # ungroup to prevent unexpected behaviour down stream
    dplyr::ungroup() -> los_numbers

spread(los_numbers, key = numbers, value = inpatient_days) -> los_wide

los_wide$mean_los=rowMeans(los_wide[,c(2:41)], na.rm=TRUE) 

los_wide%>% 
  dplyr::select(sbj_id, mean_los) -> mean_days

Merge LOS data

mean_days$sbj_id = as.numeric(mean_days$sbj_id)
use_9$sbj_id = as.numeric(use_9$sbj_id)
left_join(use_9, mean_days, by = "sbj_id") -> ut_clean1

ut_clean1$num_admits = as.factor(ut_clean1$num_admits)
ut_clean1$ed_count = as.factor(ut_clean1$ed_count)
ut_clean1$readmit_count = as.factor(ut_clean1$readmit_count)

Baseline Characteristics

ut_clean1 %>% 
  dplyr::select(ibd_3, age_yrs, gender, race_5, ethnic_3, lang_3, max_ch, IC, act_tob, num_admits, ed_count, surg_count, readmit_count, mean_los) -> ut_baseline
ut_baseline %>% tbl_summary(label = list(ibd_3 ~ "IBD Diagnosis", age_yrs ~ "Age", gender~ "Gender", race_5 ~ "Race", ethnic_3 ~ "Ethnicity", lang_3 ~ "Preferred Language"),
        statistic = list(all_continuous() ~ "{mean} ({sd})"),
        missing_text = "(Missing)")

Characteristic	N = 15,2081
IBD Diagnosis
CD	7,712 (51%)
UC	7,345 (48%)
Unspecified	151 (1.0%)
Age	49 (19)
Gender
Male	6,962 (46%)
Female	8,246 (54%)
Race
White	13,202 (87%)
Black	946 (6.2%)
Asian	373 (2.5%)
Native	56 (0.4%)
Other	631 (4.1%)
Ethnicity
NonHispanic	14,371 (98%)
UNKNOWN	0 (0%)
Hispanic	299 (2.0%)
CHOOSE NOT TO DISCLOSE	0 (0%)
(Missing)	538
Preferred Language
English	15,046 (99%)
Other	162 (1.1%)
max_ch	3.3 (4.9)
(Missing)	435
IC	7,111 (47%)
act_tob
No	12,593 (87%)
Yes	1,878 (13%)
NOT ASKED	0 (0%)
(Missing)	737
num_admits
0	11,137 (73%)
1	2,159 (14%)
2	859 (5.6%)
3	432 (2.8%)
4	217 (1.4%)
5	133 (0.9%)
6	76 (0.5%)
7	52 (0.3%)
8	32 (0.2%)
9	26 (0.2%)
10	12 (<0.1%)
11	17 (0.1%)
12	8 (<0.1%)
13	13 (<0.1%)
14	9 (<0.1%)
15	6 (<0.1%)
16	1 (<0.1%)
17	2 (<0.1%)
18	1 (<0.1%)
19	3 (<0.1%)
20	4 (<0.1%)
22	3 (<0.1%)
26	2 (<0.1%)
29	1 (<0.1%)
30	1 (<0.1%)
32	1 (<0.1%)
35	1 (<0.1%)
ed_count
0	13,231 (87%)
1	1,308 (8.6%)
2	370 (2.4%)
3	114 (0.7%)
4	62 (0.4%)
5	33 (0.2%)
6	26 (0.2%)
7	17 (0.1%)
8	7 (<0.1%)
9	11 (<0.1%)
10	8 (<0.1%)
11	1 (<0.1%)
12	3 (<0.1%)
13	2 (<0.1%)
14	5 (<0.1%)
15	2 (<0.1%)
16	3 (<0.1%)
18	3 (<0.1%)
19	1 (<0.1%)
22	1 (<0.1%)
surg_count
0	13,454 (88%)
1	1,079 (7.1%)
2	432 (2.8%)
3	178 (1.2%)
4	40 (0.3%)
5	18 (0.1%)
6	3 (<0.1%)
7	2 (<0.1%)
8	2 (<0.1%)
readmit_count
0	14,356 (94%)
1	564 (3.7%)
2	145 (1.0%)
3	71 (0.5%)
4	25 (0.2%)
5	12 (<0.1%)
6	10 (<0.1%)
7	4 (<0.1%)
8	6 (<0.1%)
9	1 (<0.1%)
10	5 (<0.1%)
11	3 (<0.1%)
12	1 (<0.1%)
14	2 (<0.1%)
16	1 (<0.1%)
18	2 (<0.1%)
mean_los	5.3 (6.1)
(Missing)	10,946
1 n (%); Mean (SD)

Bivariate analysis

Hospitalizations

Hospitalizations

Total SVI

library(MASS)
ut_clean1$num_admits = as.numeric(ut_clean1$num_admits)
admits.nb <- glm.nb(num_admits ~ ibd_3 + age_yrs + gender + race_5 + ethnic_3 + lang_3  + mstat_5 + act_tob + max_ch + IC + RPL_THEMES,
                    offset(log(pyears)),
               data = ut_clean1) 
summary(admits.nb)

Call:
glm.nb(formula = num_admits ~ ibd_3 + age_yrs + gender + race_5 + 
    ethnic_3 + lang_3 + mstat_5 + act_tob + max_ch + IC + RPL_THEMES, 
    data = ut_clean1, weights = offset(log(pyears)), init.theta = 7.223965147, 
    link = log)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-2.7926  -0.7421  -0.2962   0.0526  10.4417  

Coefficients:
                           Estimate Std. Error z value
(Intercept)               0.3917295  0.0326036  12.015
ibd_3UC                  -0.0909401  0.0138504  -6.566
ibd_3Unspecified         -0.3492108  0.1112171  -3.140
age_yrs                  -0.0056051  0.0004602 -12.180
genderFemale              0.0295595  0.0133026   2.222
race_5Black               0.1714237  0.0251101   6.827
race_5Asian              -0.1644412  0.0477163  -3.446
race_5Native             -0.2102921  0.1113673  -1.888
race_5Other              -0.0102631  0.0411783  -0.249
ethnic_3Hispanic         -0.0699396  0.0503400  -1.389
lang_3Other               0.0099557  0.0749167   0.133
mstat_5Unknown            0.1258660  0.0199318   6.315
mstat_5Unmarried         -0.0214087  0.0173712  -1.232
mstat_5DivorcedSeparated  0.1020796  0.0354365   2.881
mstat_5Widow              0.2053962  0.0448163   4.583
act_tobYes                0.0934379  0.0198512   4.707
max_ch                    0.0415718  0.0013175  31.553
IC                        0.3768615  0.0148294  25.413
RPL_THEMES                0.1974175  0.0266912   7.396
                         Pr(>|z|)    
(Intercept)               < 2e-16 ***
ibd_3UC                  5.17e-11 ***
ibd_3Unspecified         0.001690 ** 
age_yrs                   < 2e-16 ***
genderFemale             0.026277 *  
race_5Black              8.68e-12 ***
race_5Asian              0.000568 ***
race_5Native             0.058989 .  
race_5Other              0.803178    
ethnic_3Hispanic         0.164728    
lang_3Other              0.894280    
mstat_5Unknown           2.70e-10 ***
mstat_5Unmarried         0.217791    
mstat_5DivorcedSeparated 0.003969 ** 
mstat_5Widow             4.58e-06 ***
act_tobYes               2.51e-06 ***
max_ch                    < 2e-16 ***
IC                        < 2e-16 ***
RPL_THEMES               1.40e-13 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for Negative Binomial(7.224) family taken to be 1)

    Null deviance: 13419  on 10930  degrees of freedom
Residual deviance: 11124  on 10912  degrees of freedom
  (4269 observations deleted due to missingness)
AIC: 52476

Number of Fisher Scoring iterations: 1

              Theta:  7.224 
          Std. Err.:  0.283 

 2 x log-likelihood:  -52436.276 

broom::glance(admits.nb)

broom::tidy(admits.nb, exponentiate = TRUE)

model_performance(admits.nb)

# Indices of model performance

AIC       |       BIC | Nagelkerke's R2 |  RMSE | Sigma | Score_log | Score_spherical
-------------------------------------------------------------------------------------
52476.276 | 52622.263 |           0.268 | 1.708 | 1.010 |    -1.578 |           0.009

tbl_regression(admits.nb, label = list(age_yrs ~ "Age", gender~ "Gender", ethnic_3 ~ "Ethnicity", lang_3 ~ "Preferred Language", race_5 ~ "Race", RPL_THEMES ~ "Total SVI", act_tob ~"Active Tobacco Use", max_ch ~ "Charlson Comorbidity Index", IC ~ "Moderate to Severe Disease", ibd_3 ~ "IBD Type"), exponentiate = TRUE)

Characteristic	IRR1	95% CI1	p-value
IBD Type
CD	—	—
UC	0.91	0.89, 0.94	<0.001
Unspecified	0.71	0.56, 0.87	0.002
Age	0.99	0.99, 1.00	<0.001
Gender
Male	—	—
Female	1.03	1.00, 1.06	0.026
Race
White	—	—
Black	1.19	1.13, 1.25	<0.001
Asian	0.85	0.77, 0.93	<0.001
Native	0.81	0.65, 1.00	0.059
Other	0.99	0.91, 1.07	0.8
Ethnicity
NonHispanic	—	—
Hispanic	0.93	0.84, 1.03	0.2
Preferred Language
English	—	—
Other	1.01	0.87, 1.17	0.9
mstat_5
Married	—	—
Unknown	1.13	1.09, 1.18	<0.001
Unmarried	0.98	0.95, 1.01	0.2
DivorcedSeparated	1.11	1.03, 1.19	0.004
Widow	1.23	1.12, 1.34	<0.001
Active Tobacco Use
No	—	—
Yes	1.10	1.06, 1.14	<0.001
Charlson Comorbidity Index	1.04	1.04, 1.05	<0.001
Moderate to Severe Disease	1.46	1.42, 1.50	<0.001
Total SVI	1.22	1.16, 1.28	<0.001
1 IRR = Incidence Rate Ratio, CI = Confidence Interval

NA

All Themes

admits.nb2 <- glm.nb(num_admits ~ ibd_3 + age_yrs + gender + race_5 + ethnic_3 + lang_3  + mstat_5  + max_ch + IC + + act_tob + RPL_THEME1
                    + RPL_THEME2 + RPL_THEME3 + RPL_THEME4,
                    offset(log(pyears)),
               data = ut_clean1) 
summary(admits.nb2)

Call:
glm.nb(formula = num_admits ~ ibd_3 + age_yrs + gender + race_5 + 
    ethnic_3 + lang_3 + mstat_5 + max_ch + IC + +act_tob + RPL_THEME1 + 
    RPL_THEME2 + RPL_THEME3 + RPL_THEME4, data = ut_clean1, weights = offset(log(pyears)), 
    init.theta = 7.183204388, link = log)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-2.8477  -0.7385  -0.2989   0.0556  10.4045  

Coefficients:
                          Estimate Std. Error z value
(Intercept)               0.389750   0.035358  11.023
ibd_3UC                  -0.090462   0.013902  -6.507
ibd_3Unspecified         -0.350544   0.111272  -3.150
age_yrs                  -0.005596   0.000462 -12.113
genderFemale              0.028861   0.013349   2.162
race_5Black               0.166647   0.025626   6.503
race_5Asian              -0.149410   0.048403  -3.087
race_5Native             -0.215938   0.111496  -1.937
race_5Other              -0.008467   0.041454  -0.204
ethnic_3Hispanic         -0.060883   0.050507  -1.205
lang_3Other              -0.008891   0.075818  -0.117
mstat_5Unknown            0.125518   0.019991   6.279
mstat_5Unmarried         -0.021131   0.017458  -1.210
mstat_5DivorcedSeparated  0.103796   0.035506   2.923
mstat_5Widow              0.209667   0.044907   4.669
max_ch                    0.041594   0.001323  31.449
IC                        0.375128   0.014892  25.191
act_tobYes                0.087452   0.019938   4.386
RPL_THEME1                0.228140   0.038053   5.995
RPL_THEME2                0.022217   0.033291   0.667
RPL_THEME3               -0.002511   0.024186  -0.104
RPL_THEME4               -0.022670   0.027659  -0.820
                         Pr(>|z|)    
(Intercept)               < 2e-16 ***
ibd_3UC                  7.65e-11 ***
ibd_3Unspecified          0.00163 ** 
age_yrs                   < 2e-16 ***
genderFemale              0.03062 *  
race_5Black              7.87e-11 ***
race_5Asian               0.00202 ** 
race_5Native              0.05278 .  
race_5Other               0.83817    
ethnic_3Hispanic          0.22803    
lang_3Other               0.90665    
mstat_5Unknown           3.41e-10 ***
mstat_5Unmarried          0.22614    
mstat_5DivorcedSeparated  0.00346 ** 
mstat_5Widow             3.03e-06 ***
max_ch                    < 2e-16 ***
IC                        < 2e-16 ***
act_tobYes               1.15e-05 ***
RPL_THEME1               2.03e-09 ***
RPL_THEME2                0.50454    
RPL_THEME3                0.91732    
RPL_THEME4                0.41242    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for Negative Binomial(7.1832) family taken to be 1)

    Null deviance: 13351  on 10869  degrees of freedom
Residual deviance: 11053  on 10848  degrees of freedom
  (4330 observations deleted due to missingness)
AIC: 52168

Number of Fisher Scoring iterations: 1

              Theta:  7.183 
          Std. Err.:  0.281 

 2 x log-likelihood:  -52122.063 

broom::glance(admits.nb2)

broom::tidy(admits.nb2, exponentiate = TRUE)

model_performance(admits.nb2)

# Indices of model performance

AIC       |       BIC | Nagelkerke's R2 |  RMSE | Sigma | Score_log | Score_spherical
-------------------------------------------------------------------------------------
52168.063 | 52335.819 |           0.270 | 1.711 | 1.009 |    -1.579 |           0.009

tbl_regression(admits.nb2, label = list(age_yrs ~ "Age", gender~ "Gender", ethnic_3 ~ "Ethnicity", lang_3 ~ "Preferred Language", race_5 ~ "Race", RPL_THEME1 ~ "Soceioeconomic Status", RPL_THEME2 ~ "Household Composition", RPL_THEME3 ~ "Minority Status and Language", RPL_THEME4 ~ "Housing and Transportation", max_ch ~ "Charlson Comorbidity Index", IC ~ "Moderate to Severe Disease", ibd_3 ~ "IBD Type", act_tob ~ "Active Tobacco Use"), exponentiate = TRUE)

Characteristic	IRR1	95% CI1	p-value
IBD Type
CD	—	—
UC	0.91	0.89, 0.94	<0.001
Unspecified	0.70	0.56, 0.87	0.002
Age	0.99	0.99, 1.00	<0.001
Gender
Male	—	—
Female	1.03	1.00, 1.06	0.031
Race
White	—	—
Black	1.18	1.12, 1.24	<0.001
Asian	0.86	0.78, 0.95	0.002
Native	0.81	0.64, 1.00	0.053
Other	0.99	0.91, 1.07	0.8
Ethnicity
NonHispanic	—	—
Hispanic	0.94	0.85, 1.04	0.2
Preferred Language
English	—	—
Other	0.99	0.85, 1.15	>0.9
mstat_5
Married	—	—
Unknown	1.13	1.09, 1.18	<0.001
Unmarried	0.98	0.95, 1.01	0.2
DivorcedSeparated	1.11	1.03, 1.19	0.003
Widow	1.23	1.13, 1.35	<0.001
Charlson Comorbidity Index	1.04	1.04, 1.05	<0.001
Moderate to Severe Disease	1.46	1.41, 1.50	<0.001
Active Tobacco Use
No	—	—
Yes	1.09	1.05, 1.13	<0.001
Soceioeconomic Status	1.26	1.17, 1.35	<0.001
Household Composition	1.02	0.96, 1.09	0.5
Minority Status and Language	1.00	0.95, 1.05	>0.9
Housing and Transportation	0.98	0.93, 1.03	0.4
1 IRR = Incidence Rate Ratio, CI = Confidence Interval

NA

ED Visits

SVI

ut_clean1$ed_count = as.numeric(ut_clean1$ed_count)
ed.nb <- glm.nb(ed_count ~ ibd_3 + age_yrs + gender + race_5 + ethnic_3 + lang_3 + act_tob + mstat_5 + max_ch + IC + RPL_THEMES,
                    offset(log(pyears)),
               data = use_7) 
summary(ed.nb)

Call:
glm.nb(formula = ed_count ~ ibd_3 + age_yrs + gender + race_5 + 
    ethnic_3 + lang_3 + act_tob + mstat_5 + max_ch + IC + RPL_THEMES, 
    data = use_7, weights = offset(log(pyears)), init.theta = 0.1514506594, 
    link = log)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-1.0473  -0.7522  -0.6451  -0.4113   6.4678  

Coefficients:
                          Estimate Std. Error z value Pr(>|z|)    
(Intercept)              -1.582752   0.128723 -12.296  < 2e-16 ***
ibd_3UC                  -0.033155   0.054751  -0.606 0.544809    
ibd_3Unspecified          0.378116   0.346825   1.090 0.275616    
age_yrs                   0.001509   0.001828   0.826 0.408855    
genderFemale              0.184720   0.052802   3.498 0.000468 ***
race_5Black              -0.282472   0.112957  -2.501 0.012394 *  
race_5Asian               0.057952   0.169176   0.343 0.731935    
race_5Native             -0.111512   0.417718  -0.267 0.789504    
race_5Other               0.261556   0.154482   1.693 0.090433 .  
ethnic_3Hispanic         -0.124695   0.194623  -0.641 0.521716    
lang_3Other               0.427002   0.258857   1.650 0.099032 .  
act_tobYes               -0.020520   0.082313  -0.249 0.803133    
mstat_5Unknown           -0.737323   0.089768  -8.214  < 2e-16 ***
mstat_5Unmarried         -0.039180   0.067317  -0.582 0.560552    
mstat_5DivorcedSeparated  0.140240   0.136489   1.027 0.304193    
mstat_5Widow              0.143838   0.174693   0.823 0.410293    
max_ch                    0.004108   0.005527   0.743 0.457340    
IC                        0.076803   0.057309   1.340 0.180192    
RPL_THEMES                0.109823   0.107205   1.024 0.305637    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for Negative Binomial(0.1515) family taken to be 1)

    Null deviance: 6370.4  on 10930  degrees of freedom
Residual deviance: 6259.5  on 10912  degrees of freedom
  (4269 observations deleted due to missingness)
AIC: 18001

Number of Fisher Scoring iterations: 1

              Theta:  0.15145 
          Std. Err.:  0.00566 

 2 x log-likelihood:  -17960.88500 

broom::glance(ed.nb)

broom::tidy(ed.nb, exponentiate = TRUE)

model_performance(ed.nb)

# Indices of model performance

AIC       |       BIC | Nagelkerke's R2 |  RMSE | Sigma | Score_log | Score_spherical
-------------------------------------------------------------------------------------
18000.885 | 18146.872 |           0.023 | 0.927 | 0.757 |    -0.647 |           0.009

tbl_regression(ed.nb, label = list(age_yrs ~ "Age", gender~ "Gender", ethnic_3 ~ "Ethnicity", lang_3 ~ "Preferred Language", race_5 ~ "Race", RPL_THEMES ~ "Total SVI", max_ch ~ "Charlson Comorbidity Index", IC ~ "Moderate to Severe Disease", ibd_3 ~ "IBD Type"), exponentiate = TRUE)

Characteristic	IRR1	95% CI1	p-value
IBD Type
CD	—	—
UC	0.97	0.87, 1.08	0.5
Unspecified	1.46	0.75, 3.00	0.3
Age	1.00	1.00, 1.01	0.4
Gender
Male	—	—
Female	1.20	1.08, 1.33	<0.001
Race
White	—	—
Black	0.75	0.60, 0.94	0.012
Asian	1.06	0.76, 1.49	0.7
Native	0.89	0.39, 2.11	0.8
Other	1.30	0.96, 1.78	0.090
Ethnicity
NonHispanic	—	—
Hispanic	0.88	0.61, 1.30	0.5
Preferred Language
English	—	—
Other	1.53	0.93, 2.61	0.10
act_tob
No	—	—
Yes	0.98	0.83, 1.15	0.8
mstat_5
Married	—	—
Unknown	0.48	0.40, 0.57	<0.001
Unmarried	0.96	0.84, 1.10	0.6
DivorcedSeparated	1.15	0.88, 1.51	0.3
Widow	1.15	0.82, 1.64	0.4
Charlson Comorbidity Index	1.00	0.99, 1.02	0.5
Moderate to Severe Disease	1.08	0.96, 1.21	0.2
Total SVI	1.12	0.90, 1.38	0.3
1 IRR = Incidence Rate Ratio, CI = Confidence Interval

View(IBD_surgery_11_17_22)
View(use_7)
View(use_7)

Quartiles

ut_clean1$ed_count = as.numeric(ut_clean1$ed_count)
ed.nb2 <- glm.nb(ed_count ~ ibd_3 + age_yrs + gender + race_5 + ethnic_3 + lang_3 + act_tob + mstat_5 + max_ch + IC + RPL_4,
                    offset(log(pyears)),
               data = use_7) 
summary(ed.nb2)

Call:
glm.nb(formula = ed_count ~ ibd_3 + age_yrs + gender + race_5 + 
    ethnic_3 + lang_3 + act_tob + mstat_5 + max_ch + IC + RPL_4, 
    data = use_7, weights = offset(log(pyears)), init.theta = 0.1515050523, 
    link = log)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-1.0516  -0.7523  -0.6438  -0.4112   6.3669  

Coefficients:
                          Estimate Std. Error z value Pr(>|z|)    
(Intercept)              -1.573679   0.127708 -12.322  < 2e-16 ***
ibd_3UC                  -0.035929   0.054754  -0.656 0.511699    
ibd_3Unspecified          0.402735   0.345870   1.164 0.244257    
age_yrs                   0.001693   0.001827   0.927 0.354182    
genderFemale              0.185660   0.052808   3.516 0.000439 ***
race_5Black              -0.267866   0.113266  -2.365 0.018033 *  
race_5Asian               0.055434   0.169196   0.328 0.743192    
race_5Native             -0.105068   0.416800  -0.252 0.800977    
race_5Other               0.264951   0.154383   1.716 0.086127 .  
ethnic_3Hispanic         -0.115186   0.194301  -0.593 0.553299    
lang_3Other               0.429221   0.259187   1.656 0.097717 .  
act_tobYes               -0.011579   0.082258  -0.141 0.888058    
mstat_5Unknown           -0.730288   0.089764  -8.136  4.1e-16 ***
mstat_5Unmarried         -0.026653   0.067278  -0.396 0.691982    
mstat_5DivorcedSeparated  0.145395   0.136534   1.065 0.286921    
mstat_5Widow              0.151982   0.174510   0.871 0.383806    
max_ch                    0.004115   0.005528   0.744 0.456621    
IC                        0.074345   0.057299   1.297 0.194466    
RPL_4Second               0.053945   0.061623   0.875 0.381352    
RPL_4Third               -0.045668   0.073688  -0.620 0.535421    
RPL_4Fourth               0.069829   0.100254   0.697 0.486104    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for Negative Binomial(0.1515) family taken to be 1)

    Null deviance: 6371.6  on 10930  degrees of freedom
Residual deviance: 6259.5  on 10910  degrees of freedom
  (4269 observations deleted due to missingness)
AIC: 18004

Number of Fisher Scoring iterations: 1

              Theta:  0.15151 
          Std. Err.:  0.00567 

 2 x log-likelihood:  -17959.72200 

broom::glance(ed.nb2)

broom::tidy(ed.nb2, exponentiate = TRUE)

model_performance(ed.nb2)

# Indices of model performance

AIC       |       BIC | Nagelkerke's R2 |  RMSE | Sigma | Score_log | Score_spherical
-------------------------------------------------------------------------------------
18003.722 | 18164.308 |           0.023 | 0.927 | 0.757 |    -0.646 |           0.009

tbl_regression(ed.nb2, label = list(age_yrs ~ "Age", gender~ "Gender", ethnic_3 ~ "Ethnicity", lang_3 ~ "Preferred Language", race_5 ~ "Race", RPL_4 ~ "SVI Quartiles", max_ch ~ "Charlson Comorbidity Index", IC ~ "Moderate to Severe Disease", ibd_3 ~ "IBD Type"), exponentiate = TRUE)

Characteristic	IRR1	95% CI1	p-value
IBD Type
CD	—	—
UC	0.96	0.87, 1.07	0.5
Unspecified	1.50	0.77, 3.08	0.2
Age	1.00	1.00, 1.01	0.4
Gender
Male	—	—
Female	1.20	1.09, 1.34	<0.001
Race
White	—	—
Black	0.77	0.61, 0.96	0.018
Asian	1.06	0.76, 1.48	0.7
Native	0.90	0.39, 2.13	0.8
Other	1.30	0.96, 1.78	0.086
Ethnicity
NonHispanic	—	—
Hispanic	0.89	0.61, 1.31	0.6
Preferred Language
English	—	—
Other	1.54	0.93, 2.62	0.10
act_tob
No	—	—
Yes	0.99	0.84, 1.16	0.9
mstat_5
Married	—	—
Unknown	0.48	0.40, 0.57	<0.001
Unmarried	0.97	0.85, 1.11	0.7
DivorcedSeparated	1.16	0.89, 1.52	0.3
Widow	1.16	0.83, 1.65	0.4
Charlson Comorbidity Index	1.00	0.99, 1.02	0.5
Moderate to Severe Disease	1.08	0.96, 1.21	0.2
SVI Quartiles
First	—	—
Second	1.06	0.94, 1.19	0.4
Third	0.96	0.83, 1.11	0.5
Fourth	1.07	0.88, 1.31	0.5
1 IRR = Incidence Rate Ratio, CI = Confidence Interval

All themes

ut_clean1$ed_count = as.numeric(ut_clean1$ed_count)
ed.nb3 <- glm.nb(ed_count ~ ibd_3 + age_yrs + gender + race_5 + ethnic_3 + lang_3 + act_tob + mstat_5 + max_ch + IC + RPL_THEME1 + RPL_THEME2 + RPL_THEME3 + RPL_THEME4,
                    offset(log(pyears)),
               data = use_7) 
summary(ed.nb3)

Call:
glm.nb(formula = ed_count ~ ibd_3 + age_yrs + gender + race_5 + 
    ethnic_3 + lang_3 + act_tob + mstat_5 + max_ch + IC + RPL_THEME1 + 
    RPL_THEME2 + RPL_THEME3 + RPL_THEME4, data = use_7, weights = offset(log(pyears)), 
    init.theta = 0.1517000316, link = log)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-1.0634  -0.7508  -0.6428  -0.4109   6.4619  

Coefficients:
                          Estimate Std. Error z value Pr(>|z|)    
(Intercept)              -1.554546   0.139369 -11.154  < 2e-16 ***
ibd_3UC                  -0.030151   0.054938  -0.549 0.583132    
ibd_3Unspecified          0.362560   0.347374   1.044 0.296617    
age_yrs                   0.001318   0.001834   0.719 0.472419    
genderFemale              0.185748   0.052980   3.506 0.000455 ***
race_5Black              -0.316058   0.114763  -2.754 0.005887 ** 
race_5Asian               0.048257   0.172315   0.280 0.779439    
race_5Native             -0.118942   0.417861  -0.285 0.775916    
race_5Other               0.290557   0.154634   1.879 0.060245 .  
ethnic_3Hispanic         -0.139477   0.195480  -0.714 0.475530    
lang_3Other               0.326029   0.262810   1.241 0.214772    
act_tobYes               -0.032128   0.082656  -0.389 0.697507    
mstat_5Unknown           -0.734552   0.089913  -8.170 3.09e-16 ***
mstat_5Unmarried         -0.042803   0.067646  -0.633 0.526894    
mstat_5DivorcedSeparated  0.142611   0.136680   1.043 0.296767    
mstat_5Widow              0.157399   0.174933   0.900 0.368245    
max_ch                    0.004088   0.005549   0.737 0.461332    
IC                        0.063982   0.057534   1.112 0.266103    
RPL_THEME1                0.427581   0.151668   2.819 0.004815 ** 
RPL_THEME2               -0.156605   0.132466  -1.182 0.237118    
RPL_THEME3                0.024274   0.095333   0.255 0.799016    
RPL_THEME4               -0.155772   0.110143  -1.414 0.157280    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for Negative Binomial(0.1517) family taken to be 1)

    Null deviance: 6334.4  on 10869  degrees of freedom
Residual deviance: 6217.9  on 10848  degrees of freedom
  (4330 observations deleted due to missingness)
AIC: 17881

Number of Fisher Scoring iterations: 1

              Theta:  0.15170 
          Std. Err.:  0.00570 

 2 x log-likelihood:  -17834.60700 

broom::glance(ed.nb3)

broom::tidy(ed.nb3, exponentiate = TRUE)

model_performance(ed.nb3)

# Indices of model performance

AIC       |       BIC | Nagelkerke's R2 |  RMSE | Sigma | Score_log | Score_spherical
-------------------------------------------------------------------------------------
17880.607 | 18048.364 |           0.024 | 0.927 | 0.757 |    -0.647 |           0.009

tbl_regression(ed.nb3, label = list(age_yrs ~ "Age", gender~ "Gender", ethnic_3 ~ "Ethnicity", lang_3 ~ "Preferred Language", race_5 ~ "Race", max_ch ~ "Charlson Comorbidity Index", IC ~ "Moderate to Severe Disease", ibd_3 ~ "IBD Type"), exponentiate = TRUE)

Characteristic	IRR1	95% CI1	p-value
IBD Type
CD	—	—
UC	0.97	0.87, 1.08	0.6
Unspecified	1.44	0.74, 2.96	0.3
Age	1.00	1.00, 1.01	0.5
Gender
Male	—	—
Female	1.20	1.09, 1.34	<0.001
Race
White	—	—
Black	0.73	0.58, 0.91	0.006
Asian	1.05	0.75, 1.48	0.8
Native	0.89	0.39, 2.09	0.8
Other	1.34	0.98, 1.83	0.060
Ethnicity
NonHispanic	—	—
Hispanic	0.87	0.60, 1.28	0.5
Preferred Language
English	—	—
Other	1.39	0.83, 2.38	0.2
act_tob
No	—	—
Yes	0.97	0.82, 1.14	0.7
mstat_5
Married	—	—
Unknown	0.48	0.40, 0.57	<0.001
Unmarried	0.96	0.84, 1.09	0.5
DivorcedSeparated	1.15	0.88, 1.52	0.3
Widow	1.17	0.84, 1.66	0.4
Charlson Comorbidity Index	1.00	0.99, 1.02	0.5
Moderate to Severe Disease	1.07	0.95, 1.19	0.3
RPL_THEME1	1.53	1.13, 2.09	0.005
RPL_THEME2	0.86	0.65, 1.12	0.2
RPL_THEME3	1.02	0.85, 1.24	0.8
RPL_THEME4	0.86	0.68, 1.07	0.2
1 IRR = Incidence Rate Ratio, CI = Confidence Interval

Surgeries

ut_clean1$surg_count = as.numeric(ut_clean1$surg_count)
surg.nb <- glm.nb(surg_count ~ ibd_3 + age_yrs + gender + race_5 + ethnic_3 + lang_3  + mstat_5 + act_tob + max_ch + IC + RPL_THEMES,
                    offset(log(pyears)),
               data = ut_clean1) 
summary(surg.nb)

Call:
glm.nb(formula = surg_count ~ ibd_3 + age_yrs + gender + race_5 + 
    ethnic_3 + lang_3 + mstat_5 + act_tob + max_ch + IC + RPL_THEMES, 
    data = ut_clean1, weights = offset(log(pyears)), init.theta = 0.3285786911, 
    link = log)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-1.5977  -0.8019  -0.5629  -0.3522   5.0741  

Coefficients:
                          Estimate Std. Error z value Pr(>|z|)    
(Intercept)              -1.394565   0.110548 -12.615  < 2e-16 ***
ibd_3UC                  -0.140971   0.046376  -3.040  0.00237 ** 
ibd_3Unspecified          0.313783   0.319754   0.981  0.32643    
age_yrs                  -0.023594   0.001594 -14.805  < 2e-16 ***
genderFemale             -0.084468   0.043752  -1.931  0.05353 .  
race_5Black               0.119419   0.080838   1.477  0.13961    
race_5Asian              -0.877973   0.196266  -4.473 7.70e-06 ***
race_5Native             -0.509517   0.406911  -1.252  0.21051    
race_5Other              -0.647724   0.157355  -4.116 3.85e-05 ***
ethnic_3Hispanic         -0.016385   0.166006  -0.099  0.92137    
lang_3Other              -0.236229   0.301073  -0.785  0.43267    
mstat_5Unknown            0.120547   0.067414   1.788  0.07375 .  
mstat_5Unmarried         -0.005805   0.057707  -0.101  0.91987    
mstat_5DivorcedSeparated  0.333602   0.112981   2.953  0.00315 ** 
mstat_5Widow              0.187385   0.178803   1.048  0.29464    
act_tobYes                0.182114   0.062275   2.924  0.00345 ** 
max_ch                    0.061792   0.004635  13.332  < 2e-16 ***
IC                        0.989474   0.052986  18.674  < 2e-16 ***
RPL_THEMES                0.557573   0.086487   6.447 1.14e-10 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for Negative Binomial(0.3286) family taken to be 1)

    Null deviance: 9070.8  on 11250  degrees of freedom
Residual deviance: 7824.4  on 11232  degrees of freedom
  (4275 observations deleted due to missingness)
AIC: 19192

Number of Fisher Scoring iterations: 1

              Theta:  0.3286 
          Std. Err.:  0.0144 

 2 x log-likelihood:  -19152.0590 

broom::glance(surg.nb)

broom::tidy(surg.nb, exponentiate = TRUE)

model_performance(surg.nb)

# Indices of model performance

AIC       |       BIC | Nagelkerke's R2 |  RMSE | Sigma | Score_log | Score_spherical
-------------------------------------------------------------------------------------
19192.059 | 19338.623 |           0.189 | 0.667 | 0.835 |    -0.598 |           0.009

tbl_regression(surg.nb, label = list(age_yrs ~ "Age", gender~ "Gender", ethnic_3 ~ "Ethnicity", lang_3 ~ "Preferred Language", race_5 ~ "Race", RPL_THEMES ~ "Total SVI", act_tob ~"Active Tobacco Use", max_ch ~ "Charlson Comorbidity Index", IC ~ "Moderate to Severe Disease", ibd_3 ~ "IBD Type"), exponentiate = TRUE)

Characteristic	IRR1	95% CI1	p-value
IBD Type
CD	—	—
UC	0.87	0.79, 0.95	0.002
Unspecified	1.37	0.69, 2.62	0.3
Age	0.98	0.97, 0.98	<0.001
Gender
Male	—	—
Female	0.92	0.84, 1.00	0.054
Race
White	—	—
Black	1.13	0.96, 1.32	0.14
Asian	0.42	0.28, 0.60	<0.001
Native	0.60	0.27, 1.26	0.2
Other	0.52	0.38, 0.71	<0.001
Ethnicity
NonHispanic	—	—
Hispanic	0.98	0.71, 1.37	>0.9
Preferred Language
English	—	—
Other	0.79	0.44, 1.38	0.4
mstat_5
Married	—	—
Unknown	1.13	0.99, 1.29	0.074
Unmarried	0.99	0.89, 1.11	>0.9
DivorcedSeparated	1.40	1.12, 1.74	0.003
Widow	1.21	0.83, 1.72	0.3
Active Tobacco Use
No	—	—
Yes	1.20	1.06, 1.36	0.003
Charlson Comorbidity Index	1.06	1.05, 1.07	<0.001
Moderate to Severe Disease	2.69	2.42, 2.99	<0.001
Total SVI	1.75	1.47, 2.07	<0.001
1 IRR = Incidence Rate Ratio, CI = Confidence Interval

All themes

surg.nb2 <- glm.nb(surg_count ~ ibd_3 + age_yrs + gender + race_5 + ethnic_3 + lang_3  + mstat_5 + act_tob + max_ch + IC + RPL_THEME1 + RPL_THEME2 + RPL_THEME3 + RPL_THEME4,
                    offset(log(pyears)),
               data = ut_clean1) 
summary(surg.nb2)

Call:
glm.nb(formula = surg_count ~ ibd_3 + age_yrs + gender + race_5 + 
    ethnic_3 + lang_3 + mstat_5 + act_tob + max_ch + IC + RPL_THEME1 + 
    RPL_THEME2 + RPL_THEME3 + RPL_THEME4, data = ut_clean1, weights = offset(log(pyears)), 
    init.theta = 0.3308139959, link = log)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-1.6090  -0.7997  -0.5606  -0.3482   4.9820  

Coefficients:
                          Estimate Std. Error z value Pr(>|z|)    
(Intercept)              -1.359512   0.119175 -11.408  < 2e-16 ***
ibd_3UC                  -0.135858   0.046452  -2.925 0.003448 ** 
ibd_3Unspecified          0.309797   0.319725   0.969 0.332571    
age_yrs                  -0.023865   0.001597 -14.940  < 2e-16 ***
genderFemale             -0.085866   0.043806  -1.960 0.049981 *  
race_5Black               0.153340   0.082210   1.865 0.062150 .  
race_5Asian              -0.763508   0.197333  -3.869 0.000109 ***
race_5Native             -0.533707   0.408381  -1.307 0.191252    
race_5Other              -0.621963   0.158542  -3.923 8.74e-05 ***
ethnic_3Hispanic          0.026972   0.166065   0.162 0.870977    
lang_3Other              -0.216466   0.302775  -0.715 0.474647    
mstat_5Unknown            0.119926   0.067415   1.779 0.075252 .  
mstat_5Unmarried         -0.005907   0.057911  -0.102 0.918757    
mstat_5DivorcedSeparated  0.311452   0.113640   2.741 0.006131 ** 
mstat_5Widow              0.205566   0.178907   1.149 0.250551    
act_tobYes                0.171750   0.062296   2.757 0.005833 ** 
max_ch                    0.061604   0.004644  13.266  < 2e-16 ***
IC                        0.984178   0.053109  18.531  < 2e-16 ***
RPL_THEME1                0.506688   0.125040   4.052 5.07e-05 ***
RPL_THEME2                0.351011   0.108294   3.241 0.001190 ** 
RPL_THEME3               -0.122337   0.080186  -1.526 0.127091    
RPL_THEME4               -0.163951   0.091553  -1.791 0.073330 .  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for Negative Binomial(0.3308) family taken to be 1)

    Null deviance: 9052.2  on 11189  degrees of freedom
Residual deviance: 7784.8  on 11168  degrees of freedom
  (4336 observations deleted due to missingness)
AIC: 19116

Number of Fisher Scoring iterations: 1

              Theta:  0.3308 
          Std. Err.:  0.0145 

 2 x log-likelihood:  -19070.4560 

broom::glance(surg.nb2)

broom::tidy(surg.nb2, exponentiate = TRUE)

model_performance(surg.nb2)

# Indices of model performance

AIC       |       BIC | Nagelkerke's R2 |  RMSE | Sigma | Score_log | Score_spherical
-------------------------------------------------------------------------------------
19116.456 | 19284.880 |           0.193 | 0.668 | 0.835 |    -0.598 |           0.009

tbl_regression(surg.nb2, label = list(age_yrs ~ "Age", gender~ "Gender", ethnic_3 ~ "Ethnicity", lang_3 ~ "Preferred Language", race_5 ~ "Race", act_tob ~"Active Tobacco Use", max_ch ~ "Charlson Comorbidity Index", IC ~ "Moderate to Severe Disease", ibd_3 ~ "IBD Type"), exponentiate = TRUE)

Characteristic	IRR1	95% CI1	p-value
IBD Type
CD	—	—
UC	0.87	0.80, 0.96	0.003
Unspecified	1.36	0.69, 2.60	0.3
Age	0.98	0.97, 0.98	<0.001
Gender
Male	—	—
Female	0.92	0.84, 1.00	0.050
Race
White	—	—
Black	1.17	0.99, 1.37	0.062
Asian	0.47	0.31, 0.68	<0.001
Native	0.59	0.27, 1.23	0.2
Other	0.54	0.39, 0.73	<0.001
Ethnicity
NonHispanic	—	—
Hispanic	1.03	0.74, 1.43	0.9
Preferred Language
English	—	—
Other	0.81	0.45, 1.41	0.5
mstat_5
Married	—	—
Unknown	1.13	0.99, 1.29	0.075
Unmarried	0.99	0.89, 1.11	>0.9
DivorcedSeparated	1.37	1.10, 1.70	0.006
Widow	1.23	0.85, 1.75	0.3
Active Tobacco Use
No	—	—
Yes	1.19	1.05, 1.34	0.006
Charlson Comorbidity Index	1.06	1.05, 1.07	<0.001
Moderate to Severe Disease	2.68	2.41, 2.98	<0.001
RPL_THEME1	1.66	1.30, 2.13	<0.001
RPL_THEME2	1.42	1.15, 1.76	0.001
RPL_THEME3	0.88	0.76, 1.04	0.13
RPL_THEME4	0.85	0.71, 1.02	0.073
1 IRR = Incidence Rate Ratio, CI = Confidence Interval

Readmissions

SVI

ut_clean1$readmit_count = as.numeric(ut_clean1$readmit_count)
readmit.nb <- glm.nb(readmit_count ~ ibd_3 + age_yrs + gender + race_5 + ethnic_3 + lang_3  + mstat_5 + act_tob + max_ch + IC + RPL_THEMES,
                    offset(log(pyears)),
               data = ut_clean1) 

Warning: iteration limit reachedWarning: iteration limit reached

summary(readmit.nb)

Call:
glm.nb(formula = readmit_count ~ ibd_3 + age_yrs + gender + race_5 + 
    ethnic_3 + lang_3 + mstat_5 + act_tob + max_ch + IC + RPL_THEMES, 
    data = ut_clean1, weights = offset(log(pyears)), init.theta = 49510.38307, 
    link = log)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-0.9232  -0.2401  -0.1316  -0.0072  10.1003  

Coefficients:
                           Estimate Std. Error z value
(Intercept)               0.1091401  0.0362461   3.011
ibd_3UC                  -0.0154505  0.0154689  -0.999
ibd_3Unspecified         -0.0560662  0.1091470  -0.514
age_yrs                  -0.0032447  0.0005142  -6.310
genderFemale              0.0075701  0.0148541   0.510
race_5Black               0.1073172  0.0286804   3.742
race_5Asian              -0.0616009  0.0506021  -1.217
race_5Native             -0.0279398  0.1165397  -0.240
race_5Other               0.0526996  0.0446290   1.181
ethnic_3Hispanic         -0.0714176  0.0557525  -1.281
lang_3Other              -0.0476137  0.0826788  -0.576
mstat_5Unknown            0.0481555  0.0225502   2.135
mstat_5Unmarried         -0.0265728  0.0193416  -1.374
mstat_5DivorcedSeparated  0.0653761  0.0395136   1.655
mstat_5Widow              0.0658805  0.0519369   1.268
act_tobYes                0.0631236  0.0223934   2.819
max_ch                    0.0173344  0.0015171  11.426
IC                        0.1286881  0.0163188   7.886
RPL_THEMES                0.1009949  0.0300033   3.366
                         Pr(>|z|)    
(Intercept)              0.002603 ** 
ibd_3UC                  0.317887    
ibd_3Unspecified         0.607479    
age_yrs                  2.78e-10 ***
genderFemale             0.610310    
race_5Black              0.000183 ***
race_5Asian              0.223469    
race_5Native             0.810528    
race_5Other              0.237667    
ethnic_3Hispanic         0.200202    
lang_3Other              0.564691    
mstat_5Unknown           0.032722 *  
mstat_5Unmarried         0.169483    
mstat_5DivorcedSeparated 0.098022 .  
mstat_5Widow             0.204630    
act_tobYes               0.004820 ** 
max_ch                    < 2e-16 ***
IC                       3.12e-15 ***
RPL_THEMES               0.000762 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for Negative Binomial(49510.38) family taken to be 1)

    Null deviance: 3821.3  on 10931  degrees of freedom
Residual deviance: 3527.4  on 10913  degrees of freedom
  (4269 observations deleted due to missingness)
AIC: 37166

Number of Fisher Scoring iterations: 1

              Theta:  49510 
          Std. Err.:  60103 
Warning while fitting theta: iteration limit reached 

 2 x log-likelihood:  -37126.07 

broom::glance(readmit.nb)

broom::tidy(readmit.nb, exponentiate = TRUE)

model_performance(readmit.nb)

# Indices of model performance

AIC       |       BIC | Nagelkerke's R2 |  RMSE | Sigma | Score_log | Score_spherical
-------------------------------------------------------------------------------------
37166.073 | 37312.062 |           0.090 | 0.675 | 0.569 |    -1.326 |           0.009

tbl_regression(readmit.nb, label = list(age_yrs ~ "Age", gender~ "Gender", ethnic_3 ~ "Ethnicity", lang_3 ~ "Preferred Language", race_5 ~ "Race", RPL_THEMES ~ "Total SVI", act_tob ~"Active Tobacco Use", max_ch ~ "Charlson Comorbidity Index", IC ~ "Moderate to Severe Disease", ibd_3 ~ "IBD Type"), exponentiate = TRUE)

Characteristic	IRR1	95% CI1	p-value
IBD Type
CD	—	—
UC	0.98	0.96, 1.01	0.3
Unspecified	0.95	0.76, 1.16	0.6
Age	1.00	1.00, 1.00	<0.001
Gender
Male	—	—
Female	1.01	0.98, 1.04	0.6
Race
White	—	—
Black	1.11	1.05, 1.18	<0.001
Asian	0.94	0.85, 1.04	0.2
Native	0.97	0.77, 1.21	0.8
Other	1.05	0.96, 1.15	0.2
Ethnicity
NonHispanic	—	—
Hispanic	0.93	0.83, 1.04	0.2
Preferred Language
English	—	—
Other	0.95	0.81, 1.12	0.6
mstat_5
Married	—	—
Unknown	1.05	1.00, 1.10	0.033
Unmarried	0.97	0.94, 1.01	0.2
DivorcedSeparated	1.07	0.99, 1.15	0.10
Widow	1.07	0.96, 1.18	0.2
Active Tobacco Use
No	—	—
Yes	1.07	1.02, 1.11	0.005
Charlson Comorbidity Index	1.02	1.01, 1.02	<0.001
Moderate to Severe Disease	1.14	1.10, 1.17	<0.001
Total SVI	1.11	1.04, 1.17	<0.001
1 IRR = Incidence Rate Ratio, CI = Confidence Interval

All themes

readmit.nb2 <- glm.nb(readmit_count ~ ibd_3 + age_yrs + gender + race_5 + ethnic_3 + lang_3  + mstat_5 + act_tob + max_ch + IC + RPL_THEME1 + RPL_THEME2 + RPL_THEME3 + RPL_THEME4,
                    offset(log(pyears)),
               data = ut_clean1) 

Warning: iteration limit reachedWarning: iteration limit reached

summary(readmit.nb2)

Call:
glm.nb(formula = readmit_count ~ ibd_3 + age_yrs + gender + race_5 + 
    ethnic_3 + lang_3 + mstat_5 + act_tob + max_ch + IC + RPL_THEME1 + 
    RPL_THEME2 + RPL_THEME3 + RPL_THEME4, data = ut_clean1, weights = offset(log(pyears)), 
    init.theta = 49491.15318, link = log)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-0.9409  -0.2417  -0.1299  -0.0064  10.0386  

Coefficients:
                          Estimate Std. Error z value
(Intercept)               0.096163   0.039340   2.444
ibd_3UC                  -0.014901   0.015520  -0.960
ibd_3Unspecified         -0.056561   0.109156  -0.518
age_yrs                  -0.003231   0.000516  -6.261
genderFemale              0.007223   0.014902   0.485
race_5Black               0.098433   0.029240   3.366
race_5Asian              -0.062208   0.051308  -1.212
race_5Native             -0.028323   0.116557  -0.243
race_5Other               0.050844   0.044885   1.133
ethnic_3Hispanic         -0.069651   0.055960  -1.245
lang_3Other              -0.059464   0.083403  -0.713
mstat_5Unknown            0.047657   0.022610   2.108
mstat_5Unmarried         -0.027323   0.019430  -1.406
mstat_5DivorcedSeparated  0.066691   0.039570   1.685
mstat_5Widow              0.066730   0.052051   1.282
act_tobYes                0.061494   0.022481   2.735
max_ch                    0.017360   0.001523  11.402
IC                        0.128656   0.016384   7.853
RPL_THEME1                0.120746   0.042620   2.833
RPL_THEME2                0.008033   0.037248   0.216
RPL_THEME3                0.033301   0.026996   1.234
RPL_THEME4               -0.024206   0.030946  -0.782
                         Pr(>|z|)    
(Intercept)              0.014510 *  
ibd_3UC                  0.336971    
ibd_3Unspecified         0.604343    
age_yrs                  3.84e-10 ***
genderFemale             0.627892    
race_5Black              0.000762 ***
race_5Asian              0.225343    
race_5Native             0.808007    
race_5Other              0.257318    
ethnic_3Hispanic         0.213259    
lang_3Other              0.475863    
mstat_5Unknown           0.035048 *  
mstat_5Unmarried         0.159659    
mstat_5DivorcedSeparated 0.091919 .  
mstat_5Widow             0.199840    
act_tobYes               0.006230 ** 
max_ch                    < 2e-16 ***
IC                       4.08e-15 ***
RPL_THEME1               0.004611 ** 
RPL_THEME2               0.829256    
RPL_THEME3               0.217372    
RPL_THEME4               0.434096    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for Negative Binomial(49491.15) family taken to be 1)

    Null deviance: 3814.3  on 10870  degrees of freedom
Residual deviance: 3516.5  on 10849  degrees of freedom
  (4330 observations deleted due to missingness)
AIC: 36939

Number of Fisher Scoring iterations: 1

              Theta:  49491 
          Std. Err.:  60384 
Warning while fitting theta: iteration limit reached 

 2 x log-likelihood:  -36892.62 

broom::glance(readmit.nb2)

broom::tidy(readmit.nb2, exponentiate = TRUE)

model_performance(readmit.nb2)

# Indices of model performance

AIC       |       BIC | Nagelkerke's R2 |  RMSE | Sigma | Score_log | Score_spherical
-------------------------------------------------------------------------------------
36938.618 | 37106.376 |           0.091 | 0.677 | 0.569 |    -1.326 |           0.009

tbl_regression(readmit.nb2, label = list(age_yrs ~ "Age", gender~ "Gender", ethnic_3 ~ "Ethnicity", lang_3 ~ "Preferred Language", race_5 ~ "Race", act_tob ~"Active Tobacco Use", max_ch ~ "Charlson Comorbidity Index", IC ~ "Moderate to Severe Disease", ibd_3 ~ "IBD Type"), exponentiate = TRUE)

Characteristic	IRR1	95% CI1	p-value
IBD Type
CD	—	—
UC	0.99	0.96, 1.02	0.3
Unspecified	0.95	0.76, 1.16	0.6
Age	1.00	1.00, 1.00	<0.001
Gender
Male	—	—
Female	1.01	0.98, 1.04	0.6
Race
White	—	—
Black	1.10	1.04, 1.17	<0.001
Asian	0.94	0.85, 1.04	0.2
Native	0.97	0.77, 1.21	0.8
Other	1.05	0.96, 1.15	0.3
Ethnicity
NonHispanic	—	—
Hispanic	0.93	0.83, 1.04	0.2
Preferred Language
English	—	—
Other	0.94	0.80, 1.11	0.5
mstat_5
Married	—	—
Unknown	1.05	1.00, 1.10	0.035
Unmarried	0.97	0.94, 1.01	0.2
DivorcedSeparated	1.07	0.99, 1.15	0.092
Widow	1.07	0.96, 1.18	0.2
Active Tobacco Use
No	—	—
Yes	1.06	1.02, 1.11	0.006
Charlson Comorbidity Index	1.02	1.01, 1.02	<0.001
Moderate to Severe Disease	1.14	1.10, 1.17	<0.001
RPL_THEME1	1.13	1.04, 1.23	0.005
RPL_THEME2	1.01	0.94, 1.08	0.8
RPL_THEME3	1.03	0.98, 1.09	0.2
RPL_THEME4	0.98	0.92, 1.04	0.4
1 IRR = Incidence Rate Ratio, CI = Confidence Interval

LOS

SVI

LOS1 <- lm(mean_los ~ ibd_3 + age_yrs + gender + race_5 + ethnic_3 + lang_3  
           + mstat_5 + act_tob + max_ch + IC + RPL_THEMES,
               data = ut_clean1) 
summary(LOS1)

Call:
lm(formula = mean_los ~ ibd_3 + age_yrs + gender + race_5 + ethnic_3 + 
    lang_3 + mstat_5 + act_tob + max_ch + IC + RPL_THEMES, data = ut_clean1)

Residuals:
    Min      1Q  Median      3Q     Max 
 -7.241  -2.880  -1.308   1.041 190.494 

Coefficients:
                          Estimate Std. Error t value
(Intercept)               4.617986   0.473440   9.754
ibd_3UC                   0.383490   0.198973   1.927
ibd_3Unspecified         -0.330775   1.472417  -0.225
age_yrs                  -0.001101   0.006445  -0.171
genderFemale             -0.741483   0.191195  -3.878
race_5Black              -0.249550   0.351987  -0.709
race_5Asian              -0.225567   0.754134  -0.299
race_5Native              0.715623   1.746401   0.410
race_5Other              -0.493310   0.581214  -0.849
ethnic_3Hispanic          0.943411   0.667113   1.414
lang_3Other               0.179694   0.958232   0.188
mstat_5Unknown            0.655235   0.260748   2.513
mstat_5Unmarried          0.491969   0.254576   1.933
mstat_5DivorcedSeparated  0.018758   0.504688   0.037
mstat_5Widow              0.493465   0.618031   0.798
act_tobYes               -0.397432   0.287285  -1.383
max_ch                    0.106951   0.019330   5.533
IC                       -0.257667   0.215850  -1.194
RPL_THEMES                1.003162   0.375009   2.675
                         Pr(>|t|)    
(Intercept)               < 2e-16 ***
ibd_3UC                  0.054005 .  
ibd_3Unspecified         0.822265    
age_yrs                  0.864319    
genderFemale             0.000107 ***
race_5Black              0.478380    
race_5Asian              0.764874    
race_5Native             0.681996    
race_5Other              0.396066    
ethnic_3Hispanic         0.157388    
lang_3Other              0.851257    
mstat_5Unknown           0.012012 *  
mstat_5Unmarried         0.053367 .  
mstat_5DivorcedSeparated 0.970354    
mstat_5Widow             0.424658    
act_tobYes               0.166616    
max_ch                   3.34e-08 ***
IC                       0.232651    
RPL_THEMES               0.007502 ** 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 6.028 on 4106 degrees of freedom
  (11084 observations deleted due to missingness)
Multiple R-squared:  0.02022,   Adjusted R-squared:  0.01592 
F-statistic: 4.707 on 18 and 4106 DF,  p-value: 1.679e-10

broom::glance(LOS1)

broom::tidy(LOS1)

model_performance(LOS1)

# Indices of model performance

AIC       |       BIC |    R2 | R2 (adj.) |  RMSE | Sigma
---------------------------------------------------------
26548.016 | 26674.512 | 0.020 |     0.016 | 6.014 | 6.028

tbl_regression(LOS1, label = list(age_yrs ~ "Age", gender~ "Gender", ibd_3 ~ "IBD Type"))

Characteristic	Beta	95% CI1	p-value
IBD Type
CD	—	—
UC	0.38	-0.01, 0.77	0.054
Unspecified	-0.33	-3.2, 2.6	0.8
Age	0.00	-0.01, 0.01	0.9
Gender
Male	—	—
Female	-0.74	-1.1, -0.37	<0.001
race_5
White	—	—
Black	-0.25	-0.94, 0.44	0.5
Asian	-0.23	-1.7, 1.3	0.8
Native	0.72	-2.7, 4.1	0.7
Other	-0.49	-1.6, 0.65	0.4
ethnic_3
NonHispanic	—	—
Hispanic	0.94	-0.36, 2.3	0.2
lang_3
English	—	—
Other	0.18	-1.7, 2.1	0.9
mstat_5
Married	—	—
Unknown	0.66	0.14, 1.2	0.012
Unmarried	0.49	-0.01, 0.99	0.053
DivorcedSeparated	0.02	-0.97, 1.0	>0.9
Widow	0.49	-0.72, 1.7	0.4
act_tob
No	—	—
Yes	-0.40	-0.96, 0.17	0.2
max_ch	0.11	0.07, 0.14	<0.001
IC	-0.26	-0.68, 0.17	0.2
RPL_THEMES	1.0	0.27, 1.7	0.008
1 CI = Confidence Interval

All themes

LOS2 <- lm(mean_los ~ ibd_3 + age_yrs + gender + race_5 + ethnic_3 + lang_3  
           + mstat_5 + act_tob + max_ch + IC + RPL_THEME1 + RPL_THEME2 + RPL_THEME3 + RPL_THEME4,
               data = ut_clean1) 
summary(LOS2)

Call:
lm(formula = mean_los ~ ibd_3 + age_yrs + gender + race_5 + ethnic_3 + 
    lang_3 + mstat_5 + act_tob + max_ch + IC + RPL_THEME1 + RPL_THEME2 + 
    RPL_THEME3 + RPL_THEME4, data = ut_clean1)

Residuals:
    Min      1Q  Median      3Q     Max 
 -7.453  -2.857  -1.242   1.107 189.850 

Coefficients:
                           Estimate Std. Error t value
(Intercept)               4.7660360  0.5077701   9.386
ibd_3UC                   0.3669323  0.1993127   1.841
ibd_3Unspecified         -0.3511920  1.4705746  -0.239
age_yrs                  -0.0005296  0.0064560  -0.082
genderFemale             -0.7427117  0.1914089  -3.880
race_5Black               0.0035331  0.3585456   0.010
race_5Asian               0.3282376  0.7658361   0.429
race_5Native              0.7940245  1.7442224   0.455
race_5Other              -0.2975341  0.5850926  -0.509
ethnic_3Hispanic          0.9932351  0.6667261   1.490
lang_3Other               0.6100951  0.9725180   0.627
mstat_5Unknown            0.7023405  0.2614293   2.687
mstat_5Unmarried          0.5613192  0.2555348   2.197
mstat_5DivorcedSeparated -0.0090254  0.5041620  -0.018
mstat_5Widow              0.5162948  0.6200364   0.833
act_tobYes               -0.4143760  0.2874691  -1.441
max_ch                    0.1045399  0.0193927   5.391
IC                       -0.2482566  0.2161495  -1.149
RPL_THEME1                0.3676135  0.5412177   0.679
RPL_THEME2                1.0232129  0.4736323   2.160
RPL_THEME3               -1.1303273  0.3486785  -3.242
RPL_THEME4                0.3138061  0.3962756   0.792
                         Pr(>|t|)    
(Intercept)               < 2e-16 ***
ibd_3UC                  0.065696 .  
ibd_3Unspecified         0.811263    
age_yrs                  0.934629    
genderFemale             0.000106 ***
race_5Black              0.992138    
race_5Asian              0.668237    
race_5Native             0.648967    
race_5Other              0.611113    
ethnic_3Hispanic         0.136375    
lang_3Other              0.530474    
mstat_5Unknown           0.007249 ** 
mstat_5Unmarried         0.028102 *  
mstat_5DivorcedSeparated 0.985718    
mstat_5Widow             0.405071    
act_tobYes               0.149531    
max_ch                   7.42e-08 ***
IC                       0.250812    
RPL_THEME1               0.497028    
RPL_THEME2               0.030803 *  
RPL_THEME3               0.001198 ** 
RPL_THEME4               0.428472    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 6.018 on 4081 degrees of freedom
  (11106 observations deleted due to missingness)
Multiple R-squared:  0.02529,   Adjusted R-squared:  0.02028 
F-statistic: 5.043 on 21 and 4081 DF,  p-value: 4.065e-13

broom::glance(LOS2)

broom::tidy(LOS2)

model_performance(LOS2)

# Indices of model performance

AIC       |       BIC |    R2 | R2 (adj.) |  RMSE | Sigma
---------------------------------------------------------
26396.064 | 26541.411 | 0.025 |     0.020 | 6.002 | 6.018

tbl_regression(LOS2, label = list(age_yrs ~ "Age", gender~ "Gender", ibd_3 ~ "IBD Type"))

Characteristic	Beta	95% CI1	p-value
IBD Type
CD	—	—
UC	0.37	-0.02, 0.76	0.066
Unspecified	-0.35	-3.2, 2.5	0.8
Age	0.00	-0.01, 0.01	>0.9
Gender
Male	—	—
Female	-0.74	-1.1, -0.37	<0.001
race_5
White	—	—
Black	0.00	-0.70, 0.71	>0.9
Asian	0.33	-1.2, 1.8	0.7
Native	0.79	-2.6, 4.2	0.6
Other	-0.30	-1.4, 0.85	0.6
ethnic_3
NonHispanic	—	—
Hispanic	0.99	-0.31, 2.3	0.14
lang_3
English	—	—
Other	0.61	-1.3, 2.5	0.5
mstat_5
Married	—	—
Unknown	0.70	0.19, 1.2	0.007
Unmarried	0.56	0.06, 1.1	0.028
DivorcedSeparated	-0.01	-1.0, 0.98	>0.9
Widow	0.52	-0.70, 1.7	0.4
act_tob
No	—	—
Yes	-0.41	-0.98, 0.15	0.15
max_ch	0.10	0.07, 0.14	<0.001
IC	-0.25	-0.67, 0.18	0.3
RPL_THEME1	0.37	-0.69, 1.4	0.5
RPL_THEME2	1.0	0.09, 2.0	0.031
RPL_THEME3	-1.1	-1.8, -0.45	0.001
RPL_THEME4	0.31	-0.46, 1.1	0.4
1 CI = Confidence Interval

save(ut_clean1, file = "ut_clean1.rda")