はじめに
生存時間分析において、時間共変量型のモデルを組むときに、時間変化する変数を縦持ちにするのにsurvival::tmergeが便利です。
その中で、可変長引数に当たる部分の変数名を動的に変更したいと思いました。
しかし、なかなか苦戦したため、ここに残しておきます。
使用したコードはmako_sharkさんの記事にあったコードを使用させていただきました。
問題
今回利用するデータを生成します
test <- data.frame(id = 1:400,
baseline= 0,
smoke = c(rbinom(400, 1, prob=0.3)),
treat = c(rep(0, 200), rep(1, 200)),
censor = sample(10:180, 400, replace = T),
outcome = c(rbinom(200, 1, prob=0.1), rbinom(200, 1, prob=0.2)))
test1 <- test %>%
mutate(visit1 = 0,
smoke1 = c(rbinom(400, 1, prob=0.3)),
exercise1 = c(rpois(400,2))) %>%
mutate(visit2 = visit1 + sample(21:42, 400, replace = T),
smoke2 = c(rbinom(400, 1, prob=0.3)),
exercise2 = c(rpois(400,2))) %>%
mutate(visit3 = visit2 + sample(21:42, 400, replace = T),
smoke3 = c(rbinom(400, 1, prob=0.3)),
exercise3 = c(rpois(400,2))) %>%
mutate(visit1 = ifelse(visit1<censor, visit1, NA),
smoke1 = ifelse(visit1<censor, smoke1, NA),
exercise1 = ifelse(visit1<censor, exercise1, NA),
visit2 = ifelse(visit2<censor, visit2, NA),
smoke2 = ifelse(visit2<censor, smoke2, NA),
exercise2 = ifelse(visit2<censor, exercise2, NA),
visit3 = ifelse(visit2<censor, visit3, NA),
smoke3 = ifelse(visit2<censor, smoke3, NA),
exercise3 = ifelse(visit2<censor, exercise3, NA))
head(test1)
# id baseline smoke treat censor outcome visit1 smoke1 exercise1 visit2 smoke2 exercise2 visit3 smoke3 exercise3
#1 1 0 0 0 90 0 0 1 3 41 0 7 80 0 3
#2 2 0 1 0 141 0 0 1 0 24 0 1 47 1 2
#3 3 0 0 0 138 0 0 0 5 36 0 3 66 0 1
#4 4 0 1 0 89 1 0 0 3 38 0 1 65 0 2
#5 5 0 0 0 148 0 0 0 4 26 1 2 51 0 4
#6 6 0 0 0 18 0 0 0 3 NA NA NA NA NA NA
survival::tmergeの使い方は、次のような形です。
詳しくはリファレンスを。
library(survival)
test2 = tmerge(data1 = test, data2=test1, id=id, tstop=censor, event = event(censor, outcome))
test2 = tmerge(data1 = test2, data2=test1, id=id, smoke=tdc(visit1, smoke1))
test2 = tmerge(data1 = test2, data2=test1, id=id, smoke=tdc(visit2, smoke2))
test2 = tmerge(data1 = test2, data2=test1, id=id, smoke=tdc(visit3, smoke3))
head(test2)
# id baseline smoke treat censor outcome tstart tstop event
#1 1 0 1 0 90 0 0 41 0
#2 1 0 0 0 90 0 41 80 0
#3 1 0 0 0 90 0 80 90 0
#4 2 0 1 0 141 0 0 24 0
#5 2 0 0 0 141 0 24 47 0
#6 2 0 1 0 141 0 47 141 0
縦持ちにした変数名=tdc(期間, 変数)のように引数に記述していくことで、この期間ごとに別々のデータとしてい縦持ちにすることができます。
この数名は可変長引数であり、複数指定することができます。
この時、カラム名(変数名)をlist等で保持して、ループを回して変換したいと思いました。
これができると、多数のカラムの追加や共変量の時間変化が多数ある場合にも簡単に対応できると思いました。
今までの経験から動的に変数名を定義するような記述しましたが、上手くできませんでした。
解決策
いろいろ頑張りましたが、最終的には私のお師匠さんがr-wakarangに質問を投げて解決しました。
#library(survival)
library(rlang)
feature_list <- list("smoke","exercise")
test2 = tmerge(data1 = test1 %>% select(id,treat,censor,outcome),
data2 = test1,
id = id,
tstop = censor,
event = event(censor, outcome))
for(feature in feature_list){
for(i in 1:3){
# expressionのリストを組み立てる
extra_args <- exprs(!! feature := tdc(get(paste0("visit",i)), get(paste0(feature,i))))
# コード文字列を生成
tmerge_quo <- quo(tmerge(data1 = test2, data2=test1, id = id, !!! extra_args))
test2 <- eval_tidy(tmerge_quo)
}
}
head(test2,9)
# id treat censor outcome tstart tstop event smoke exercise
#1 1 0 90 0 0 26 0 1 3
#2 1 0 90 0 26 60 0 1 3
#3 1 0 90 0 60 90 0 0 1
#4 2 0 141 0 0 23 0 1 2
#5 2 0 141 0 23 46 0 0 4
#6 2 0 141 0 46 141 0 0 2
#7 3 0 138 0 0 33 0 0 1
#8 3 0 138 0 33 56 0 0 0
#9 3 0 138 0 56 138 0 0 2
私一人だとこの回答にはたどり着きませんでした。
わからない時には賢人に頼るものですね。