This is a report focusing on end of life care preferences and their implementation by the clinical care team. In particular, we will focus on the National Quality Forum Measure #1626:
Percentage of vulnerable adults admitted to ICU who survive at least 48 hours who have their care preferences documented within 48 hours OR documentation as to why this was not done.
For our purposeses:
Care preference documentation within 48hrs of admission, age
>74.
Care provider/Patient interactions from the MIMIC-III Intensive Care Unit Database are analyzed within the framework of a hospital admission utilizing a series of Mixed Effects Models. The patient is the random effect, so we can explain the difference in variance explained by patients compared to that explained by other effects (when the patient is introduced to other demographic and clinical variables) in the model.
The use of fixed and random effects is deliberate. A fixed effect is used when a categorical variable contains all levels of interest. A random effect is used when a categorical variable contains only a sample of all levels of interest, the samples which do not exist are unobserved variables and, since their distributions do not exist, they will be estimated. The levels of the random effect are assumed to be independent of each other; this assumption requires a priori knowledge of the phenomenon being modeled.
We will use the Sequential Organ Failure Assessment (SOFA) from Illness Severity Scores (github) to determine the severity of illness at admission.
attending_check <- function(dat){
## Temporary data frame
tmp_frame <- data.frame()
## Results frame
res <- data.frame()
## For each hospital admission
for (name in unique(dat$HADM_ID)){
## Subset admission
tmp_frame <- dat[dat$HADM_ID == name, ]
## If any care providers are "Attending"
if (any("Attending" %in% tmp_frame$CG_DESCRIPTION)){
## add hospital admission to results
res <- rbind(res, tmp_frame)
}
}
## Return control to outer level
return(res)
}
provider_caremeasure_check <- function(dat){
## Create Caremeasure Variable
dat$NQF <- rep(0, each = nrow(dat))
## Temporary hospital admission data frame
tmp_frame <- data.frame()
## Temporary careprovider data frame
cg_frame <- data.frame()
## Results frame
res <- data.frame()
## For each hospital admission
for (name in unique(dat$HADM_ID)){
## Subset admission
tmp_frame <- dat[dat$HADM_ID == name, ]
## For each care provider
for (id in unique(tmp_frame$CGID)){
## Subset provider
cg_frame <- tmp_frame[tmp_frame$CGID == id, ]
## if care measure were implemented during the admission by the provider
if (any(cg_frame$CIM.machine == 1)){
## Apply credit for entire admission
cg_frame$NQF <- rep(1, each = nrow(cg_frame))
}
## Bind
res <- rbind(res, cg_frame)
}
}
## Return control to outer level
return(res)
}
team_caremeasure_check <- function(dat){
## Create Caremeasure Variable
dat$NQF <- rep(0, each = nrow(dat))
## Temporary hospital admission data frame
tmp_frame <- data.frame()
## Results frame
res <- data.frame()
## For each hospital admission
for (name in unique(dat$HADM_ID)){
## Subset admission
tmp_frame <- dat[dat$HADM_ID == name, ]
## If ANY care provider documented care preferences
if (any(tmp_frame$CIM.machine == 1)){
## Apply credit for entire admission
tmp_frame$NQF <- rep(1, each = nrow(tmp_frame))
}
## Bind
res <- rbind(res, tmp_frame)
}
## Return control to outer level
return(res)
}
plotDat <- function(dat, column, x_col, bs, mn, xl, yl){
tmp <- as.matrix(table(dat[[column]], dat[[x_col]]))
prop <- prop.table(tmp, margin = 2)#2 for column-wise proportions
par(mar = c(5.0, 4.0, 4.0, 15), xpd = TRUE)
barplot(prop, col = cm.colors(length(rownames(prop))), beside = bs, width = 2, main = mn, xlab = xl, ylab = yl)
legend("topright", inset = c(-0.5,0), fill = cm.colors(length(rownames(prop))), legend=rownames(prop))
}
expire_check <- function(dat){
## Create Caremeasure Variable
dat$HOSP_DEATH <- rep(0, each = nrow(dat))
## Temporary data frame
tmp_frame <- data.frame()
## Results frame
res <- data.frame()
## For each hospital admission
for (name in unique(dat$HADM_ID)){
## Subset admission
tmp_frame <- dat[dat$HADM_ID == name, ]
## If any admission has a death
## Check if pt expired in hospital
if (any(tmp_frame$HOSPITAL_EXPIRE_FLAG == 1)){
tmp_frame$HOSP_DEATH <- rep(1, each = nrow(tmp_frame))
}
## add hospital admission to results
res <- rbind(res, tmp_frame)
}
## Return control to outer level
return(res)
}
rocplot <- function(pred, truth, ...) {
predob = prediction(pred, truth)
perf = performance(predob, "tpr", "fpr")
plot(perf, ...)
area <- auc(truth, pred)
area <- format(round(area, 4), nsmall = 4)
text(x=0.7, y=0.1, labels = paste("AUC =", area))
# the reference x=y line
segments(x0=0, y0=0, x1=1, y1=1, col="gray", lty=2)
}
model_info <- function(fit){
#Summary info
model_sum <- summary(fit)
#Odds ratio, confidence interval
odds_ratio <- cbind(OR = exp(fit$coef), exp(confint(fit)))
#Create list for return
my_list <- list(model_sum, odds_ratio)
#names
names(my_list) <- c("Model Summary","OR Summary")
return(my_list)
}
## Latest Dataset of NeuroNER Predictions
dat <- read.csv("~/nqf_caregivers/data/20180607_EOL_data_ICU.csv", header = T, stringsAsFactors = F)
dim(dat)
## [1] 10250 57
## Note: Notes had been logged by multiple Care providers, we will reintroduce those annotations
## Load Labeled Note Data for NQF Caremeasure Cohort (From NOTEEVENTS table)
tmp <- read.csv("~/nqf_caregivers/data/note_labels_over75.csv", header = T, stringsAsFactors = F)
dim(tmp)
## [1] 11575 25
## Keep only TEXT and ROW_ID from tmp
tmp <- tmp[ ,c("ROW_ID", "TEXT")]
## Inner join
dat <- merge(tmp, dat, by = "TEXT")
## Clean tmp
rm(tmp)
## Check column names
colnames(dat)
## [1] "TEXT" "ROW_ID.x" "SUBJECT_ID"
## [4] "HADM_ID" "ROW_ID.y" "CHARTDATE"
## [7] "CHARTTIME" "STORETIME" "CATEGORY"
## [10] "DESCRIPTION" "CGID" "ISERROR"
## [13] "ADMITTIME" "DISCHTIME" "DEATHTIME"
## [16] "ADMISSION_TYPE" "ADMISSION_LOCATION" "DISCHARGE_LOCATION"
## [19] "INSURANCE" "LANGUAGE" "RELIGION"
## [22] "MARITAL_STATUS" "ETHNICITY" "EDREGTIME"
## [25] "EDOUTTIME" "DIAGNOSIS" "HOSPITAL_EXPIRE_FLAG"
## [28] "HAS_CHARTEVENTS_DATA" "GENDER" "DOB"
## [31] "DOD" "DOD_HOSP" "DOD_SSN"
## [34] "EXPIRE_FLAG" "ICUSTAY_ID" "DBSOURCE"
## [37] "FIRST_CAREUNIT" "LAST_CAREUNIT" "FIRST_WARDID"
## [40] "LAST_WARDID" "INTIME" "OUTTIME"
## [43] "LOS" "AGE" "ADMISSION_NUMBER"
## [46] "DAYS_UNTIL_DEATH" "TIME_SINCE_ADMIT" "CGID.1"
## [49] "HADM_ID.1" "FAM.machine" "CIM.machine"
## [52] "LIM.machine" "CAR.machine" "COD.machine"
## [55] "check.CGID" "check.dadm_id" "CIM.or.FAM"
## [58] "Died.in.Hospital"
## What is HADM_ID.1?
head(table(dat$HADM_ID.1))
##
## #N/A 100102 100153 100347 100391 100525
## 32 8 3 12 15 2
## What is HADM_ID?
head(table(dat$HADM_ID))
##
## 100102 100153 100347 100391 100525 100575
## 8 3 12 15 2 11
## #N/A? Clean HADM_ID.1
dat$HADM_ID.1 <- NULL
## What is CGID.1?
head(table(dat$CGID.1))
##
## #N/A 14010 14022 14037 14045 14056
## 32 22 1 93 37 14
## What is CGID
head(table(dat$CGID))
##
## 14010 14022 14037 14045 14056 14080
## 22 1 93 37 14 6
## #N/A? Clean CGID.1
dat$CGID.1 <- NULL
## What is check.CGID
head(table(dat$check.CGID))
## 0
## 11575
## Clean it
dat$check.CGID <- NULL
## What is check.dadm_id?
head(table(dat$check.dadm_id))
## 0
## 11575
## Clean it
dat$check.dadm_id <- NULL
## Clean column names
dat$ROW_ID.y <- NULL
colnames(dat)[which(colnames(dat) == "ROW_ID.x")] <- "ROW_ID"
## Load CAREGIVERS Table for join on CGID
cg <- read.csv("~/nqf_caregivers/data/mimic/CAREGIVERS.csv",
header = T, stringsAsFactors = F)
## Change column name of "NOTEEVENTS.DESCRIPTION" to explicitly mention that it describes the note
colnames(dat)[which(colnames(dat) == "DESCRIPTION")] <- "NOTE_DESCRIPTION"
## Change column name of "CAREGIVERS. DESCRIPTION" to explicitly mention that it describes the careprovider
colnames(cg)[which(colnames(cg) == "DESCRIPTION")] <- "CG_DESCRIPTION"
## Remove ROW_ID from CG
cg$ROW_ID <- NULL
## Remove TEXT
dat$TEXT <- NULL
## Merge to caregivers
dat <- merge(dat, cg, by = "CGID")
dim(dat)
## [1] 11575 54
## Clean CG
rm(cg)
sofa <- read.csv("~/nqf_caregivers/data/sofa.csv", header = T, stringsAsFactors = F)
#oasis <- read.csv("~/nqf_caregivers/data/oasis.csv", header = T, stringsAsFactors = F)
#saps <- read.csv("~/nqf_caregivers/data/saps.csv", header = T, stringsAsFactors = F)
colnames(sofa) <- toupper(colnames(sofa))
dat <- merge(dat, sofa, by = c("SUBJECT_ID", "HADM_ID", "ICUSTAY_ID"))
dim(dat)
## [1] 11575 61
## Clean environment
rm(sofa)
## Clean ethnicity to Black/White/Other
dat[(grepl("WHITE|PORTUGUESE", dat$ETHNICITY)),]$ETHNICITY <- "WHITE"
dat[(grepl("ASIAN", dat$ETHNICITY)),]$ETHNICITY <- "OTHER"
dat[(grepl("BLACK", dat$ETHNICITY)),]$ETHNICITY <- "BLACK"
dat[(grepl("HISPANIC", dat$ETHNICITY)),]$ETHNICITY <- "OTHER"
dat[(grepl("MIDDLE|NATIVE|MULTI|DECLINED|UNABLE|OTHER|NOT", dat$ETHNICITY)),]$ETHNICITY <- "OTHER"
## Clean Marital Status to Married, Single, Widowed, Unknown
dat$MARITAL_STATUS[dat$MARITAL_STATUS == ""] <- "UNKNOWN (DEFAULT)"
dat$MARITAL_STATUS[dat$MARITAL_STATUS == "UNKNOWN (DEFAULT)"] <- "UNKNOWN"
dat$MARITAL_STATUS[dat$MARITAL_STATUS == "SEPARATED"] <- "SINGLE"
dat$MARITAL_STATUS[dat$MARITAL_STATUS == "DIVORCED"] <- "SINGLE"
## Clean admission data
dat$ADMISSION_TYPE[dat$ADMISSION_TYPE == "URGENT"] <- "EMERGENCY"
At this level, we will see if providers documented patient care preferences at least once in the first 48hrs of the hospital admission. Here, the clinician/patient interaction is the unit of analysis.
temp <- provider_caremeasure_check(attending_check(dat))
temp <- expire_check(temp)
plotDat(temp, "NQF", "GENDER", F, "Gender (Provider-Level)", "Gender", "Frequency")
test <- table(temp$GENDER, temp$NQF)
test
##
## 0 1
## F 2187 3344
## M 2768 2860
chisq.test(test)
##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: test
## X-squared = 104.66, df = 1, p-value < 2.2e-16
pairwiseNominalIndependence(
as.matrix(test),
fisher = F, gtest = F, chisq = T, method = "fdr")
## Comparison p.Chisq p.adj.Chisq
## 1 F : M 1.45e-24 1.45e-24
plotDat(temp, "NQF", "ADMISSION_TYPE", F, "Admission Type (Provider-Level)", "Admission Type", "Frequency")
test <- table(temp$ADMISSION_TYPE, temp$NQF)
test
##
## 0 1
## ELECTIVE 294 28
## EMERGENCY 4661 6176
chisq.test(test)
##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: test
## X-squared = 293.48, df = 1, p-value < 2.2e-16
pairwiseNominalIndependence(
as.matrix(test),
fisher = F, gtest = F, chisq = T, method = "fdr")
## Comparison p.Chisq p.adj.Chisq
## 1 ELECTIVE : EMERGENCY 8.65e-66 8.65e-66
plotDat(temp, "NQF","ETHNICITY", F, "Ethnicity (Provider-Level)", "Ethnicity", "Frequency")
test <- table(temp$ETHNICITY, temp$NQF)
test
##
## 0 1
## BLACK 405 619
## OTHER 475 636
## WHITE 4075 4949
chisq.test(test)
##
## Pearson's Chi-squared test
##
## data: test
## X-squared = 13.069, df = 2, p-value = 0.001452
pairwiseNominalIndependence(
as.matrix(test),
fisher = F, gtest = F, chisq = T, method = "fdr")
## Comparison p.Chisq p.adj.Chisq
## 1 BLACK : OTHER 0.145000 0.14500
## 2 BLACK : WHITE 0.000706 0.00212
## 3 OTHER : WHITE 0.137000 0.14500
plotDat(temp, "NQF", "MARITAL_STATUS", F, "Marital Status (Provider-Level)", "Marital Status", "Frequency")
test <- table(temp$MARITAL_STATUS, temp$NQF)
test
##
## 0 1
## MARRIED 2315 2343
## SINGLE 884 1469
## UNKNOWN 227 187
## WIDOWED 1529 2205
chisq.test(test)
##
## Pearson's Chi-squared test
##
## data: test
## X-squared = 133.74, df = 3, p-value < 2.2e-16
pairwiseNominalIndependence(
as.matrix(test),
fisher = F, gtest = F, chisq = T, method = "fdr")
## Comparison p.Chisq p.adj.Chisq
## 1 MARRIED : SINGLE 7.70e-22 4.62e-21
## 2 MARRIED : UNKNOWN 5.12e-02 5.12e-02
## 3 MARRIED : WIDOWED 1.53e-15 4.59e-15
## 4 SINGLE : UNKNOWN 5.65e-11 1.13e-10
## 5 SINGLE : WIDOWED 9.39e-03 1.13e-02
## 6 UNKNOWN : WIDOWED 7.79e-08 1.17e-07
x <- barplot(table(mtcars$cyl), xaxt="n")
labs <- paste(names(table(mtcars$cyl)), "cylinders")
text(cex=1, x=x-.25, y=-1.25, labs, xpd=TRUE, srt=45)
plotDat(temp, "NQF", "FIRST_CAREUNIT", F, "First Careunit (Provider-Level)", "First Careunit", "Frequency")
test <- table(temp$FIRST_CAREUNIT, temp$NQF)
test
##
## 0 1
## CCU 696 905
## CSRU 273 111
## MICU 3031 4312
## SICU 593 546
## TSICU 362 330
chisq.test(test)
##
## Pearson's Chi-squared test
##
## data: test
## X-squared = 185.04, df = 4, p-value < 2.2e-16
pairwiseNominalIndependence(
as.matrix(test),
fisher = F, gtest = F, chisq = T, method = "fdr")
## Comparison p.Chisq p.adj.Chisq
## 1 CCU : CSRU 4.13e-22 2.06e-21
## 2 CCU : MICU 1.13e-01 1.26e-01
## 3 CCU : SICU 1.08e-05 1.54e-05
## 4 CCU : TSICU 1.17e-04 1.46e-04
## 5 CSRU : MICU 2.10e-30 2.10e-29
## 6 CSRU : SICU 1.10e-10 2.75e-10
## 7 CSRU : TSICU 2.91e-09 5.82e-09
## 8 MICU : SICU 9.46e-12 3.15e-11
## 9 MICU : TSICU 2.43e-08 4.05e-08
## 10 SICU : TSICU 9.56e-01 9.56e-01
boxplot(temp$AGE ~ temp$NQF,
main = "Caremeasure Implementation\n by Age (Provider-Level)",
xlab = "Implementation (1 == Yes)",
ylab = "Age (Years)")
t.test(temp$AGE ~ temp$NQF)
##
## Welch Two Sample t-test
##
## data: temp$AGE by temp$NQF
## t = -12.58, df = 10825, p-value < 2.2e-16
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -1.425367 -1.041058
## sample estimates:
## mean in group 0 mean in group 1
## 83.40638 84.63959
boxplot(temp$SOFA ~ temp$NQF,
main = "Caremeasure Implementation\n by SOFA Score (Provider-Level)",
xlab = "Implementation (1 == Yes)",
ylab = "SOFA Score")
t.test(temp$SOFA ~ temp$NQF)
##
## Welch Two Sample t-test
##
## data: temp$SOFA by temp$NQF
## t = 0.11386, df = 10337, p-value = 0.9094
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -0.1062057 0.1193046
## sample estimates:
## mean in group 0 mean in group 1
## 5.037336 5.030787
Continuous data, Age and SOFA, will be standardized in the form:
\[z = \frac{x - \mu}{\sigma}\]
temp$AGE <- (temp$AGE - mean(temp$AGE))/sd(temp$AGE)
temp$SOFA <- (temp$SOFA - mean(temp$SOFA)/sd(temp$SOFA))
gm_initial <- glmer(NQF ~ (1 | CGID),
data = temp,
family = binomial,
control = glmerControl(optimizer = "bobyqa"),
nAGQ = 10)
## View
summary(gm_initial)
## Generalized linear mixed model fit by maximum likelihood (Adaptive
## Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
## Family: binomial ( logit )
## Formula: NQF ~ (1 | CGID)
## Data: temp
## Control: glmerControl(optimizer = "bobyqa")
##
## AIC BIC logLik deviance df.resid
## 13715.3 13730.0 -6855.7 13711.3 11157
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.2922 -0.8253 0.3867 0.7561 3.6611
##
## Random effects:
## Groups Name Variance Std.Dev.
## CGID (Intercept) 2.138 1.462
## Number of obs: 11159, groups: CGID, 493
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.16423 0.07617 2.156 0.0311 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
sjstats::icc(gm_initial)
##
## Generalized linear mixed model
## Family: binomial (logit)
## Formula: NQF ~ (1 | CGID)
##
## ICC (CGID): 0.393846
gm_a <- glmer(NQF ~ GENDER +
AGE +
ETHNICITY +
MARITAL_STATUS +
(1 | CGID),
data = temp,
family = binomial,
control = glmerControl(optimizer = "bobyqa"),
nAGQ = 10)
## View
summary(gm_a)
## Generalized linear mixed model fit by maximum likelihood (Adaptive
## Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
## Family: binomial ( logit )
## Formula: NQF ~ GENDER + AGE + ETHNICITY + MARITAL_STATUS + (1 | CGID)
## Data: temp
## Control: glmerControl(optimizer = "bobyqa")
##
## AIC BIC logLik deviance df.resid
## 13472.3 13538.1 -6727.1 13454.3 11150
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.4389 -0.7877 0.3476 0.7307 4.1859
##
## Random effects:
## Groups Name Variance Std.Dev.
## CGID (Intercept) 2.112 1.453
## Number of obs: 11159, groups: CGID, 493
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.22100 0.11608 1.904 0.05692 .
## GENDERM -0.33510 0.04934 -6.792 1.11e-11 ***
## AGE 0.25396 0.02348 10.817 < 2e-16 ***
## ETHNICITYOTHER 0.03976 0.10762 0.369 0.71182
## ETHNICITYWHITE -0.02416 0.08113 -0.298 0.76585
## MARITAL_STATUSSINGLE 0.38954 0.06319 6.164 7.08e-10 ***
## MARITAL_STATUSUNKNOWN -0.14376 0.12465 -1.153 0.24878
## MARITAL_STATUSWIDOWED 0.17954 0.05818 3.086 0.00203 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) GENDER AGE ETHNICITYO ETHNICITYW MARITAL_STATUSS
## GENDERM -0.266
## AGE 0.089 -0.054
## ETHNICITYOT -0.495 -0.037 -0.051
## ETHNICITYWH -0.647 -0.064 -0.064 0.675
## MARITAL_STATUSS -0.302 0.248 0.020 0.093 0.095
## MARITAL_STATUSU -0.130 0.066 -0.077 -0.068 0.048 0.161
## MARITAL_STATUSW -0.369 0.368 -0.177 0.101 0.096 0.432
## MARITAL_STATUSU
## GENDERM
## AGE
## ETHNICITYOT
## ETHNICITYWH
## MARITAL_STATUSS
## MARITAL_STATUSU
## MARITAL_STATUSW 0.205
sjt.glmer(gm_a)
| NQF | ||||
| Odds Ratio | CI | p | ||
| Fixed Parts | ||||
| (Intercept) | 1.25 | 0.99 – 1.57 | .057 | |
| GENDER (M) | 0.72 | 0.65 – 0.79 | <.001 | |
| AGE | 1.29 | 1.23 – 1.35 | <.001 | |
| ETHNICITY (OTHER) | 1.04 | 0.84 – 1.28 | .712 | |
| ETHNICITY (WHITE) | 0.98 | 0.83 – 1.14 | .766 | |
| MARITAL_STATUS (SINGLE) | 1.48 | 1.30 – 1.67 | <.001 | |
| MARITAL_STATUS (UNKNOWN) | 0.87 | 0.68 – 1.11 | .249 | |
| MARITAL_STATUS (WIDOWED) | 1.20 | 1.07 – 1.34 | .002 | |
| Random Parts | ||||
| τ00, CGID | 2.112 | |||
| NCGID | 493 | |||
| ICCCGID | 0.391 | |||
| Observations | 11159 | |||
| Deviance | 12244.947 | |||
gm_b <- glmer(NQF ~ SOFA +
(1 | CGID),
data = temp,
family = binomial,
control = glmerControl(optimizer = "bobyqa"),
nAGQ = 10)
## View
summary(gm_b)
## Generalized linear mixed model fit by maximum likelihood (Adaptive
## Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
## Family: binomial ( logit )
## Formula: NQF ~ SOFA + (1 | CGID)
## Data: temp
## Control: glmerControl(optimizer = "bobyqa")
##
## AIC BIC logLik deviance df.resid
## 13717.3 13739.3 -6855.7 13711.3 11156
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.2906 -0.8251 0.3866 0.7564 3.6627
##
## Random effects:
## Groups Name Variance Std.Dev.
## CGID (Intercept) 2.138 1.462
## Number of obs: 11159, groups: CGID, 493
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.1633086 0.0801375 2.038 0.0416 *
## SOFA 0.0002795 0.0075712 0.037 0.9706
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## SOFA -0.311
sjt.glmer(gm_b)
| NQF | ||||
| Odds Ratio | CI | p | ||
| Fixed Parts | ||||
| (Intercept) | 1.18 | 1.01 – 1.38 | .042 | |
| SOFA | 1.00 | 0.99 – 1.02 | .971 | |
| Random Parts | ||||
| τ00, CGID | 2.138 | |||
| NCGID | 493 | |||
| ICCCGID | 0.394 | |||
| Observations | 11159 | |||
| Deviance | 12489.152 | |||
gm_2a <- glmer(NQF ~ GENDER +
AGE +
ETHNICITY +
MARITAL_STATUS +
FIRST_CAREUNIT +
(1 | CGID),
data = temp,
family = binomial,
control = glmerControl(optimizer = "bobyqa"),
nAGQ = 10)
## View
summary(gm_2a)
## Generalized linear mixed model fit by maximum likelihood (Adaptive
## Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
## Family: binomial ( logit )
## Formula:
## NQF ~ GENDER + AGE + ETHNICITY + MARITAL_STATUS + FIRST_CAREUNIT +
## (1 | CGID)
## Data: temp
## Control: glmerControl(optimizer = "bobyqa")
##
## AIC BIC logLik deviance df.resid
## 13425.1 13520.2 -6699.5 13399.1 11146
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.3984 -0.7883 0.3440 0.7329 4.7920
##
## Random effects:
## Groups Name Variance Std.Dev.
## CGID (Intercept) 2.184 1.478
## Number of obs: 11159, groups: CGID, 493
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.089391 0.136182 -0.656 0.511562
## GENDERM -0.317199 0.049573 -6.399 1.57e-10 ***
## AGE 0.259135 0.023635 10.964 < 2e-16 ***
## ETHNICITYOTHER 0.084372 0.108248 0.779 0.435723
## ETHNICITYWHITE 0.006944 0.081463 0.085 0.932066
## MARITAL_STATUSSINGLE 0.402595 0.063371 6.353 2.11e-10 ***
## MARITAL_STATUSUNKNOWN -0.189578 0.125578 -1.510 0.131135
## MARITAL_STATUSWIDOWED 0.199230 0.058503 3.405 0.000660 ***
## FIRST_CAREUNITCSRU -0.453809 0.171432 -2.647 0.008117 **
## FIRST_CAREUNITMICU 0.301187 0.078861 3.819 0.000134 ***
## FIRST_CAREUNITSICU 0.347727 0.116786 2.977 0.002906 **
## FIRST_CAREUNITTSICU 0.764958 0.142851 5.355 8.56e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) GENDER AGE ETHNICITYO ETHNICITYW MARITAL_STATUSS
## GENDERM -0.242
## AGE 0.044 -0.050
## ETHNICITYOT -0.455 -0.033 -0.045
## ETHNICITYWH -0.577 -0.061 -0.062 0.673
## MARITAL_STATUSS -0.272 0.246 0.022 0.096 0.097
## MARITAL_STATUSU -0.083 0.064 -0.078 -0.069 0.043 0.158
## MARITAL_STATUSW -0.310 0.367 -0.176 0.103 0.098 0.431
## FIRST_CAREUNITC -0.205 -0.031 0.002 0.011 0.014 0.000
## FIRST_CAREUNITM -0.475 0.031 0.083 0.070 0.035 0.027
## FIRST_CAREUNITS -0.386 0.033 0.014 0.026 0.043 0.015
## FIRST_CAREUNITT -0.359 0.009 0.011 0.045 0.066 0.027
## MARITAL_STATUSU MARITAL_STATUSW FIRST_CAREUNITC
## GENDERM
## AGE
## ETHNICITYOT
## ETHNICITYWH
## MARITAL_STATUSS
## MARITAL_STATUSU
## MARITAL_STATUSW 0.202
## FIRST_CAREUNITC -0.031 -0.049
## FIRST_CAREUNITM -0.033 -0.010 0.300
## FIRST_CAREUNITS -0.041 -0.018 0.290
## FIRST_CAREUNITT -0.076 0.021 0.266
## FIRST_CAREUNITM FIRST_CAREUNITS
## GENDERM
## AGE
## ETHNICITYOT
## ETHNICITYWH
## MARITAL_STATUSS
## MARITAL_STATUSU
## MARITAL_STATUSW
## FIRST_CAREUNITC
## FIRST_CAREUNITM
## FIRST_CAREUNITS 0.540
## FIRST_CAREUNITT 0.441 0.556
sjt.glmer(gm_2a)
| NQF | ||||
| Odds Ratio | CI | p | ||
| Fixed Parts | ||||
| (Intercept) | 0.91 | 0.70 – 1.19 | .512 | |
| GENDER (M) | 0.73 | 0.66 – 0.80 | <.001 | |
| AGE | 1.30 | 1.24 – 1.36 | <.001 | |
| ETHNICITY (OTHER) | 1.09 | 0.88 – 1.35 | .436 | |
| ETHNICITY (WHITE) | 1.01 | 0.86 – 1.18 | .932 | |
| MARITAL_STATUS (SINGLE) | 1.50 | 1.32 – 1.69 | <.001 | |
| MARITAL_STATUS (UNKNOWN) | 0.83 | 0.65 – 1.06 | .131 | |
| MARITAL_STATUS (WIDOWED) | 1.22 | 1.09 – 1.37 | <.001 | |
| FIRST_CAREUNIT (CSRU) | 0.64 | 0.45 – 0.89 | .008 | |
| FIRST_CAREUNIT (MICU) | 1.35 | 1.16 – 1.58 | <.001 | |
| FIRST_CAREUNIT (SICU) | 1.42 | 1.13 – 1.78 | .003 | |
| FIRST_CAREUNIT (TSICU) | 2.15 | 1.62 – 2.84 | <.001 | |
| Random Parts | ||||
| τ00, CGID | 2.184 | |||
| NCGID | 493 | |||
| ICCCGID | 0.399 | |||
| Observations | 11159 | |||
| Deviance | 12180.260 | |||
gm_2b <- glmer(NQF ~ SOFA +
FIRST_CAREUNIT +
(1 | CGID),
data = temp,
family = binomial,
control = glmerControl(optimizer = "bobyqa"),
nAGQ = 10)
## View
summary(gm_2b)
## Generalized linear mixed model fit by maximum likelihood (Adaptive
## Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
## Family: binomial ( logit )
## Formula: NQF ~ SOFA + FIRST_CAREUNIT + (1 | CGID)
## Data: temp
## Control: glmerControl(optimizer = "bobyqa")
##
## AIC BIC logLik deviance df.resid
## 13673.1 13724.3 -6829.5 13659.1 11152
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.3161 -0.8161 0.3780 0.7567 4.2014
##
## Random effects:
## Groups Name Variance Std.Dev.
## CGID (Intercept) 2.243 1.498
## Number of obs: 11159, groups: CGID, 493
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.077708 0.103960 -0.747 0.454771
## SOFA 0.003185 0.007625 0.418 0.676141
## FIRST_CAREUNITCSRU -0.466875 0.168642 -2.768 0.005633 **
## FIRST_CAREUNITMICU 0.241646 0.077318 3.125 0.001776 **
## FIRST_CAREUNITSICU 0.382926 0.115541 3.314 0.000919 ***
## FIRST_CAREUNITTSICU 0.734990 0.140761 5.222 1.77e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) SOFA FIRST_CAREUNITC FIRST_CAREUNITM
## SOFA -0.245
## FIRST_CAREUNITC -0.270 -0.020
## FIRST_CAREUNITM -0.570 -0.018 0.302
## FIRST_CAREUNITS -0.474 0.051 0.290 0.536
## FIRST_CAREUNITT -0.421 0.047 0.268 0.438
## FIRST_CAREUNITS
## SOFA
## FIRST_CAREUNITC
## FIRST_CAREUNITM
## FIRST_CAREUNITS
## FIRST_CAREUNITT 0.557
sjt.glmer(gm_2b)
| NQF | ||||
| Odds Ratio | CI | p | ||
| Fixed Parts | ||||
| (Intercept) | 0.93 | 0.75 – 1.13 | .455 | |
| SOFA | 1.00 | 0.99 – 1.02 | .676 | |
| FIRST_CAREUNIT (CSRU) | 0.63 | 0.45 – 0.87 | .006 | |
| FIRST_CAREUNIT (MICU) | 1.27 | 1.09 – 1.48 | .002 | |
| FIRST_CAREUNIT (SICU) | 1.47 | 1.17 – 1.84 | <.001 | |
| FIRST_CAREUNIT (TSICU) | 2.09 | 1.58 – 2.75 | <.001 | |
| Random Parts | ||||
| τ00, CGID | 2.243 | |||
| NCGID | 493 | |||
| ICCCGID | 0.405 | |||
| Observations | 11159 | |||
| Deviance | 12421.574 | |||
gm_icu <- glmer(NQF ~
GENDER +
AGE +
ETHNICITY +
MARITAL_STATUS +
SOFA +
FIRST_CAREUNIT +
(1 | CGID),
data = temp,
family = binomial,
control = glmerControl(optimizer = "bobyqa"),
nAGQ = 10)
## View
summary(gm_icu)
## Generalized linear mixed model fit by maximum likelihood (Adaptive
## Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
## Family: binomial ( logit )
## Formula:
## NQF ~ GENDER + AGE + ETHNICITY + MARITAL_STATUS + SOFA + FIRST_CAREUNIT +
## (1 | CGID)
## Data: temp
## Control: glmerControl(optimizer = "bobyqa")
##
## AIC BIC logLik deviance df.resid
## 13420.5 13523.0 -6696.3 13392.5 11145
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.2806 -0.7926 0.3433 0.7313 4.9852
##
## Random effects:
## Groups Name Variance Std.Dev.
## CGID (Intercept) 2.194 1.481
## Number of obs: 11159, groups: CGID, 493
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.158075 0.138990 -1.137 0.255410
## GENDERM -0.331264 0.049905 -6.638 3.18e-11 ***
## AGE 0.263075 0.023705 11.098 < 2e-16 ***
## ETHNICITYOTHER 0.080285 0.108411 0.741 0.458960
## ETHNICITYWHITE 0.014362 0.081614 0.176 0.860314
## MARITAL_STATUSSINGLE 0.408120 0.063430 6.434 1.24e-10 ***
## MARITAL_STATUSUNKNOWN -0.179696 0.125682 -1.430 0.152785
## MARITAL_STATUSWIDOWED 0.203522 0.058556 3.476 0.000510 ***
## SOFA 0.020113 0.007848 2.563 0.010383 *
## FIRST_CAREUNITCSRU -0.462606 0.171437 -2.698 0.006967 **
## FIRST_CAREUNITMICU 0.297664 0.078919 3.772 0.000162 ***
## FIRST_CAREUNITSICU 0.362566 0.117035 3.098 0.001949 **
## FIRST_CAREUNITTSICU 0.782931 0.143072 5.472 4.44e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation matrix not shown by default, as p = 13 > 12.
## Use print(x, correlation=TRUE) or
## vcov(x) if you need it
sjt.glmer(gm_icu)
| NQF | ||||
| Odds Ratio | CI | p | ||
| Fixed Parts | ||||
| (Intercept) | 0.85 | 0.65 – 1.12 | .255 | |
| GENDER (M) | 0.72 | 0.65 – 0.79 | <.001 | |
| AGE | 1.30 | 1.24 – 1.36 | <.001 | |
| ETHNICITY (OTHER) | 1.08 | 0.88 – 1.34 | .459 | |
| ETHNICITY (WHITE) | 1.01 | 0.86 – 1.19 | .860 | |
| MARITAL_STATUS (SINGLE) | 1.50 | 1.33 – 1.70 | <.001 | |
| MARITAL_STATUS (UNKNOWN) | 0.84 | 0.65 – 1.07 | .153 | |
| MARITAL_STATUS (WIDOWED) | 1.23 | 1.09 – 1.37 | <.001 | |
| SOFA | 1.02 | 1.00 – 1.04 | .010 | |
| FIRST_CAREUNIT (CSRU) | 0.63 | 0.45 – 0.88 | .007 | |
| FIRST_CAREUNIT (MICU) | 1.35 | 1.15 – 1.57 | <.001 | |
| FIRST_CAREUNIT (SICU) | 1.44 | 1.14 – 1.81 | .002 | |
| FIRST_CAREUNIT (TSICU) | 2.19 | 1.65 – 2.90 | <.001 | |
| Random Parts | ||||
| τ00, CGID | 2.194 | |||
| NCGID | 493 | |||
| ICCCGID | 0.400 | |||
| Observations | 11159 | |||
| Deviance | 12172.161 | |||
par(mfrow=c(2,3))
rocplot(as.numeric(predict(gm_initial, type="response")), temp$NQF, col="blue", main = "Initial Model")
rocplot(as.numeric(predict(gm_a, type="response")), temp$NQF, col="blue", main = "Model 1a (Demographics)")
rocplot(as.numeric(predict(gm_b, type="response")), temp$NQF, col="blue", main = "Model 1b (Clinical Variables)")
rocplot(as.numeric(predict(gm_2a, type="response")), temp$NQF, col="blue", main = "Model 2a (Demographics + ICU)")
rocplot(as.numeric(predict(gm_2b, type="response")), temp$NQF, col="blue", main = "Model 2b (Clinical Variables + ICU)")
rocplot(as.numeric(predict(gm_icu, type="response")), temp$NQF, col="blue", main = "Model 3 (All Covariates)")
Not meaningful.
gm_a <- glm(NQF ~ GENDER +
AGE +
ETHNICITY +
MARITAL_STATUS,
data = temp,
family = binomial(link = "logit"))
## View
model_info(gm_a)
## Waiting for profiling to be done...
## $`Model Summary`
##
## Call:
## glm(formula = NQF ~ GENDER + AGE + ETHNICITY + MARITAL_STATUS,
## family = binomial(link = "logit"), data = temp)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.6955 -1.2224 0.8935 1.0738 1.5704
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.41451 0.07676 5.400 6.67e-08 ***
## GENDERM -0.31305 0.04191 -7.470 8.04e-14 ***
## AGE 0.25069 0.02009 12.476 < 2e-16 ***
## ETHNICITYOTHER -0.00460 0.09117 -0.050 0.959760
## ETHNICITYWHITE -0.16072 0.06906 -2.327 0.019959 *
## MARITAL_STATUSSINGLE 0.40458 0.05420 7.464 8.37e-14 ***
## MARITAL_STATUSUNKNOWN -0.39663 0.10611 -3.738 0.000185 ***
## MARITAL_STATUSWIDOWED 0.10679 0.04917 2.172 0.029868 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 15330 on 11158 degrees of freedom
## Residual deviance: 14980 on 11151 degrees of freedom
## AIC: 14996
##
## Number of Fisher Scoring iterations: 4
##
##
## $`OR Summary`
## OR 2.5 % 97.5 %
## (Intercept) 1.5136245 1.3026783 1.7601383
## GENDERM 0.7312139 0.6735309 0.7937934
## AGE 1.2849074 1.2353570 1.3366010
## ETHNICITYOTHER 0.9954107 0.8324493 1.1901204
## ETHNICITYWHITE 0.8515325 0.7433614 0.9745490
## MARITAL_STATUSSINGLE 1.4986755 1.3478250 1.6669219
## MARITAL_STATUSUNKNOWN 0.6725828 0.5458913 0.8276527
## MARITAL_STATUSWIDOWED 1.1127017 1.0104664 1.2252801
gm_b <- glm(NQF ~ SOFA,
data = temp,
family = binomial(link = "logit"))
## View
model_info(gm_b)
## Waiting for profiling to be done...
## $`Model Summary`
##
## Call:
## glm(formula = NQF ~ SOFA, family = binomial(link = "logit"),
## data = temp)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.276 -1.274 1.083 1.084 1.087
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.2272396 0.0285897 7.948 1.89e-15 ***
## SOFA -0.0007278 0.0063508 -0.115 0.909
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 15330 on 11158 degrees of freedom
## Residual deviance: 15330 on 11157 degrees of freedom
## AIC: 15334
##
## Number of Fisher Scoring iterations: 3
##
##
## $`OR Summary`
## OR 2.5 % 97.5 %
## (Intercept) 1.2551306 1.1867773 1.327527
## SOFA 0.9992725 0.9869175 1.011798
gm_2a <- glm(NQF ~ GENDER +
AGE +
ETHNICITY +
MARITAL_STATUS +
FIRST_CAREUNIT,
data = temp,
family = binomial(link = "logit"))
## View
model_info(gm_2a)
## Waiting for profiling to be done...
## $`Model Summary`
##
## Call:
## glm(formula = NQF ~ GENDER + AGE + ETHNICITY + MARITAL_STATUS +
## FIRST_CAREUNIT, family = binomial(link = "logit"), data = temp)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.7632 -1.2201 0.8364 1.0560 1.7705
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.36198 0.09268 3.906 9.39e-05 ***
## GENDERM -0.29432 0.04240 -6.941 3.90e-12 ***
## AGE 0.25485 0.02026 12.577 < 2e-16 ***
## ETHNICITYOTHER 0.07258 0.09240 0.785 0.43217
## ETHNICITYWHITE -0.13145 0.06950 -1.891 0.05858 .
## MARITAL_STATUSSINGLE 0.39380 0.05451 7.225 5.01e-13 ***
## MARITAL_STATUSUNKNOWN -0.29419 0.10809 -2.722 0.00649 **
## MARITAL_STATUSWIDOWED 0.15119 0.04982 3.034 0.00241 **
## FIRST_CAREUNITCSRU -1.08414 0.12547 -8.641 < 2e-16 ***
## FIRST_CAREUNITMICU 0.13121 0.05697 2.303 0.02127 *
## FIRST_CAREUNITSICU -0.36046 0.07956 -4.531 5.88e-06 ***
## FIRST_CAREUNITTSICU -0.28322 0.09372 -3.022 0.00251 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 15330 on 11158 degrees of freedom
## Residual deviance: 14804 on 11147 degrees of freedom
## AIC: 14828
##
## Number of Fisher Scoring iterations: 4
##
##
## $`OR Summary`
## OR 2.5 % 97.5 %
## (Intercept) 1.4361771 1.1980433 1.7229352
## GENDERM 0.7450390 0.6856020 0.8095905
## AGE 1.2902742 1.2401053 1.3426299
## ETHNICITYOTHER 1.0752767 0.8970970 1.2887367
## ETHNICITYWHITE 0.8768212 0.7647764 1.0043495
## MARITAL_STATUSSINGLE 1.4826089 1.3325846 1.6500387
## MARITAL_STATUSUNKNOWN 0.7451372 0.6024862 0.9205865
## MARITAL_STATUSWIDOWED 1.1632135 1.0550064 1.2825623
## FIRST_CAREUNITCSRU 0.3381933 0.2636760 0.4313495
## FIRST_CAREUNITMICU 1.1402075 1.0195868 1.2747458
## FIRST_CAREUNITSICU 0.6973538 0.5965642 0.8149293
## FIRST_CAREUNITTSICU 0.7533563 0.6268141 0.9051693
gm_2b <- glm(NQF ~ SOFA +
FIRST_CAREUNIT,
data = temp,
family = binomial(link = "logit"))
## View
model_info(gm_2b)
## Waiting for profiling to be done...
## $`Model Summary`
##
## Call:
## glm(formula = NQF ~ SOFA + FIRST_CAREUNIT, family = binomial(link = "logit"),
## data = temp)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.338 -1.323 1.027 1.036 1.579
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.273813 0.054575 5.017 5.24e-07 ***
## SOFA -0.003461 0.006437 -0.538 0.591
## FIRST_CAREUNITCSRU -1.159929 0.123438 -9.397 < 2e-16 ***
## FIRST_CAREUNITMICU 0.090876 0.055739 1.630 0.103
## FIRST_CAREUNITSICU -0.346994 0.077920 -4.453 8.46e-06 ***
## FIRST_CAREUNITTSICU -0.357304 0.091384 -3.910 9.23e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 15330 on 11158 degrees of freedom
## Residual deviance: 15143 on 11153 degrees of freedom
## AIC: 15155
##
## Number of Fisher Scoring iterations: 4
##
##
## $`OR Summary`
## OR 2.5 % 97.5 %
## (Intercept) 1.3149687 1.1817847 1.4637357
## SOFA 0.9965445 0.9840557 1.0092069
## FIRST_CAREUNITCSRU 0.3135084 0.2453797 0.3982432
## FIRST_CAREUNITMICU 1.0951331 0.9816228 1.2213724
## FIRST_CAREUNITSICU 0.7068092 0.6066043 0.8233301
## FIRST_CAREUNITTSICU 0.6995600 0.5847083 0.8366620
gm_icu <- glm(NQF ~
GENDER +
AGE +
ETHNICITY +
MARITAL_STATUS +
SOFA +
FIRST_CAREUNIT,
data = temp,
family = binomial(link = "logit"))
## View
model_info(gm_icu)
## Waiting for profiling to be done...
## $`Model Summary`
##
## Call:
## glm(formula = NQF ~ GENDER + AGE + ETHNICITY + MARITAL_STATUS +
## SOFA + FIRST_CAREUNIT, family = binomial(link = "logit"),
## data = temp)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.7742 -1.2190 0.8377 1.0578 1.7829
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.323744 0.095265 3.398 0.000678 ***
## GENDERM -0.302932 0.042703 -7.094 1.30e-12 ***
## AGE 0.257538 0.020330 12.668 < 2e-16 ***
## ETHNICITYOTHER 0.070420 0.092457 0.762 0.446265
## ETHNICITYWHITE -0.126631 0.069594 -1.820 0.068826 .
## MARITAL_STATUSSINGLE 0.396581 0.054539 7.272 3.55e-13 ***
## MARITAL_STATUSUNKNOWN -0.289381 0.108153 -2.676 0.007458 **
## MARITAL_STATUSWIDOWED 0.153350 0.049840 3.077 0.002092 **
## SOFA 0.011589 0.006655 1.741 0.081608 .
## FIRST_CAREUNITCSRU -1.091802 0.125550 -8.696 < 2e-16 ***
## FIRST_CAREUNITMICU 0.128070 0.057006 2.247 0.024664 *
## FIRST_CAREUNITSICU -0.355293 0.079631 -4.462 8.13e-06 ***
## FIRST_CAREUNITTSICU -0.276171 0.093780 -2.945 0.003231 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 15330 on 11158 degrees of freedom
## Residual deviance: 14801 on 11146 degrees of freedom
## AIC: 14827
##
## Number of Fisher Scoring iterations: 4
##
##
## $`OR Summary`
## OR 2.5 % 97.5 %
## (Intercept) 1.3822940 1.1472187 1.6666589
## GENDERM 0.7386491 0.6793220 0.8031131
## AGE 1.2937404 1.2432756 1.3464150
## ETHNICITYOTHER 1.0729589 0.8950591 1.2861032
## ETHNICITYWHITE 0.8810590 0.7683379 1.0093910
## MARITAL_STATUSSINGLE 1.4867326 1.3362079 1.6547356
## MARITAL_STATUSUNKNOWN 0.7487268 0.6053150 0.9251427
## MARITAL_STATUSWIDOWED 1.1657328 1.0572550 1.2853853
## SOFA 1.0116560 0.9985574 1.0249514
## FIRST_CAREUNITCSRU 0.3356111 0.2616212 0.4281265
## FIRST_CAREUNITMICU 1.1366331 1.0163182 1.2708378
## FIRST_CAREUNITSICU 0.7009680 0.5995752 0.8192669
## FIRST_CAREUNITTSICU 0.7586832 0.6311746 0.9116708
par(mfrow=c(2,3))
rocplot(as.numeric(predict(gm_a, type="response")), temp$NQF, col="blue", main = "Model 1a (Demographics)")
rocplot(as.numeric(predict(gm_b, type="response")), temp$NQF, col="blue", main = "Model 1b (Clinical Variables)")
rocplot(as.numeric(predict(gm_2a, type="response")), temp$NQF, col="blue", main = "Model 2a (Demographics + ICU)")
rocplot(as.numeric(predict(gm_2b, type="response")), temp$NQF, col="blue", main = "Model 2b (Clinical Variables + ICU)")
rocplot(as.numeric(predict(gm_icu, type="response")), temp$NQF, col="blue", main = "Model 3 (All Covariates)")
For this analysis, we will look at the clinical care team. If any care provider documented the patient care preferences during the admission, then the entire care team will get credit for the documentation.
temp <- team_caremeasure_check(attending_check(dat))
temp <- expire_check(temp)
plotDat(temp, "NQF", "GENDER", F, "Gender (Care Team Level)", "Gender", "Frequency")
test <- table(temp$GENDER, temp$NQF)
test
##
## 0 1
## F 1338 4193
## M 1812 3816
chisq.test(test)
##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: test
## X-squared = 87.841, df = 1, p-value < 2.2e-16
pairwiseNominalIndependence(
as.matrix(test),
fisher = F, gtest = F, chisq = T, method = "fdr")
## Comparison p.Chisq p.adj.Chisq
## 1 F : M 7.09e-21 7.09e-21
plotDat(temp, "NQF", "ADMISSION_TYPE", F, "Admission Type (Care Team Level)", "Admission Type", "Frequency")
test <- table(temp$ADMISSION_TYPE, temp$NQF)
test
##
## 0 1
## ELECTIVE 275 47
## EMERGENCY 2875 7962
chisq.test(test)
##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: test
## X-squared = 532.1, df = 1, p-value < 2.2e-16
pairwiseNominalIndependence(
as.matrix(test),
fisher = F, gtest = F, chisq = T, method = "fdr")
## Comparison p.Chisq p.adj.Chisq
## 1 ELECTIVE : EMERGENCY 9.89e-118 9.89e-118
plotDat(temp, "NQF","ETHNICITY", F, "Ethnicity (Care Team Level)", "Ethnicity", "Frequency")
test <- table(temp$ETHNICITY, temp$NQF)
test
##
## 0 1
## BLACK 250 774
## OTHER 255 856
## WHITE 2645 6379
chisq.test(test)
##
## Pearson's Chi-squared test
##
## data: test
## X-squared = 27.836, df = 2, p-value = 9.024e-07
pairwiseNominalIndependence(
as.matrix(test),
fisher = F, gtest = F, chisq = T, method = "fdr")
## Comparison p.Chisq p.adj.Chisq
## 1 BLACK : OTHER 4.57e-01 4.57e-01
## 2 BLACK : WHITE 1.19e-03 1.78e-03
## 3 OTHER : WHITE 1.14e-05 3.42e-05
plotDat(temp, "NQF", "MARITAL_STATUS", F, "Marital Status (Care Team Level)", "Marital Status", "Frequency")
test <- table(temp$MARITAL_STATUS, temp$NQF)
test
##
## 0 1
## MARRIED 1509 3149
## SINGLE 490 1863
## UNKNOWN 175 239
## WIDOWED 976 2758
chisq.test(test)
##
## Pearson's Chi-squared test
##
## data: test
## X-squared = 151.94, df = 3, p-value < 2.2e-16
pairwiseNominalIndependence(
as.matrix(test),
fisher = F, gtest = F, chisq = T, method = "fdr")
## Comparison p.Chisq p.adj.Chisq
## 1 MARRIED : SINGLE 5.20e-24 3.12e-23
## 2 MARRIED : UNKNOWN 5.48e-05 5.48e-05
## 3 MARRIED : WIDOWED 5.10e-10 7.65e-10
## 4 SINGLE : UNKNOWN 8.35e-21 2.50e-20
## 5 SINGLE : WIDOWED 2.72e-06 3.26e-06
## 6 UNKNOWN : WIDOWED 5.29e-12 1.06e-11
plotDat(temp, "NQF", "FIRST_CAREUNIT", F, "First Careunit (Care Team Level)", "First Careunit", "Frequency")
test <- table(temp$FIRST_CAREUNIT, temp$NQF)
test
##
## 0 1
## CCU 406 1195
## CSRU 211 173
## MICU 1749 5594
## SICU 494 645
## TSICU 290 402
chisq.test(test)
##
## Pearson's Chi-squared test
##
## data: test
## X-squared = 405.13, df = 4, p-value < 2.2e-16
pairwiseNominalIndependence(
as.matrix(test),
fisher = F, gtest = F, chisq = T, method = "fdr")
## Comparison p.Chisq p.adj.Chisq
## 1 CCU : CSRU 4.60e-29 1.53e-28
## 2 CCU : MICU 2.03e-01 2.26e-01
## 3 CCU : SICU 6.67e-23 1.33e-22
## 4 CCU : TSICU 3.79e-15 6.32e-15
## 5 CSRU : MICU 3.64e-42 1.82e-41
## 6 CSRU : SICU 1.06e-04 1.32e-04
## 7 CSRU : TSICU 5.24e-05 7.49e-05
## 8 MICU : SICU 7.77e-44 7.77e-43
## 9 MICU : TSICU 2.29e-25 5.72e-25
## 10 SICU : TSICU 5.72e-01 5.72e-01
boxplot(temp$AGE ~ temp$NQF,
main = "Caremeasure Implementation\n by Age (Care Team Level)",
xlab = "Implementation (1 == Yes)",
ylab = "Age (Years)")
t.test(temp$AGE ~ temp$NQF)
##
## Welch Two Sample t-test
##
## data: temp$AGE by temp$NQF
## t = -11.503, df = 6320.5, p-value < 2.2e-16
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -1.4056701 -0.9963286
## sample estimates:
## mean in group 0 mean in group 1
## 83.23003 84.43103
boxplot(temp$SOFA ~ temp$NQF,
main = "Caremeasure Implementation\n by SOFA Score (Care Team Level)",
xlab = "Implementation (1 == Yes)",
ylab = "SOFA Score")
t.test(temp$SOFA ~ temp$NQF)
##
## Welch Two Sample t-test
##
## data: temp$SOFA by temp$NQF
## t = -4.5397, df = 5777.5, p-value = 5.748e-06
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -0.4092173 -0.1623825
## sample estimates:
## mean in group 0 mean in group 1
## 4.828571 5.114371
Continuous data, Age and SOFA, will be standardized in the form:
\[z = \frac{x - \mu}{\sigma}\]
temp$AGE <- (temp$AGE - mean(temp$AGE))/sd(temp$AGE)
temp$SOFA <- (temp$SOFA - mean(temp$SOFA)/sd(temp$SOFA))
gm_initial <- glmer(NQF ~ (1 | CGID),
data = temp,
family = binomial,
control = glmerControl(optimizer = "bobyqa"),
nAGQ = 10)
## View
summary(gm_initial)
## Generalized linear mixed model fit by maximum likelihood (Adaptive
## Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
## Family: binomial ( logit )
## Formula: NQF ~ (1 | CGID)
## Data: temp
## Control: glmerControl(optimizer = "bobyqa")
##
## AIC BIC logLik deviance df.resid
## 11907.8 11922.4 -5951.9 11903.8 11157
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.6384 -0.7102 0.4009 0.5699 2.2754
##
## Random effects:
## Groups Name Variance Std.Dev.
## CGID (Intercept) 2.011 1.418
## Number of obs: 11159, groups: CGID, 493
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.05472 0.07638 13.81 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
sjstats::icc(gm_initial)
##
## Generalized linear mixed model
## Family: binomial (logit)
## Formula: NQF ~ (1 | CGID)
##
## ICC (CGID): 0.379369
gm_a <- glmer(NQF ~ GENDER +
AGE +
ETHNICITY +
MARITAL_STATUS +
(1 | CGID),
data = temp,
family = binomial,
control = glmerControl(optimizer = "bobyqa"),
nAGQ = 10)
## View
summary(gm_a)
## Generalized linear mixed model fit by maximum likelihood (Adaptive
## Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
## Family: binomial ( logit )
## Formula: NQF ~ GENDER + AGE + ETHNICITY + MARITAL_STATUS + (1 | CGID)
## Data: temp
## Control: glmerControl(optimizer = "bobyqa")
##
## AIC BIC logLik deviance df.resid
## 11663.7 11729.6 -5822.9 11645.7 11150
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -6.3725 -0.6454 0.3750 0.5813 2.7621
##
## Random effects:
## Groups Name Variance Std.Dev.
## CGID (Intercept) 1.937 1.392
## Number of obs: 11159, groups: CGID, 493
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.06543 0.12256 8.693 < 2e-16 ***
## GENDERM -0.36430 0.05429 -6.710 1.95e-11 ***
## AGE 0.25305 0.02607 9.706 < 2e-16 ***
## ETHNICITYOTHER 0.48539 0.12337 3.934 8.34e-05 ***
## ETHNICITYWHITE 0.01746 0.09065 0.193 0.84727
## MARITAL_STATUSSINGLE 0.50199 0.07124 7.047 1.83e-12 ***
## MARITAL_STATUSUNKNOWN -0.38472 0.12847 -2.995 0.00275 **
## MARITAL_STATUSWIDOWED 0.14126 0.06340 2.228 0.02588 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) GENDER AGE ETHNICITYO ETHNICITYW MARITAL_STATUSS
## GENDERM -0.275
## AGE 0.108 -0.061
## ETHNICITYOT -0.507 -0.050 -0.030
## ETHNICITYWH -0.676 -0.078 -0.067 0.665
## MARITAL_STATUSS -0.289 0.226 0.034 0.093 0.087
## MARITAL_STATUSU -0.133 0.068 -0.078 -0.081 0.053 0.150
## MARITAL_STATUSW -0.373 0.361 -0.176 0.102 0.095 0.390
## MARITAL_STATUSU
## GENDERM
## AGE
## ETHNICITYOT
## ETHNICITYWH
## MARITAL_STATUSS
## MARITAL_STATUSU
## MARITAL_STATUSW 0.202
sjt.glmer(gm_a)
| NQF | ||||
| Odds Ratio | CI | p | ||
| Fixed Parts | ||||
| (Intercept) | 2.90 | 2.28 – 3.69 | <.001 | |
| GENDER (M) | 0.69 | 0.62 – 0.77 | <.001 | |
| AGE | 1.29 | 1.22 – 1.36 | <.001 | |
| ETHNICITY (OTHER) | 1.62 | 1.28 – 2.07 | <.001 | |
| ETHNICITY (WHITE) | 1.02 | 0.85 – 1.22 | .847 | |
| MARITAL_STATUS (SINGLE) | 1.65 | 1.44 – 1.90 | <.001 | |
| MARITAL_STATUS (UNKNOWN) | 0.68 | 0.53 – 0.88 | .003 | |
| MARITAL_STATUS (WIDOWED) | 1.15 | 1.02 – 1.30 | .026 | |
| Random Parts | ||||
| τ00, CGID | 1.937 | |||
| NCGID | 493 | |||
| ICCCGID | 0.371 | |||
| Observations | 11159 | |||
| Deviance | 10536.307 | |||
gm_b <- glmer(NQF ~ SOFA +
(1 | CGID),
data = temp,
family = binomial,
control = glmerControl(optimizer = "bobyqa"),
nAGQ = 10)
## View
summary(gm_b)
## Generalized linear mixed model fit by maximum likelihood (Adaptive
## Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
## Family: binomial ( logit )
## Formula: NQF ~ SOFA + (1 | CGID)
## Data: temp
## Control: glmerControl(optimizer = "bobyqa")
##
## AIC BIC logLik deviance df.resid
## 11898.9 11920.9 -5946.5 11892.9 11156
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.4794 -0.7124 0.3938 0.5783 2.4456
##
## Random effects:
## Groups Name Variance Std.Dev.
## CGID (Intercept) 2.006 1.416
## Number of obs: 11159, groups: CGID, 493
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.964591 0.080949 11.92 < 2e-16 ***
## SOFA 0.027701 0.008446 3.28 0.00104 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## SOFA -0.334
sjt.glmer(gm_b)
| NQF | ||||
| Odds Ratio | CI | p | ||
| Fixed Parts | ||||
| (Intercept) | 2.62 | 2.24 – 3.07 | <.001 | |
| SOFA | 1.03 | 1.01 – 1.05 | .001 | |
| Random Parts | ||||
| τ00, CGID | 2.006 | |||
| NCGID | 493 | |||
| ICCCGID | 0.379 | |||
| Observations | 11159 | |||
| Deviance | 10761.908 | |||
gm_2a <- glmer(NQF ~ GENDER +
AGE +
ETHNICITY +
MARITAL_STATUS +
FIRST_CAREUNIT +
(1 | CGID),
data = temp,
family = binomial,
control = glmerControl(optimizer = "bobyqa"),
nAGQ = 10)
## View
summary(gm_2a)
## Generalized linear mixed model fit by maximum likelihood (Adaptive
## Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
## Family: binomial ( logit )
## Formula:
## NQF ~ GENDER + AGE + ETHNICITY + MARITAL_STATUS + FIRST_CAREUNIT +
## (1 | CGID)
## Data: temp
## Control: glmerControl(optimizer = "bobyqa")
##
## AIC BIC logLik deviance df.resid
## 11620.3 11715.5 -5797.1 11594.3 11146
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -6.3780 -0.6342 0.3752 0.5795 3.1486
##
## Random effects:
## Groups Name Variance Std.Dev.
## CGID (Intercept) 1.8 1.342
## Number of obs: 11159, groups: CGID, 493
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.97904 0.14253 6.869 6.47e-12 ***
## GENDERM -0.34575 0.05451 -6.343 2.25e-10 ***
## AGE 0.25758 0.02615 9.849 < 2e-16 ***
## ETHNICITYOTHER 0.53361 0.12383 4.309 1.64e-05 ***
## ETHNICITYWHITE 0.03406 0.09052 0.376 0.70675
## MARITAL_STATUSSINGLE 0.51937 0.07140 7.275 3.48e-13 ***
## MARITAL_STATUSUNKNOWN -0.40106 0.12901 -3.109 0.00188 **
## MARITAL_STATUSWIDOWED 0.17425 0.06371 2.735 0.00624 **
## FIRST_CAREUNITCSRU -0.76667 0.17108 -4.481 7.42e-06 ***
## FIRST_CAREUNITMICU 0.13239 0.08679 1.525 0.12718
## FIRST_CAREUNITSICU -0.21185 0.12372 -1.712 0.08683 .
## FIRST_CAREUNITTSICU 0.32709 0.14833 2.205 0.02745 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) GENDER AGE ETHNICITYO ETHNICITYW MARITAL_STATUSS
## GENDERM -0.254
## AGE 0.071 -0.058
## ETHNICITYOT -0.466 -0.044 -0.020
## ETHNICITYWH -0.603 -0.076 -0.066 0.661
## MARITAL_STATUSS -0.262 0.226 0.036 0.099 0.091
## MARITAL_STATUSU -0.086 0.067 -0.078 -0.081 0.048 0.147
## MARITAL_STATUSW -0.315 0.361 -0.176 0.105 0.098 0.390
## FIRST_CAREUNITC -0.231 -0.043 -0.003 0.002 0.016 -0.010
## FIRST_CAREUNITM -0.500 0.037 0.056 0.068 0.033 0.023
## FIRST_CAREUNITS -0.391 0.034 -0.006 0.021 0.038 0.007
## FIRST_CAREUNITT -0.364 0.009 -0.001 0.052 0.068 0.032
## MARITAL_STATUSU MARITAL_STATUSW FIRST_CAREUNITC
## GENDERM
## AGE
## ETHNICITYOT
## ETHNICITYWH
## MARITAL_STATUSS
## MARITAL_STATUSU
## MARITAL_STATUSW 0.200
## FIRST_CAREUNITC -0.037 -0.061
## FIRST_CAREUNITM -0.035 -0.009 0.351
## FIRST_CAREUNITS -0.047 -0.025 0.350
## FIRST_CAREUNITT -0.083 0.020 0.320
## FIRST_CAREUNITM FIRST_CAREUNITS
## GENDERM
## AGE
## ETHNICITYOT
## ETHNICITYWH
## MARITAL_STATUSS
## MARITAL_STATUSU
## MARITAL_STATUSW
## FIRST_CAREUNITC
## FIRST_CAREUNITM
## FIRST_CAREUNITS 0.559
## FIRST_CAREUNITT 0.467 0.592
sjt.glmer(gm_2a)
| NQF | ||||
| Odds Ratio | CI | p | ||
| Fixed Parts | ||||
| (Intercept) | 2.66 | 2.01 – 3.52 | <.001 | |
| GENDER (M) | 0.71 | 0.64 – 0.79 | <.001 | |
| AGE | 1.29 | 1.23 – 1.36 | <.001 | |
| ETHNICITY (OTHER) | 1.71 | 1.34 – 2.17 | <.001 | |
| ETHNICITY (WHITE) | 1.03 | 0.87 – 1.24 | .707 | |
| MARITAL_STATUS (SINGLE) | 1.68 | 1.46 – 1.93 | <.001 | |
| MARITAL_STATUS (UNKNOWN) | 0.67 | 0.52 – 0.86 | .002 | |
| MARITAL_STATUS (WIDOWED) | 1.19 | 1.05 – 1.35 | .006 | |
| FIRST_CAREUNIT (CSRU) | 0.46 | 0.33 – 0.65 | <.001 | |
| FIRST_CAREUNIT (MICU) | 1.14 | 0.96 – 1.35 | .127 | |
| FIRST_CAREUNIT (SICU) | 0.81 | 0.63 – 1.03 | .087 | |
| FIRST_CAREUNIT (TSICU) | 1.39 | 1.04 – 1.85 | .027 | |
| Random Parts | ||||
| τ00, CGID | 1.800 | |||
| NCGID | 493 | |||
| ICCCGID | 0.354 | |||
| Observations | 11159 | |||
| Deviance | 10513.905 | |||
gm_2b <- glmer(NQF ~ SOFA +
FIRST_CAREUNIT +
(1 | CGID),
data = temp,
family = binomial,
control = glmerControl(optimizer = "bobyqa"),
nAGQ = 10)
## View
summary(gm_2b)
## Generalized linear mixed model fit by maximum likelihood (Adaptive
## Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
## Family: binomial ( logit )
## Formula: NQF ~ SOFA + FIRST_CAREUNIT + (1 | CGID)
## Data: temp
## Control: glmerControl(optimizer = "bobyqa")
##
## AIC BIC logLik deviance df.resid
## 11861.9 11913.1 -5923.9 11847.9 11152
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.5217 -0.6889 0.3914 0.5803 2.7499
##
## Random effects:
## Groups Name Variance Std.Dev.
## CGID (Intercept) 1.913 1.383
## Number of obs: 11159, groups: CGID, 493
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.947851 0.105952 8.946 < 2e-16 ***
## SOFA 0.028912 0.008497 3.403 0.000667 ***
## FIRST_CAREUNITCSRU -0.800718 0.168384 -4.755 1.98e-06 ***
## FIRST_CAREUNITMICU 0.068061 0.085470 0.796 0.425849
## FIRST_CAREUNITSICU -0.136505 0.122717 -1.112 0.265985
## FIRST_CAREUNITTSICU 0.308234 0.145968 2.112 0.034716 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) SOFA FIRST_CAREUNITC FIRST_CAREUNITM
## SOFA -0.266
## FIRST_CAREUNITC -0.321 -0.017
## FIRST_CAREUNITM -0.623 -0.006 0.358
## FIRST_CAREUNITS -0.501 0.060 0.352 0.557
## FIRST_CAREUNITT -0.438 0.051 0.324 0.467
## FIRST_CAREUNITS
## SOFA
## FIRST_CAREUNITC
## FIRST_CAREUNITM
## FIRST_CAREUNITS
## FIRST_CAREUNITT 0.596
sjt.glmer(gm_2b)
| NQF | ||||
| Odds Ratio | CI | p | ||
| Fixed Parts | ||||
| (Intercept) | 2.58 | 2.10 – 3.18 | <.001 | |
| SOFA | 1.03 | 1.01 – 1.05 | <.001 | |
| FIRST_CAREUNIT (CSRU) | 0.45 | 0.32 – 0.62 | <.001 | |
| FIRST_CAREUNIT (MICU) | 1.07 | 0.91 – 1.27 | .426 | |
| FIRST_CAREUNIT (SICU) | 0.87 | 0.69 – 1.11 | .266 | |
| FIRST_CAREUNIT (TSICU) | 1.36 | 1.02 – 1.81 | .035 | |
| Random Parts | ||||
| τ00, CGID | 1.913 | |||
| NCGID | 493 | |||
| ICCCGID | 0.368 | |||
| Observations | 11159 | |||
| Deviance | 10736.210 | |||
gm_icu <- glmer(NQF ~
GENDER +
AGE +
ETHNICITY +
MARITAL_STATUS +
SOFA +
FIRST_CAREUNIT +
(1 | CGID),
data = temp,
family = binomial,
control = glmerControl(optimizer = "bobyqa"),
nAGQ = 10)
## View
summary(gm_icu)
## Generalized linear mixed model fit by maximum likelihood (Adaptive
## Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
## Family: binomial ( logit )
## Formula:
## NQF ~ GENDER + AGE + ETHNICITY + MARITAL_STATUS + SOFA + FIRST_CAREUNIT +
## (1 | CGID)
## Data: temp
## Control: glmerControl(optimizer = "bobyqa")
##
## AIC BIC logLik deviance df.resid
## 11595.5 11698.0 -5783.7 11567.5 11145
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -6.4330 -0.6264 0.3722 0.5762 3.4517
##
## Random effects:
## Groups Name Variance Std.Dev.
## CGID (Intercept) 1.812 1.346
## Number of obs: 11159, groups: CGID, 493
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.825834 0.145765 5.666 1.47e-08 ***
## GENDERM -0.374802 0.054877 -6.830 8.50e-12 ***
## AGE 0.265987 0.026273 10.124 < 2e-16 ***
## ETHNICITYOTHER 0.531797 0.124263 4.280 1.87e-05 ***
## ETHNICITYWHITE 0.048033 0.090726 0.529 0.59651
## MARITAL_STATUSSINGLE 0.531716 0.071543 7.432 1.07e-13 ***
## MARITAL_STATUSUNKNOWN -0.377833 0.129057 -2.928 0.00342 **
## MARITAL_STATUSWIDOWED 0.184446 0.063817 2.890 0.00385 **
## SOFA 0.045084 0.008769 5.141 2.73e-07 ***
## FIRST_CAREUNITCSRU -0.782858 0.171328 -4.569 4.89e-06 ***
## FIRST_CAREUNITMICU 0.128886 0.086892 1.483 0.13800
## FIRST_CAREUNITSICU -0.179204 0.124222 -1.443 0.14913
## FIRST_CAREUNITTSICU 0.361549 0.148659 2.432 0.01501 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation matrix not shown by default, as p = 13 > 12.
## Use print(x, correlation=TRUE) or
## vcov(x) if you need it
sjt.glmer(gm_icu)
| NQF | ||||
| Odds Ratio | CI | p | ||
| Fixed Parts | ||||
| (Intercept) | 2.28 | 1.72 – 3.04 | <.001 | |
| GENDER (M) | 0.69 | 0.62 – 0.77 | <.001 | |
| AGE | 1.30 | 1.24 – 1.37 | <.001 | |
| ETHNICITY (OTHER) | 1.70 | 1.33 – 2.17 | <.001 | |
| ETHNICITY (WHITE) | 1.05 | 0.88 – 1.25 | .597 | |
| MARITAL_STATUS (SINGLE) | 1.70 | 1.48 – 1.96 | <.001 | |
| MARITAL_STATUS (UNKNOWN) | 0.69 | 0.53 – 0.88 | .003 | |
| MARITAL_STATUS (WIDOWED) | 1.20 | 1.06 – 1.36 | .004 | |
| SOFA | 1.05 | 1.03 – 1.06 | <.001 | |
| FIRST_CAREUNIT (CSRU) | 0.46 | 0.33 – 0.64 | <.001 | |
| FIRST_CAREUNIT (MICU) | 1.14 | 0.96 – 1.35 | .138 | |
| FIRST_CAREUNIT (SICU) | 0.84 | 0.66 – 1.07 | .149 | |
| FIRST_CAREUNIT (TSICU) | 1.44 | 1.07 – 1.92 | .015 | |
| Random Parts | ||||
| τ00, CGID | 1.812 | |||
| NCGID | 493 | |||
| ICCCGID | 0.355 | |||
| Observations | 11159 | |||
| Deviance | 10485.611 | |||
par(mfrow=c(2,3))
rocplot(as.numeric(predict(gm_initial, type="response")), temp$NQF, col="blue", main = "Initial Model")
rocplot(as.numeric(predict(gm_a, type="response")), temp$NQF, col="blue", main = "Model 1a (Demographics)")
rocplot(as.numeric(predict(gm_b, type="response")), temp$NQF, col="blue", main = "Model 1b (Clinical Variables)")
rocplot(as.numeric(predict(gm_2a, type="response")), temp$NQF, col="blue", main = "Model 2a (Demographics + ICU)")
rocplot(as.numeric(predict(gm_2b, type="response")), temp$NQF, col="blue", main = "Model 2b (Clinical Variables + ICU)")
rocplot(as.numeric(predict(gm_icu, type="response")), temp$NQF, col="blue", main = "Model 3 (All Covariates)")
Not meaningful.
gm_a <- glm(NQF ~ GENDER +
AGE +
ETHNICITY +
MARITAL_STATUS,
data = temp,
family = binomial(link = "logit"))
## View
model_info(gm_a)
## Waiting for profiling to be done...
## $`Model Summary`
##
## Call:
## glm(formula = NQF ~ GENDER + AGE + ETHNICITY + MARITAL_STATUS,
## family = binomial(link = "logit"), data = temp)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.0524 -1.3558 0.7131 0.8514 1.3900
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.14273 0.08652 13.208 < 2e-16 ***
## GENDERM -0.32540 0.04645 -7.006 2.46e-12 ***
## AGE 0.25791 0.02235 11.542 < 2e-16 ***
## ETHNICITYOTHER 0.27144 0.10569 2.568 0.0102 *
## ETHNICITYWHITE -0.17328 0.07843 -2.209 0.0272 *
## MARITAL_STATUSSINGLE 0.52593 0.06214 8.463 < 2e-16 ***
## MARITAL_STATUSUNKNOWN -0.68581 0.10836 -6.329 2.47e-10 ***
## MARITAL_STATUSWIDOWED 0.05847 0.05395 1.084 0.2785
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 13281 on 11158 degrees of freedom
## Residual deviance: 12914 on 11151 degrees of freedom
## AIC: 12930
##
## Number of Fisher Scoring iterations: 4
##
##
## $`OR Summary`
## OR 2.5 % 97.5 %
## (Intercept) 3.1353184 2.6496285 3.7197749
## GENDERM 0.7222378 0.6593390 0.7910238
## AGE 1.2942246 1.2388622 1.3522800
## ETHNICITYOTHER 1.3118546 1.0664176 1.6140752
## ETHNICITYWHITE 0.8409056 0.7199808 0.9792482
## MARITAL_STATUSSINGLE 1.6920322 1.4988336 1.9123353
## MARITAL_STATUSUNKNOWN 0.5036813 0.4075010 0.6233374
## MARITAL_STATUSWIDOWED 1.0602117 0.9538886 1.1785708
gm_b <- glm(NQF ~ SOFA,
data = temp,
family = binomial(link = "logit"))
## View
model_info(gm_b)
## Waiting for profiling to be done...
## $`Model Summary`
##
## Call:
## glm(formula = NQF ~ SOFA, family = binomial(link = "logit"),
## data = temp)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.7808 -1.5470 0.8029 0.8254 0.8717
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.826230 0.031302 26.395 < 2e-16 ***
## SOFA 0.032490 0.007181 4.525 6.05e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 13281 on 11158 degrees of freedom
## Residual deviance: 13261 on 11157 degrees of freedom
## AIC: 13265
##
## Number of Fisher Scoring iterations: 4
##
##
## $`OR Summary`
## OR 2.5 % 97.5 %
## (Intercept) 2.284689 2.149072 2.429646
## SOFA 1.033024 1.018642 1.047725
gm_2a <- glm(NQF ~ GENDER +
AGE +
ETHNICITY +
MARITAL_STATUS +
FIRST_CAREUNIT,
data = temp,
family = binomial(link = "logit"))
## View
model_info(gm_2a)
## Waiting for profiling to be done...
## $`Model Summary`
##
## Call:
## glm(formula = NQF ~ GENDER + AGE + ETHNICITY + MARITAL_STATUS +
## FIRST_CAREUNIT, family = binomial(link = "logit"), data = temp)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.1567 -1.1270 0.6710 0.8155 1.5558
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.17362 0.10565 11.108 < 2e-16 ***
## GENDERM -0.31896 0.04750 -6.716 1.87e-11 ***
## AGE 0.26826 0.02277 11.779 < 2e-16 ***
## ETHNICITYOTHER 0.41483 0.10855 3.822 0.000133 ***
## ETHNICITYWHITE -0.13291 0.07975 -1.667 0.095600 .
## MARITAL_STATUSSINGLE 0.50508 0.06306 8.010 1.15e-15 ***
## MARITAL_STATUSUNKNOWN -0.50087 0.11204 -4.470 7.81e-06 ***
## MARITAL_STATUSWIDOWED 0.11838 0.05528 2.141 0.032235 *
## FIRST_CAREUNITCSRU -1.18681 0.12034 -9.863 < 2e-16 ***
## FIRST_CAREUNITMICU 0.13976 0.06504 2.149 0.031650 *
## FIRST_CAREUNITSICU -0.84086 0.08503 -9.889 < 2e-16 ***
## FIRST_CAREUNITTSICU -0.66598 0.09906 -6.723 1.78e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 13281 on 11158 degrees of freedom
## Residual deviance: 12545 on 11147 degrees of freedom
## AIC: 12569
##
## Number of Fisher Scoring iterations: 4
##
##
## $`OR Summary`
## OR 2.5 % 97.5 %
## (Intercept) 3.2336696 2.6320563 3.9828340
## GENDERM 0.7269040 0.6622378 0.7977673
## AGE 1.3076887 1.2507043 1.3675004
## ETHNICITYOTHER 1.5141176 1.2240603 1.8735154
## ETHNICITYWHITE 0.8755423 0.7477099 1.0222409
## MARITAL_STATUSSINGLE 1.6571136 1.4652722 1.8762073
## MARITAL_STATUSUNKNOWN 0.6060026 0.4868062 0.7554407
## MARITAL_STATUSWIDOWED 1.1256749 1.0101810 1.2546353
## FIRST_CAREUNITCSRU 0.3051929 0.2408667 0.3861290
## FIRST_CAREUNITMICU 1.1500031 1.0115715 1.3054206
## FIRST_CAREUNITSICU 0.4313385 0.3649816 0.5093978
## FIRST_CAREUNITTSICU 0.5137711 0.4231037 0.6239222
gm_2b <- glm(NQF ~ SOFA +
FIRST_CAREUNIT,
data = temp,
family = binomial(link = "logit"))
## View
model_info(gm_2b)
## Waiting for profiling to be done...
## $`Model Summary`
##
## Call:
## glm(formula = NQF ~ SOFA + FIRST_CAREUNIT, family = binomial(link = "logit"),
## data = temp)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.8415 -1.2707 0.7305 0.7639 1.3155
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.997956 0.061864 16.131 < 2e-16 ***
## SOFA 0.025667 0.007342 3.496 0.000473 ***
## FIRST_CAREUNITCSRU -1.299375 0.117808 -11.030 < 2e-16 ***
## FIRST_CAREUNITMICU 0.076010 0.063712 1.193 0.232859
## FIRST_CAREUNITSICU -0.800513 0.083025 -9.642 < 2e-16 ***
## FIRST_CAREUNITTSICU -0.738141 0.096241 -7.670 1.72e-14 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 13281 on 11158 degrees of freedom
## Residual deviance: 12892 on 11153 degrees of freedom
## AIC: 12904
##
## Number of Fisher Scoring iterations: 4
##
##
## $`OR Summary`
## OR 2.5 % 97.5 %
## (Intercept) 2.7127318 2.4051091 3.0653430
## SOFA 1.0259989 1.0113926 1.0409289
## FIRST_CAREUNITCSRU 0.2727022 0.2162857 0.3433094
## FIRST_CAREUNITMICU 1.0789731 0.9515195 1.2215327
## FIRST_CAREUNITSICU 0.4490984 0.3815144 0.5282975
## FIRST_CAREUNITTSICU 0.4780018 0.3958215 0.5772770
gm_icu <- glm(NQF ~
GENDER +
AGE +
ETHNICITY +
MARITAL_STATUS +
SOFA +
FIRST_CAREUNIT,
data = temp,
family = binomial(link = "logit"))
## View
model_info(gm_icu)
## Waiting for profiling to be done...
## $`Model Summary`
##
## Call:
## glm(formula = NQF ~ GENDER + AGE + ETHNICITY + MARITAL_STATUS +
## SOFA + FIRST_CAREUNIT, family = binomial(link = "logit"),
## data = temp)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.2506 -1.1409 0.6676 0.8177 1.5492
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.037661 0.108429 9.570 < 2e-16 ***
## GENDERM -0.348734 0.047881 -7.283 3.26e-13 ***
## AGE 0.278813 0.022917 12.166 < 2e-16 ***
## ETHNICITYOTHER 0.411646 0.108942 3.779 0.000158 ***
## ETHNICITYWHITE -0.118724 0.079960 -1.485 0.137599
## MARITAL_STATUSSINGLE 0.517150 0.063203 8.182 2.78e-16 ***
## MARITAL_STATUSUNKNOWN -0.482863 0.112254 -4.302 1.70e-05 ***
## MARITAL_STATUSWIDOWED 0.127479 0.055333 2.304 0.021231 *
## SOFA 0.042354 0.007607 5.568 2.58e-08 ***
## FIRST_CAREUNITCSRU -1.218994 0.120658 -10.103 < 2e-16 ***
## FIRST_CAREUNITMICU 0.129467 0.065185 1.986 0.047016 *
## FIRST_CAREUNITSICU -0.824054 0.085209 -9.671 < 2e-16 ***
## FIRST_CAREUNITTSICU -0.644605 0.099133 -6.502 7.90e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 13281 on 11158 degrees of freedom
## Residual deviance: 12514 on 11146 degrees of freedom
## AIC: 12540
##
## Number of Fisher Scoring iterations: 4
##
##
## $`OR Summary`
## OR 2.5 % 97.5 %
## (Intercept) 2.8226063 2.2847918 3.4951818
## GENDERM 0.7055805 0.6423184 0.7749432
## AGE 1.3215604 1.2636228 1.3824017
## ETHNICITYOTHER 1.5093001 1.2192309 1.8690038
## ETHNICITYWHITE 0.8880526 0.7580906 1.0372807
## MARITAL_STATUSSINGLE 1.6772412 1.4826441 1.8995524
## MARITAL_STATUSUNKNOWN 0.6170144 0.4954494 0.7694974
## MARITAL_STATUSWIDOWED 1.1359613 1.0193093 1.2662317
## SOFA 1.0432640 1.0278801 1.0589950
## FIRST_CAREUNITCSRU 0.2955274 0.2330900 0.3741354
## FIRST_CAREUNITMICU 1.1382211 1.0009305 1.2924047
## FIRST_CAREUNITSICU 0.4386497 0.3710425 0.5182143
## FIRST_CAREUNITTSICU 0.5248698 0.4321841 0.6374953
par(mfrow=c(2,3))
rocplot(as.numeric(predict(gm_a, type="response")), temp$NQF, col="blue", main = "Model 1a (Demographics)")
rocplot(as.numeric(predict(gm_b, type="response")), temp$NQF, col="blue", main = "Model 1b (Clinical Variables)")
rocplot(as.numeric(predict(gm_2a, type="response")), temp$NQF, col="blue", main = "Model 2a (Demographics + ICU)")
rocplot(as.numeric(predict(gm_2b, type="response")), temp$NQF, col="blue", main = "Model 2b (Clinical Variables + ICU)")
rocplot(as.numeric(predict(gm_icu, type="response")), temp$NQF, col="blue", main = "Model 3 (All Covariates)")
Vincent JL, de Mendonça A, Cantraine F, Moreno R, Takala J, Suter PM, Sprung CL, Colardyn F, Blecher S. Use of the SOFA score to assess the incidence of organ dysfunction/failure in intensive care units: results of a multicenter, prospective study. Working group on “sepsis-related problems” of the European Society of Intensive Care Medicine. Crit Care Med 1998 Nov;26(11):1793-800. PMID 9824069.
Moreno R, Vincent JL, Matos R, Mendonça A, Cantraine F, Thijs L, Takala J, Sprung C, Antonelli M, Bruining H, Willatts S. The use of maximum SOFA score to quantify organ dysfunction/failure in intensive care. Results of a prospective, multicentre study. Working Group on Sepsis related Problems of the ESICM. Intensive Care Med 1999 Jul;25(7):686-96. PMID 10470572.
Benjamini, Y., and Hochberg, Y. (1995). Controlling the false discovery rate: a practical and powerful approach to multiple testing. Journal of the Royal Statistical Society Series B 57, 289–300.