varix main analysis
市田 親正
2023/4/20
必要なlibraryの読み込み
library(tidyverse)
library(readxl)
library(tableone)
library(lubridate)
library(skimr)
library(svglite)
library(haven)
library(summarytools)## Warning in system2("/usr/bin/otool", c("-L", shQuote(DSO)), stdout = TRUE):
## 命令 ''/usr/bin/otool' -L
## '/Library/Frameworks/R.framework/Resources/library/tcltk/libs//tcltk.so''
## の実行は状態 1 を持ちましたlibrary(naniar)
library(devtools)
library(reader)
library(stringr)
library(missForest)
library(mice)
library(geepack)
library(ggplot2)
library(cobalt)
library(tableone)
library(sandwich)
library(survey)
library(WeightIt)
library(MatchIt)
library(twang)
library(naniar)
pacman::p_load(norm2, miceadds, lmtest, car, ROCR, pROC, knitr)データ読み込み
data <- read_excel("varix_icu_ok.xlsx")データの確認
str(data)## tibble [790 × 50] (S3: tbl_df/tbl/data.frame)
## $ hosp_id : num [1:790] 1001 1001 1001 1001 1001 ...
## $ pt_id : num [1:790] 1 2 4 5 6 9 12 14 16 17 ...
## $ ad_num : num [1:790] 1 1 1 1 1 1 1 1 1 1 ...
## $ year : num [1:790] 2012 2011 2011 2010 2017 ...
## $ antibio_exposure : num [1:790] 0 0 0 0 0 0 0 0 1 0 ...
## $ antibio_type_a : chr [1:790] "99" "99" "99" "99" ...
## $ antibio_type_b : chr [1:790] "99" "99" "99" "99" ...
## $ antibio_term : num [1:790] 0 0 0 0 0 0 0 0 5 0 ...
## $ age_num : num [1:790] 50 80 44 67 47 73 65 41 38 62 ...
## $ age_cate : num [1:790] 0 2 0 1 0 1 1 0 0 0 ...
## $ sex : chr [1:790] "M" "M" "M" "F" ...
## $ bmi_num : num [1:790] NA 25.3 14.5 NA 21 ...
## $ smoking : num [1:790] 0 0 240 0 270 1000 0 210 0 0 ...
## $ brthel_cate : num [1:790] NA 2 1 NA NA 0 1 1 1 2 ...
## $ child_numl : num [1:790] 11 6 8 NA 11 9 6 8 NA 9 ...
## $ child_score : chr [1:790] "c" "a" "b" NA ...
## $ jcs_all : num [1:790] 0 0 0 0 0 0 0 0 0 3 ...
## $ cci_num : num [1:790] 4 4 4 4 4 4 5 5 4 4 ...
## $ malignancy : num [1:790] 0 0 0 0 0 0 0 0 0 0 ...
## $ hd : num [1:790] 0 0 0 0 0 0 0 0 0 0 ...
## $ hcc : num [1:790] 0 0 0 0 0 0 0 0 0 0 ...
## $ alocohol : num [1:790] 1 0 0 0 1 1 0 1 1 0 ...
## $ past_varix_rup : num [1:790] 0 0 0 0 0 0 0 0 0 0 ...
## $ antiplate_user_01: num [1:790] 0 0 0 0 0 0 0 0 0 0 ...
## $ anticoag_user01 : num [1:790] 0 0 0 0 0 0 0 0 0 0 ...
## $ nsaids_user01 : num [1:790] 0 0 0 0 0 0 0 0 0 0 ...
## $ steroid_user01 : num [1:790] 0 0 0 0 0 0 0 0 0 0 ...
## $ ppi_user01 : num [1:790] 1 1 1 1 1 1 1 1 1 1 ...
## $ vaso : num [1:790] 0 0 0 0 0 0 0 0 0 0 ...
## $ map_day0 : num [1:790] 0 2 6 2 0 4 0 6 0 8 ...
## $ shock : num [1:790] 1 0 1 1 0 1 0 1 0 NA ...
## $ t_bil : num [1:790] 2.2 1.2 3.4 1.2 7.5 2.2 0.4 1.5 1.8 0.5 ...
## $ ast : num [1:790] 217 31 129 52 112 55 55 56 169 56 ...
## $ alt : num [1:790] 63 22 46 36 53 20 67 12 36 30 ...
## $ wbc : num [1:790] 7400 5000 9100 3900 9800 11100 2900 8700 9300 14600 ...
## $ hb : num [1:790] 6.9 10.8 10.7 6.3 12.2 5.6 5.7 4.9 10.9 3.7 ...
## $ plt : num [1:790] 115 77 162 63 69 124 61 129 90 437 ...
## $ alb : num [1:790] 2.2 3.2 2.9 2.8 2.6 2.3 4 2.2 3.2 2.1 ...
## $ eGFR : num [1:790] 58 57.7 112.4 123.6 89.3 ...
## $ eGFR30 : num [1:790] 0 0 0 0 0 0 1 0 0 0 ...
## $ crp : num [1:790] 0.92 0.29 0.29 0.29 0.1 NA 0.04 0.44 0.09 0.17 ...
## $ pt : num [1:790] 37.8 55 37.8 74.6 28.5 45.9 71.5 43.6 66.7 51.3 ...
## $ aptt : num [1:790] 29.4 29 35.3 30.1 37.8 27.5 25 24.8 20.8 30.9 ...
## $ outcome : num [1:790] 0 0 0 0 0 0 0 0 0 0 ...
## $ death_day2 : num [1:790] 0 0 0 0 0 0 0 0 0 0 ...
## $ rebleed_end_hemo : num [1:790] 0 0 0 0 0 0 0 0 0 0 ...
## $ sbp : num [1:790] 0 0 0 0 0 0 0 0 0 0 ...
## $ infection : num [1:790] 0 0 0 0 0 1 0 0 0 0 ...
## $ infection_all : num [1:790] 0 0 0 0 0 1 0 0 0 0 ...
## $ los : num [1:790] 12 7 10 3 6 8 5 8 5 2 ...データの整理
df <-
data|> #miceパッケージが順序ありの因子型でないと読めないためこれを指定。
mutate(
hosp_id=as.integer(hosp_id),
year=as.integer(year),
pt_id=as.integer(pt_id),
ad_num=as.integer(ad_num),
antibio_exposure=as.integer(antibio_exposure), #geeglmを回すにはfactorでは回らない
antibio_type_a= factor(antibio_type_a, levels = c("a", "b","c","d","e","f","g","h","i","99")),
antibio_type_b= factor(antibio_type_b, levels = c("c","d","f","h","i","99")),
antibio_term= as.integer(antibio_term ),
age_num=as.integer(age_num),
age_cate=factor(age_cate,levels = c("0", "1", "2", "3")),
sex= factor(sex, levels = c("M", "F")),
smoking= as.integer(smoking),
brthel_cate= factor(brthel_cate, levels = c("0", "1", "2")),
child_numl= as.integer(child_numl),
child_score=factor(child_score, levels = c("a", "b", "c")),
jcs_all=factor(jcs_all, levels = c("0", "1", "2","3","10","20","30","100","200","300")),
cci_num=as.integer(cci_num),
malignancy=factor(malignancy),
hd=factor(hd),
hcc=factor(hcc),
alocohol=factor(alocohol),
past_varix_rup=factor(past_varix_rup),
antiplate_user_01=factor(antiplate_user_01),
anticoag_user01=factor(anticoag_user01),
nsaids_user01=factor(nsaids_user01),
steroid_user01=factor(steroid_user01),
ppi_user01=factor(ppi_user01),
vaso=factor(vaso),
map_day0= as.integer(map_day0),
shock=factor(shock),
eGFR30=factor(eGFR30),
outcome=factor(outcome),
death_day2=factor(death_day2),
rebleed_end_hemo=factor(rebleed_end_hemo),
sbp=factor(sbp),
infection=factor(infection),
infection_all=factor(infection_all),
los=as.integer(los)
)整理後の確認
str(df)## tibble [790 × 50] (S3: tbl_df/tbl/data.frame)
## $ hosp_id : int [1:790] 1001 1001 1001 1001 1001 1001 1001 1001 1001 1001 ...
## $ pt_id : int [1:790] 1 2 4 5 6 9 12 14 16 17 ...
## $ ad_num : int [1:790] 1 1 1 1 1 1 1 1 1 1 ...
## $ year : int [1:790] 2012 2011 2011 2010 2017 2010 2013 2014 2016 2019 ...
## $ antibio_exposure : int [1:790] 0 0 0 0 0 0 0 0 1 0 ...
## $ antibio_type_a : Factor w/ 10 levels "a","b","c","d",..: 10 10 10 10 10 10 10 10 4 10 ...
## $ antibio_type_b : Factor w/ 6 levels "c","d","f","h",..: 6 6 6 6 6 6 6 6 6 6 ...
## $ antibio_term : int [1:790] 0 0 0 0 0 0 0 0 5 0 ...
## $ age_num : int [1:790] 50 80 44 67 47 73 65 41 38 62 ...
## $ age_cate : Factor w/ 4 levels "0","1","2","3": 1 3 1 2 1 2 2 1 1 1 ...
## $ sex : Factor w/ 2 levels "M","F": 1 1 1 2 1 1 1 2 1 2 ...
## $ bmi_num : num [1:790] NA 25.3 14.5 NA 21 ...
## $ smoking : int [1:790] 0 0 240 0 270 1000 0 210 0 0 ...
## $ brthel_cate : Factor w/ 3 levels "0","1","2": NA 3 2 NA NA 1 2 2 2 3 ...
## $ child_numl : int [1:790] 11 6 8 NA 11 9 6 8 NA 9 ...
## $ child_score : Factor w/ 3 levels "a","b","c": 3 1 2 NA 3 2 1 2 NA 2 ...
## $ jcs_all : Factor w/ 10 levels "0","1","2","3",..: 1 1 1 1 1 1 1 1 1 4 ...
## $ cci_num : int [1:790] 4 4 4 4 4 4 5 5 4 4 ...
## $ malignancy : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
## $ hd : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
## $ hcc : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
## $ alocohol : Factor w/ 2 levels "0","1": 2 1 1 1 2 2 1 2 2 1 ...
## $ past_varix_rup : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
## $ antiplate_user_01: Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
## $ anticoag_user01 : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
## $ nsaids_user01 : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
## $ steroid_user01 : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
## $ ppi_user01 : Factor w/ 2 levels "0","1": 2 2 2 2 2 2 2 2 2 2 ...
## $ vaso : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
## $ map_day0 : int [1:790] 0 2 6 2 0 4 0 6 0 8 ...
## $ shock : Factor w/ 2 levels "0","1": 2 1 2 2 1 2 1 2 1 NA ...
## $ t_bil : num [1:790] 2.2 1.2 3.4 1.2 7.5 2.2 0.4 1.5 1.8 0.5 ...
## $ ast : num [1:790] 217 31 129 52 112 55 55 56 169 56 ...
## $ alt : num [1:790] 63 22 46 36 53 20 67 12 36 30 ...
## $ wbc : num [1:790] 7400 5000 9100 3900 9800 11100 2900 8700 9300 14600 ...
## $ hb : num [1:790] 6.9 10.8 10.7 6.3 12.2 5.6 5.7 4.9 10.9 3.7 ...
## $ plt : num [1:790] 115 77 162 63 69 124 61 129 90 437 ...
## $ alb : num [1:790] 2.2 3.2 2.9 2.8 2.6 2.3 4 2.2 3.2 2.1 ...
## $ eGFR : num [1:790] 58 57.7 112.4 123.6 89.3 ...
## $ eGFR30 : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 2 1 1 1 ...
## $ crp : num [1:790] 0.92 0.29 0.29 0.29 0.1 NA 0.04 0.44 0.09 0.17 ...
## $ pt : num [1:790] 37.8 55 37.8 74.6 28.5 45.9 71.5 43.6 66.7 51.3 ...
## $ aptt : num [1:790] 29.4 29 35.3 30.1 37.8 27.5 25 24.8 20.8 30.9 ...
## $ outcome : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
## $ death_day2 : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
## $ rebleed_end_hemo : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
## $ sbp : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
## $ infection : Factor w/ 2 levels "0","1": 1 1 1 1 1 2 1 1 1 1 ...
## $ infection_all : Factor w/ 2 levels "0","1": 1 1 1 1 1 2 1 1 1 1 ...
## $ los : int [1:790] 12 7 10 3 6 8 5 8 5 2 ...確認その2
dfSummary(df)## Data Frame Summary
## df
## Dimensions: 790 x 50
## Duplicates: 0
##
## ------------------------------------------------------------------------------------------------------------------------
## No Variable Stats / Values Freqs (% of Valid) Graph Valid Missing
## ---- ------------------- ----------------------------- --------------------- ---------------------- ---------- ---------
## 1 hosp_id Mean (sd) : 1021.6 (21.9) 44 distinct values : 790 0
## [integer] min < med < max: : (100.0%) (0.0%)
## 1001 < 1017 < 1075 :
## IQR (CV) : 20 (0) : : :
## : : : : : :
##
## 2 pt_id Mean (sd) : 416.6 (243.3) 679 distinct values : . 790 0
## [integer] min < med < max: : : : . : : : : (100.0%) (0.0%)
## 1 < 425.5 < 839 : : : : : : : :
## IQR (CV) : 431.5 (0.6) : : : : : : : : .
## : : : : : : : : :
##
## 3 ad_num Mean (sd) : 1.2 (0.6) 1 : 672 (85.1%) IIIIIIIIIIIIIIIII 790 0
## [integer] min < med < max: 2 : 79 (10.0%) II (100.0%) (0.0%)
## 1 < 1 < 6 3 : 26 ( 3.3%)
## IQR (CV) : 0 (0.5) 4 : 9 ( 1.1%)
## 5 : 2 ( 0.3%)
## 6 : 2 ( 0.3%)
##
## 4 year Mean (sd) : 2015.9 (3.8) 13 distinct values : . : 790 0
## [integer] min < med < max: : : : (100.0%) (0.0%)
## 2010 < 2016 < 2022 : : . : . :
## IQR (CV) : 7 (0) : : : : : : : : : :
## : : : : : : : : : :
##
## 5 antibio_exposure Min : 0 0 : 558 (70.6%) IIIIIIIIIIIIII 790 0
## [integer] Mean : 0.3 1 : 232 (29.4%) IIIII (100.0%) (0.0%)
## Max : 1
##
## 6 antibio_type_a 1. a 4 ( 0.5%) 790 0
## [factor] 2. b 32 ( 4.1%) (100.0%) (0.0%)
## 3. c 51 ( 6.5%) I
## 4. d 104 (13.2%) II
## 5. e 2 ( 0.3%)
## 6. f 12 ( 1.5%)
## 7. g 2 ( 0.3%)
## 8. h 22 ( 2.8%)
## 9. i 3 ( 0.4%)
## 10. 99 558 (70.6%) IIIIIIIIIIIIII
##
## 7 antibio_type_b 1. c 1 ( 0.1%) 790 0
## [factor] 2. d 5 ( 0.6%) (100.0%) (0.0%)
## 3. f 2 ( 0.3%)
## 4. h 4 ( 0.5%)
## 5. i 2 ( 0.3%)
## 6. 99 776 (98.2%) IIIIIIIIIIIIIIIIIII
##
## 8 antibio_term Mean (sd) : 1.3 (3.3) 23 distinct values : 790 0
## [integer] min < med < max: : (100.0%) (0.0%)
## 0 < 0 < 35 :
## IQR (CV) : 1 (2.5) :
## : .
##
## 9 age_num Mean (sd) : 61.1 (13) 65 distinct values : 790 0
## [integer] min < med < max: . : : : (100.0%) (0.0%)
## 24 < 62 < 93 : : : : .
## IQR (CV) : 19 (0.2) : : : : : : :
## . : : : : : : : .
##
## 10 age_cate 1. 0 465 (58.9%) IIIIIIIIIII 790 0
## [factor] 2. 1 190 (24.1%) IIII (100.0%) (0.0%)
## 3. 2 133 (16.8%) III
## 4. 3 2 ( 0.3%)
##
## 11 sex 1. M 598 (75.7%) IIIIIIIIIIIIIII 790 0
## [factor] 2. F 192 (24.3%) IIII (100.0%) (0.0%)
##
## 12 bmi_num Mean (sd) : 23.8 (12.6) 554 distinct values : 706 84
## [numeric] min < med < max: : (89.4%) (10.6%)
## 13.7 < 22.8 < 339.8 :
## IQR (CV) : 5 (0.5) :
## :
##
## 13 smoking Mean (sd) : 257.5 (390) 116 distinct values : 694 96
## [integer] min < med < max: : (87.8%) (12.2%)
## 0 < 0 < 2200 :
## IQR (CV) : 410 (1.5) :
## : : . .
##
## 14 brthel_cate 1. 0 269 (41.2%) IIIIIIII 653 137
## [factor] 2. 1 215 (32.9%) IIIIII (82.7%) (17.3%)
## 3. 2 169 (25.9%) IIIII
##
## 15 child_numl Mean (sd) : 8.3 (2) 11 distinct values : . : 689 101
## [integer] min < med < max: : : : : (87.2%) (12.8%)
## 5 < 8 < 15 : : : : :
## IQR (CV) : 3 (0.2) : : : : : . .
## : : : : : : : .
##
## 16 child_score 1. a 135 (19.0%) III 711 79
## [factor] 2. b 376 (52.9%) IIIIIIIIII (90.0%) (10.0%)
## 3. c 200 (28.1%) IIIII
##
## 17 jcs_all 1. 0 628 (79.5%) IIIIIIIIIIIIIII 790 0
## [factor] 2. 1 93 (11.8%) II (100.0%) (0.0%)
## 3. 2 20 ( 2.5%)
## 4. 3 17 ( 2.2%)
## 5. 10 14 ( 1.8%)
## 6. 20 2 ( 0.3%)
## 7. 30 6 ( 0.8%)
## 8. 100 6 ( 0.8%)
## 9. 200 2 ( 0.3%)
## 10. 300 2 ( 0.3%)
##
## 18 cci_num Mean (sd) : 4.5 (1.3) 11 distinct values : 790 0
## [integer] min < med < max: : (100.0%) (0.0%)
## 3 < 4 < 13 :
## IQR (CV) : 1 (0.3) : :
## : : : .
##
## 19 malignancy 1. 0 696 (88.1%) IIIIIIIIIIIIIIIII 790 0
## [factor] 2. 1 94 (11.9%) II (100.0%) (0.0%)
##
## 20 hd 1. 0 779 (98.6%) IIIIIIIIIIIIIIIIIII 790 0
## [factor] 2. 1 11 ( 1.4%) (100.0%) (0.0%)
##
## 21 hcc 1. 0 640 (81.0%) IIIIIIIIIIIIIIII 790 0
## [factor] 2. 1 150 (19.0%) III (100.0%) (0.0%)
##
## 22 alocohol 1. 0 417 (52.8%) IIIIIIIIII 790 0
## [factor] 2. 1 373 (47.2%) IIIIIIIII (100.0%) (0.0%)
##
## 23 past_varix_rup 1. 0 600 (75.9%) IIIIIIIIIIIIIII 790 0
## [factor] 2. 1 190 (24.1%) IIII (100.0%) (0.0%)
##
## 24 antiplate_user_01 1. 0 779 (98.6%) IIIIIIIIIIIIIIIIIII 790 0
## [factor] 2. 1 11 ( 1.4%) (100.0%) (0.0%)
##
## 25 anticoag_user01 1. 0 777 (98.4%) IIIIIIIIIIIIIIIIIII 790 0
## [factor] 2. 1 13 ( 1.6%) (100.0%) (0.0%)
##
## 26 nsaids_user01 1. 0 772 (97.7%) IIIIIIIIIIIIIIIIIII 790 0
## [factor] 2. 1 18 ( 2.3%) (100.0%) (0.0%)
##
## 27 steroid_user01 1. 0 788 (99.7%) IIIIIIIIIIIIIIIIIII 790 0
## [factor] 2. 1 2 ( 0.3%) (100.0%) (0.0%)
##
## 28 ppi_user01 1. 0 91 (11.5%) II 790 0
## [factor] 2. 1 699 (88.5%) IIIIIIIIIIIIIIIII (100.0%) (0.0%)
##
## 29 vaso 1. 0 764 (96.7%) IIIIIIIIIIIIIIIIIII 790 0
## [factor] 2. 1 26 ( 3.3%) (100.0%) (0.0%)
##
## 30 map_day0 Mean (sd) : 2.9 (2.9) 13 distinct values : 790 0
## [integer] min < med < max: : (100.0%) (0.0%)
## 0 < 4 < 26 : :
## IQR (CV) : 4 (1) : : .
## : : : .
##
## 31 shock 1. 0 475 (62.0%) IIIIIIIIIIII 766 24
## [factor] 2. 1 291 (38.0%) IIIIIII (97.0%) (3.0%)
##
## 32 t_bil Mean (sd) : 2.1 (1.9) 254 distinct values : 765 25
## [numeric] min < med < max: : (96.8%) (3.2%)
## 0.2 < 1.4 < 13.7 :
## IQR (CV) : 1.6 (0.9) : :
## : : : .
##
## 33 ast Mean (sd) : 77.4 (90.1) 195 distinct values : 773 17
## [numeric] min < med < max: : (97.8%) (2.2%)
## 13 < 49 < 984 :
## IQR (CV) : 56 (1.2) :
## : .
##
## 34 alt Mean (sd) : 39.9 (44.5) 124 distinct values : 773 17
## [numeric] min < med < max: : (97.8%) (2.2%)
## 7 < 28 < 562 :
## IQR (CV) : 25 (1.1) :
## : .
##
## 35 wbc Mean (sd) : 8417.8 (4639.6) 345 distinct values : 776 14
## [numeric] min < med < max: : : (98.2%) (1.8%)
## 1400 < 7400 < 58400 : :
## IQR (CV) : 5027.5 (0.6) : :
## : : :
##
## 36 hb Mean (sd) : 8.8 (2.4) 119 distinct values . . : 776 14
## [numeric] min < med < max: : : : (98.2%) (1.8%)
## 2.1 < 8.6 < 16.5 : : : :
## IQR (CV) : 3.3 (0.3) : : : : : .
## . : : : : : : . .
##
## 37 plt Mean (sd) : 116.5 (69.4) 217 distinct values : 776 14
## [numeric] min < med < max: : (98.2%) (1.8%)
## 18 < 102 < 1073 : .
## IQR (CV) : 68 (0.6) : :
## : : .
##
## 38 alb Mean (sd) : 2.9 (0.6) 35 distinct values : 754 36
## [numeric] min < med < max: : : (95.4%) (4.6%)
## 1.1 < 2.9 < 4.9 : : :
## IQR (CV) : 0.8 (0.2) : : : :
## : : : : : .
##
## 39 eGFR Mean (sd) : 71.6 (30.8) 737 distinct values : : 775 15
## [numeric] min < med < max: : : . (98.1%) (1.9%)
## 4.3 < 69.2 < 196 : : :
## IQR (CV) : 40.6 (0.4) : : : : .
## : : : : : : . .
##
## 40 eGFR30 1. 0 716 (92.4%) IIIIIIIIIIIIIIIIII 775 15
## [factor] 2. 1 59 ( 7.6%) I (98.1%) (1.9%)
##
## 41 crp Mean (sd) : 0.8 (1.4) 360 distinct values : 753 37
## [numeric] min < med < max: : (95.3%) (4.7%)
## 0 < 0.3 < 18.6 :
## IQR (CV) : 0.7 (1.9) :
## : .
##
## 42 pt Mean (sd) : 57.1 (16.7) 383 distinct values . . : 744 46
## [numeric] min < med < max: : : : (94.2%) (5.8%)
## 13.4 < 57.9 < 107.9 : : : : :
## IQR (CV) : 23.8 (0.3) . : : : : : .
## . : : : : : : : . .
##
## 43 aptt Mean (sd) : 32.7 (14.3) 234 distinct values : 700 90
## [numeric] min < med < max: : (88.6%) (11.4%)
## 17.6 < 30.2 < 200 :
## IQR (CV) : 7.6 (0.4) :
## : :
##
## 44 outcome 1. 0 724 (91.6%) IIIIIIIIIIIIIIIIII 790 0
## [factor] 2. 1 66 ( 8.4%) I (100.0%) (0.0%)
##
## 45 death_day2 1. 0 742 (93.9%) IIIIIIIIIIIIIIIIII 790 0
## [factor] 2. 1 48 ( 6.1%) I (100.0%) (0.0%)
##
## 46 rebleed_end_hemo 1. 0 775 (98.1%) IIIIIIIIIIIIIIIIIII 790 0
## [factor] 2. 1 15 ( 1.9%) (100.0%) (0.0%)
##
## 47 sbp 1. 0 775 (98.1%) IIIIIIIIIIIIIIIIIII 790 0
## [factor] 2. 1 15 ( 1.9%) (100.0%) (0.0%)
##
## 48 infection 1. 0 753 (95.3%) IIIIIIIIIIIIIIIIIII 790 0
## [factor] 2. 1 37 ( 4.7%) (100.0%) (0.0%)
##
## 49 infection_all 1. 0 744 (94.2%) IIIIIIIIIIIIIIIIII 790 0
## [factor] 2. 1 46 ( 5.8%) I (100.0%) (0.0%)
##
## 50 los Mean (sd) : 13 (15) 58 distinct values : 790 0
## [integer] min < med < max: : (100.0%) (0.0%)
## 2 < 9 < 217 :
## IQR (CV) : 9 (1.2) :
## : .
## ------------------------------------------------------------------------------------------------------------------------確認その3
dfSummary(df) %>% view()## Switching method to 'browser'## Output file written: /var/folders/n9/tf_wmwpn3gl2t7l1cz4tqk000000gn/T//RtmptELf1k/filec4176a972eac.html欠測評価
gg_miss_var(df,show_pct = T)## Warning: The `guide` argument in `scale_*()` cannot be `FALSE`. This was deprecated in
## ggplot2 3.3.4.
## ℹ Please use "none" instead.
## ℹ The deprecated feature was likely used in the naniar package.
## Please report the issue at <]8;;https://github.com/njtierney/naniar/issueshttps://github.com/njtierney/naniar/issues]8;;>.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
欠測評価2
vis_miss(df,cluster = T)## Warning: `gather_()` was deprecated in tidyr 1.2.0.
## ℹ Please use `gather()` instead.
## ℹ The deprecated feature was likely used in the visdat package.
## Please report the issue at <]8;;https://github.com/ropensci/visdat/issueshttps://github.com/ropensci/visdat/issues]8;;>.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
tableoneアウトカムも含めて作成
col_cont <-
c("antibio_term","child_numl","cci_num", "map_day0", "t_bil", "ast", "alt", "wbc", "hb", "plt", "alb", "crp", "pt", "aptt", "los")
col_fact <-
c("antibio_type_a", "antibio_type_b", "age_cate","sex", "brthel_cate", "malignancy", "child_score", "hd", "hcc", "alocohol", "past_varix_rup", "antiplate_user_01", "anticoag_user01", "nsaids_user01" ,"steroid_user01", "ppi_user01", "shock", "eGFR30", "outcome", "death_day2", "rebleed_end_hemo", "sbp", "infection", "infection_all")
df|> select(c(col_cont, col_fact, "antibio_exposure"))|>
CreateTableOne(vars = c(col_cont, col_fact), strata = "antibio_exposure", factorVars = col_fact) -> tableone2## Warning: Using an external vector in selections was deprecated in tidyselect 1.1.0.
## ℹ Please use `all_of()` or `any_of()` instead.
## # Was:
## data %>% select(col_cont)
##
## # Now:
## data %>% select(all_of(col_cont))
##
## See <https://tidyselect.r-lib.org/reference/faq-external-vector.html>.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.## Warning: Using an external vector in selections was deprecated in tidyselect 1.1.0.
## ℹ Please use `all_of()` or `any_of()` instead.
## # Was:
## data %>% select(col_fact)
##
## # Now:
## data %>% select(all_of(col_fact))
##
## See <https://tidyselect.r-lib.org/reference/faq-external-vector.html>.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.print(tableone2, smd = T, missing = T, test = T, explain = T) ## Stratified by antibio_exposure
## 0 1 p test
## n 558 232
## antibio_term (mean (SD)) 0.00 (0.00) 4.59 (4.83) <0.001
## child_numl (mean (SD)) 8.30 (1.96) 8.41 (2.03) 0.508
## cci_num (mean (SD)) 4.51 (1.32) 4.52 (1.33) 0.917
## map_day0 (mean (SD)) 2.86 (2.92) 3.13 (2.98) 0.231
## t_bil (mean (SD)) 2.01 (1.88) 2.26 (2.03) 0.098
## ast (mean (SD)) 74.97 (82.14) 83.16 (106.55) 0.248
## alt (mean (SD)) 38.83 (42.49) 42.57 (49.00) 0.285
## wbc (mean (SD)) 8194.22 (4242.15) 8945.32 (5438.02) 0.039
## hb (mean (SD)) 8.69 (2.47) 8.98 (2.39) 0.133
## plt (mean (SD)) 116.49 (61.73) 116.46 (84.84) 0.996
## alb (mean (SD)) 2.91 (0.57) 2.92 (0.64) 0.789
## crp (mean (SD)) 0.76 (1.39) 0.72 (1.57) 0.699
## pt (mean (SD)) 57.00 (16.56) 57.31 (17.01) 0.816
## aptt (mean (SD)) 31.95 (12.69) 34.43 (17.33) 0.035
## los (mean (SD)) 12.92 (13.73) 13.06 (17.84) 0.905
## antibio_type_a (%) <0.001
## a 0 ( 0.0) 4 ( 1.7)
## b 0 ( 0.0) 32 (13.8)
## c 0 ( 0.0) 51 (22.0)
## d 0 ( 0.0) 104 (44.8)
## e 0 ( 0.0) 2 ( 0.9)
## f 0 ( 0.0) 12 ( 5.2)
## g 0 ( 0.0) 2 ( 0.9)
## h 0 ( 0.0) 22 ( 9.5)
## i 0 ( 0.0) 3 ( 1.3)
## 99 558 (100.0) 0 ( 0.0)
## antibio_type_b (%) <0.001
## c 0 ( 0.0) 1 ( 0.4)
## d 0 ( 0.0) 5 ( 2.2)
## f 0 ( 0.0) 2 ( 0.9)
## h 0 ( 0.0) 4 ( 1.7)
## i 0 ( 0.0) 2 ( 0.9)
## 99 558 (100.0) 218 (94.0)
## age_cate (%) 0.235
## 0 322 ( 57.7) 143 (61.6)
## 1 144 ( 25.8) 46 (19.8)
## 2 90 ( 16.1) 43 (18.5)
## 3 2 ( 0.4) 0 ( 0.0)
## sex = F (%) 141 ( 25.3) 51 (22.0) 0.374
## brthel_cate (%) 0.486
## 0 186 ( 41.5) 83 (40.5)
## 1 152 ( 33.9) 63 (30.7)
## 2 110 ( 24.6) 59 (28.8)
## malignancy = 1 (%) 65 ( 11.6) 29 (12.5) 0.829
## child_score (%) 0.776
## a 93 ( 18.8) 42 (19.4)
## b 266 ( 53.7) 110 (50.9)
## c 136 ( 27.5) 64 (29.6)
## hd = 1 (%) 8 ( 1.4) 3 ( 1.3) 1.000
## hcc = 1 (%) 112 ( 20.1) 38 (16.4) 0.269
## alocohol = 1 (%) 246 ( 44.1) 127 (54.7) 0.008
## past_varix_rup = 1 (%) 127 ( 22.8) 63 (27.2) 0.221
## antiplate_user_01 = 1 (%) 8 ( 1.4) 3 ( 1.3) 1.000
## anticoag_user01 = 1 (%) 8 ( 1.4) 5 ( 2.2) 0.675
## nsaids_user01 = 1 (%) 13 ( 2.3) 5 ( 2.2) 1.000
## steroid_user01 = 1 (%) 2 ( 0.4) 0 ( 0.0) 0.892
## ppi_user01 = 1 (%) 486 ( 87.1) 213 (91.8) 0.077
## shock = 1 (%) 197 ( 36.5) 94 (41.4) 0.236
## eGFR30 = 1 (%) 40 ( 7.3) 19 ( 8.3) 0.769
## outcome = 1 (%) 44 ( 7.9) 22 ( 9.5) 0.550
## death_day2 = 1 (%) 34 ( 6.1) 14 ( 6.0) 1.000
## rebleed_end_hemo = 1 (%) 9 ( 1.6) 6 ( 2.6) 0.531
## sbp = 1 (%) 10 ( 1.8) 5 ( 2.2) 0.957
## infection = 1 (%) 33 ( 5.9) 4 ( 1.7) 0.019
## infection_all = 1 (%) 38 ( 6.8) 8 ( 3.4) 0.095
## Stratified by antibio_exposure
## SMD Missing
## n
## antibio_term (mean (SD)) 1.345 0.0
## child_numl (mean (SD)) 0.054 12.8
## cci_num (mean (SD)) 0.008 0.0
## map_day0 (mean (SD)) 0.093 0.0
## t_bil (mean (SD)) 0.129 3.2
## ast (mean (SD)) 0.086 2.2
## alt (mean (SD)) 0.082 2.2
## wbc (mean (SD)) 0.154 1.8
## hb (mean (SD)) 0.119 1.8
## plt (mean (SD)) <0.001 1.8
## alb (mean (SD)) 0.021 4.6
## crp (mean (SD)) 0.030 4.7
## pt (mean (SD)) 0.018 5.8
## aptt (mean (SD)) 0.163 11.4
## los (mean (SD)) 0.009 0.0
## antibio_type_a (%) 10.677 0.0
## a
## b
## c
## d
## e
## f
## g
## h
## i
## 99
## antibio_type_b (%) 0.358 0.0
## c
## d
## f
## h
## i
## 99
## age_cate (%) 0.170 0.0
## 0
## 1
## 2
## 3
## sex = F (%) 0.077 0.0
## brthel_cate (%) 0.101 17.3
## 0
## 1
## 2
## malignancy = 1 (%) 0.026 0.0
## child_score (%) 0.058 10.0
## a
## b
## c
## hd = 1 (%) 0.012 0.0
## hcc = 1 (%) 0.096 0.0
## alocohol = 1 (%) 0.214 0.0
## past_varix_rup = 1 (%) 0.102 0.0
## antiplate_user_01 = 1 (%) 0.012 0.0
## anticoag_user01 = 1 (%) 0.054 0.0
## nsaids_user01 = 1 (%) 0.012 0.0
## steroid_user01 = 1 (%) 0.085 0.0
## ppi_user01 = 1 (%) 0.154 0.0
## shock = 1 (%) 0.100 3.0
## eGFR30 = 1 (%) 0.034 1.9
## outcome = 1 (%) 0.057 0.0
## death_day2 = 1 (%) 0.002 0.0
## rebleed_end_hemo = 1 (%) 0.068 0.0
## sbp = 1 (%) 0.026 0.0
## infection = 1 (%) 0.220 0.0
## infection_all = 1 (%) 0.153 0.0連続量をアウトカムも含め可視化
df |> #全体
select(col_cont) |>
pivot_longer(cols = col_cont, names_to = "name", values_to = "value") |>
ggplot()+
geom_histogram(aes(x = value), color = "black")+
facet_wrap(~ name, scales = "free", ncol = 5) +
theme_bw()+
theme(text = element_text(size = 12))## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.## Warning: Removed 411 rows containing non-finite values (`stat_bin()`).
素データでのtable1の作成
col_continuous <-
c("child_numl","cci_num", "map_day0", "t_bil", "ast", "alt", "wbc", "hb", "plt", "alb", "crp", "pt", "aptt")
col_factors <-
c("age_cate","sex", "brthel_cate", "malignancy", "child_score", "hd", "hcc", "alocohol", "past_varix_rup", "antiplate_user_01", "anticoag_user01", "nsaids_user01" ,"steroid_user01", "ppi_user01","vaso", "shock", "eGFR30")
df|> select(c(col_continuous, col_factors, "antibio_exposure"))|>
CreateTableOne(vars = c(col_continuous, col_factors), strata = "antibio_exposure", factorVars = col_factors) -> tableone## Warning: Using an external vector in selections was deprecated in tidyselect 1.1.0.
## ℹ Please use `all_of()` or `any_of()` instead.
## # Was:
## data %>% select(col_continuous)
##
## # Now:
## data %>% select(all_of(col_continuous))
##
## See <https://tidyselect.r-lib.org/reference/faq-external-vector.html>.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.## Warning: Using an external vector in selections was deprecated in tidyselect 1.1.0.
## ℹ Please use `all_of()` or `any_of()` instead.
## # Was:
## data %>% select(col_factors)
##
## # Now:
## data %>% select(all_of(col_factors))
##
## See <https://tidyselect.r-lib.org/reference/faq-external-vector.html>.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.print(tableone, smd = T, missing = T, test = T, explain = T) ## Stratified by antibio_exposure
## 0 1 p test
## n 558 232
## child_numl (mean (SD)) 8.30 (1.96) 8.41 (2.03) 0.508
## cci_num (mean (SD)) 4.51 (1.32) 4.52 (1.33) 0.917
## map_day0 (mean (SD)) 2.86 (2.92) 3.13 (2.98) 0.231
## t_bil (mean (SD)) 2.01 (1.88) 2.26 (2.03) 0.098
## ast (mean (SD)) 74.97 (82.14) 83.16 (106.55) 0.248
## alt (mean (SD)) 38.83 (42.49) 42.57 (49.00) 0.285
## wbc (mean (SD)) 8194.22 (4242.15) 8945.32 (5438.02) 0.039
## hb (mean (SD)) 8.69 (2.47) 8.98 (2.39) 0.133
## plt (mean (SD)) 116.49 (61.73) 116.46 (84.84) 0.996
## alb (mean (SD)) 2.91 (0.57) 2.92 (0.64) 0.789
## crp (mean (SD)) 0.76 (1.39) 0.72 (1.57) 0.699
## pt (mean (SD)) 57.00 (16.56) 57.31 (17.01) 0.816
## aptt (mean (SD)) 31.95 (12.69) 34.43 (17.33) 0.035
## age_cate (%) 0.235
## 0 322 (57.7) 143 (61.6)
## 1 144 (25.8) 46 (19.8)
## 2 90 (16.1) 43 (18.5)
## 3 2 ( 0.4) 0 ( 0.0)
## sex = F (%) 141 (25.3) 51 (22.0) 0.374
## brthel_cate (%) 0.486
## 0 186 (41.5) 83 (40.5)
## 1 152 (33.9) 63 (30.7)
## 2 110 (24.6) 59 (28.8)
## malignancy = 1 (%) 65 (11.6) 29 (12.5) 0.829
## child_score (%) 0.776
## a 93 (18.8) 42 (19.4)
## b 266 (53.7) 110 (50.9)
## c 136 (27.5) 64 (29.6)
## hd = 1 (%) 8 ( 1.4) 3 ( 1.3) 1.000
## hcc = 1 (%) 112 (20.1) 38 (16.4) 0.269
## alocohol = 1 (%) 246 (44.1) 127 (54.7) 0.008
## past_varix_rup = 1 (%) 127 (22.8) 63 (27.2) 0.221
## antiplate_user_01 = 1 (%) 8 ( 1.4) 3 ( 1.3) 1.000
## anticoag_user01 = 1 (%) 8 ( 1.4) 5 ( 2.2) 0.675
## nsaids_user01 = 1 (%) 13 ( 2.3) 5 ( 2.2) 1.000
## steroid_user01 = 1 (%) 2 ( 0.4) 0 ( 0.0) 0.892
## ppi_user01 = 1 (%) 486 (87.1) 213 (91.8) 0.077
## vaso = 1 (%) 19 ( 3.4) 7 ( 3.0) 0.953
## shock = 1 (%) 197 (36.5) 94 (41.4) 0.236
## eGFR30 = 1 (%) 40 ( 7.3) 19 ( 8.3) 0.769
## Stratified by antibio_exposure
## SMD Missing
## n
## child_numl (mean (SD)) 0.054 12.8
## cci_num (mean (SD)) 0.008 0.0
## map_day0 (mean (SD)) 0.093 0.0
## t_bil (mean (SD)) 0.129 3.2
## ast (mean (SD)) 0.086 2.2
## alt (mean (SD)) 0.082 2.2
## wbc (mean (SD)) 0.154 1.8
## hb (mean (SD)) 0.119 1.8
## plt (mean (SD)) <0.001 1.8
## alb (mean (SD)) 0.021 4.6
## crp (mean (SD)) 0.030 4.7
## pt (mean (SD)) 0.018 5.8
## aptt (mean (SD)) 0.163 11.4
## age_cate (%) 0.170 0.0
## 0
## 1
## 2
## 3
## sex = F (%) 0.077 0.0
## brthel_cate (%) 0.101 17.3
## 0
## 1
## 2
## malignancy = 1 (%) 0.026 0.0
## child_score (%) 0.058 10.0
## a
## b
## c
## hd = 1 (%) 0.012 0.0
## hcc = 1 (%) 0.096 0.0
## alocohol = 1 (%) 0.214 0.0
## past_varix_rup = 1 (%) 0.102 0.0
## antiplate_user_01 = 1 (%) 0.012 0.0
## anticoag_user01 = 1 (%) 0.054 0.0
## nsaids_user01 = 1 (%) 0.012 0.0
## steroid_user01 = 1 (%) 0.085 0.0
## ppi_user01 = 1 (%) 0.154 0.0
## vaso = 1 (%) 0.022 0.0
## shock = 1 (%) 0.100 3.0
## eGFR30 = 1 (%) 0.034 1.9連続量の可視化
df |> #全体
select(col_continuous) |>
pivot_longer(cols = col_continuous, names_to = "name", values_to = "value") |>
ggplot()+
geom_histogram(aes(x = value), color = "black")+
facet_wrap(~ name, scales = "free", ncol = 5) +
theme_bw()+
theme(text = element_text(size = 12))## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.## Warning: Removed 411 rows containing non-finite values (`stat_bin()`).
losのmedianとIQRを求める
library(dplyr)
# antibio_exposureの0.1を基準にデータフレームを群別
grouped_df <- df %>%
mutate(group = ifelse(antibio_exposure >= 0.1, "Group1", "Group0"))
# 各群のlosのmedianとIQRを計算
results <- grouped_df %>%
group_by(group) %>%
summarise(median_los = median(los),
lower_quartile = quantile(los, 0.25),
upper_quartile = quantile(los, 0.75)) %>%
mutate(IQR_los = paste0(lower_quartile, " - ", upper_quartile))
# 不要な列を削除
results <- results[, c("group", "median_los", "IQR_los")]
print(results)## # A tibble: 2 × 3
## group median_los IQR_los
## <chr> <dbl> <chr>
## 1 Group0 9 6 - 15
## 2 Group1 8 5 - 15単変量の解析をみてみる
# Load necessary libraries
library(tidyverse)
library(MASS)##
## 次のパッケージを付け加えます: 'MASS'## 以下のオブジェクトは 'package:dplyr' からマスクされています:
##
## selectlibrary(broom)
# Univariate logistic regression for outcome, death_day2, rebleed_end_hemo, and sepsis
outcome_model <- glm(outcome ~ antibio_exposure, data = df, family = binomial(link = "logit"))
death_day2_model <- glm(death_day2 ~ antibio_exposure, data = df, family = binomial(link = "logit"))
rebleed_end_hemo_model <- glm(rebleed_end_hemo ~ antibio_exposure, data = df, family = binomial(link = "logit"))
sbp_model <- glm(sbp ~ antibio_exposure, data = df, family = binomial(link = "logit"))
# Univariate negative binomial regression for los
los_model <- glm.nb(los ~ antibio_exposure, data = df)
# Calculate odds ratios
outcome_or <- coef(outcome_model) %>% exp()
death_day2_or <- coef(death_day2_model) %>% exp()
rebleed_end_hemo_or <- coef(rebleed_end_hemo_model) %>% exp()
sbp_or <- coef(sbp_model) %>% exp()
# Calculate rate ratio for LOS
los_rr <- coef(los_model) %>% exp()
# Get odds ratios and 95% confidence intervals for logistic regression models
outcome_ci <- confint(outcome_model, method = "Wald") %>% as.data.frame() %>%
rownames_to_column("term") %>%
mutate(estimate = exp(`2.5 %`), conf.low = exp(`2.5 %`), conf.high = exp(`97.5 %`))## Waiting for profiling to be done...death_day2_ci <- confint(death_day2_model, method = "Wald") %>% as.data.frame() %>%
rownames_to_column("term") %>%
mutate(estimate = exp(`2.5 %`), conf.low = exp(`2.5 %`), conf.high = exp(`97.5 %`))## Waiting for profiling to be done...rebleed_end_hemo_ci <- confint(rebleed_end_hemo_model, method = "Wald") %>% as.data.frame() %>%
rownames_to_column("term") %>%
mutate(estimate = exp(`2.5 %`), conf.low = exp(`2.5 %`), conf.high = exp(`97.5 %`))## Waiting for profiling to be done...sbp_ci <- confint(sbp_model, method = "Wald") %>% as.data.frame() %>%
rownames_to_column("term") %>%
mutate(estimate = exp(`2.5 %`), conf.low = exp(`2.5 %`), conf.high = exp(`97.5 %`))## Waiting for profiling to be done...# Get rate ratios and 95% confidence intervals for negative binomial regression model
los_ci <- confint(los_model, method = "Wald") %>% as.data.frame() %>%
rownames_to_column("term") %>%
mutate(estimate = exp(`2.5 %`), conf.low = exp(`2.5 %`), conf.high = exp(`97.5 %`))## Waiting for profiling to be done...# Print results with rounding
print("Outcome:")## [1] "Outcome:"print(paste("OR:", round(outcome_or["antibio_exposure"], 4), "95% CI:", round(outcome_ci$conf.low[2], 4), "to", round(outcome_ci$conf.high[2], 4)))## [1] "OR: 1.2238 95% CI: 0.7045 to 2.0704"print("Death day 2:")## [1] "Death day 2:"print(paste("OR:", round(death_day2_or["antibio_exposure"], 4), "95% CI:", round(death_day2_ci$conf.low[2], 4), "to", round(death_day2_ci$conf.high[2], 4)))## [1] "OR: 0.9897 95% CI: 0.5052 to 1.8429"print("Rebleed end hemo:")## [1] "Rebleed end hemo:"print(paste("OR:", round(rebleed_end_hemo_or["antibio_exposure"], 4), "95% CI:", round(rebleed_end_hemo_ci$conf.low[2], 4), "to", round(rebleed_end_hemo_ci$conf.high[2], 4)))## [1] "OR: 1.6195 95% CI: 0.5374 to 4.5435"print("SBP:")## [1] "SBP:"print(paste("OR:", round(sbp_or["antibio_exposure"], 4), "95% CI:", round(sbp_ci$conf.low[2], 4), "to", round(sbp_ci$conf.high[2], 4)))## [1] "OR: 1.207 95% CI: 0.3725 to 3.4367"print("LOS:")## [1] "LOS:"print(paste("RR:", round(los_rr["antibio_exposure"], 4), "95% CI:", round(los_ci$conf.low[2], 4), "to", round(los_ci$conf.high[2], 4)))## [1] "RR: 1.0109 95% CI: 0.9 to 1.137"多重代入
imp1 <- mice(df,
m=100, #mは代入の結果として作成するデータセットの数
maxit =200, #maxitは反復回数
method="pmm",
printFlag = FALSE,
seed = 10
)## Warning: Number of logged events: 360000plot(imp1)
ATTでIPTW
# Function to perform IPTW analysis on different outcomes
perform_iptw_analysis <- function(outcome_name) {
gee_results <- list()
summary_tables <- list()
conf_ints <- list()
weighted_data <- list()
ps_data <- list()
for (i in 1:100) {
imputed_data <- complete(imp1, i)
# Fit GEE propensity score model
propensity_model <- geeglm(antibio_exposure ~ age_cate + sex + brthel_cate + child_numl + cci_num + malignancy + hd + hcc + alocohol + past_varix_rup + antiplate_user_01 + anticoag_user01 + nsaids_user01 + steroid_user01 + ppi_user01 + vaso + map_day0 + shock + t_bil + ast + alt + wbc + hb + plt + alb + eGFR30 + crp + pt + aptt,
data = imputed_data,
id = hosp_id,
family = binomial("logit"),
corstr = "exchangeable")
# Calculate propensity scores
ps <- predict(propensity_model, type = "response", newdata = imputed_data)
imputed_data$ps <- ps
# Calculate IPTW and fit the logistic regression model by ATT
stabilized_weights <- with(imputed_data, antibio_exposure + (1 - antibio_exposure) * ps / (1 - ps))
iptw_model <- glm(as.formula(paste0(outcome_name, " ~ antibio_exposure")), data = imputed_data, weights = stabilized_weights, family = binomial("logit"))
# Store the results
gee_results[[i]] <- iptw_model
summary_tables[[i]] <- summary(iptw_model)
conf_ints[[i]] <- confint(iptw_model, level = 0.95, method = "Wald")
# Store weighted data
imputed_data$weight <- stabilized_weights
weighted_data[[i]] <- imputed_data
ps_data[[i]] <- imputed_data[, c("ps", "antibio_exposure")]
}
# Save the results from each imputed dataset
params <- list()
for (i in 1:100) {
params[[i]] <- data.frame(Estimate = coef(gee_results[[i]])[2],
Std.Error = coef(summary_tables[[i]])[2, "Std. Error"],
df = Inf)
}
#Calculate the pooled estimates
Qbar <- mean(sapply(params, function(x) x$Estimate))
se_Qbar <- sqrt(sum(sapply(params, function(x) x$Std.Error^2))/length(params))
#Calculate the pooled odds ratio and confidence intervals
pooled_or <- exp(Qbar)
pooled_lower <- exp(Qbar - 1.96 * se_Qbar)
pooled_upper <- exp(Qbar + 1.96 * se_Qbar)
#Calculate the pooled p-value (two-tailed)
pooled_p <- 2 * (1 - pnorm(abs(Qbar / se_Qbar)))
return(list(OR = pooled_or, Lower_CI = pooled_lower, Upper_CI = pooled_upper, P_value = pooled_p, PS_data = ps_data))
}
#Perform IPTW analysis for each outcome
outcomes <- c("outcome", "death_day2", "rebleed_end_hemo", "sbp")
results <- lapply(outcomes, perform_iptw_analysis)## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...#Combine results into a data frame for plotting
forest_plot_data <- data.frame(
Outcome = outcomes,
OR = sapply(results, function(x) x$OR),
Lower_CI = sapply(results, function(x) x$Lower_CI),
Upper_CI = sapply(results, function(x) x$Upper_CI),
P_value = sapply(results, function(x) x$P_value)
)
#Display results as a table
print(forest_plot_data)## Outcome OR Lower_CI Upper_CI P_value
## 1 outcome 1.0678232 0.5415712 2.105441 0.8497386
## 2 death_day2 0.8305911 0.3781029 1.824587 0.6438717
## 3 rebleed_end_hemo 1.3650093 0.3537599 5.266992 0.6515137
## 4 sbp 1.3361396 0.3081420 5.793656 0.6986266#Create forest plot with the specified order and adjusted p-value position
forest_plot_ordered <- ggplot(forest_plot_data, aes(x = Outcome, y = OR, ymin = Lower_CI, ymax = Upper_CI)) +
geom_pointrange(size = 0.5) +
geom_hline(yintercept = 1, linetype = "dashed") +
scale_y_log10() +
coord_flip() +
theme_bw() +
labs(
x = "Outcome",
y = "Odds Ratio (95% Confidence Interval)"
) +
geom_text(aes(label = sprintf("p = %.3f", P_value)), vjust = -1.0, size = 3) +
scale_x_discrete(limits = c("sbp", "rebleed_end_hemo", "death_day2","outcome"), labels = c("SBP", "Rebleeding", "In hospital mortality", "Composite outcome"))
#Display forest plot
forest_plot_ordered
#Calculate and display c-statistics for each imputed dataset
c_statistics <- lapply(results[[1]]$PS_data, function(data) {
roc_obj <- roc(data$antibio_exposure, data$ps)
auc(roc_obj)
})## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases## Setting levels: control = 0, case = 1## Setting direction: controls < cases#Print c-statistics
print(c_statistics)## [[1]]
## Area under the curve: 0.6193
##
## [[2]]
## Area under the curve: 0.6236
##
## [[3]]
## Area under the curve: 0.6255
##
## [[4]]
## Area under the curve: 0.6228
##
## [[5]]
## Area under the curve: 0.6236
##
## [[6]]
## Area under the curve: 0.6266
##
## [[7]]
## Area under the curve: 0.6233
##
## [[8]]
## Area under the curve: 0.6223
##
## [[9]]
## Area under the curve: 0.624
##
## [[10]]
## Area under the curve: 0.621
##
## [[11]]
## Area under the curve: 0.6155
##
## [[12]]
## Area under the curve: 0.6333
##
## [[13]]
## Area under the curve: 0.624
##
## [[14]]
## Area under the curve: 0.6323
##
## [[15]]
## Area under the curve: 0.626
##
## [[16]]
## Area under the curve: 0.6229
##
## [[17]]
## Area under the curve: 0.6258
##
## [[18]]
## Area under the curve: 0.6145
##
## [[19]]
## Area under the curve: 0.6191
##
## [[20]]
## Area under the curve: 0.6244
##
## [[21]]
## Area under the curve: 0.6168
##
## [[22]]
## Area under the curve: 0.6247
##
## [[23]]
## Area under the curve: 0.6323
##
## [[24]]
## Area under the curve: 0.6258
##
## [[25]]
## Area under the curve: 0.6297
##
## [[26]]
## Area under the curve: 0.6247
##
## [[27]]
## Area under the curve: 0.6233
##
## [[28]]
## Area under the curve: 0.6237
##
## [[29]]
## Area under the curve: 0.619
##
## [[30]]
## Area under the curve: 0.6279
##
## [[31]]
## Area under the curve: 0.6252
##
## [[32]]
## Area under the curve: 0.6315
##
## [[33]]
## Area under the curve: 0.62
##
## [[34]]
## Area under the curve: 0.631
##
## [[35]]
## Area under the curve: 0.6296
##
## [[36]]
## Area under the curve: 0.624
##
## [[37]]
## Area under the curve: 0.6291
##
## [[38]]
## Area under the curve: 0.6315
##
## [[39]]
## Area under the curve: 0.6203
##
## [[40]]
## Area under the curve: 0.6382
##
## [[41]]
## Area under the curve: 0.6234
##
## [[42]]
## Area under the curve: 0.6212
##
## [[43]]
## Area under the curve: 0.6174
##
## [[44]]
## Area under the curve: 0.6316
##
## [[45]]
## Area under the curve: 0.6213
##
## [[46]]
## Area under the curve: 0.6205
##
## [[47]]
## Area under the curve: 0.6194
##
## [[48]]
## Area under the curve: 0.6236
##
## [[49]]
## Area under the curve: 0.6171
##
## [[50]]
## Area under the curve: 0.6228
##
## [[51]]
## Area under the curve: 0.6176
##
## [[52]]
## Area under the curve: 0.6215
##
## [[53]]
## Area under the curve: 0.6156
##
## [[54]]
## Area under the curve: 0.6195
##
## [[55]]
## Area under the curve: 0.6166
##
## [[56]]
## Area under the curve: 0.6134
##
## [[57]]
## Area under the curve: 0.6175
##
## [[58]]
## Area under the curve: 0.625
##
## [[59]]
## Area under the curve: 0.6164
##
## [[60]]
## Area under the curve: 0.6261
##
## [[61]]
## Area under the curve: 0.6302
##
## [[62]]
## Area under the curve: 0.627
##
## [[63]]
## Area under the curve: 0.6234
##
## [[64]]
## Area under the curve: 0.6239
##
## [[65]]
## Area under the curve: 0.6217
##
## [[66]]
## Area under the curve: 0.6184
##
## [[67]]
## Area under the curve: 0.6221
##
## [[68]]
## Area under the curve: 0.6223
##
## [[69]]
## Area under the curve: 0.6147
##
## [[70]]
## Area under the curve: 0.6218
##
## [[71]]
## Area under the curve: 0.6113
##
## [[72]]
## Area under the curve: 0.6231
##
## [[73]]
## Area under the curve: 0.6193
##
## [[74]]
## Area under the curve: 0.6226
##
## [[75]]
## Area under the curve: 0.6283
##
## [[76]]
## Area under the curve: 0.6173
##
## [[77]]
## Area under the curve: 0.6367
##
## [[78]]
## Area under the curve: 0.6359
##
## [[79]]
## Area under the curve: 0.63
##
## [[80]]
## Area under the curve: 0.624
##
## [[81]]
## Area under the curve: 0.6186
##
## [[82]]
## Area under the curve: 0.6174
##
## [[83]]
## Area under the curve: 0.6193
##
## [[84]]
## Area under the curve: 0.6241
##
## [[85]]
## Area under the curve: 0.6241
##
## [[86]]
## Area under the curve: 0.6275
##
## [[87]]
## Area under the curve: 0.6267
##
## [[88]]
## Area under the curve: 0.6276
##
## [[89]]
## Area under the curve: 0.621
##
## [[90]]
## Area under the curve: 0.6228
##
## [[91]]
## Area under the curve: 0.6322
##
## [[92]]
## Area under the curve: 0.6179
##
## [[93]]
## Area under the curve: 0.6163
##
## [[94]]
## Area under the curve: 0.6188
##
## [[95]]
## Area under the curve: 0.6209
##
## [[96]]
## Area under the curve: 0.6204
##
## [[97]]
## Area under the curve: 0.6144
##
## [[98]]
## Area under the curve: 0.6226
##
## [[99]]
## Area under the curve: 0.6232
##
## [[100]]
## Area under the curve: 0.6332#Plot the propensity score overlap for each imputed dataset
for (i in 1:5) {
ggplot(results[[1]]$PS_data[[i]], aes(x = ps, fill = as.factor(antibio_exposure))) +
geom_histogram(position = "identity", alpha = 0.5, bins = 30) +
labs(
title = paste0("Imputed dataset ", i),
x = "Propensity Score",
y = "Count",
fill = "Antibio Exposure"
) +
theme_bw() +
scale_fill_discrete(labels = c("No Exposure", "Exposure")) +
xlim(0, 1) +
guides(fill = guide_legend(reverse = TRUE)) -> p
print(p)
}
ATTで負の二項分布で解析してみる
# Install and load the MASS package
if (!requireNamespace("MASS", quietly = TRUE)) {
install.packages("MASS")
}
library(MASS)
# LOS
# 本解析: IPTW models on imputed data with Negative Binomial distribution
gee_results <- list()
summary_tables <- list()
conf_ints <- list()
weighted_data <- list()
for (i in 1:100) {
imputed_data <- complete(imp1, i)
# Fit GEE propensity score model
propensity_model <- geeglm(antibio_exposure ~ age_cate + sex + brthel_cate + child_numl + cci_num + malignancy + hd + hcc + alocohol + past_varix_rup + antiplate_user_01 + anticoag_user01 + nsaids_user01 + steroid_user01 + ppi_user01 + vaso + map_day0 + shock + t_bil + ast + alt + wbc + hb + plt + alb + eGFR30 + crp + pt + aptt,
data = imputed_data,
id = hosp_id,
family = binomial("logit"),
corstr = "exchangeable")
# Calculate propensity scores
ps <- predict(propensity_model, type = "response", newdata = imputed_data)
imputed_data$ps <- ps
# Calculate IPTW and fit the negative binomial regression model by ATT
stabilized_weights <- with(imputed_data, antibio_exposure + (1 - antibio_exposure) * ps / (1 - ps))
iptw_model <- glm.nb(los ~ antibio_exposure, data = imputed_data, weights = stabilized_weights)
# Store the results
gee_results[[i]] <- iptw_model
summary_tables[[i]] <- summary(iptw_model)
conf_ints[[i]] <- confint(iptw_model, level = 0.95, method = "Wald")
# Store weighted data
imputed_data$weight <- stabilized_weights
weighted_data[[i]] <- imputed_data
}## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...
## Waiting for profiling to be done...# Save the results from each imputed dataset
params <- list()
for (i in 1:100) {
params[[i]] <- data.frame(estimate = coef(gee_results[[i]])[2],
std.error = coef(summary_tables[[i]])[2, "Std. Error"],
df = Inf) # Add infinite degrees of freedom
}
# Manually pool the results
Qbar <- mean(sapply(params, function(x) x$estimate))
se_Qbar <- sqrt(sum(sapply(params, function(x) x$std.error^2))/length(params))
# Calculate the pooled rate ratio and confidence intervals
pooled_rr <- exp(Qbar)
pooled_lower <- exp(Qbar - 1.96 * se_Qbar)
pooled_upper <- exp(Qbar + 1.96 * se_Qbar)
# Calculate the pooled z-value
pooled_z <- Qbar / se_Qbar
# Calculate the pooled p-value (two-tailed)
pooled_p <- 2 * (1 - pnorm(abs(pooled_z)))
cat("\nPooled results (IPTW with Negative Binomial distribution)\n")##
## Pooled results (IPTW with Negative Binomial distribution)cat("\nRate Ratio:", pooled_rr, "\n")##
## Rate Ratio: 1.008119cat("\n95% Confidence Interval: (", pooled_lower, ", ", pooled_upper, ")\n")##
## 95% Confidence Interval: ( 0.8641974 , 1.176009 )cat("\nP-value:", pooled_p, "\n")##
## P-value: 0.9180501サマリー表の作成 補完後IPTW前
# 補完データの1つを選択 (例: 最初の補完データ)
selected_data <- complete(imp1, 1)
# IPTW前のデータのTable 1
unweighted_tableone <- CreateTableOne(vars = c(col_continuous, col_factors), strata = "antibio_exposure",data = selected_data, factorVars = col_factors)
print(unweighted_tableone, smd = T, missing = T, test = T, explain = T)## Stratified by antibio_exposure
## 0 1 p test
## n 558 232
## child_numl (mean (SD)) 8.34 (1.98) 8.53 (2.03) 0.230
## cci_num (mean (SD)) 4.51 (1.32) 4.52 (1.33) 0.917
## map_day0 (mean (SD)) 2.86 (2.92) 3.13 (2.98) 0.231
## t_bil (mean (SD)) 2.00 (1.86) 2.26 (2.02) 0.073
## ast (mean (SD)) 75.18 (82.66) 82.89 (106.15) 0.274
## alt (mean (SD)) 38.84 (42.18) 42.47 (48.81) 0.294
## wbc (mean (SD)) 8163.49 (4242.94) 8930.91 (5430.68) 0.034
## hb (mean (SD)) 8.74 (2.49) 8.97 (2.39) 0.225
## plt (mean (SD)) 116.16 (61.53) 116.85 (84.86) 0.899
## alb (mean (SD)) 2.92 (0.58) 2.94 (0.65) 0.752
## crp (mean (SD)) 0.77 (1.40) 0.71 (1.55) 0.614
## pt (mean (SD)) 57.68 (16.70) 57.57 (17.14) 0.936
## aptt (mean (SD)) 32.05 (13.50) 34.20 (16.68) 0.058
## age_cate (%) 0.235
## 0 322 (57.7) 143 (61.6)
## 1 144 (25.8) 46 (19.8)
## 2 90 (16.1) 43 (18.5)
## 3 2 ( 0.4) 0 ( 0.0)
## sex = F (%) 141 (25.3) 51 (22.0) 0.374
## brthel_cate (%) 0.200
## 0 236 (42.3) 94 (40.5)
## 1 189 (33.9) 69 (29.7)
## 2 133 (23.8) 69 (29.7)
## malignancy = 1 (%) 65 (11.6) 29 (12.5) 0.829
## child_score (%) 0.600
## a 104 (18.6) 43 (18.5)
## b 307 (55.0) 120 (51.7)
## c 147 (26.3) 69 (29.7)
## hd = 1 (%) 8 ( 1.4) 3 ( 1.3) 1.000
## hcc = 1 (%) 112 (20.1) 38 (16.4) 0.269
## alocohol = 1 (%) 246 (44.1) 127 (54.7) 0.008
## past_varix_rup = 1 (%) 127 (22.8) 63 (27.2) 0.221
## antiplate_user_01 = 1 (%) 8 ( 1.4) 3 ( 1.3) 1.000
## anticoag_user01 = 1 (%) 8 ( 1.4) 5 ( 2.2) 0.675
## nsaids_user01 = 1 (%) 13 ( 2.3) 5 ( 2.2) 1.000
## steroid_user01 = 1 (%) 2 ( 0.4) 0 ( 0.0) 0.892
## ppi_user01 = 1 (%) 486 (87.1) 213 (91.8) 0.077
## vaso = 1 (%) 19 ( 3.4) 7 ( 3.0) 0.953
## shock = 1 (%) 204 (36.6) 94 (40.5) 0.335
## eGFR30 = 1 (%) 41 ( 7.3) 19 ( 8.2) 0.795
## Stratified by antibio_exposure
## SMD Missing
## n
## child_numl (mean (SD)) 0.093 0.0
## cci_num (mean (SD)) 0.008 0.0
## map_day0 (mean (SD)) 0.093 0.0
## t_bil (mean (SD)) 0.138 0.0
## ast (mean (SD)) 0.081 0.0
## alt (mean (SD)) 0.080 0.0
## wbc (mean (SD)) 0.157 0.0
## hb (mean (SD)) 0.096 0.0
## plt (mean (SD)) 0.009 0.0
## alb (mean (SD)) 0.024 0.0
## crp (mean (SD)) 0.039 0.0
## pt (mean (SD)) 0.006 0.0
## aptt (mean (SD)) 0.142 0.0
## age_cate (%) 0.170 0.0
## 0
## 1
## 2
## 3
## sex = F (%) 0.077 0.0
## brthel_cate (%) 0.139 0.0
## 0
## 1
## 2
## malignancy = 1 (%) 0.026 0.0
## child_score (%) 0.079 0.0
## a
## b
## c
## hd = 1 (%) 0.012 0.0
## hcc = 1 (%) 0.096 0.0
## alocohol = 1 (%) 0.214 0.0
## past_varix_rup = 1 (%) 0.102 0.0
## antiplate_user_01 = 1 (%) 0.012 0.0
## anticoag_user01 = 1 (%) 0.054 0.0
## nsaids_user01 = 1 (%) 0.012 0.0
## steroid_user01 = 1 (%) 0.085 0.0
## ppi_user01 = 1 (%) 0.154 0.0
## vaso = 1 (%) 0.022 0.0
## shock = 1 (%) 0.081 0.0
## eGFR30 = 1 (%) 0.031 0.0補完後IPTW後のSMD算出
# 1. Calculate SMD for all covariates
calculate_smd <- function(data, weights, col_continuous, col_factors) {
smd_continuous <- sapply(col_continuous, function(x) {
cov_unexposed <- data[data$antibio_exposure == 0, x]
cov_exposed <- data[data$antibio_exposure == 1, x]
w_unexposed <- weights[data$antibio_exposure == 0]
w_exposed <- weights[data$antibio_exposure == 1]
mean_unexposed <- sum(w_unexposed * cov_unexposed) / sum(w_unexposed)
mean_exposed <- sum(w_exposed * cov_exposed) / sum(w_exposed)
sd_unexposed <- sqrt(sum(w_unexposed * (cov_unexposed - mean_unexposed)^2) / sum(w_unexposed))
sd_exposed <- sqrt(sum(w_exposed * (cov_exposed - mean_exposed)^2) / sum(w_exposed))
smd <- (mean_exposed - mean_unexposed) / sqrt((sd_exposed^2 + sd_unexposed^2) / 2)
return(smd)
})
smd_factors <- sapply(col_factors, function(x) {
prop_unexposed <- tapply(weights[data$antibio_exposure == 0], data[data$antibio_exposure == 0, x], sum) / sum(weights[data$antibio_exposure == 0])
prop_exposed <- tapply(weights[data$antibio_exposure == 1], data[data$antibio_exposure == 1, x], sum) / sum(weights[data$antibio_exposure == 1])
common_levels <- intersect(names(prop_unexposed), names(prop_exposed))
prop_unexposed <- prop_unexposed[common_levels]
prop_exposed <- prop_exposed[common_levels]
smd <- sum((prop_exposed - prop_unexposed) / sqrt((prop_exposed * (1 - prop_exposed) + prop_unexposed * (1 - prop_unexposed)) / 2), na.rm = TRUE)
return(smd)
})
return(c(smd_continuous, smd_factors))
}
# 2. Calculate SMD for each imputed dataset
smd_results <- lapply(weighted_data, function(data) {
calculate_smd(data, data$weight, col_continuous, col_factors)
})
# 3. Create a data frame with the mean SMD for each covariate
mean_smd_dataframe <- function(smd_results, col_continuous, col_factors) {
mean_smd_continuous <- sapply(1:length(col_continuous), function(i) {
mean(sapply(smd_results, function(x) x[i]))
})
mean_smd_factors <- sapply((length(col_continuous) + 1):(length(col_continuous) + length(col_factors)), function(i) {
mean(sapply(smd_results, function(x) x[i]))
})
mean_smd_continuous_rounded <- round(mean_smd_continuous, 4)
mean_smd_factors_rounded <- round(mean_smd_factors, 4)
smd_df_continuous <- data.frame(Covariate = col_continuous, Mean_SMD = mean_smd_continuous_rounded)
smd_df_factors <- data.frame(Covariate = col_factors, Mean_SMD = mean_smd_factors_rounded)
smd_df <- rbind(smd_df_continuous, smd_df_factors)
smd_df <- smd_df[order(abs(smd_df$Mean_SMD), decreasing = TRUE),]
return(smd_df)
}
# 4. Display the mean SMD dataframe
mean_smd_df <- mean_smd_dataframe(smd_results, col_continuous, col_factors)
cat("Mean SMD for all covariates (rounded to 4 decimal places):\n")## Mean SMD for all covariates (rounded to 4 decimal places):print(mean_smd_df)## Covariate Mean_SMD
## 13 aptt -0.1693
## 10 alb 0.0842
## 11 crp -0.0780
## 8 hb 0.0562
## 12 pt -0.0562
## 4 t_bil 0.0508
## 9 plt -0.0401
## 1 child_numl 0.0327
## 3 map_day0 0.0309
## 14 age_cate -0.0199
## 2 cci_num -0.0189
## 6 alt -0.0120
## 5 ast -0.0115
## 7 wbc -0.0077
## 18 child_score 0.0060
## 16 brthel_cate -0.0031
## 26 steroid_user01 0.0005
## 15 sex 0.0000
## 17 malignancy 0.0000
## 19 hd 0.0000
## 20 hcc 0.0000
## 21 alocohol 0.0000
## 22 past_varix_rup 0.0000
## 23 antiplate_user_01 0.0000
## 24 anticoag_user01 0.0000
## 25 nsaids_user01 0.0000
## 27 ppi_user01 0.0000
## 28 vaso 0.0000
## 29 shock 0.0000
## 30 eGFR30 0.0000love plotで比較
# Calculate pre-imputation SMD
pre_smd_continuous <- sapply(col_continuous, function(x) {
cov_unexposed <- data[data$antibio_exposure == 0, x, drop = FALSE]
cov_exposed <- data[data$antibio_exposure == 1, x, drop = FALSE]
mean_unexposed <- mean(cov_unexposed[[x]], na.rm = TRUE)
mean_exposed <- mean(cov_exposed[[x]], na.rm = TRUE)
sd_unexposed <- sd(cov_unexposed[[x]], na.rm = TRUE)
sd_exposed <- sd(cov_exposed[[x]], na.rm = TRUE)
smd <- (mean_exposed - mean_unexposed) / sqrt((sd_exposed^2 + sd_unexposed^2) / 2)
return(smd)
})
pre_smd_factors <- calculate_smd(data, rep(1, nrow(data)), c(), col_factors)
pre_smd <- c(pre_smd_continuous, pre_smd_factors)
# Calculate the absolute SMD
abs_pre_smd_continuous <- abs(pre_smd_continuous)
abs_pre_smd_factors <- abs(unlist(pre_smd_factors))
abs_mean_smd <- abs(mean_smd_df$Mean_SMD)
# Combine pre- and post-imputation SMDs
smd_df_pre_post <- data.frame(
Covariate = rep(c(col_continuous, col_factors), 2),
SMD = c(abs_pre_smd_continuous, abs_pre_smd_factors, abs_mean_smd),
Imputation = rep(c("Pre IPTW", "Post IPTW"), each = length(col_continuous) + length(col_factors))
)plot作成のためのpreデータ作成
# パッケージの読み込みは不要
# データの読み込み
pre_table <- read.csv("before.csv", col.names = c("covariate", "non_user", "user", "smd", "missing"))
# NAを削除
pre_table <- pre_table[!is.na(pre_table$smd),]
# Imputation列を追加
pre_table$Imputation <- "Pre IPTW"
# 新しい行名を指定
new_row_names <- c("child_numl", "cci_num", "map_day0", "t_bil", "ast", "alt", "wbc", "hb", "plt", "alb", "crp", "pt", "aptt", "age_cate", "sex", "brthel_cate", "malignancy", "child_score", "hd", "hcc", "alocohol", "past_varix_rup", "antiplate_user_01", "anticoag_user01", "nsaids_user01", "steroid_user01", "ppi_user01", "vaso", "shock", "eGFR30")
# 新しい行名のマッピングデータフレームを作成
mapping <- data.frame(covariate = pre_table$covariate, new_covariate = new_row_names)
# covariate列の名前を変更
pre_table <- merge(pre_table, mapping, by = "covariate", all.x = TRUE)
pre_table$covariate <- NULL
colnames(pre_table)[colnames(pre_table) == "new_covariate"] <- "covariate"
# covariate, smd, Imputation列だけを選択
pre_table <- pre_table[, c("covariate", "smd", "Imputation")]
pre_table## covariate smd Imputation
## 1 age_cate 0.170 Pre IPTW
## 2 alb 0.021 Pre IPTW
## 3 alocohol 0.214 Pre IPTW
## 4 alt 0.082 Pre IPTW
## 5 anticoag_user01 0.054 Pre IPTW
## 6 antiplate_user_01 0.012 Pre IPTW
## 7 aptt 0.163 Pre IPTW
## 8 ast 0.086 Pre IPTW
## 9 brthel_cate 0.101 Pre IPTW
## 10 cci_num 0.008 Pre IPTW
## 11 child_numl 0.054 Pre IPTW
## 12 child_score 0.058 Pre IPTW
## 13 crp 0.030 Pre IPTW
## 14 eGFR30 0.034 Pre IPTW
## 15 hb 0.119 Pre IPTW
## 16 hcc 0.096 Pre IPTW
## 17 hd 0.012 Pre IPTW
## 18 malignancy 0.026 Pre IPTW
## 19 map_day0 0.093 Pre IPTW
## 20 nsaids_user01 0.012 Pre IPTW
## 21 past_varix_rup 0.102 Pre IPTW
## 22 plt 0.001 Pre IPTW
## 23 ppi_user01 0.154 Pre IPTW
## 24 pt 0.018 Pre IPTW
## 25 sex 0.077 Pre IPTW
## 26 shock 0.100 Pre IPTW
## 27 steroid_user01 0.085 Pre IPTW
## 28 t_bil 0.129 Pre IPTW
## 29 vaso 0.022 Pre IPTW
## 30 wbc 0.154 Pre IPTWpostSMDも計算する。
mean_smd_df %>% mutate(Imputation = "Post IPTW") ->mean_smd_df
mean_smd_df <- mean_smd_df[, c("Covariate", "Mean_SMD", "Imputation")]
colnames(mean_smd_df) <- c("covariate", "smd", "Imputation")
# smd列の値を絶対値に変換
mean_smd_df$smd <- abs(mean_smd_df$smd)
mean_smd_df## covariate smd Imputation
## 13 aptt 0.1693 Post IPTW
## 10 alb 0.0842 Post IPTW
## 11 crp 0.0780 Post IPTW
## 8 hb 0.0562 Post IPTW
## 12 pt 0.0562 Post IPTW
## 4 t_bil 0.0508 Post IPTW
## 9 plt 0.0401 Post IPTW
## 1 child_numl 0.0327 Post IPTW
## 3 map_day0 0.0309 Post IPTW
## 14 age_cate 0.0199 Post IPTW
## 2 cci_num 0.0189 Post IPTW
## 6 alt 0.0120 Post IPTW
## 5 ast 0.0115 Post IPTW
## 7 wbc 0.0077 Post IPTW
## 18 child_score 0.0060 Post IPTW
## 16 brthel_cate 0.0031 Post IPTW
## 26 steroid_user01 0.0005 Post IPTW
## 15 sex 0.0000 Post IPTW
## 17 malignancy 0.0000 Post IPTW
## 19 hd 0.0000 Post IPTW
## 20 hcc 0.0000 Post IPTW
## 21 alocohol 0.0000 Post IPTW
## 22 past_varix_rup 0.0000 Post IPTW
## 23 antiplate_user_01 0.0000 Post IPTW
## 24 anticoag_user01 0.0000 Post IPTW
## 25 nsaids_user01 0.0000 Post IPTW
## 27 ppi_user01 0.0000 Post IPTW
## 28 vaso 0.0000 Post IPTW
## 29 shock 0.0000 Post IPTW
## 30 eGFR30 0.0000 Post IPTWpre,postの結合
smd_df_pre_post_main <- rbind(pre_table, mean_smd_df)
smd_df_pre_post_main %>% arrange(smd)## covariate smd Imputation
## 151 sex 0.0000 Post IPTW
## 171 malignancy 0.0000 Post IPTW
## 191 hd 0.0000 Post IPTW
## 201 hcc 0.0000 Post IPTW
## 211 alocohol 0.0000 Post IPTW
## 221 past_varix_rup 0.0000 Post IPTW
## 231 antiplate_user_01 0.0000 Post IPTW
## 241 anticoag_user01 0.0000 Post IPTW
## 251 nsaids_user01 0.0000 Post IPTW
## 271 ppi_user01 0.0000 Post IPTW
## 281 vaso 0.0000 Post IPTW
## 291 shock 0.0000 Post IPTW
## 301 eGFR30 0.0000 Post IPTW
## 261 steroid_user01 0.0005 Post IPTW
## 22 plt 0.0010 Pre IPTW
## 161 brthel_cate 0.0031 Post IPTW
## 181 child_score 0.0060 Post IPTW
## 71 wbc 0.0077 Post IPTW
## 10 cci_num 0.0080 Pre IPTW
## 51 ast 0.0115 Post IPTW
## 6 antiplate_user_01 0.0120 Pre IPTW
## 17 hd 0.0120 Pre IPTW
## 20 nsaids_user01 0.0120 Pre IPTW
## 61 alt 0.0120 Post IPTW
## 24 pt 0.0180 Pre IPTW
## 210 cci_num 0.0189 Post IPTW
## 141 age_cate 0.0199 Post IPTW
## 2 alb 0.0210 Pre IPTW
## 29 vaso 0.0220 Pre IPTW
## 18 malignancy 0.0260 Pre IPTW
## 13 crp 0.0300 Pre IPTW
## 31 map_day0 0.0309 Post IPTW
## 110 child_numl 0.0327 Post IPTW
## 14 eGFR30 0.0340 Pre IPTW
## 91 plt 0.0401 Post IPTW
## 41 t_bil 0.0508 Post IPTW
## 5 anticoag_user01 0.0540 Pre IPTW
## 11 child_numl 0.0540 Pre IPTW
## 81 hb 0.0562 Post IPTW
## 121 pt 0.0562 Post IPTW
## 12 child_score 0.0580 Pre IPTW
## 25 sex 0.0770 Pre IPTW
## 111 crp 0.0780 Post IPTW
## 4 alt 0.0820 Pre IPTW
## 101 alb 0.0842 Post IPTW
## 27 steroid_user01 0.0850 Pre IPTW
## 8 ast 0.0860 Pre IPTW
## 19 map_day0 0.0930 Pre IPTW
## 16 hcc 0.0960 Pre IPTW
## 26 shock 0.1000 Pre IPTW
## 9 brthel_cate 0.1010 Pre IPTW
## 21 past_varix_rup 0.1020 Pre IPTW
## 15 hb 0.1190 Pre IPTW
## 28 t_bil 0.1290 Pre IPTW
## 23 ppi_user01 0.1540 Pre IPTW
## 30 wbc 0.1540 Pre IPTW
## 7 aptt 0.1630 Pre IPTW
## 131 aptt 0.1693 Post IPTW
## 1 age_cate 0.1700 Pre IPTW
## 3 alocohol 0.2140 Pre IPTWlove plot作成
# Assuming smd_df_pre_post_main is your dataframe with columns:
# covariate, smd, Imputation (with values "Pre IPTW" and "Post IPTW")
# Create a mapping of covariate names
covariate_name_mapping <- c(
"age_cate"="Age",
"child_numl" = "Child-Pugh score",
"cci_num" = "Charlson Comorbidity Index",
"map_day0" = "RBC transfusion",
"t_bil" = "Total bilirubin",
"ast" = "Aspartate aminotransferase",
"alt" = "Alanine aminotransferase",
"wbc" = "White blood cell",
"hb" = "Hemoglobin",
"plt" = "Platelet",
"alb" = "Albumin",
"crp" = "C-reactive protein",
"pt" = "Prothrombin time",
"aptt" = "Activated partial thromboplastin time",
"sex" = "Sex",
"brthel_cate" = "Barthel index",
"malignancy" = "Malignant tomor history",
"child_score" = "Child-Pugh classification",
"hd" = "Maintenance hemodialysis",
"hcc" = "Hepatic cancer",
"alocohol" = "Alcohol related disease",
"past_varix_rup" = "Past varix rupture history",
"antiplate_user_01" = "Antiplatelet use",
"anticoag_user01" = "Anticoagulant use",
"nsaids_user01" = "NSAIDs use",
"steroid_user01" = "Corticosteroid use",
"ppi_user01" = "Acid blocker use",
"vaso" = "Vasopressor",
"shock" = "Shock index > 1",
"eGFR30" = "eGFR < 30"
)
# Replace covariate names with the new names
smd_df_pre_post_main$covariate <- covariate_name_mapping[smd_df_pre_post_main$covariate]
# Sort the dataframe by absolute SMD of Pre IPTW in descending order
smd_df_pre_post_main_sorted <- smd_df_pre_post_main %>%
filter(Imputation == "Pre IPTW") %>%
arrange(abs(smd)) %>%
mutate(covariate = factor(covariate, levels = covariate))
# Merge sorted Pre IPTW data with Post IPTW data
smd_df_pre_post_main_sorted <- smd_df_pre_post_main %>%
filter(Imputation == "Post IPTW") %>%
mutate(covariate = factor(covariate, levels = smd_df_pre_post_main_sorted$covariate)) %>%
bind_rows(smd_df_pre_post_main_sorted)
# Create Love plot
love_plot_main <- ggplot(smd_df_pre_post_main_sorted, aes(x = covariate, y = smd, color = Imputation)) +
geom_point() +
geom_hline(yintercept = 0, linetype = "dashed") +
coord_flip() +
theme_bw() +
labs(
x = "Covariate",
y = "Absolute Standardized Mean Difference"
) +
guides(color = guide_legend(title = NULL))
# Show the love plot
print(love_plot_main)
素データでの背景サマリーの作成
# Rename columns
col_continuous <- c("Child-Pugh score", "Charlson Comorbidity Index", "RBC transfusion", "Total bilirubin", "Aspartate aminotransferase", "Alanine aminotransferase", "White blood cell", "Hemoglobin", "Platelet", "Albumin", "C-reactive protein", "Prothrombin time", "Activated partial thromboplastin time")
col_factors <- c("Age", "Sex", "Barthel index", "Child-Pugh classification", "Maintenance hemodialysis", "Hepatic cancer", "Malignant tomor history", "Alcohol related disease", "Past varix rupture history", "Antiplatelet use", "Anticoagulant use", "NSAIDs use", "Corticosteroid use", "Acid blocker use", "Vasopressor","Shock index > 1", "eGFR < 30")
# Update column names in the dataframe
colnames(df)[colnames(df) %in% c("child_numl", "cci_num", "map_day0", "t_bil", "ast", "alt", "wbc", "hb", "plt", "alb", "crp", "pt", "aptt")] <- col_continuous
colnames(df)[colnames(df) %in% c("age_cate", "sex", "brthel_cate", "child_score", "hd", "hcc", "malignancy", "alocohol", "past_varix_rup", "antiplate_user_01", "anticoag_user01", "nsaids_user01", "steroid_user01", "ppi_user01", "vaso", "shock", "eGFR30")] <- col_factors
# Change order of columns
ordered_colnames <- c("Age", "Sex", "Barthel index", "Child-Pugh classification", "Child-Pugh score", "Charlson Comorbidity Index", "Maintenance hemodialysis", "Hepatic cancer", "Malignant tomor history", "Alcohol related disease", "Past varix rupture history", "Antiplatelet use", "Anticoagulant use", "NSAIDs use", "Corticosteroid use", "Acid blocker use","Vasopressor", "RBC transfusion", "Total bilirubin", "Aspartate aminotransferase", "Alanine aminotransferase", "Albumin", "White blood cell", "Hemoglobin", "Platelet", "eGFR < 30", "C-reactive protein", "Prothrombin time", "Activated partial thromboplastin time", "Shock index > 1")
# Create table
tableone <- CreateTableOne(vars = c(ordered_colnames, "antibio_exposure"), strata = "antibio_exposure", data = df, factorVars = col_factors)
# Print table
print(tableone, smd = T, missing = T, test = T, explain = T)## Stratified by antibio_exposure
## 0
## n 558
## Age (%)
## 0 322 (57.7)
## 1 144 (25.8)
## 2 90 (16.1)
## 3 2 ( 0.4)
## Sex = F (%) 141 (25.3)
## Barthel index (%)
## 0 186 (41.5)
## 1 152 (33.9)
## 2 110 (24.6)
## Child-Pugh classification (%)
## a 93 (18.8)
## b 266 (53.7)
## c 136 (27.5)
## Child-Pugh score (mean (SD)) 8.30 (1.96)
## Charlson Comorbidity Index (mean (SD)) 4.51 (1.32)
## Maintenance hemodialysis = 1 (%) 65 (11.6)
## Hepatic cancer = 1 (%) 8 ( 1.4)
## Malignant tomor history = 1 (%) 112 (20.1)
## Alcohol related disease = 1 (%) 246 (44.1)
## Past varix rupture history = 1 (%) 127 (22.8)
## Antiplatelet use = 1 (%) 8 ( 1.4)
## Anticoagulant use = 1 (%) 8 ( 1.4)
## NSAIDs use = 1 (%) 13 ( 2.3)
## Corticosteroid use = 1 (%) 2 ( 0.4)
## Acid blocker use = 1 (%) 486 (87.1)
## Vasopressor = 1 (%) 19 ( 3.4)
## RBC transfusion (mean (SD)) 2.86 (2.92)
## Total bilirubin (mean (SD)) 2.01 (1.88)
## Aspartate aminotransferase (mean (SD)) 74.97 (82.14)
## Alanine aminotransferase (mean (SD)) 38.83 (42.49)
## Albumin (mean (SD)) 2.91 (0.57)
## White blood cell (mean (SD)) 8194.22 (4242.15)
## Hemoglobin (mean (SD)) 8.69 (2.47)
## Platelet (mean (SD)) 116.49 (61.73)
## eGFR < 30 = 1 (%) 40 ( 7.3)
## C-reactive protein (mean (SD)) 0.76 (1.39)
## Prothrombin time (mean (SD)) 57.00 (16.56)
## Activated partial thromboplastin time (mean (SD)) 31.95 (12.69)
## Shock index > 1 = 1 (%) 197 (36.5)
## antibio_exposure (mean (SD)) 0.00 (0.00)
## Stratified by antibio_exposure
## 1 p
## n 232
## Age (%) 0.235
## 0 143 (61.6)
## 1 46 (19.8)
## 2 43 (18.5)
## 3 0 ( 0.0)
## Sex = F (%) 51 (22.0) 0.374
## Barthel index (%) 0.486
## 0 83 (40.5)
## 1 63 (30.7)
## 2 59 (28.8)
## Child-Pugh classification (%) 0.776
## a 42 (19.4)
## b 110 (50.9)
## c 64 (29.6)
## Child-Pugh score (mean (SD)) 8.41 (2.03) 0.508
## Charlson Comorbidity Index (mean (SD)) 4.52 (1.33) 0.917
## Maintenance hemodialysis = 1 (%) 29 (12.5) 0.829
## Hepatic cancer = 1 (%) 3 ( 1.3) 1.000
## Malignant tomor history = 1 (%) 38 (16.4) 0.269
## Alcohol related disease = 1 (%) 127 (54.7) 0.008
## Past varix rupture history = 1 (%) 63 (27.2) 0.221
## Antiplatelet use = 1 (%) 3 ( 1.3) 1.000
## Anticoagulant use = 1 (%) 5 ( 2.2) 0.675
## NSAIDs use = 1 (%) 5 ( 2.2) 1.000
## Corticosteroid use = 1 (%) 0 ( 0.0) 0.892
## Acid blocker use = 1 (%) 213 (91.8) 0.077
## Vasopressor = 1 (%) 7 ( 3.0) 0.953
## RBC transfusion (mean (SD)) 3.13 (2.98) 0.231
## Total bilirubin (mean (SD)) 2.26 (2.03) 0.098
## Aspartate aminotransferase (mean (SD)) 83.16 (106.55) 0.248
## Alanine aminotransferase (mean (SD)) 42.57 (49.00) 0.285
## Albumin (mean (SD)) 2.92 (0.64) 0.789
## White blood cell (mean (SD)) 8945.32 (5438.02) 0.039
## Hemoglobin (mean (SD)) 8.98 (2.39) 0.133
## Platelet (mean (SD)) 116.46 (84.84) 0.996
## eGFR < 30 = 1 (%) 19 ( 8.3) 0.769
## C-reactive protein (mean (SD)) 0.72 (1.57) 0.699
## Prothrombin time (mean (SD)) 57.31 (17.01) 0.816
## Activated partial thromboplastin time (mean (SD)) 34.43 (17.33) 0.035
## Shock index > 1 = 1 (%) 94 (41.4) 0.236
## antibio_exposure (mean (SD)) 1.00 (0.00) <0.001
## Stratified by antibio_exposure
## test SMD Missing
## n
## Age (%) 0.170 0.0
## 0
## 1
## 2
## 3
## Sex = F (%) 0.077 0.0
## Barthel index (%) 0.101 17.3
## 0
## 1
## 2
## Child-Pugh classification (%) 0.058 10.0
## a
## b
## c
## Child-Pugh score (mean (SD)) 0.054 12.8
## Charlson Comorbidity Index (mean (SD)) 0.008 0.0
## Maintenance hemodialysis = 1 (%) 0.026 0.0
## Hepatic cancer = 1 (%) 0.012 0.0
## Malignant tomor history = 1 (%) 0.096 0.0
## Alcohol related disease = 1 (%) 0.214 0.0
## Past varix rupture history = 1 (%) 0.102 0.0
## Antiplatelet use = 1 (%) 0.012 0.0
## Anticoagulant use = 1 (%) 0.054 0.0
## NSAIDs use = 1 (%) 0.012 0.0
## Corticosteroid use = 1 (%) 0.085 0.0
## Acid blocker use = 1 (%) 0.154 0.0
## Vasopressor = 1 (%) 0.022 0.0
## RBC transfusion (mean (SD)) 0.093 0.0
## Total bilirubin (mean (SD)) 0.129 3.2
## Aspartate aminotransferase (mean (SD)) 0.086 2.2
## Alanine aminotransferase (mean (SD)) 0.082 2.2
## Albumin (mean (SD)) 0.021 4.6
## White blood cell (mean (SD)) 0.154 1.8
## Hemoglobin (mean (SD)) 0.119 1.8
## Platelet (mean (SD)) <0.001 1.8
## eGFR < 30 = 1 (%) 0.034 1.9
## C-reactive protein (mean (SD)) 0.030 4.7
## Prothrombin time (mean (SD)) 0.018 5.8
## Activated partial thromboplastin time (mean (SD)) 0.163 11.4
## Shock index > 1 = 1 (%) 0.100 3.0
## antibio_exposure (mean (SD)) Inf 0.0tableobe作成の詳細
library(dplyr)
library(tidyr)
# Continuous variables for which we want median and IQR
continuous_vars <- c("Child-Pugh score", "Charlson Comorbidity Index", "Total bilirubin", "Aspartate aminotransferase", "Alanine aminotransferase", "Albumin", "White blood cell", "Hemoglobin", "Platelet", "C-reactive protein", "Prothrombin time", "Activated partial thromboplastin time")
# Function to calculate median and IQR
calc_median_iqr <- function(x) {
q <- quantile(x, c(0.25, 0.5, 0.75), na.rm = TRUE)
paste0("(", round(q[2], 1), ", ", round(q[1], 1), " - ", round(q[3], 1), ")")
}
# Calculate median and IQR for continuous variables
continuous_summary <- df %>%
group_by(antibio_exposure) %>%
summarise(across(all_of(continuous_vars), calc_median_iqr), .groups = "drop") %>%
pivot_longer(cols = continuous_vars, names_to = "Variable", values_to = "Value")
# Calculate percentages for categorical variables
categorical_summary <- df %>%
group_by(antibio_exposure) %>%
summarise(across(all_of(col_factors), ~mean(. == "M", na.rm = TRUE)), .groups = "drop") %>%
pivot_longer(cols = col_factors, names_to = "Variable", values_to = "Value") %>%
mutate(Value = paste0(round(Value * 100, 1), "%"))
# Combine continuous and categorical summaries
tableone <- bind_rows(continuous_summary, categorical_summary) %>%
spread(key = antibio_exposure, value = Value)
# Print table
kable(tableone, row.names = FALSE)| Variable | 0 | 1 |
|---|---|---|
| Acid blocker use | 0% | 0% |
| Activated partial thromboplastin time | (29.6, 26.4 - 33.8) | (31.1, 28.4 - 35.4) |
| Age | 0% | 0% |
| Alanine aminotransferase | (27, 19 - 42) | (30.5, 20 - 47) |
| Albumin | (2.9, 2.5 - 3.3) | (2.9, 2.5 - 3.3) |
| Alcohol related disease | 0% | 0% |
| Anticoagulant use | 0% | 0% |
| Antiplatelet use | 0% | 0% |
| Aspartate aminotransferase | (47, 31 - 83) | (54.5, 32.2 - 94.8) |
| Barthel index | 0% | 0% |
| C-reactive protein | (0.3, 0.1 - 0.8) | (0.3, 0.1 - 0.7) |
| Charlson Comorbidity Index | (4, 4 - 5) | (4, 4 - 5) |
| Child-Pugh classification | 0% | 0% |
| Child-Pugh score | (8, 7 - 10) | (8, 7 - 10) |
| Corticosteroid use | 0% | 0% |
| eGFR < 30 | 0% | 0% |
| Hemoglobin | (8.5, 6.9 - 10.2) | (9, 7.3 - 10.5) |
| Hepatic cancer | 0% | 0% |
| Maintenance hemodialysis | 0% | 0% |
| Malignant tomor history | 0% | 0% |
| NSAIDs use | 0% | 0% |
| Past varix rupture history | 0% | 0% |
| Platelet | (103, 75 - 144) | (99, 72 - 139) |
| Prothrombin time | (57, 44.8 - 68) | (59, 44.2 - 69.1) |
| Sex | 74.7% | 78% |
| Shock index > 1 | 0% | 0% |
| Total bilirubin | (1.4, 0.9 - 2.4) | (1.6, 1 - 2.9) |
| Vasopressor | 0% | 0% |
| White blood cell | (7300, 5200 - 10400) | (7720, 5900 - 10700) |
tableone作成の詳細
# Continuous variables for which we want median and IQR
continuous_vars <- c("Child-Pugh score", "Charlson Comorbidity Index", "Total bilirubin", "Aspartate aminotransferase", "Alanine aminotransferase", "Albumin", "White blood cell", "Hemoglobin", "Platelet", "C-reactive protein", "Prothrombin time", "Activated partial thromboplastin time", "RBC transfusion")
# Function to calculate median and IQR
calc_median_iqr <- function(x) {
q <- quantile(x, c(0.25, 0.5, 0.75), na.rm = TRUE)
paste0("(", round(q[2], 1), ", ", round(q[1], 1), " - ", round(q[3], 1), ")")
}
# Calculate median and IQR for continuous variables
continuous_summary <- df %>%
group_by(antibio_exposure) %>%
summarise(across(all_of(continuous_vars), calc_median_iqr)) %>%
pivot_longer(cols = continuous_vars, names_to = "Variable", values_to = "Value")
# Calculate percentages for categorical variables
categorical_summary <- df %>%
group_by(antibio_exposure) %>%
summarise(across(all_of(col_factors), ~mean(. == "M", na.rm = TRUE))) %>%
pivot_longer(cols = col_factors, names_to = "Variable", values_to = "Value") %>%
mutate(Value = paste0(round(Value * 100, 1), "%"))
# Combine continuous and categorical summaries
tableone <- bind_rows(continuous_summary, categorical_summary) %>%
spread(key = antibio_exposure, value = Value)
# Print table
kable(tableone, row.names = FALSE)| Variable | 0 | 1 |
|---|---|---|
| Acid blocker use | 0% | 0% |
| Activated partial thromboplastin time | (29.6, 26.4 - 33.8) | (31.1, 28.4 - 35.4) |
| Age | 0% | 0% |
| Alanine aminotransferase | (27, 19 - 42) | (30.5, 20 - 47) |
| Albumin | (2.9, 2.5 - 3.3) | (2.9, 2.5 - 3.3) |
| Alcohol related disease | 0% | 0% |
| Anticoagulant use | 0% | 0% |
| Antiplatelet use | 0% | 0% |
| Aspartate aminotransferase | (47, 31 - 83) | (54.5, 32.2 - 94.8) |
| Barthel index | 0% | 0% |
| C-reactive protein | (0.3, 0.1 - 0.8) | (0.3, 0.1 - 0.7) |
| Charlson Comorbidity Index | (4, 4 - 5) | (4, 4 - 5) |
| Child-Pugh classification | 0% | 0% |
| Child-Pugh score | (8, 7 - 10) | (8, 7 - 10) |
| Corticosteroid use | 0% | 0% |
| eGFR < 30 | 0% | 0% |
| Hemoglobin | (8.5, 6.9 - 10.2) | (9, 7.3 - 10.5) |
| Hepatic cancer | 0% | 0% |
| Maintenance hemodialysis | 0% | 0% |
| Malignant tomor history | 0% | 0% |
| NSAIDs use | 0% | 0% |
| Past varix rupture history | 0% | 0% |
| Platelet | (103, 75 - 144) | (99, 72 - 139) |
| Prothrombin time | (57, 44.8 - 68) | (59, 44.2 - 69.1) |
| RBC transfusion | (2.5, 0 - 4) | (4, 0 - 4) |
| Sex | 74.7% | 78% |
| Shock index > 1 | 0% | 0% |
| Total bilirubin | (1.4, 0.9 - 2.4) | (1.6, 1 - 2.9) |
| Vasopressor | 0% | 0% |
| White blood cell | (7300, 5200 - 10400) | (7720, 5900 - 10700) |