, Woojoo Lee
Department of Public Health Science, Graduate School of Public Health, Seoul National University, Seoul, Korea
Copyright © 2022 The Korean Society for Preventive Medicine
This is an Open Access article distributed under the terms of the Creative Commons Attribution Non-Commercial License (https://creativecommons.org/licenses/by-nc/4.0/) which permits unrestricted non-commercial use, distribution, and reproduction in any medium, provided the original work is properly cited.
# Load libraries
library (readxl)
library (dplyr)
library (boot)
library (stdReg)
# Import liver synthetic data
liver_syn_data1 <- read.csv (‘~/liver_syn_data.csv’)
# [1] Build function for standardization
standardization <- function( data, indices ) {
liver_syn_data0 <- data[indices, ]
# [1]-1. data expansion
# original
liver_syn_data0$data.class <- ‘ori’
# 2nd copy
liver_syn_data2 <- liver_syn_data0 %>%
mutate (data.class=‘T0’, Treatment=0, death=NA)
# 3rd copy
liver_syn_data3 <- liver_syn_data0 %>%
mutate (data.class=‘T1’, Treatment=1, death=NA)
# Combine all data
onesample <- rbind ( liver_syn_data0, liver_syn_
data2, liver_syn_data3)
# [1]-2. Fit the model
fit1 <-glm (death ~ Treatment * i_bclc+age+Liver_
Cancer_Cause+MELD+cpc_cat+platelet_cat+
Sodium_level+AFP_level+Ascites_status, family=
‘binomial’,data = onesample)
# Confounders in the dataset (i_bclc: BCLC stage; age: Age;
Liver_Cancer_Cause: cause of liver cancer; MELD: MELD
score; cpc_cat: child pugh classification; platelet_cat:
Platelet count; Sodium_level: sodium level, AFP_level:
alpha-fetoprotein; Ascites_status: ascites status)
# [1]-3. Predict the outcome Y
onesample$predicted.meanY 1 <- predict.glm (fit1,
onesample, type=“response”)
Y1T1=mean( onesample$predicted.meanY1
[onesample$data.class==‘T1’])
Y1T0=mean( onesample$predicted.meanY1
[onesample$data.class==‘T0’])
# [1]-4. Calculate the causal estimands
ATE=Y1T1 - Y1T0 #Risk difference
RR=Y1T1/Y1T0 #Relative ratio
OR=(Y1T1/(1-Y1T1))/(Y1T0/(1-Y1T0)) #Odds ratio
return(c(ATE, RR, OR))}
# [2] Generate confidence intervals
# [2]-1. Calculate the 95% confidence interval
set.seed(1234)
results <- boot(data=liver_syn_data1, statistic=
standardization, R=1,000, parallel=“multicore”)
se <- c(sd(results$t[, 1]), sd(results$t[, 2]),
sd(results$t[, 3]))
mean <- results$t0
# 95% normal confidence interval using se
ll1 <- mean - qnorm(0.975) * se
ul1 <- mean + qnorm(0.975) * se
# 95% percentile confidence interval
ll2 <- c (quantile (results$t[,1], 0.025), quantile
(results$t[,2], 0.025), quantile (results$t[,3], 0.025))
ul2 <- c (quantile (results$t[,1], 0.975), quantile
(results$t[,2], 0.975), quantile (results$t[,3], 0.975))
# [2]-2. Present the result
bootstrap <-data.frame(cbind(c(“ATE”, “RR”, “OR”), round
(mean, 4), round (se, 4), round (ll1, 4), round (ul1, 4),
round (ll2, 4), round (ul2, 4)), row.names=NULL)
colnames (bootstrap) <- c(“Estimand”, “mean”, “se”,
“Lower1”, “Upper1”, “Lower2”, “Upper2”)
# [3] Standardization using the stdReg package
fit1 <-glm (death ~ Treatment * i_bclc+age+Liver_
Cancer_Cause+MELD+cpc_cat+platelet_cat+
Sodium_level+AFP_level+Ascites_status, family=
‘binomial’, data=liver_syn_data1)
fit.std <- stdGlm (fit=fit1, data=liver_syn_data1,
X=“Treatment”, x= seq (0,1,1))
summary (fit.std, contrast=‘difference’, reference=0)
#Risk difference
summary (fit.std, contrast=‘ratio’, reference=0)
#Relative ratio
summary (fit.std, transform=‘odds’, contrast=“ratio”,
reference=0) #Odds ratio
CONFLICT OF INTEREST
The authors have no conflicts of interest associated with the material presented in this paper.
AUTHOR CONTRIBUTIONS
Both authors contributed equally to conceiving the study, analyzing the data, and writing this paper.
FUNDING
This study was supported by the National Research Foundation of Korea (NRF) grant funded by the Korea government (MSIT) (No.2021R1A2C1014409).
| Estimand | Model 11 | Model 22 | Model 33 |
|---|---|---|---|
| Risk difference | 0.02 (−0.01, 0.04) | 0.02 (−0.01, 0.04) | 0.02 (−0.01, 0.04) |
| Relative risk | 1.09 (0.94, 1.25) | 1.10 (0.94, 1.25) | 1.10 (0.94, 1.25) |
| Odds ratio | 1.12 (0.92, 1.32) | 1.12 (0.92, 1.32) | 1.12 (0.92, 1.32) |
Values are presented as estimate (95% confidence interval).
BCLC, Barcelona Clinic Liver Cancer.
1 Interaction term between treatment, BCLC stage, and alpha-fetoprotein level and other confounders are included in the model.
2 Interaction term between treatment, BCLC stage, and cause of liver cancer and other confounders are included in the model.
3 Interaction term between treatment and Child-Pugh classification and other confounders are included in the model.
| Name | Measures of causation | Measures of association |
|---|---|---|
| Risk difference | Pr[Y1=1]–Pr[Y0=1] | Pr[Y=1|T=1]–Pr[Y=1|T=0] |
| Risk ratio | ||
| Odds ratio |
| Name | Causal estimand | Standardization method |
|---|---|---|
| Risk difference | Pr[Y1=1]–Pr[Y0=1] | ∑c Pr[Y=1|T=1,C=c] Pr[C=c]– ∑c Pr[Y=1|T=0,C=c] Pr[C=c] |
| Risk ratio | ||
| Odds ratio |
| Name | Causal estimand for subgroups | Standardization for subgroups |
|---|---|---|
| Risk difference | Pr[Y1=1|S=s]–Pr[Y0=1|S=s] | ∑c Pr[Y=1|T=1,C=c, S=s] Pr[C=c|S=s]–∑c Pr[Y=1|T=0,C=c, S=s] Pr[C=c|S=s] |
| Risk ratio | ||
| Odds ratio |
| Estimand | Standardization function | stdGlm function |
|---|---|---|
| Risk difference | 0.02 (−0.01, 0.05) | 0.02 (−0.01, 0.05) |
| Relative risk | 1.11 (0.94, 1.27) | 1.11 (0.95, 1.27) |
| Odds ratio | 1.13 (0.93, 1.34) | 1.13 (0.93, 1.33) |
| Estimand | Model 1 |
Model 2 |
Model 3 |
|---|---|---|---|
| Risk difference | 0.02 (−0.01, 0.04) | 0.02 (−0.01, 0.04) | 0.02 (−0.01, 0.04) |
| Relative risk | 1.09 (0.94, 1.25) | 1.10 (0.94, 1.25) | 1.10 (0.94, 1.25) |
| Odds ratio | 1.12 (0.92, 1.32) | 1.12 (0.92, 1.32) | 1.12 (0.92, 1.32) |
Values are presented as estimate (95% confidence interval).
Values are presented as estimate (95% confidence interval). BCLC, Barcelona Clinic Liver Cancer. Interaction term between treatment, BCLC stage, and alpha-fetoprotein level and other confounders are included in the model. Interaction term between treatment, BCLC stage, and cause of liver cancer and other confounders are included in the model. Interaction term between treatment and Child-Pugh classification and other confounders are included in the model.