BMJ在2020年發(fā)表了一篇關(guān)于預(yù)測模型樣本量計算的文章（10.1136/bmj.m441，這篇文章是免費下載的，記得把補(bǔ)充文件也下載下來），算是目前樣本量計算的指南性文件。

并且作者還提供了一個R包(pmsampsize)用于計算樣本量，使用起來也非常方便，有點類似于pwr包計算樣本量

但是要注意該文獻(xiàn)的樣本量計算是針對開發(fā)臨床預(yù)測模型需要的樣本量，也就是訓(xùn)練集的樣本量，不是外部驗證集的樣本量。如果是驗證集的樣本量，作者專門又寫了3篇文章，分別針對回歸、二分類、生存數(shù)據(jù)：

• Riley R D, Collins G S, Ensor J, et al. Minimum sample size calculations for external validation of a clinical prediction model with a time-to-event outcome[J]. Statistics in Medicine, 2022, 41(7): 1280–1295. DOI:10.1002/sim.9275.
• Riley R D, Debray T P A, Collins G S, et al. Minimum sample size for external validation of a clinical prediction model with a binary outcome[J]. Statistics in Medicine, 2021, 40(19): 4230–4251. DOI:10.1002/sim.9025.
• Archer L, Snell K I E, Ensor J, et al. Minimum sample size for external validation of a clinical prediction model with a continuous outcome[J]. Statistics in Medicine, 2021, 40(1): 133–146. DOI:10.1002/sim.8766.

另外，該文獻(xiàn)的樣本量計算方法是針對3大回歸模型的：線性回歸、邏輯回歸、cox回歸。如果是一些機(jī)器學(xué)習(xí)方法（比如隨機(jī)森林、支持向量機(jī)等）則需要更多的樣本。

經(jīng)典方法

Harrell老爺子在他的書《Regression Modeling Strategies》中介紹的開發(fā)模型的樣本量計算方法是：

在開發(fā)數(shù)據(jù)集(也就是訓(xùn)練集)中，連續(xù)型結(jié)果的有效樣本量由研究參與者的總數(shù)決定(有多少用多少)。對于二分類結(jié)果，有效樣本量通常被認(rèn)為大約等于事件（有結(jié)果的事件）和非事件（沒有結(jié)果的事件）的最小值; time-to-event數(shù)據(jù)中，樣本量可以粗略等于陽性事件的數(shù)量。

在為二分類或time-to-event數(shù)據(jù)開發(fā)預(yù)測模型時，所需要的樣本量常用的計算方法是10EPV法，即陽性事件的數(shù)量至少是預(yù)測變量個數(shù)的10倍（10 events per variable，10EPV）。

但是“variable”一詞具有誤導(dǎo)性，因為在模型中一個預(yù)測變量可能有多個β（即回歸系數(shù)），例如，具有三個類別的分類型預(yù)測變量就會有兩個β（例如腫瘤等級1、2、3，那么就會有β(2和1比)、β(3和1比)，因為分類變量在回歸分析中需要進(jìn)行啞變量編碼）。還有就是在建模過程中使用了多項式轉(zhuǎn)換和樣條變換等也會使得同一個變量有多個β，如果變量之間有交互項也會產(chǎn)生同樣的結(jié)果。

由于預(yù)測模型的參數(shù)（也就是回歸系數(shù)β）通常多于實際的預(yù)測變量個數(shù)，所以最好使用10EPP（10 events per candidate predictor parameter）法，即陽性事件的數(shù)量至少是“候選預(yù)測變量的參數(shù)”的10倍?！昂蜻x”一詞很重要，因為模型過擬合的程度與預(yù)測變量參數(shù)的數(shù)量有關(guān)，而不是最終模型方程中的參數(shù)數(shù)量。

但是10EPP原則目前也有一些爭議，也有大佬建議5EPP或者15、20、50EPP。這些數(shù)量的使用都是和具體的情況有關(guān)的，也沒個金標(biāo)準(zhǔn)，不僅取決于相對于候選預(yù)測變量參數(shù)數(shù)量的事件數(shù)量，還取決于參與者總數(shù)、研究人群中的結(jié)果比例（發(fā)生率）以及模型的預(yù)期預(yù)測性能等。

4步法

Van Smeden等和Riley等人最近的工作描述了如何計算預(yù)測模型開發(fā)所需的樣本量，使用條件是用戶指定目標(biāo)人群中的總體結(jié)果風(fēng)險或平均結(jié)果值、候選預(yù)測變量參數(shù)的數(shù)量以及總體模型擬合方面的預(yù)期模型性能。具體實施起來總結(jié)為4個步驟：

1. 第一步是確定需要多少樣本才能準(zhǔn)確估計平均風(fēng)險（也就是平均概率，對應(yīng)二分類或者生存數(shù)據(jù)）或者平均值（對應(yīng)連續(xù)型數(shù)據(jù)）；
2. 第二步是確定需要多少樣本才能讓模型的預(yù)測值和真實值之間的誤差更?。?/section>
3. 第三步是確定有足夠的的樣本量以防止過擬合；
4. 第四步是確定有足夠的的樣本量使模型的表面性能和真實性能的誤差更小。

第一步

樣本大小必須讓預(yù)測模型的截距能被精確估計，以確保開發(fā)的模型可以準(zhǔn)確預(yù)測平均結(jié)果值(對應(yīng)連續(xù)型數(shù)據(jù))或總體結(jié)果比例(對應(yīng)二分類或者生存數(shù)據(jù))。一個簡單的方法是計算：能夠準(zhǔn)確估計“沒有預(yù)測變量的空模型(null model)的截距”所需要的樣本量。

這里涉及一個簡單的數(shù)學(xué)知識，就是線性模型的截距反映了模型預(yù)測的平均值。

這個“準(zhǔn)確估計”一般要求誤差在0.05以內(nèi)，也就是說預(yù)測值最好不要超過均值的95%可信區(qū)間。

下圖是計算公式和一個例子。假如一個二分類數(shù)據(jù)，它的陽性事件比例是0.5，為了控制誤差在0.05以內(nèi)，根據(jù)以下公式計算，需要的樣本量最少是385個。

第二步

預(yù)測值和真實值之間的誤差可以用很多指標(biāo)衡量，比如平均絕對百分比誤差（Mean Absolute Percentage Error MAPE），這個指標(biāo)其實是衡量回歸模型的常用指標(biāo)，對于二分類數(shù)據(jù)如果使用的是概率的話也能用這個指標(biāo)衡量。

下面是計算公式和一個例子。假如一個二分類數(shù)據(jù)，它的陽性事件比例是0.3，預(yù)測變量有10個，為了控制誤差（MAPE）在0.05以內(nèi)，根據(jù)以下公式計算，需要的樣本量最少是461個。

第三步

樣本量越少且預(yù)測變量數(shù)量越多，則越容易過擬合，因此需要足夠的樣本量防止過擬合。

建模過程中通常會使用收縮法（Shrinkage，或者被稱為懲罰(penalisation)或正則化(regularisation)）來降低過擬合的風(fēng)險。Riley等人建議使用一個較小的收縮值（≤10%），并計算此時所需要的樣本量。并且還需要指定候選預(yù)測變量參數(shù)的個數(shù)以及一個模型性能指標(biāo)，比如Cox-Snell R²(記為CS-R方，屬于偽R方的一種)。CS-R方可以反應(yīng)信噪比（signal:noise），從而反應(yīng)模型是否過擬合。

對于連續(xù)型數(shù)據(jù)來說，CS-R方就是決定系數(shù)，反應(yīng)模型所能解釋的方差（或者叫變異）百分比，范圍是0到1之間，越接近1越好，說明模型能夠準(zhǔn)確識別數(shù)據(jù)內(nèi)部的模式，不會被噪聲（誤差）干擾，如果CS-R方接近0則說明模型很有可能過擬合。

對于二分類數(shù)據(jù)和生存數(shù)據(jù)來說，CS-R方的范圍是0到max(CS-R方)。對于邏輯回歸模型來說，如果陽性事件發(fā)生率為0.5，0.4，0.3，0.2，0.1，0.05，0.01，那么對應(yīng)的max(CS-R方)分別是0.75，0.74，0.71，0.63，0.48，0.33，0.11。所以即使模型的預(yù)期性能非常好，這個CS-R方的值也一般會選擇比較小的值。

以下是二分類和生存數(shù)據(jù)的樣本量計算公式和一個示例。對于一個邏輯回歸模型，如果有20個候選預(yù)測變量參數(shù)（EPP），CS-R方選擇0.1，那么為了使收縮值保持在10%，最少的樣本量是1698。

Cox-Snell R²(也就是CS-R方)的選擇有多種方法，以下是作者比較推薦的幾種：

• 直接使用別人文獻(xiàn)里報道的值
• 使用其他指標(biāo)近似，比如使用C指數(shù)、AUC值、其他偽R方等
• 根據(jù)max(CS-R方)計算

本文獻(xiàn)的附件5提供了詳細(xì)的公式用于計算max(CS-R方)，感興趣的自己查看一下吧。

第四步

應(yīng)該有足夠多的樣本量保證模型的表面指標(biāo)和真實指標(biāo)之間的差異足夠小。

表面指標(biāo)（apparent values），假如我們用訓(xùn)練集開發(fā)了一個模型，然后讓這個模型對訓(xùn)練集進(jìn)行預(yù)測，這樣得到的指標(biāo)就是表面指標(biāo)，這種計算模型指標(biāo)的方法叫做重代入法（resubstitution）。真實指標(biāo)是指模型在其他數(shù)據(jù)中（就是模型開發(fā)時沒用過的數(shù)據(jù)）得到的更真實、更接近模型真實性能的指標(biāo)。

本篇文獻(xiàn)中采用的指標(biāo)是另外一種調(diào)整的R方，即：Nagelkerke-R方（也是偽R方的一種），Nagelkerke-R方=CS-R方/max(CS-R方)。

下面是二分類和生存數(shù)據(jù)的樣本量計算公式。對于一個邏輯回歸模型，假設(shè)陽性事件的比例是0.05（此時對應(yīng)的max(CS-R方)是0.33），指定CS-R方為0.2，那么為了使真實的Nagelkerke-R方和表面Nagelkerke-R方的差異保持在0.05，至少需要的樣本量是1079。

總結(jié)

下面是一個總結(jié)，對于連續(xù)型數(shù)據(jù)，推薦使用4步法(C1-C4)，對于二分類數(shù)據(jù)推薦使用4步法(B1-B4)，對于生存數(shù)據(jù)推薦使用3步法(T1-T4)。

除此之外作者專門寫了一個R包用于計算臨床預(yù)測模型的樣本量：pmsampsize，這個R包可以計算以上每一個步驟（除了B2這個步驟不行，這一步是通過網(wǎng)頁計算的，網(wǎng)址是：https://mvansmeden./BeyondEPV/）所需要的樣本量，并選擇最大的一個作為開發(fā)模型所需要的最少樣本量。

R包使用方法

下面用3個實例演示這個R包的使用方法。

library(pmsampsize)

二分類數(shù)據(jù)

假如我們要根據(jù)妊娠15周時測定的各種指標(biāo)預(yù)測孕婦發(fā)生子癇前期的風(fēng)險，這是一個二分類數(shù)據(jù)，結(jié)果變量是發(fā)生子癇/不發(fā)生子癇。

假設(shè)該數(shù)據(jù)中，發(fā)生子癇的比例是0.05（陽性事件的比例），候選預(yù)測變量的參數(shù)數(shù)量是30（EPP是30），max(CS-R方)是0.33。如果我們預(yù)期模型能夠解釋15%的變異，根據(jù)第3步中介紹的CS-R方的計算方法，可以得到CS-R方=0.15*0.33=0.05。

有了這幾個數(shù)據(jù)，就可以計算樣本量了：

pmsampsize(
  type = "b",         # 二分類數(shù)據(jù)
  csrsquared = 0.05,  # CS-R方
  parameters = 30,    # EPP的數(shù)量
  prevalence = 0.05,  # 陽性事件發(fā)生率
  seed = 123          # 設(shè)置隨機(jī)數(shù)種子
)

    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 

對于二分類數(shù)據(jù)，使用4步法計算樣本量，其中B2這一步不能通過這個包計算，所以這個包給出了其他3個步驟所需要的樣本量，B2這個步驟算出來是需要544例，因為要同時滿足4個步驟的要求，所以最終需要的樣本量是5249例。

生存數(shù)據(jù)

假如我們要預(yù)測治療停止一段時間后，靜脈血栓栓塞復(fù)發(fā)的風(fēng)險。這是一個time-to-event類型的數(shù)據(jù)，結(jié)局是復(fù)發(fā)/不復(fù)發(fā)，時間就是治療停止后的時長。

假設(shè)該數(shù)據(jù)中，C指數(shù)是0.69，CS-R方是0.051，EPP=30，平均隨訪時間是2.07年，陽性事件發(fā)生比例是0.065，需要進(jìn)行預(yù)測的時間點選擇2年，那么樣本量計算如下：

pmsampsize(
  type = "s",         # 生存數(shù)據(jù)
  csrsquared = 0.051, # CS-R方
  parameters = 30,    # EPP的數(shù)量
  rate = 0.065,       # 陽性事件發(fā)生率
  timepoint = 2,      # 指定要預(yù)測的時間點
  meanfup = 2.07      # 平均隨訪時間
)
    NB: Assuming 0.05 acceptable difference in apparent & adjusted R-squared 
    NB: Assuming 0.05 margin of error in estimation of overall risk at time point = 2  
    NB: Events per Predictor Parameter (EPP) assumes overall event rate = 0.065  
     
                 Samp_size Shrinkage Parameter CS_Rsq Max_Rsq Nag_Rsq   EPP
    Criteria 1        5143     0.900        30  0.051   0.555   0.092 23.07
    Criteria 2        1039     0.648        30  0.051   0.555   0.092  4.66
    Criteria 3 *      5143     0.900        30  0.051   0.555   0.092 23.07
    Final SS          5143     0.900        30  0.051   0.555   0.092 23.07
     
     Minimum sample size required for new model development based on user inputs = 5143, 
     corresponding to 10646 person-time** of follow-up, with 692 outcome events 
     assuming an overall event rate = 0.065 and therefore an EPP = 23.07  
     
     * 95% CI for overall risk = (0.113, 0.13), for true value of 0.122 and sample size n = 5143 
     **where time is in the units mean follow-up time was specified in

生存數(shù)據(jù)的樣本量計算遵循4步法，所以結(jié)果中給出了4個步驟每一步驟所需要的樣本量，最終需要的樣本量是5143例。

連續(xù)型數(shù)據(jù)

假如我們要預(yù)測青少年的無脂肪體重，該任務(wù)很明顯是一個回歸任務(wù)，結(jié)果變量是數(shù)值型的。

假設(shè)該數(shù)據(jù)中，CS-R方為0.9，EPP=20，總體的平均無脂肪體重是26.7kg（截距的值），總體的無脂肪體重的標(biāo)準(zhǔn)差是8.7kg，那么計算樣本量的代碼為：

pmsampsize(
  type = "c",       # 連續(xù)型數(shù)據(jù)
  rsquared = 0.9,   # 連續(xù)型數(shù)據(jù)的CS-R方=R方
  parameters = 20,  # EPP的數(shù)量
  intercept = 26.7, # 截距，也就是均值
  sd = 8.7          # 總體的標(biāo)準(zhǔn)差
)

    NB: Assuming 0.05 acceptable difference in apparent & adjusted R-squared 
    NB: Assuming MMOE <= 1.1 in estimation of intercept & residual standard deviation 
    SPP - Subjects per Predictor Parameter 
     
                Samp_size Shrinkage Parameter Rsq   SPP
    Criteria 1         68     0.900        20 0.9  3.40
    Criteria 2         41     0.853        20 0.9  2.05
    Criteria 3        254     0.970        20 0.9 12.70
    Criteria 4*       254     0.970        20 0.9 12.70
    Final             254     0.970        20 0.9 12.70
     
     Minimum sample size required for new model development based on user inputs = 254  
     
     * 95% CI for intercept = (26.36, 27.04), for sample size n = 254

連續(xù)型數(shù)據(jù)的樣本量計算也遵循4步法，最終所需要的樣本量是254。