SG_SIRI_analysis
Stephen Murata
10/24/2021
1 Overall summary statistics
options(width=190)
compareGroups::descrTable(~. , data = SG_df_new, hide.no = '0')##
## --------Summary descriptives table ---------
##
## __________________________
## [ALL] N
## N=79
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯
## Pt_Group: 79
## HC 32 (40.5%)
## Patient 47 (59.5%)
## Exp_Group: 47
## ESC_PBO 21 (44.7%)
## ESC_CBX 26 (55.3%)
## SII_BL 489 (230) 79
## SIRI_BL 1.03 (0.60) 79
## SII_WK8 462 (219) 47
## SIRI_WK8 1.07 (0.61) 47
## Ethnicity 1.52 (0.86) 46
## BMI 31.6 (5.78) 45
## Age 42.0 (12.7) 46
## Sex: 46
## 1 16 (34.8%)
## 2 30 (65.2%)
## HAMD17_BL 22.5 (6.34) 46
## HAMD17_WK8 10.0 (6.32) 43
## Remission 0.42 (0.50) 43
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯2 Summary of sample by patient status
createTable(compareGroups(Pt_Group ~ ., data = SG_df_new, method = NA), hide.no = '0', show.p.mul= T)##
## --------Summary descriptives table by 'Pt_Group'---------
##
## _______________________________________________________
## HC Patient p.overall
## N=32 N=47
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯
## Exp_Group: .
## ESC_PBO 0 (.%) 21 (44.7%)
## ESC_CBX 0 (.%) 26 (55.3%)
## SII_BL 491 [305;611] 477 [320;640] 0.834
## SIRI_BL 0.76 [0.58;1.21] 0.98 [0.64;1.39] 0.201
## SII_WK8 . [.;.] 408 [309;579] .
## SIRI_WK8 . [.;.] 1.02 [0.62;1.38] .
## Ethnicity: .
## 1 0 (.%) 30 (65.2%)
## 2 0 (.%) 10 (21.7%)
## 3 0 (.%) 5 (10.9%)
## 5 0 (.%) 1 (2.17%)
## BMI . (.) 31.6 (5.78) .
## Age . (.) 42.0 (12.7) .
## Sex: .
## 1 0 (.%) 16 (34.8%)
## 2 0 (.%) 30 (65.2%)
## HAMD17_BL . (.) 22.5 (6.34) .
## HAMD17_WK8 . [.;.] 9.00 [6.00;14.0] .
## Remission 0 (.%) 18 (41.9%) .
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯3 Summary of sample by remission
createTable(compareGroups(Remission ~ ., data = SG_df_new, method = NA), hide.no = '0', show.p.mul= T)##
## --------Summary descriptives table by 'Remission'---------
##
## _____________________________________________________________
## 0 1 p.overall
## N=25 N=18
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯
## Pt_Group: Patient 25 (100%) 18 (100%) .
## Exp_Group: <0.001
## ESC_PBO 17 (68.0%) 1 (5.56%)
## ESC_CBX 8 (32.0%) 17 (94.4%)
## SII_BL 494 [336;646] 436 [296;612] 0.649
## SIRI_BL 1.01 [0.78;1.38] 0.90 [0.48;1.17] 0.445
## SII_WK8 451 [307;639] 407 [302;539] 0.453
## SIRI_WK8 1.12 [0.69;1.69] 0.94 [0.52;1.23] 0.143
## Ethnicity: 0.057
## 1 15 (62.5%) 12 (66.7%)
## 2 3 (12.5%) 6 (33.3%)
## 3 5 (20.8%) 0 (0.00%)
## 5 1 (4.17%) 0 (0.00%)
## BMI 32.0 (4.85) 30.7 (7.00) 0.503
## Age 44.0 [34.8;53.8] 37.5 [31.8;43.2] 0.186
## Sex: 1.000
## 1 9 (37.5%) 7 (38.9%)
## 2 15 (62.5%) 11 (61.1%)
## HAMD17_BL 21.8 (5.86) 22.8 (6.00) 0.613
## HAMD17_WK8 13.0 [10.0;16.0] 5.50 [2.25;7.00] <0.001
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯4 Summary of sample by treatment
##
## --------Summary descriptives table by 'Exp_Group'---------
##
## _____________________________________________________________
## ESC_PBO ESC_CBX p.overall
## N=21 N=26
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯
## Pt_Group: Patient 21 (100%) 26 (100%) .
## SII_BL 523 [357;759] 379 [296;618] 0.168
## SIRI_BL 1.06 [0.78;1.53] 0.94 [0.59;1.25] 0.294
## SII_WK8 399 [320;639] 415 [297;518] 0.386
## SIRI_WK8 1.04 [0.73;1.69] 0.94 [0.57;1.22] 0.118
## Ethnicity: 0.658
## 1 12 (60.0%) 18 (69.2%)
## 2 4 (20.0%) 6 (23.1%)
## 3 3 (15.0%) 2 (7.69%)
## 5 1 (5.00%) 0 (0.00%)
## BMI 32.2 (4.89) 31.1 (6.40) 0.515
## Age 46.7 (13.0) 38.3 (11.3) 0.028
## Sex: 0.776
## 1 6 (30.0%) 10 (38.5%)
## 2 14 (70.0%) 16 (61.5%)
## HAMD17_BL 21.6 (7.31) 23.1 (5.56) 0.472
## HAMD17_WK8 11.5 [9.00;16.0] 7.00 [5.00;11.0] 0.001
## Remission 1 (5.56%) 17 (68.0%) <0.001
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯5 Visual inspection of outcome variables: histograms
p<-ggplot(SG_df_new, aes(x=HAMD17_BL)) +
geom_histogram(aes(y=..density..), colour="black", fill="white")+
geom_density(alpha=.2, fill="#FF6666")
p+geom_vline(aes(xintercept=mean(HAMD17_BL)), color="blue", linetype="dashed", size=1)## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
p<-ggplot(SG_df_new, aes(x=HAMD17_WK8)) +
geom_histogram(aes(y=..density..), colour="black", fill="white")+
geom_density(alpha=.2, fill="#FF6666")
p+geom_vline(aes(xintercept=mean(HAMD17_WK8)), color="blue", linetype="dashed", size=1)## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
6 Visual inspection of biomarker distributions: histograms
# Histogram with density plotz
p<-ggplot(SG_df_new, aes(x=SII_BL)) +
geom_histogram(aes(y=..density..), colour="black", fill="white")+
geom_density(alpha=.2, fill="#FF6666")
p+geom_vline(aes(xintercept=mean(SII_BL)), color="blue", linetype="dashed", size=1)## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
p<-ggplot(SG_df_new, aes(x=SII_WK8)) +
geom_histogram(aes(y=..density..), colour="black", fill="white")+
geom_density(alpha=.2, fill="#FF6666")
p+geom_vline(aes(xintercept=mean(SII_WK8)), color="blue", linetype="dashed", size=1)## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
p<-ggplot(SG_df_new, aes(x=SIRI_BL)) +
geom_histogram(aes(y=..density..), colour="black", fill="white")+
geom_density(alpha=.2, fill="#FF6666")
p+geom_vline(aes(xintercept=mean(SIRI_BL)), color="blue", linetype="dashed", size=1)## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
p<-ggplot(SG_df_new, aes(x=SIRI_WK8)) +
geom_histogram(aes(y=..density..), colour="black", fill="white")+
geom_density(alpha=.2, fill="#FF6666")
p+geom_vline(aes(xintercept=mean(SIRI_WK8)), color="blue", linetype="dashed", size=1)## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
7 Variable transformations
sg_bx_vars<-SG_df_new %>% dplyr::select(contains("SII"), contains("SIRI")) %>% names()
SG_df_new[paste(sg_bx_vars,sep="_", "log")] <- log(SG_df_new[sg_bx_vars])
SG_df_new[paste(sg_bx_vars,sep="_", "inverse")] <- 1/(SG_df_new[sg_bx_vars])
SG_df_new[paste(sg_bx_vars,sep="_", "sqrt")] <- sqrt(SG_df_new[sg_bx_vars])8 Visual inspection of biomarker residuals (modelled by hypothesis)
Rcoef_list <- list()
beta_list<-list()
se_list<-list()
t_list<-list()
p_list<-list()
var_list<-list()
sg_bx_trans_vars<-SG_df_new %>% dplyr::select(contains("SII"), contains("SIRI")) %>% names() %>% sort()
for (x in sg_bx_trans_vars) {
LM1 <- lm(substitute(HAMD17_WK8 ~ Sex+Age+BMI+Ethnicity+Exp_Group+HAMD17_BL+i, list(i = as.name(x))), data = SG_df_new)
Rcoef_list[[x]]<- summary(LM1)$r.squared[1]
beta_list[[x]]<-summary(LM1)$coefficients[8,1]
se_list[[x]]<-summary(LM1)$coefficients[8,2]
t_list[[x]]<-summary(LM1)$coefficients[8,3]
p_list[[x]]<-summary(LM1)$coefficients[8,4]
var_list[[x]]<-x
plot(LM1, which=2, main=x)
}
SG_bx_vars_selected<-c("SII_BL", "SII_WK8", "SIRI_BL", "SIRI_WK8_inverse") Note: Generally raw variables are preferable here, except for SIRI_WK8 possibly
9 Tabulating biomarker coefficients (models above)
do.call(rbind, Map(data.frame,
Var=var_list,
Beta=beta_list,
t_test=t_list,
SE=se_list,
p_value=p_list,
R_coeff=Rcoef_list)) %>%
knitr::kable(digits=2) | Var | Beta | t_test | SE | p_value | R_coeff | |
|---|---|---|---|---|---|---|
| SII_BL | SII_BL | 0.00 | 0.09 | 0.00 | 0.93 | 0.35 |
| SII_BL_inverse | SII_BL_inverse | -353.93 | -0.48 | 735.46 | 0.63 | 0.35 |
| SII_BL_log | SII_BL_log | 0.51 | 0.25 | 2.01 | 0.80 | 0.35 |
| SII_BL_sqrt | SII_BL_sqrt | 0.03 | 0.16 | 0.19 | 0.87 | 0.35 |
| SII_WK8 | SII_WK8 | 0.01 | 1.54 | 0.00 | 0.13 | 0.39 |
| SII_WK8_inverse | SII_WK8_inverse | -503.11 | -0.70 | 714.63 | 0.49 | 0.35 |
| SII_WK8_log | SII_WK8_log | 2.06 | 1.07 | 1.92 | 0.29 | 0.37 |
| SII_WK8_sqrt | SII_WK8_sqrt | 0.24 | 1.30 | 0.18 | 0.20 | 0.38 |
| SIRI_BL | SIRI_BL | 0.00 | 0.00 | 1.56 | 1.00 | 0.35 |
| SIRI_BL_inverse | SIRI_BL_inverse | -1.78 | -1.53 | 1.16 | 0.13 | 0.39 |
| SIRI_BL_log | SIRI_BL_log | 1.48 | 0.87 | 1.69 | 0.39 | 0.36 |
| SIRI_BL_sqrt | SIRI_BL_sqrt | 1.50 | 0.44 | 3.45 | 0.67 | 0.35 |
| SIRI_WK8 | SIRI_WK8 | 2.10 | 1.38 | 1.52 | 0.18 | 0.38 |
| SIRI_WK8_inverse | SIRI_WK8_inverse | -1.43 | -1.44 | 0.99 | 0.16 | 0.38 |
| SIRI_WK8_log | SIRI_WK8_log | 2.18 | 1.42 | 1.53 | 0.16 | 0.38 |
| SIRI_WK8_sqrt | SIRI_WK8_sqrt | 4.54 | 1.41 | 3.21 | 0.17 | 0.38 |
10 Multivariate linear model of HAMD17_WK8 by SIRI_BL
SG_model<-lm(HAMD17_WK8~Sex+Age+BMI+Ethnicity+Exp_Group+HAMD17_BL+SIRI_BL, data=SG_df_new)
car::Anova(SG_model, type = "II")## Anova Table (Type II tests)
##
## Response: HAMD17_WK8
## Sum Sq Df F value Pr(>F)
## Sex 47.00 1 1.4531 0.2363
## Age 72.57 1 2.2436 0.1434
## BMI 5.02 1 0.1552 0.6961
## Ethnicity 40.13 1 1.2408 0.2731
## Exp_Group 205.39 1 6.3500 0.0166 *
## HAMD17_BL 5.43 1 0.1680 0.6845
## SIRI_BL 0.00 1 0.0000 0.9989
## Residuals 1099.71 34
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1summary(SG_model)##
## Call:
## lm(formula = HAMD17_WK8 ~ Sex + Age + BMI + Ethnicity + Exp_Group +
## HAMD17_BL + SIRI_BL, data = SG_df_new)
##
## Residuals:
## Min 1Q Median 3Q Max
## -9.0195 -3.9494 -0.8909 2.9092 11.5150
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 8.19247 8.15691 1.004 0.3223
## Sex2 -2.42556 2.01214 -1.205 0.2363
## Age 0.11565 0.07721 1.498 0.1434
## BMI -0.06772 0.17189 -0.394 0.6961
## Ethnicity 1.27559 1.14515 1.114 0.2731
## Exp_GroupESC_CBX -5.08426 2.01763 -2.520 0.0166 *
## HAMD17_BL 0.07209 0.17589 0.410 0.6845
## SIRI_BL 0.00223 1.56425 0.001 0.9989
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 5.687 on 34 degrees of freedom
## (37 observations deleted due to missingness)
## Multiple R-squared: 0.345, Adjusted R-squared: 0.2102
## F-statistic: 2.558 on 7 and 34 DF, p-value: 0.03136Note: Here we are asking the question whether HAMD17_WK8 is predicted by baseline factors, independently of treatment arm. Notice that we are controlling for baseline HAMD17_BL, in addition to demographics. As expected, ESC+CBX arm significantly associates with lower HAMD17_WK8. There are no other significant factors, however.
11 Multivariate linear model (same as above) including interactions
SG_model<-lm(HAMD17_WK8~
Sex
+Age
+BMI
+Ethnicity
+Exp_Group
+HAMD17_BL*SIRI_BL
+SIRI_BL
, data=SG_df_new, na.action=na.exclude)
car::Anova(SG_model, type = "II")## Anova Table (Type II tests)
##
## Response: HAMD17_WK8
## Sum Sq Df F value Pr(>F)
## Sex 123.71 1 5.5026 0.0251482 *
## Age 17.17 1 0.7638 0.3884604
## BMI 1.73 1 0.0769 0.7833386
## Ethnicity 83.26 1 3.7036 0.0629548 .
## Exp_Group 137.58 1 6.1195 0.0186870 *
## HAMD17_BL 5.43 1 0.2417 0.6262430
## SIRI_BL 0.00 1 0.0000 0.9986461
## HAMD17_BL:SIRI_BL 357.80 1 15.9152 0.0003467 ***
## Residuals 741.91 33
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1summary(SG_model)##
## Call:
## lm(formula = HAMD17_WK8 ~ Sex + Age + BMI + Ethnicity + Exp_Group +
## HAMD17_BL * SIRI_BL + SIRI_BL, data = SG_df_new, na.action = na.exclude)
##
## Residuals:
## Min 1Q Median 3Q Max
## -9.9081 -2.8551 0.0042 2.7696 8.3298
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 31.83697 9.02081 3.529 0.001251 **
## Sex2 -4.04929 1.72622 -2.346 0.025148 *
## Age 0.05767 0.06599 0.874 0.388460
## BMI -0.03978 0.14348 -0.277 0.783339
## Ethnicity 1.85875 0.96586 1.924 0.062955 .
## Exp_GroupESC_CBX -4.19738 1.69675 -2.474 0.018687 *
## HAMD17_BL -1.10587 0.32968 -3.354 0.002011 **
## SIRI_BL -22.35393 5.75366 -3.885 0.000466 ***
## HAMD17_BL:SIRI_BL 1.22482 0.30702 3.989 0.000347 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4.742 on 33 degrees of freedom
## (37 observations deleted due to missingness)
## Multiple R-squared: 0.5581, Adjusted R-squared: 0.451
## F-statistic: 5.21 on 8 and 33 DF, p-value: 0.0002949reg_cohend(SG_model)## (Intercept) Sex2 Age BMI Ethnicity Exp_GroupESC_CBX HAMD17_BL SIRI_BL HAMD17_BL:SIRI_BL
## 0.545 -0.362 0.135 -0.043 0.297 -0.382 -0.518 -0.599 0.616 interact_plot(SG_model, pred = SIRI_BL, modx = HAMD17_BL, jitter=0.1, plot.points = TRUE, main.title = "The effect of SIRI_BL on HAMD17_WK8 depends on HAMD17_BL")
Note: In this model, we are screening for “interaction effects” amongst the independent variables. When evaluating for interactions between SIRI_BL and HAMD17_BL, something interesting happens: (1) elevated baseline SIRI is significantly associated with lower HAMD17_WK8 (main effect) and (2) elevated SIRI_BL is associated with elevated HAMD17_WK8 in patients with elevated HAMD17_BL (interaction effect). As expected, treatment arm (CBX) is an independent predictor of lower HAMD17_WK8, but its effect does not depend on SIRI_BL. Sex1 is independently associated with higher HAMD17_WK8 as well.
12 Multivariate logistic model of remission by SIRI_BL
SG_model_logit<-glm(Remission~
Sex
+Age
+BMI
+Ethnicity
+Exp_Group
+HAMD17_BL*SIRI_BL
+SIRI_BL
, family="binomial", na.action=na.exclude, data=SG_df_new)
car::Anova(SG_model_logit, type = "II")## Analysis of Deviance Table (Type II tests)
##
## Response: Remission
## LR Chisq Df Pr(>Chisq)
## Sex 3.8000 1 0.0512526 .
## Age 0.0449 1 0.8321389
## BMI 0.0013 1 0.9715982
## Ethnicity 4.0392 1 0.0444545 *
## Exp_Group 20.0354 1 7.602e-06 ***
## HAMD17_BL 0.3465 1 0.5561241
## SIRI_BL 1.1485 1 0.2838712
## HAMD17_BL:SIRI_BL 12.0966 1 0.0005051 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1summary(SG_model_logit)##
## Call:
## glm(formula = Remission ~ Sex + Age + BMI + Ethnicity + Exp_Group +
## HAMD17_BL * SIRI_BL + SIRI_BL, family = "binomial", data = SG_df_new,
## na.action = na.exclude)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.40675 -0.23258 -0.01554 0.43021 1.40361
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -25.711800 10.383014 -2.476 0.0133 *
## Sex2 3.165269 1.888376 1.676 0.0937 .
## Age 0.012026 0.057343 0.210 0.8339
## BMI 0.003513 0.098741 0.036 0.9716
## Ethnicity -2.114256 1.193857 -1.771 0.0766 .
## Exp_GroupESC_CBX 6.706381 2.905740 2.308 0.0210 *
## HAMD17_BL 1.049533 0.414098 2.535 0.0113 *
## SIRI_BL 21.400254 8.391651 2.550 0.0108 *
## HAMD17_BL:SIRI_BL -1.131213 0.454905 -2.487 0.0129 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 57.364 on 41 degrees of freedom
## Residual deviance: 24.757 on 33 degrees of freedom
## (37 observations deleted due to missingness)
## AIC: 42.757
##
## Number of Fisher Scoring iterations: 7reg_cohend(SG_model_logit)## (Intercept) Sex2 Age BMI Ethnicity Exp_GroupESC_CBX HAMD17_BL SIRI_BL HAMD17_BL:SIRI_BL
## -0.382 0.259 0.032 0.005 -0.273 0.356 0.391 0.394 -0.384Note: Here we are reproducing the results from the above interaction model, except this time we are using a dichotomous outcome (remission vs non-remission). The estimate signs are reversed because we’re talking about remission, not HAMD17 score.This confirms the prior model. Seeing that there is something special about SIRI_BL, We can now start to ask the question whether changes in SIRI are related to treatment response. That would involve another analysis (“mediation analysis”, more involved) looking at whether CBX “causes” changes in SIRI that are linked to treatment response. Out of scope here, but for that we would need to establish (1) treatment arm is associated with SIRI at point X, where X temporally precedes outcome, and that (2) SIRI at point X is associated with treatment response.