多中心RCT样本量计算：基于R模拟数据的混合效应模型疑问咨询

我正在开展一项治疗组与安慰剂组1:1对照的多中心RCT，主要结局为妊娠延长天数，最小临床相关差异（MCRD）为5天。根据既往研究，该结局呈偏态分布，假定为对数正态（log-normal）分布。中心对妊娠延长的态度会影响主要结局，但该影响程度难以估计。基于反馈，我决定利用可比单中心试验的汇总统计量生成模拟数据，采用以中心为随机效应、治疗为固定效应的对数正态混合效应模型，计算试验所需样本量（检验效能80%，α=0.05）。

既往试验的统计信息如下：

• 安慰剂组对数转换后妊娠延长数据均值为2.25，标准差为1.23；
• 未转换数据均值为17.6天，标准差为18天；
• 总样本量180例（每组90例）。

基于MCRD=5天，我编写了如下R代码：

# Define parameters
## I used the same sample size as the other trial, 180 patients (90 per arm)
n_centers <- 9  # Number of centers
n_per_center_per_group <- 10  # Number of participants per treatment group per center
total_n <- n_centers * n_per_center_per_group * 2  # Total number of participants (with 2 treatment arms)
n <- 90

# Generate treatment variable for placebo
treatment_placebo <- rep("Placebo", n_per_center_per_group * n_centers)

# Generate treatment variable for medicine
treatment_medicine <- rep("Medicine", n_per_center_per_group * n_centers)

treatment <- factor(c(treatment_placebo, treatment_medicine))

# Generate center variable for each treatment group
center_placebo <- rep(1:n_centers, each = n_per_center_per_group)
center_medicine <- rep(1:n_centers, each = n_per_center_per_group)

center <- factor(c(center_placebo, center_medicine))

# Combine placebo and medicine data
simulated_data <- data.frame(treatment = factor(c(treatment_placebo, treatment_medicine)),
                             center = factor(c(center_placebo, center_medicine)))


# Generate log-normal prolongation of gestation with different means for treatment and placebo
mean_placebo <- 2.25 # mean of log normal transformed data in previous trial
MCRD <- log (16/11) #MCRD of 5 days on the untransformed scale: 11 days estimated as a mean prolongation in placebo and 16 days in treatment group, transformed to put it on log-scale. 
mean_treatment <- mean_placebo + MCRD 
sd_log <- 1.23  # Standard deviation on the log scale (edited, incorrectly defined as 0.7 in previous version)

log_prolongation_placebo <- rlnorm(n, meanlog = mean_placebo, sdlog = sd_log)
log_prolongation_treatment <- rlnorm(n, meanlog = mean_treatment, sdlog = sd_log)

log_prolongation_placebo <- log(log_prolongation_placebo)
log_prolongation_treatment <- log(log_prolongation_treatment)

log_prolongation <- c(log_prolongation_placebo, log_prolongation_treatment)

# Create data frame
simulated_data <- data.frame(treatment = factor(c(treatment_placebo, treatment_metformin)),
                             log_prolongation = log_prolongation,
                             center = factor(c(center_placebo, center_metformin)))

#RUN SIMULATION WITH SIMULATED DATA
library(lme4)
lmer_model <- lmer(log_prolongation ~ treatment + (1 | center), data = simulated_data)
summary(lmer_model)

powerSim(lmer_model, nsim=1000, alpha = 0.05)

由于我对R和统计学经验不足，担心模型输入存在问题（尤其是反复进行对数转换的部分），特提出以下疑问：

1. 代码中MCRD的计算是否正确？
2. 对数正态尺度的编码是否正确？是否存在尺度混叠问题？
3. 该样本量计算方法是否合理？

解答 #

1. MCRD计算的问题 #

你当前用log(16/11)对应5天差异的逻辑存在矛盾：这个计算基于安慰剂组均值11天的假设，但你提供的既往试验中安慰剂组未转换数据均值是17.6天，两者不匹配。

正确的转换逻辑应基于对数正态分布的尺度转换公式：原始尺度的均值 = exp(对数均值 + 对数标准差²/2)。已知安慰剂组对数均值=2.25、对数标准差=1.23，对应的原始尺度理论均值为exp(2.25 + 1.23²/2)≈20.2，和既往试验的17.6接近（属于样本均值与理论均值的正常差异）。

如果临床需求是原始尺度均值差5天，治疗组原始均值应为17.6+5=22.6，对应的对数均值需反推：对数均值 = log(22.6) - 1.23²/2≈3.12-0.756≈2.364，因此对数尺度的MCRD应为2.364-2.25=0.114，而非log(16/11)。

2. 对数正态尺度编码的问题 #

代码存在明显的尺度混叠和变量错误：

• 冗余的对数转换：rlnorm()生成的是原始尺度的对数正态数据，你随后又对其做log()转换，这完全多余。直接用rnorm()生成对数尺度的正态数据即可，和log(rlnorm())结果一致，但更高效：

  log_prolongation_placebo <- rnorm(n, mean = mean_placebo, sd = sd_log)
  log_prolongation_treatment <- rnorm(n, mean = mean_treatment, sd = sd_log)
  

• 变量名错误：数据框创建时用了treatment_metformin和center_metformin，但之前定义的变量是treatment_medicine和center_medicine，会导致代码运行报错。
• 缺失中心随机效应：当前模拟数据仅设置了中心分组变量，但未加入中心间的变异，会导致后续混合效应模型的中心随机效应方差估计为0，完全不符合多中心试验的实际情况。需要先生成中心随机效应再生成结局数据：

  # 假设中心间标准差为0.2（可根据类似试验调整或做敏感性分析）
  sigma_center <- 0.2
  center_effect <- rnorm(n_centers, mean = 0, sd = sigma_center)
  
  log_prolongation_placebo <- rep(center_effect, each = n_per_center_per_group) + rnorm(n, mean = mean_placebo, sd = sd_log)
  log_prolongation_treatment <- rep(center_effect, each = n_per_center_per_group) + rnorm(n, mean = mean_treatment, sd = sd_log)
  

3. 样本量计算方法的合理性 #

用模拟法计算多中心混合效应模型的样本量是合理的，但当前实现存在关键缺陷：

• 未加入中心随机效应：忽略中心间变异会低估所需样本量，因为中心间差异会增加总变异度，降低检验效能。
• 固定样本量模拟：你当前直接用180例样本模拟，而非迭代寻找达到80%效能的最小样本量。应设置不同样本量（比如从100到300，按中心数和每组每中心例数调整），分别计算效能，找到满足要求的最小样本量。
• 依赖simr包的powerSim()：需确保提前加载simr包，且模拟过程中模型拟合稳定，避免收敛问题。

优化后的核心步骤：

1. 确定合理的中心间变异度（参考同类试验或做敏感性分析，比如尝试0.1/0.2/0.3）。
2. 基于原始尺度5天的绝对差，正确转换得到对数尺度的MCRD。
3. 编写循环测试不同样本量，计算对应效能，找到达标最小样本量。

内容的提问来源于stack exchange，提问作者Catherina Muscat