1. Table 1 만들기

아래처럼 데이터 읽은후, Q_로 시작하는 변수를 범주형으로 바꾼 후, 연도별 기술통계량을 Table로 나태내어라 (jstable::CreateTableOneJS 이용)

a <- fread("https://raw.githubusercontent.com/jinseob2kim/lecture-snuhlab/master/data/example_g1e.csv")
vars.factor <- grep("Q_", names(a), value = T)

#for (v in vars.factor){
#  a[[v]] <- factor(a[[v]])
#}

a[, (vars.factor) := lapply(.SD, factor), .SDcols = vars.factor]

vars.tb1 <- names(a)[-c(1:3)]

tb1 <- CreateTableOneJS(vars = vars.tb1, strata = "EXMD_BZ_YYYY", data = a)
#knitr::kable(tb1$table, caption = tb1$caption)


#install.packages("DT")
DT::datatable(tb1$table, caption = tb1$caption)

회귀분석

선형회귀분석

time ~ rx + age + sex 선형회귀 실행 후 Table로 나타내어라(jstable::glmshow.display 이용)

library(survival)
#colon
res.reg <- glm(time ~ rx + age + sex, data = colon)

tb.reg <- glmshow.display(res.reg)
knitr::kable(tb.reg$table, caption = tb.reg$first.line)
Linear regression predicting time
crude coeff.(95%CI) crude P value adj. coeff.(95%CI) adj. P value
rx: ref.=Obs NA NA NA NA
Lev 24.66 (-79.49,128.82) 0.643 22.98 (-81.3,127.27) 0.666
Lev+5FU 271.07 (166.41,375.74) < 0.001 273.05 (168.19,377.9) < 0.001
age 0.38 (-3.22,3.99) 0.835 0.37 (-3.21,3.95) 0.84
sex 13.93 (-72.26,100.12) 0.751 32.67 (-53.21,118.55) 0.456

로지스틱

status ~ rx + age + sex 로지스틱회귀 실행 후 Table로 나타내어라(jstable::glmshow.display 이용)

library(survival)
#colon
res.logistic <- glm(status ~ rx + age + sex, data = colon, family = binomial)

tb.logistic <- glmshow.display(res.logistic)
knitr::kable(tb.logistic$table, caption = tb.logistic$first.line)
Logistic regression predicting status
crude OR.(95%CI) crude P value adj. OR.(95%CI) adj. P value
rx: ref.=Obs NA NA NA NA
Lev 0.96 (0.77,1.2) 0.709 0.96 (0.77,1.2) 0.747
Lev+5FU 0.55 (0.44,0.68) < 0.001 0.54 (0.43,0.68) < 0.001
age 1 (0.99,1) 0.296 1 (0.99,1) 0.294
sex 0.97 (0.81,1.17) 0.758 0.93 (0.77,1.12) 0.454

생존분석

Cox

Surv(time, status) ~ rx + age + sex 실행 후 Table로 나태내어라(jstable::cox2.display 이용)

#library(survival)
res.cox <- coxph(Surv(time, status) ~ rx + age + sex, data = colon, model = T)
tb.cox <- cox2.display(res.cox)
knitr::kable(tb.cox$table, caption = tb.cox$caption)
Cox model on time (‘time’) to event (‘status’)
crude HR(95%CI) crude P value adj. HR(95%CI) adj. P value
rx: ref.=Obs NA NA NA NA
Lev 0.98 (0.84,1.14) 0.786 0.98 (0.84,1.14) 0.811
Lev+5FU 0.64 (0.55,0.76) < 0.001 0.64 (0.55,0.76) < 0.001
age 1 (0.99,1) 0.382 1 (0.99,1) 0.468
sex 0.97 (0.85,1.1) 0.61 0.95 (0.84,1.08) 0.446

kaplan-meier

Surv(time, status) ~ rx 실행 후 그림으로 나타내어라(jskm::jskm 이용)

#install.packages("jskm")

res.km <- survfit(Surv(time, status) ~ rx, data = colon)

jskm::jskm(res.km, table = T, pval = T, label.nrisk = "No. at risk", size.label.nrisk = 8, 
     xlabs = "Time(Day)", ylabs = "Survival", ystratalabs = c("Obs", "Lev", "Lev + 5FU"), ystrataname = "rx",
     marks = F, timeby = 365, xlims = c(0, 3000), ylims = c(0.25, 1), showpercent = T)