R でログランク検定に必要なサンプルサイズを計算する方法

2018年8月25日2024年10月13日

ログランク検定のサンプルサイズ計算を R で行う方法

＞＞もう統計で悩むのは終わりにしませんか？

↑期間・数量限定で無料プレゼント中！

ログランク検定はどうやるか？

ログランク検定は、時間経過とともにイベントが起きていくデータの群間比較をする方法。

患者さんの死亡をイベントとしたデータ、病気が再発することをイベントとしたデータ、病気が発症することをイベントとしたデータ、などが扱える。

病気ばかりではなく、時間経過とともに起きてくる事柄なら、なんでも扱える。

サンプルサイズ計算の前に、ログランク検定はどうやるか？を見てみる。

R の survival パッケージの survdiff() を使う。

survival パッケージはインストール済みなので、library()で呼び出すだけでOK。

ovarian という卵巣がんの治療後生存時間の解析データ。

futime が生存時間、fustat が死亡発生のデータ。

rx が治療群のデータ。

library(survival)
survdiff(Surv(futime, fustat) ~ rx,data=ovarian)

結果はp=0.3で統計学的有意に異ならなかった。

> survdiff(Surv(futime, fustat) ~ rx,data=ovarian)
Call:
survdiff(formula = Surv(futime, fustat) ~ rx, data = ovarian)
N Observed Expected (O-E)^2/E (O-E)^2/V
rx=1 13        7     5.23     0.596      1.06
rx=2 13        5     6.77     0.461      1.06
Chisq= 1.1  on 1 degrees of freedom, p= 0.3

もう一つの例は、MASSパッケージの中の、gehanというデータ。

白血病の患者さんの寛解から再燃までのデータ。

timeが再燃までの時間、censが再燃したかどうか、treatが6－メルカプトプリン（6-MP）で治療したか、対照だったか。

library(MASS)
survdiff(Surv(time, cens) ~ treat, data = gehan)

結果は、統計学的有意に6-MP治療群のほうが再燃が少なかった。

> library(MASS)
> survdiff(Surv(time, cens) ~ treat, data = gehan)
Call:
survdiff(formula = Surv(time, cens) ~ treat, data = gehan)
N Observed Expected (O-E)^2/E (O-E)^2/V
treat=6-MP    21        9     19.3      5.46      16.8
treat=control 21       21     10.7      9.77      16.8
Chisq= 16.8  on 1 degrees of freedom, p= 4e-05

ログランク検定のサンプルサイズ計算はどうやるか？

ログランク検定のサンプルサイズ計算は、実はCoxの比例ハザードモデルのサンプルサイズ計算を使う。

あわせて読みたい

R と EZR で生存時間解析のサンプルサイズ計算を行う方法生存時間解析のサンプルサイズ計算の方法 Cox の比例ハザードモデルを使う前提の計算生存時間解析のサンプルサイズ計算その 1 R でスクリプトを書いてみた。 S0がコン…

ログランク検定もCoxも生存時間解析と言いながら、イベントが起きる順番しか使っていない。

あとはイベントが起きる比を使う。

サンプルサイズ計算はどちらにも流用できるということだ。

sample.size.cox <- function(S1,S0,alternative="two.sided",power=.8,
sig.level=.05){
alternative <- match.arg(alternative)
tside <- switch(alternative, one.sided=1, two.sided=2)
beta <- log(log(S1)/log(S0))
Za <- qnorm(sig.level/tside, lower.tail=FALSE)
Zb <- qnorm(power)
d <- ((Za+Zb)*(1+exp(beta))/(1-exp(beta)))^2
n <- d/(((1-S1)+(1-S0))/2)
NOTE <- "n is total number of subjects in *both* of groups"
METHOD <- "Sample Size Calculation of Cox Propotional Hazard Model"
structure(
list("Total number of death" = d,
"Total number of subjects" = n,
"Surv. rate of trt. arm" = S1,
"Surv. rate of ctl. arm" = S0,
"Hazard ratio" = exp(beta),
sig.level = sig.level,
power = power,
alternative = alternative,
note = NOTE,
method = METHOD),
class = "power.htest")
}

先ほどのovarianの結果をもとにサンプルサイズ計算をすると、両群合わせて151例の死亡、326例の患者さんが必要だった。

各群13例の試験では統計学的有意にならなくて当然だ。

> sample.size.cox(S1=8/13, S0=6/13)
Sample Size Calculation of Cox Propotional Hazard Model
Total number of death = 150.2536
Total number of subjects = 325.5495
Surv. rate of trt. arm = 0.6153846
Surv. rate of ctl. arm = 0.4615385
Hazard ratio = 0.6279283
sig.level = 0.05
power = 0.8
alternative = two.sided
NOTE: n is total number of subjects in *both* of groups

gehanデータは、全部で8例の再燃症例があればよく、必要なサンプルサイズは二群全体で11例と計算された。

> sample.size.cox(S1=(21-9)/21, S0=(21-21)/21)
Sample Size Calculation of Cox Propotional Hazard Model
Total number of death = 7.84888
Total number of subjects = 10.98843
Surv. rate of trt. arm = 0.5714286
Surv. rate of ctl. arm = 0
Hazard ratio = 0
sig.level = 0.05
power = 0.8
alternative = two.sided
NOTE: n is total number of subjects in *both* of groups