はじめに

毎年、何百もの臨床予測モデルが発表されているが、その多くはあまりにも小さなデータセットを用いて開発されている。このため、予測は不正確となり、役に立たないモデルで文献が溢れかえることになる。二値またはイベントまでの時間に関する予測モデルを開発する場合、経験則として、変数ごとに少なくとも10イベントを確保することが言われてきた。しかし、それでも不足している場合やときには多すぎる場合もあり、適切なサンプルサイズ設計が求められてきた。

Calculating the sample size required for developing a clinical prediction model BMJ 2020; 368 doi: https://doi.org/10.1136/bmj.m441

はサンプルサイズを計算する、より高度な方法を提案している。このマークダウンでは主に二値アウトカムのClinical prediction model開発のサンプルサイズ計算についてまとめる。

サンプルサイズに必要なもの

サンプルサイズ計算には以下のものが必要となる
アウトカム平均値を誤差幅を小さく推定するために必要なもの (B1 method)

全体的なアウトカムリスクまたは平均アウトカム値
誤差幅 0.05以下

全体で予測誤差が小さい予測値を得るために必要なもの(B2 method)

予測変数の候補の数
MAPE (平均予測誤差：新しい予測モデルが、平均して、標的集団における推定結果確率の予測誤差)

Over fittingの問題を最小化するために必要なもの (B3 method)

Shrinkage (収縮（ペナルティ法、正則化法とも呼ばれる）は、極端な予測（例えば、0または1に近い予測確率）を全体の平均に引き戻すように、開発したモデルの予測値の変動を減らすことでオーバーフィッティングの問題に対処するもの。0.9に設定)
　　
Cox-Snell R二乗統計量（\(R^2_{cs}\)）: 予測されるモデルの当てはまりの良さの（保守的な）値. 予測されるR2csは、複数のパラメータの推定と過剰適合の可能性に影響を与える信号-ノイズ比を反映している。信号-ノイズ比が高いと予想される場合（つまり、ノイズが少なく予測モデルでR2csが1に近い場合）.真のパターンが検出しやすくなるため、オーバーフィットの懸念が少なく、より多くの予測パラメータを推定することができる。しかし、信号-ノイズ比が低い場合（つまり、ノイズが多くR2csが0に近いと予想される）、真のパターンは識別しにくく、オーバーフィッティングの可能性が高くなり、少ない予測パラメータしか信頼性をもって推定できなくなる.

必要なサンプルサイズとEPPは、\(R^2_{cs}\)の選択に敏感で、\(R^2_{cs}\)の予想値が低いと必要なサンプルサイズが大きくなる。したがって、\(R^2_{cs}\)は保守的に選択することが推奨される。二値アウトカムのモデルではだいたい0.1~0.2.

\(R^2_{cs}\)の求め方として、

既存モデルの報告から引用
既存のモデル他の統計量から計算
期待される\(R^2_cs\)を計算の方法

がある.3は以下の方法で計算する.

\(maxR^2_{cs}\): \(R^2_{cs}\)の最大値：計算して求める
期待される\(R^2_{cs}\): \(R^2Nagelkerke * maxR^2_{cs}\)で計算.

Optimismを小さくするために必要なもの (B4 methods)

\(R^2Nagelkerke\): Nagelkerke R2乗： Cox-Snell R2乗の最大値が1にならない点を修正したもので、期待されるモデル適合度(\(R^2_{cs}/maxR^2_{cs}\)).これの差が0.05以下を目指す。つまり楽観性が小さいことを目指す。

実際に計算してみる。

一度にB1からB4を計算するため以下の手順を取る。

BMJのStepとは一致していない。

Step 1 \(maxR^2cs\) を求める

\(R^2_{cs}\)を計算するにあたり、既存報告から引用できないことを前提に算出する方法（上記の3）を用いる。

B1-B4の方法でサンプルサイズを一括で算出したいので、最初に\(maxR^2cs\) を算出しとく。

E: イベントありの総数
n: 任意に選んだサンプルサイズ

例えば、結果の割合が50％ならば任意のサンプルサイズ𝑛=100したがって𝐸=50イベントとする。

これは、n=10,E=5にしても一緒である。

Function calcurate \(max R^2_{cs}\)

## Create function
maxR2cs_fun<-function(E,n){
  lnL<-E*log(E/n)+(n-E)*log(1-E/n)
  R2cs<-1-exp(2*lnL/n)
  return(R2cs)
}
## See data supplement S5 BMJ 2022

例えば50％のイベント割合であれば、

maxR2cs_fun(50,100)

## [1] 0.75

これは割合が同じであれば何でも良い。

maxR2cs_fun(5,10)

## [1] 0.75

Step 2 \(R^2Nagelkerke\)を決める

\(R^2Nagelkerke\)は、モデルが分散の合計を何％説明するかで決定する。
医療だと10-20%程度になることが多い。
モデルにベースラインのアウトカムを含むような場合は0.5になることもある。
例えば、ベースのCreから1週間後のCreを予測する場合などである。

BMJ2020ではを保守的に15%と設定していた。

Step 3 期待される\(R^2_{cs}\)を求める

\(R^2Nagelkerke=R^2_{cs}/maxR^2_{cs}\)
なので

期待される\(R^2_{cs}\)は

\(期待されるR^2_{cs}=R^2Nagelkerke * maxR^2_{cs}\)で求まる。

ここからは、イベント割合は0.05 (E=5, n=100)として進めていく。

# anticipated R2cs
0.15*maxR2cs_fun(5,100) #0.05

## [1] 0.04915297

だいたい0.05であった。

Step 4 pmsampsizeでサンプルサイズを求める

B1, B3, B4 methodsのみ算出してくれる。B2は後に計算する。

library(pmsampsize)

では、以下の設定でサンプルサイズ計算を行う。

\(R^2_{cs}\)=0.05 (Step3で求めたもの)
候補となる予測因子の数=30 (最終的に残った予測因子でなく変数選択前の候補を含む)
イベント割合=0.05 (既報から自分で設定したもの)

#B1,B3,B4 methods
pmsampsize(type="b", #binary outcome
           rsquared = 0.05, # anticipated R^2cs (R^2Nagelkerke*maxR^2cs)
           parameters = 30, # Predictor
           prevalence = 0.05 )#Prevarence

## NB: Assuming 0.05 acceptable difference in apparent & adjusted R-squared 
## NB: Assuming 0.05 margin of error in estimation of intercept 
## NB: Events per Predictor Parameter (EPP) assumes prevalence = 0.05  
##  
##            Samp_size Shrinkage Parameter CS_Rsq Max_Rsq Nag_Rsq  EPP
## Criteria 1      5249     0.900        30   0.05   0.328   0.153 8.75
## Criteria 2      1770     0.753        30   0.05   0.328   0.153 2.95
## Criteria 3        73     0.900        30   0.05   0.328   0.153 0.12
## Final           5249     0.900        30   0.05   0.328   0.153 8.75
##  
##  Minimum sample size required for new model development based on user inputs = 5249, 
##  with 263 events (assuming an outcome prevalence = 0.05) and an EPP = 8.75 
##  
##

Events per Predictor Parameter (EPP) =Events per Predictor Vliable (EPV)

Step 5 B2 methods

B2 methodは自分で計算する。

Pr: 期待されるアウトカム割合
P: 予測因子として設定(30未満まで)
MAPE: 平均予測誤差、基本0.05で設定

B2_samplesize<-function(Pr,P,MAPE){
  n<-exp({-0.508+0.259*log(Pr)+0.504*log(P)-log(MAPE)}/0.544)
  return(n)
}

B2_samplesize(0.05,30,0.05)

## [1] 543.413

上記のB1-B4 methodsのうち、最大のサンプルサイズであり、5249人がモデルの開発に必要な人数である（EPP=8.75）。

これは、EPP=10より少なく済んでいる。

外部検証のサンプルサイズ

これらのサンプルサイズ計算はCPR開発のサンプルサイズ計算であり、外部検証のサンプルサイズ計算は別のものを参照すること。

Collins GS, Ogundimu EO, Altman DG. Sample size considerations for the external validation of a multivariable prognostic model: a resampling study. Stat Med2016;35:214-26. doi:10.1002/sim.6787 pmid:26553135
Van Calster B, Nieboer D, Vergouwe Y, De Cock B, Pencina MJ, Steyerberg EW. A calibration hierarchy for risk models was defined: from utopia to empirical data. J Clin Epidemiol2016;74:167-76. doi:10.1016/j.jclinepi.2015.12.005 pmid:26772608
Jinks RC, Royston P, Parmar MK. Discrimination-based sample size calculations for multivariable prognostic models for time-to-event data. BMC Med Res Methodol2015;15:82. doi:10.1186/s12874-015-0078-y pmid:26459415
Vergouwe Y, Steyerberg EW, Eijkemans MJ, Habbema JD. Substantial effective sample sizes were required for external validation studies of predictive logistic regression models. J Clin Epidemiol2005;58:475-83.

Reference

Calculating the sample size required for developing a clinical prediction model BMJ 2020; 368 doi: https://doi.org/10.1136/bmj.m441 (Published 18 March 2020) Cite this as: BMJ 2020;368:m441

CPR samplesize

Ryohei Yamamoto

2022/08/06