生存時間解析のサンプルサイズ設計 ― 必要イベント数・Schoenfeld式・脱落を考慮した症例数をRで実装する ―
この記事でわかること
・生存時間解析では検出力を決めるのが「症例数」ではなく「イベント数（死亡や増悪などの発生件数）」である理由
・Schoenfeld式（ショーンフェルド式）を使って、想定ハザード比から必要イベント数を求める手順
・必要イベント数を必要症例数へ変換する考え方 ― 集積期間・追跡期間・脱落（追跡不能）の織り込み方
・gsDesign パッケージの
nSurv() と powerSurvEpi パッケージによるR実装、そして手計算との突き合わせ・組み上げた設計をシミュレーションで検算し、狙った検出力が本当に出ているかを確かめる方法
はじめに
第III相のがん試験や心血管アウトカム試験のように、主要評価項目が「全生存期間（OS）」「無増悪生存期間（PFS）」といった時間依存エンドポイントである試験では、サンプルサイズ設計の考え方が大きく変わります。連続量や二値の比較では「1群何例必要か」を直接計算しますが、生存時間解析ではまず「何イベント必要か」を決め、そのイベント数を集めるために何例組み入れるかを逆算するという、発想の転換が必要になります。
この違いを押さえずに、平均値の差や割合の差と同じ感覚で症例数だけを議論すると、「予定どおり症例は集まったのに、解析時点でイベントが足りず検出力が出ない」という事態が起こります。実際、臨床試験の中間解析や解析時期（データカットオフ）の議論が「何月時点で何イベント貯まるか」を軸に進むのは、この構造があるからです。連続量・二値の場合の基本的な考え方については、検出力分析（Power Analysis）入門 ― サンプルサイズ設計に欠かせない基礎知識をSAS・Rで実装 ―で整理していますので、本記事と読み比べていただくと違いが際立つはずです。
本記事は、製薬企業・CROの生物統計担当者や、統計解析を学ぶ社会人・大学生の方を対象に、検定・p値・検出力の基礎は知っているが生存時間解析の症例数設計は初めて、という読者に向けて書いています。Schoenfeld式による必要イベント数の算出から、イベント数を症例数へ変換する計算、gsDesign・powerSurvEpi によるR実装、そしてシミュレーションによる検算までを、実際にRを実行した出力とともに一気通貫でたどります。
生存時間解析では「症例数」ではなく「イベント数」が検出力を決める
情報量はイベント数にほぼ比例する
ログランク検定やCox比例ハザードモデルの部分尤度では、ハザード比の推定精度（統計学でいうFisher情報量）が、観察されたイベント数にほぼ比例するという性質があります。直感的に言えば、部分尤度はイベントが起きた各時点で「そのとき生き残っていた人（リスク集合）の中で、誰がイベントを起こしたか」を比較して情報を積み上げていく仕組みです。したがって情報を生むのはイベントが起きた瞬間だけであり、最後まで打ち切り（censoring）で終わった被験者は、リスク集合に貢献こそすれ、それ自体はほとんど情報を持ちません。
このため「1,000例組み入れたのにイベントが50件しか発生していない」という試験は、症例数がどれだけ多くても検出力が出ません。逆に、予後の悪い集団を対象にしてイベントが早く貯まる試験であれば、少ない症例数でも必要な情報量に到達できます。ログランク検定そのものの仕組みはログランク(log-rank)検定を徹底解説を、部分尤度とハザード比の解釈はCox比例ハザードモデル入門〜数式から実務応用まで〜をあわせてご確認ください。
Schoenfeld式で必要イベント数を求める
必要イベント数の算出に最もよく使われるのが、次のSchoenfeld式です。
\[ d = \frac{(z_{1-\alpha/2} + z_{1-\beta})^2}{p(1-p)\,(\log \mathrm{HR})^2} \]
ここで \( d \) は必要イベント数（両群合計）、\( p \) は実験群への割付割合、\( \mathrm{HR} \) は想定するハザード比、\( z_{1-\alpha/2} \) は両側有意水準 \( \alpha \) に対応する標準正規分布の分位点、\( z_{1-\beta} \) は検出力 \( 1-\beta \) に対応する分位点です。1:1割付なら \( p = 0.5 \) なので分母の \( p(1-p) = 0.25 \) となり、教科書でよく見る \( d = 4(z_{1-\alpha/2} + z_{1-\beta})^2 / (\log \mathrm{HR})^2 \) の形に一致します。式に症例数 \( n \) が一切登場しないことが、生存時間解析の設計を象徴しています。
この式をそのままRの関数にしてみます。
schoenfeld_events <- function(hr, alpha = 0.05, power = 0.8, p1 = 0.5) {
za <- qnorm(1 - alpha / 2)
zb <- qnorm(power)
(za + zb)^2 / (p1 * (1 - p1) * (log(hr))^2)
}
HR = 0.7、両側 \( \alpha = 0.05 \) で、検出力80%・90%それぞれの必要イベント数を求めてみます。
cat("power=80%:", schoenfeld_events(0.7, power = 0.8), "\n")
cat("power=90%:", schoenfeld_events(0.7, power = 0.9), "\n")
cat("ceiling 80%:", ceiling(schoenfeld_events(0.7, power = 0.8)), "\n")
cat("ceiling 90%:", ceiling(schoenfeld_events(0.7, power = 0.9)), "\n")
power=80%: 246.7871
power=90%: 330.3779
ceiling 80%: 247
ceiling 90%: 331
HR = 0.7（ハザードが3割低下）を検出するには、両側5%・検出力80%で 247イベント、検出力90%なら 331イベント が必要です。検出力を80%から90%へ引き上げるだけで必要イベント数は約1.34倍に増えます。イベント数は切り上げて扱うため、246.7871 は247、330.3779 は331とします。実務では、この「247イベント」がプロトコルに書かれる目標イベント数となり、データカットオフのタイミングを規定する基準にもなります。
HRの仮定が必要イベント数を左右する
同じ計算をHRごとに繰り返すと、想定HRのわずかな違いが必要イベント数を大きく動かすことがわかります。
| 想定ハザード比 HR | 必要イベント数（検出力80%） | 必要イベント数（検出力90%） |
|---|---|---|
| 0.60 | 121 | 162 |
| 0.65 | 170 | 227 |
| 0.70 | 247 | 331 |
| 0.75 | 380 | 508 |
| 0.80 | 631 | 845 |
検出力80%で比べると、HR = 0.60 なら121イベントで済むのに対し、HR = 0.80 では631イベントと5倍以上に膨らみます。Schoenfeld式の分母が \( (\log \mathrm{HR})^2 \) であるため、HRが1に近づくほど \( \log \mathrm{HR} \) がゼロに近づき、必要イベント数が急激に発散するのです。
設計段階で「効きすぎるHR」を仮定すると、試験は高い確率で失敗します。HR = 0.60 を前提に121イベントで設計した試験で、実際の効果がHR = 0.75 程度だった場合、本来必要な380イベントの3分の1しか情報が集まらず、真に有効な薬でも有意差を示せません。想定HRは「期待値」ではなく「これを下回っても検出したい臨床的に意味のある最小の効果」として、既存の類似試験やメタ解析の実績値を根拠に、やや保守的に設定してください。
Freedman式との比較
必要イベント数の公式にはSchoenfeld式のほかにFreedman式もあり、こちらは \( (\mathrm{HR}+1)^2/(\mathrm{HR}-1)^2 \) を用いる点が異なります。同じ条件で両者を計算すると次のようになります。
freedman_events <- function(hr, alpha = 0.05, power = 0.8) {
za <- qnorm(1 - alpha / 2)
zb <- qnorm(power)
((1 + hr) / (1 - hr))^2 * (za + zb)^2
}
cat("Freedman power=80%:", freedman_events(0.7), "->", ceiling(freedman_events(0.7)), "\n")
cat("Freedman power=90%:", freedman_events(0.7, power = 0.9), "->",
ceiling(freedman_events(0.7, power = 0.9)), "\n")
Freedman power=80%: 252.0362 -> 253
Freedman power=90%: 337.405 -> 338
| 想定ハザード比 HR | Schoenfeld式（検出力80%） | Freedman式（検出力80%） |
|---|---|---|
| 0.60 | 121 | 126 |
| 0.65 | 170 | 175 |
| 0.70 | 247 | 253 |
| 0.75 | 380 | 385 |
| 0.80 | 631 | 636 |
HR = 0.7 で比べると、Schoenfeld式が247イベント、Freedman式が253イベントと、実務上はSchoenfeld式のほうがやや小さめの値になります。差はHR = 0.60 で121対126、HR = 0.80 で631対636と、いずれも数イベント程度にとどまり、どちらを採用しても設計の骨格は変わりません。ただしプロトコルや統計解析計画書（SAP）にはどちらの式を用いたかを明記し、規制当局からの照会や再計算に備えて再現できる形で残しておくことが重要です。保守的に設計したい場合は、大きめに出るFreedman式を採用するのも一案です。
必要イベント数から必要症例数を求める ― 集積・追跡・脱落の考慮
Schoenfeld式で必要イベント数 \( d \) が求まっても、それだけでは「何人組み入れればよいか」は決まりません。イベント数は「組み入れた人数」に「その人が試験期間中にイベントを起こして観察される確率」を掛けたものだからです。したがって必要症例数 \( n \) は、必要イベント数を試験期間中のイベント観察確率 \( P(\text{event}) \) で割って求めます。
\[ n = \frac{d}{P(\text{event})} \]
つまり、イベントが起きにくい集団ほど（\( P(\text{event}) \) が小さいほど）、同じイベント数を貯めるために多くの症例が必要になります。
イベント確率をどう計算するか
被験者が試験期間中にイベントを起こす確率は、その人がいつ組み入れられたかによって変わります。試験開始直後に入った人は長く追跡されますが、集積期間の終わり近くに入った人は追加追跡期間しか観察されません。そこで実務では「集積期間 \( A \) の間に一様に組み入れられる（一様集積）」「生存時間は指数分布に従う」と仮定し、次の式でイベント確率を求めます。
\[ P(\text{event}) = 1 – \frac{e^{-\lambda f} – e^{-\lambda (A+f)}}{\lambda A} \]
ここで \( \lambda \) はハザード（単位時間あたりのイベント発生率）、\( A \) は集積期間、\( f \) は最後の被験者を組み入れてからの追加追跡期間で、試験全体の期間は \( A+f \) となります。指数分布では中央値とハザードが \( \lambda = \log(2)/\text{中央値} \) という関係で結ばれるため、「対照群の生存期間中央値」という臨床的に議論しやすい量から \( \lambda \) を決められるのが利点です。指数分布の性質やWeibull分布への拡張はパラメトリック生存時間解析とは ― 指数・Weibull・AFTモデルの理論からR/SAS実装まで徹底解説 ―で解説しています。
Rで症例数まで計算する
対照群の生存期間中央値12か月、HR=0.7、集積期間24か月、追加追跡12か月（試験全体36か月）という設定で、必要症例数を求めてみます。
p_event <- function(lambda, A, f) {
1 - (exp(-lambda * f) - exp(-lambda * (A + f))) / (lambda * A)
}
lambda_C <- log(2) / 12
lambda_E <- lambda_C * 0.7
A <- 24; f <- 12
pC <- p_event(lambda_C, A, f)
pE <- p_event(lambda_E, A, f)
p_bar <- (pC + pE) / 2
d <- 247
N <- d / p_bar
cat("lambda_C:", lambda_C, " lambda_E:", lambda_E, "\n")
cat("median experimental (months):", log(2) / lambda_E, "\n")
cat("P(event) control:", pC, "\n")
cat("P(event) experimental:", pE, "\n")
cat("average P(event):", p_bar, "\n")
cat("required total N:", N, "-> ceiling:", ceiling(N), " per arm:", ceiling(N / 2), "\n")
cat("with 10% dropout: total N =", ceiling(N / 0.9), " per arm:", ceiling(N / 0.9 / 2), "\n")
lambda_C: 0.05776227 lambda_E: 0.04043359
median experimental (months): 17.14286
P(event) control: 0.7294947
P(event) experimental: 0.6060268
average P(event): 0.6677607
required total N: 369.893 -> ceiling: 370 per arm: 185
with 10% dropout: total N = 411 per arm: 206
対照群のハザードは 0.05776227、HR=0.7 を掛けた実験群は 0.04043359 で、実験群の生存期間中央値は 17.14286 か月に延びます。試験期間中にイベントが観察される確率は対照群 0.7294947、実験群 0.6060268、両群平均で 0.6677607 でした。必要イベント数247をこの平均で割ると 369.893、切り上げて総症例数370例（1群185例）となります。
脱落（追跡不能）がある場合はさらに割り増しが必要です。脱落した被験者はイベントが観察されないまま試験から抜けるため、実質的に \( P(\text{event}) \) が小さくなるのと同じ効果を持ち、\( n = d / P(\text{event}) \) の分母を押し下げて必要症例数を増やします。脱落率10%を見込むなら、割り増し前の 369.893 を 0.9 で割って切り上げ、411例（1群206例）まで積み増します。プロトコルに書く症例数はこの割り増し後の数字です。
対照群の中央値の仮定はどれくらい効くか
必要症例数は「対照群のイベントがどれくらい起きるか」という仮定に強く依存します。対照群の生存期間中央値を12・15・18か月と変えたときの結果は次のとおりです。
| 対照群の生存期間中央値 | 平均イベント確率 | 必要総症例数 |
|---|---|---|
| 12か月 | 0.6678 | 370例 |
| 15か月 | 0.5906 | 419例 |
| 18か月 | 0.528 | 468例 |
対照群のイベントが起きにくいほど（＝生存期間中央値が長いほど）、必要症例数は増えます。中央値12か月なら370例で済むところ、15か月では419例、18か月では468例が必要になります。中央値の仮定を3か月動かすごとに約50例、12か月と18か月を比べれば約100例の差が生じる計算です。過去試験のデータや標準治療の進歩を踏まえて対照群の中央値を保守的（長め）に見積もり、必ず感度分析の表をプロトコルや統計解析計画書に残しておきましょう。
症例数と試験期間はトレードオフ
もう一つ効くのが追加追跡期間です。追跡を延ばせば1人あたりのイベント観察確率が上がるため、必要症例数は減ります。
| 追加追跡期間 | 平均イベント確率 | 必要総症例数 | 試験全体の期間 |
|---|---|---|---|
| 0か月 | 0.4095 | 604例 | 24か月 |
| 6か月 | 0.5577 | 443例 | 30か月 |
| 12か月 | 0.6678 | 370例 | 36か月 |
| 24か月 | 0.8111 | 305例 | 48か月 |
症例数と試験期間はトレードオフです。追加追跡なしで組み入れ終了と同時に解析するなら604例が必要ですが、12か月追跡すれば370例、24か月追跡すれば305例まで減らせます。ただし試験期間は24か月から48か月へと倍になります。「症例数を減らしたい」という要望には常に「では試験を何か月延ばせますか」という問いが対になっており、施設のリクルート能力・開発計画・競合状況を踏まえた総合判断が必要です。試験デザイン全体の選択肢は臨床試験でよく使われる試験デザイン完全ガイドも参考にしてください。
Rによるサンプルサイズ設計の実装
ここまでは式の意味を理解するために手計算に近い形で進めてきましたが、実務では専用パッケージを使うのが安全です。生存時間解析のサンプルサイズ設計では gsDesign::nSurv() と powerSurvEpi::ssizeCT.default() がよく使われます。
gsDesign::nSurv による設計
gsDesign パッケージの nSurv() は、集積期間・追跡期間・脱落（追跡不能）を陽に指定してLachin-Foulkes法で必要症例数とイベント数を同時に返します。年5%の追跡不能を仮定して、先ほどと同じ設計を計算します。
library(gsDesign)
eta <- -log(1 - 0.05) / 12 # 年5%の追跡不能
cat("dropout hazard eta (per month):", eta, "\n")
nSurv(lambdaC = log(2) / 12, hr = 0.7, eta = eta,
T = 36, minfup = 12, ratio = 1,
alpha = 0.025, beta = 0.2, sided = 1)
dropout hazard eta (per month): 0.004274441
nSurv fixed-design summary (method=LachinFoulkes; target=Accrual rate)
HR=0.700 vs HR0=1.000 | alpha=0.025 (sided=1) | power=80.0%
N=383.8 subjects | D=245.9 events | T=36.0 study duration | accrual=24.0 Accrual duration | minfup=12.0 minimum follow-up | ratio=1 randomization ratio (experimental/control)
Key inputs (names preserved):
desc item value input
Accrual rate(s) gamma 15.992 1
Accrual rate duration(s) R 24 12
Control hazard rate(s) lambdaC 0.058 0.058
Control dropout rate(s) eta 0.004 0.004
Experimental dropout rate(s) etaE 0.004 etaE
Event and dropout rate duration(s) S NULL S
必要症例数は N=383.8、必要イベント数は D=245.9 です。手計算の370例／247イベントと非常に近い値になっており、式の理解が正しいことが確認できます。イベント数はほぼ一致する一方で症例数が約14例大きく出ているのは、
nSurv() が脱落ハザード eta=0.004274441（月あたり）を設計に組み込んでいるためです。手計算では最後に「370 / 0.9」と一括で割り増しましたが、nSurv() は脱落を時間の関数として扱うので、より現実に近い症例数が得られます。nSurv() の引数で最も間違えやすいのが有意水準の指定です。上のコードの alpha = 0.025, sided = 1 は「両側5%」に対応します。片側指定なので0.05の半分を渡す必要があり、ここに alpha = 0.05, sided = 1 と書くと片側5%（両側10%）の設計になり、症例数が過小になります。逆に sided = 2 を使う場合は alpha = 0.05 です。計算結果を採用する前に、必ず出力の alpha=0.025 (sided=1) の行を目視で確認してください。powerSurvEpi による設計
powerSurvEpi の ssizeCT.default() は、各群のイベント確率を直接与えて症例数を求めるタイプの関数です。先ほど計算した \( P(\text{event}) \) をそのまま渡せます。引数は power が検出力、k が対照群に対する実験群の割付比（1なら1:1割付）、pE・pC が実験群・対照群の試験期間全体でのイベント（failure）確率、RR が想定ハザード比、alpha が有意水準です。
library(powerSurvEpi)
ss <- ssizeCT.default(power = 0.8, k = 1, pE = 0.6060268, pC = 0.7294947,
RR = 0.7, alpha = 0.05)
ss
cat("total:", sum(ss), "\n")
nE nC
189 189
total: 378
実験群189例・対照群189例、合計378例という結果です。この関数は集積や追跡のモデルを内部に持たず、イベント確率 pE・pC を利用者が与える前提なので、その値をどう計算したかが結果の妥当性をすべて決めます。逆に言えば、Kaplan-Meier曲線から読み取ったイベント割合など、指数分布を仮定しにくい場面でも柔軟に使えるのが利点です。
alpha = 0.05 がここでは両側の指定である点に注意してください。3手法の結果を比較する
| 手法 | 必要総症例数 | 必要イベント数 | 脱落の扱い |
|---|---|---|---|
| 手計算（Schoenfeld＋イベント確率） | 370例 | 247イベント | 後から一括で割り増し |
| gsDesign::nSurv | 383.8例 | 245.9イベント | 脱落ハザードとして設計に内蔵 |
| powerSurvEpi::ssizeCT.default | 378例 | ― | イベント確率に織り込む前提 |
3手法は370例・383.8例・378例と数十例の範囲でばらついていますが、これは計算の誤りではありません。脱落を後から割り増すのか設計に内蔵するのか、イベント確率をどう与えるのかといった前提の違いが素直に反映された結果です。重要なのは「どの数字が正しいか」を争うことではなく、採用した手法と置いた仮定（対照群の中央値、集積期間、追跡期間、脱落率、有意水準の側数）をプロトコルと統計解析計画書に漏れなく明記し、後から誰でも再現できる状態にしておくことです。中間解析を伴う群逐次デザインへ拡張する場合の実装はGroup Sequential Design R実装比較 ― rpact / gsDesign / SAS PROC SEQDESIGN ―で比較しています。
シミュレーションで設計を検算する
公式で求めた設計が本当に狙った検出力を持っているかどうかは、仮想的な試験データを何度も生成して実際に検定を繰り返せば確かめられます。プロトコル作成時に「自分の計算が合っているか」を自己検算する手段として、実務でも非常に有効です。ここでは指数分布で生存時間を発生させ、一様集積と管理打ち切り（administrative censoring）を再現したうえで、ログランク(log-rank)検定を2,000回繰り返し、帰無仮説が棄却された割合を経験的検出力（empirical power）として求めます。
まずシミュレーション関数を定義します。survival パッケージの survdiff() でログランク検定を行い、棄却の有無と観測イベント数を毎回記録する構造です。
library(survival)
sim_power <- function(n_arm, lam_c, hr, A, f, eta = 0, nsim = 2000, seed = 20260721) {
set.seed(seed)
Tt <- A + f
lam_e <- lam_c * hr
rej <- logical(nsim); nev <- numeric(nsim)
for (i in seq_len(nsim)) {
entry <- runif(2 * n_arm, 0, A)
arm <- rep(c(0, 1), each = n_arm)
tev <- c(rexp(n_arm, lam_c), rexp(n_arm, lam_e))
tdrop <- if (eta > 0) rexp(2 * n_arm, eta) else rep(Inf, 2 * n_arm)
adm <- Tt - entry
obs <- pmin(tev, tdrop, adm)
ev <- as.numeric(tev <= pmin(tdrop, adm))
fit <- survdiff(Surv(obs, ev) ~ arm)
rej[i] <- pchisq(fit$chisq, df = 1, lower.tail = FALSE) < 0.05
nev[i] <- sum(ev)
}
list(power = mean(rej), mean_events = mean(nev))
}
まず脱落を一切考慮しない、1群185例（合計370例）の設計を検算します。
r1 <- sim_power(n_arm = 185, lam_c = log(2) / 12, hr = 0.7, A = 24, f = 12, eta = 0)
cat("[185/arm, no dropout ] power:", r1$power, " mean events:", r1$mean_events, "\n")
[185/arm, no dropout ] power: 0.8005 mean events: 247.2675
理論上の必要イベント数247イベント・検出力80%という設計に対し、シミュレーションの結果は経験的検出力0.8005、平均観測イベント数247.2675でした。小数点以下までほぼ完全に一致しており、Schoenfeld式と指数分布に基づくイベント確率の計算が数値的にも妥当であることが確認できます。公式を「信じる」のではなく「確かめる」ことができるのがシミュレーションの価値です。
次に、同じ症例数のまま年5%の追跡不能（脱落）が起きる状況を再現します。さらに、脱落を見込んで10%上乗せした1群206例の設計も検算します。
eta <- -log(1 - 0.05) / 12 # 年5%の追跡不能
r2 <- sim_power(n_arm = 185, lam_c = log(2) / 12, hr = 0.7, A = 24, f = 12, eta = eta)
r3 <- sim_power(n_arm = 206, lam_c = log(2) / 12, hr = 0.7, A = 24, f = 12, eta = eta)
cat("[185/arm, 5%/yr drop ] power:", r2$power, " mean events:", r2$mean_events, "\n")
cat("[206/arm, 5%/yr drop ] power:", r3$power, " mean events:", r3$mean_events, "\n")
[185/arm, 5%/yr drop ] power: 0.7805 mean events: 236.6805
[206/arm, 5%/yr drop ] power: 0.824 mean events: 263.86
脱落を無視して1群185例のまま試験を実施すると、平均イベント数は247.2675から236.6805へ約10イベント減り、検出力は80%から0.7805（78.05%）へ実際に低下します。一方、10%上乗せして1群206例とすると平均イベント数は263.86まで回復し、経験的検出力は0.824（82.4%）となり目標の80%を確保できます。脱落率の見込みは机上の保険ではなく、検出力に直結する設計パラメータであることが数値で示されています。
生存曲線そのものの挙動を確認しておきたい場合は、生存時間解析の基礎：カプラン–マイヤー法と SAS・R による実装を徹底解説も併せてご覧ください。
実務でのポイントと落とし穴
計算そのものは公式に数値を入れるだけですが、実務で設計を誤らせるのは公式の外側にある前提です。ここでは特に頻出する論点を整理します。
比例ハザードの仮定 — Schoenfeld式は、試験期間を通じてハザード比が一定であることを前提にしています。免疫チェックポイント阻害薬のように効果の立ち上がりが遅れる（遅発効果、delayed effect）薬剤では、序盤の生存曲線が重なるため平均的なハザード比が1に近づき、想定した検出力を得られません。比例ハザード性の診断と破綻時の対処はCox比例ハザード性のチェックとRMSTへの切り替え判断に、代替指標としてのRMSTはRMST（制限付き平均生存時間）解析ガイドにまとめています。
中央値ではなくハザード比で考える — 「中央値を12か月から17か月に延ばしたい」という臨床的な目標が先にあることは多いのですが、中央値の比がそのままハザード比になるのは指数分布を仮定したときだけです。ワイブル分布などを想定する場合、中央値比とハザード比は一致しません。
競合リスク — 目的イベント以外の死亡などが起こると、目的イベントの発生確率の見積もりそのものが変わり、必要症例数が過小になります。詳しくは競合リスク（Competing Risks）とはをご参照ください。
estimandとの整合 — 何を「イベント」として数えるか、治療中止や後治療といった中間事象（intercurrent event）をどう扱うかが決まらなければ、イベント確率も必要症例数も定まりません。設計段階でICH E9(R1) Estimandフレームワーク徹底解説の考え方と揃えておく必要があります。
中間解析があるならαを消費する — 群逐次デザインでは中間時点で有意水準の一部を使うため、同じ検出力を保つには最終解析までの必要イベント数が固定デザインより増えます。実装はGroup Sequential Design R実装比較で扱っています。
サンプルサイズ設計で仮定した検定と、実際の主解析の手法は必ず揃えてください。設計を単純なログランク検定で行いながら主解析を層別ログランク検定や共変量調整済みCox回帰にすると、検出力は設計時の想定からずれます。層別因子を増やしすぎて各層のイベント数が乏しくなると、かえって検出力を失うこともあります。シミュレーションで検算する際も、主解析と同じモデル・同じ層別構造で検定を回すことが、設計どおりの検出力を担保する最も確実な方法です。
イベント駆動型試験の運用 — 生存時間解析の試験では、登録が終わった時点ではなく、必要イベント数に到達した時点でデータカットオフを設定する運用（event-driven design）が標準です。したがって試験期間中は、ブラインドを維持したまま両群合計のイベント数だけをモニタリングし、解析実施時期を予測することになります。群別のイベント数を見てしまうと盲検性が損なわれるため、モニタリングの手順も設計と同時に決めておく必要があります。
サンプルサイズ設計は「数値を1つ出す作業」ではなく、比例ハザード性・イベント定義・脱落率・中間解析の有無といった前提を明示し、それぞれが外れたときの影響を示す作業です。プロトコルには最終的な症例数だけでなく、必要イベント数・仮定した中央値とハザード比・脱落率・感度分析の結果まで記載しておくと、照会対応が格段に楽になります。
📚 この記事をより深く理解するための参考書籍
生存時間解析とサンプルサイズ設計をさらに深く学びたい方に、おすすめの書籍をご紹介します。
まとめ
本記事では、生存時間解析におけるサンプルサイズ設計を、必要イベント数の算出から実際の症例数、そしてシミュレーションによる検算まで一貫して紹介しました。検出力を決めるのは症例数ではなくイベント数であること、Schoenfeld式でハザード比と検出力から必要イベント数が直接求まること、そこから集積期間と追跡期間を踏まえたイベント発生確率で割り戻して必要症例数を得ること、さらに脱落を見込んで割り増しを行うこと、という流れをRの実装とともに追いました。最後にシミュレーションで検算し、理論値とほぼ一致することを数値で確認しています。
この一連の流れを自分の手で再現できるようになると、統計解析計画書の症例数根拠を人の計算に頼らず組み立てられるようになります。前提を変えたときに症例数がどう動くかを即座に示せることは、臨床チームとの議論や規制当局からの照会対応において大きな武器になります。
関連するテーマとして、サンプルサイズ設計全般の基礎は検出力分析（Power Analysis）入門 ― サンプルサイズ設計に欠かせない基礎知識をSAS・Rで実装 ―で、ハザード比を推定する解析手法そのものはCox比例ハザードモデル入門〜数式から実務応用まで〜で、比例ハザードが成り立たない場合の代替アプローチはRMST（制限付き平均生存時間）解析ガイド ― Cox比例ハザードが使えない時の判断フローとR・SAS実装 ―で詳しく解説しています。あわせてお読みいただき、設計から解析までを一本の線として理解していただければと思います。
