はじめに

formula関数とsurv_fit関数の関係をまとめる。

`surv_fit`関数

カプランマイヤー曲線を計算して可視化する際、一般的にはsurvivalパッケージのsurvfit関数を利用し、survminerパッケージのggsurvplot関数で可視化することが多いはずで、surv_fit関数は使うことはない。下記のようにすれば問題がなく計算できる。

library("survival")
library("survminer")

fit <- survfit(Surv(time, status) ~ sex, data = colon)
ggsurvplot(fit)

文字列で作成した式をformula型に変換すれば、文字列から計算もできると思ったが、これはできない。Error in x$formula : object of type 'symbol' is not subsettableと表示さてしまい可視化できない。ggsurvplot関数に渡される場合に、symbolで渡されるのがよろしくない模様。

# f <- as.formula('Surv(time, status) ~ sex')
# fit2 <- survfit(f, data = colon)
# ggsurvplot(fit2)

# Error in x$formula : object of type 'symbol' is not subsettable

このようなケースでsurv_fit関数を利用する。この関数であれば問題なく可視化できる。surv_fit関数標準のsurvfit関数のラッパー関数で、こちらであれば機能する。おそらくドキュメントに記載されている、このあたりが関わっていると思われる。

If the formula names are not available, the variables in the formulas are extracted and used to build the name of survfit object.

ggsurvplot関数側の細かい詳細を確認してないので、正確かどうかはわからない。

f <- as.formula('Surv(time, status) ~ sex')
fit2 <- surv_fit(f, data = colon)
ggsurvplot(fit2)

このように考えて、for-loopで複数のカプランマイヤー曲線を計算して可視化することもできる。

treats <- c('rx', 'sex')
res <- vector(mode = 'list', length = length(treats))
for (i in seq_along(treats)) {
  treat <- treats[i]
  survie <- as.formula(paste0('Surv(time, status) ~ ', treat))
  fit <- surv_fit(survie, data = colon)
  res[[i]] <- ggsurvplot(fit, colon)
}
res[1]

## [[1]]

res[2]

## [[1]]

surv_fit関数は他にも便利な機能があり、可視化だけでなく、複数の式を一度に計算できる。これはsurvfit関数ではできない。

formulas <- list(
 sex = Surv(time, status) ~ sex,
 rx = Surv(time, status) ~ rx
)

# Fit survival curves for each formula
fit <- surv_fit(formulas, data = colon)
surv_pvalue(fit)

## $`colon::sex`
##   variable      pval   method pval.txt
## 1      sex 0.6107936 Log-rank p = 0.61
## 
## $`colon::rx`
##   variable         pval   method   pval.txt
## 1       rx 4.990735e-08 Log-rank p < 0.0001

使いみちはわからないが、データ側を変えることができる。

fit <- surv_fit(Surv(time, status) ~ sex,
               data = list(colon, lung))
surv_pvalue(fit)

## $`colon::sex`
##   variable      pval   method pval.txt
## 1      sex 0.6107936 Log-rank p = 0.61
## 
## $`lung::sex`
##   variable        pval   method   pval.txt
## 1      sex 0.001311165 Log-rank p = 0.0013

他にもグループ化に対応している。ここではrxでグループ化している。

fit <- surv_fit(Surv(time, status) ~ sex,
               data = colon, group.by = "rx")
ggsurvplot(fit, colon)

## $`rx.Obs::sex`

## 
## $`rx.Lev::sex`

## 
## $`rx.Lev+5FU::sex`

## 
## attr(,"class")
## [1] "list"            "ggsurvplot_list"

formula関数とsurv_fit関数の話

はじめに

surv_fit関数

`surv_fit`関数