Posted by on | Featured

在你计划进行的为期30天的实验中，第14天时，你所测试的AI产品实验结果已经具备了统计学上的显著性。这个实验旨在验证这样一个因果关系：基于大语言模型的功能是否确实改善了实验结果？在场的每一位产品经理都希望尽快将这一功能推向市场。但你的统计专家建议再等待整整30天，因为只有这样，得出的p值才会具有可靠性。

你遵从了专家的建议，等待到了第30天，发现该功能的效果确实存在。然而，你却浪费了16天的时间去测试一个你已经确定有95%置信度能够有效果的功能，这使得后续的实验不得不推迟，从而造成了机会成本的损失。

从技术角度来看，如果使用的是传统的固定样本检验方法，统计专家的说法确实是正确的。在标准的t检验中，p值只有在事先确定了样本规模并且只进行一次测试的情况下才具有有效性。但如果你在测试过程中提前查看结果并在p值小于0.05时就停止测试，那么假阳性率就会上升到30%左右。

p值的设定本来就是用于那种只进行一次测试的情况的——它适用于那些具有固定终点、且测试过程是一次性的实验。而当你需要在一个可以随时查看结果的动态实验环境中使用p值时，就需要采用完全不同的数学方法来进行分析。

正是为了应对这种动态实验环境，人们才设计了序贯检验方法。其中，混合序贯概率比检验（mSPRT）（Johari等人提出）使用一种被称为“e值”的数学指标，能够确保每次测试结果都具有可靠性。你可以每天查看实验结果，一旦获得足够的证据就可以停止测试，而且这种方法的假阳性率会保持在5%左右。

Netflix已经记录了他们在实际生产中如何使用这种序贯检验方法（Lindon等人介绍），而这种方法的理论基础可以追溯到Wald在1945年关于序贯分析的研究以及Ville在1939年提出的不等式。

本教程会明确解释这些概念。你将通过模拟实验来直观地了解错误率会上升的原因，还会从零开始用Python实现mSPRT方法，并将其应用到那个共享的合成大语言模型产品数据集上，从而清楚地理解在什么情况下序贯检验方法会失效。

配套笔记本：本文中提到的每一段代码都可以在配套仓库中的msprt_demo.ipynb文件中端到端运行。

目录

为什么在经典检验中采用可选停止机制

如果频繁查看p值，你的假阳性率会升高到30%左右。这个数字确实应该引起你的注意，而你在下面的第一步中也会亲自验证这一现象。

在经典的假设检验中，p值用于回答一个具体问题：在原假设为真的前提下，当你严格按照既定的计划和样本量进行实验时，观察到如此极端的数据结果的概率是多少？

问题就出在“严格按照计划进行实验”这个要求上。如果你在第5天、第10天或第14天检查结果，然后因为p值小于0.05而停止实验，那么你实际上并没有按照最初的计划来进行实验。你其实进行了14次不同的实验，只选择了其中结果符合标准的那一次作为最终结论。但p值的计算公式并不了解这一点。

我们可以这样理解：在原假设下（即没有实际效果），p值会在0到1之间随机波动，它不会一直保持在0.5这个数值上。如果进行30天的实验，在某个时刻p值低于0.05的概率是相当高的。如果你每天都在检查结果，并且一看到p值小于0.05就立即停止实验，那么你几乎肯定会遇到这种情况，从而得出一个“有效果”的结论。但实际上这种效果并不存在。

减少检查频率同样会引发同样的问题。你需要频繁地进行检测：因为产品的发展速度很快，如果实验时间比必要的时间长16天，就会造成实际的损失，推迟产品的上市时间，还会浪费各种机会。因此，你需要一种无论你在何时停止实验，其结果都仍然有效的检验方法。

序列检验究竟能起到什么作用

序列检验正是为了解决可选停止的问题而设计的，它用另一种统计量——e值来代替p值。

与p值不同，e值是非负的，并且随着时间的推移，由一系列e值构成的序列在原假设下满足“超级鞅性质”：根据已有的数据结果来看，下一个e值的期望值最多只会等于当前的值。

正是这种“超级鞅性质”使得在实验过程中选择任意时间停止都不会带来问题。虽然每个步骤中e值的平均值低于1是必要的条件，但这个条件并不充分；“超级鞅性质”的要求更为严格，它能够确保在任何停止时间点下结果都是可靠的。

原因如下：如果e值序列满足非负超级鞅的性质，并且在原假设H0下E[e_t] ≤ 1，那么根据Ville不等式的结论，该序列的运行最大值超过1/α的概率最多为α。以α = 0.05、停止阈值1/α = 20为例，那么e值序列达到20这个值的概率最多只有5%。

无论你在何时停止实验或检查多少次，这个错误率的上限都是不变的。这种保证具有时间一致性，它适用于所有可能的停止时间点。

而传统的p值所具有的保证仅适用于事先确定的样本量。如果反复进行检测，这个上限就会失效，因此传统方法并不具备时间一致性的优势。

mSPRT通过计算贝叶斯因子来得到e值：这个贝叶斯因子表示在备择假设下观测数据出现的概率与在零假设下出现的概率之比。

“混合”这一概念意味着，在H1假设下，你并不指定一个具体的效应大小。而是需要根据预先设定的效应大小分布，来计算这些概率比的均值。

对于伯努利型结果（比如任务是否完成，只能是“是”或“否”），如果为每个实验组的完成率设定一个Beta(1,1)的先验分布，那么就可以利用对数贝塔函数将贝叶斯因子表示为封闭形式的数学表达式。其实这种计算方法并不复杂：正如第二步所展示的，整个计算过程只需要调用四次betaln函数即可完成。

实际应用中，你可以每天收集数据并计算当时的e值，当这个数值超过20时就可以停止实验。但如果在最大样本量范围内，e值始终低于20，那就说明你无法拒绝零假设。无论你是每天、每小时还是每分钟检查一次数据，I型错误的概率都会保持在5%的水平。

识别假设

mSPRT的有效性取决于四个条件。如果其中任何一个条件不满足，mSPRT的统计功效就会受到影响。下文中会将各种可能导致的失效情况与相应的条件进行对应说明。

  1. 在H0假设下，数据生成过程必须具有非负超鞅性质。也就是说，在H0假设下，E[e_{t+1} | e_1, …, e_t] 必须小于或等于e_t。对于这里使用的Beta-Binomial贝叶斯因子来说，只要先验分布是合理的（Beta(1,1)符合这个要求），且各实验组内的观测数据是独立同分布的，这一条件就成立。

  2. 数据生成过程必须具有稳定性。在整个实验过程中，数据生成机制必须保持不变。如果由于某些外部因素（比如模型更新、营销活动的变化，或者不同工作日的影响）导致实验组的完成率发生变化，那么e值就会受到这些干扰因素的影响，从而使得你无法准确区分处理效应与这些干扰因素带来的影响。

  3. 每个实验组内的观测数据必须是独立的。

    1. 这意味着不同用户的实验结果不能相互影响。如果存在网络效应、共享工作空间，或者推荐系统产生的间接影响，就会违反这一条件。

    2. 必须明确指定先验分布。

      1. 这里使用的Beta(1,1)先验分布其实是一个建模假设。mSPRT的统计功效取决于这个先验分布是否能够合理地反映真实的效应大小。如果先验分布设定不当，虽然不会影响I型错误的概率，但可能会导致e值增长得非常缓慢，从而使得你在样本量耗尽之前始终无法达到判断阈值。

      先决条件

      • Python 3.11或更高版本

      • pandas 2.x版本（通过pip install pandas安装）

      • numpy 1.26或更高版本（通过pip install numpy安装）

      • scipy 1.12或更高版本（通过pip install scipy安装）

      • matplotlib 3.8或更高版本（通过pip install matplotlib安装）

      如果你需要获取合成数据集，可以克隆以下仓库：

      git clone https://github.com/RudrenduPaul/product-experimentation-causal-inference-genai-llm.git
cd product-experimentation-causal-inference-genai-llm
python data/generate_data.py --seed 42 --n-users 50000 --out data/synthetic_llm_logs.csv

      具体操作过程如下：首先会克隆包含该系列所有13份配套笔记本的数据仓库，然后生成一个可供50,000名用户使用的合成数据集，并将其保存到文件《data/synthetic_llm_logs.csv》中。该系列中的每篇文章都会使用这个相同的CSV数据进行分析，因此各种方法的结果可以直接进行比较。此外，数据生成工具还会在数据中人为设置一个+5个百分点的效应，以便用于评估第一组用户使用新功能后任务完成率的变化情况。

      设置实验环境

      这个合成数据集模拟了一个拥有50,000名用户的SaaS人工智能辅助产品。其中，《task_completed》这一列记录了人工智能是否成功完成了用户的任务（值为1）；而《wave》这一列则将用户分为不同的组：第一组用户会率先使用新功能，而第二组用户则作为对照组。

      422306d0-efbc-44d6-a0a1-f8415f2d5e6d

      图1：概念性的e值变化轨迹。蓝色曲线代表真实效果，其数值会逐渐上升，并在绿色虚线处超过临界值；紫色曲线代表的效果较弱，在30天内始终未超过临界值；灰色曲线表示零效应，其数值几乎一直保持在1附近；红色虚线则是临界边界，对应于1/α=20这个条件。请将此图与下方的图2进行对比，后者展示了在实际数据集上e值的变化情况。

      import pandas as pd
import numpy as np

df = pd.read_csv("data/synthetic_llm_logs.csv")

treated = df[df["wave"] == 1]["task_completed"].values
control = df[df["wave"] == 2]["task_completed"].values

print(f"实验组：样本数量={len(treated):,}, 平均值={treated.mean():.4f}")
print(f>对照组：样本数量={len(control):,}, 平均值={control.mean():.4f}")
print(f>观察到的效果提升幅度：{treated.mean() - control.mean():.4f})

      预期输出结果：

      实验组：样本数量=24,937, 平均值=0.6202
对照组：样本数量=25,063, 平均值=0.5718
观察到的效果提升幅度：0.0485

      具体操作过程如下：首先会加载这个包含50,000条记录的数据集，并将其按照用户所属组别进行划分。实验组共有24,937名用户，他们的任务完成率为62.0%；对照组也有25,063名用户，他们的任务完成率为57.2%。观察到的4.85个百分点的效果提升幅度与数据生成工具中预设的5个百分点数值较为接近，这一细微差距可能是由于抽样误差造成的。接下来，这些数据将会被依次用于后续的统计分析中，具体步骤如下所述。

      第一步：模拟“窥探效应”问题

      “窥探效应”是一个真实存在且可以被测量的现象：如果进行30天的连续监测，虚假正例的比例会从4.2%上升到30.2%，下面的模拟实验可以验证这一结论。

      在这个模拟实验中，共进行了1,000次零效应实验（即处理措施对结果没有任何影响），并且每天都会检查当前的p值是否已经低于0.05。实验设置为每组每天有60名用户参与，整个实验持续30天，因此每组总共会有1,800条观测数据，这个规模对于一个中等规模的SaaS产品来说是比较现实的。

      from scipy import stats
import numpy as np

np.random.seed(42)

N_SIMS = 1000
N_days = 30
USERS_PER_ARM_PER_DAY = 60
NULL_RATE = 0.60

false_positives_peeking = 0
false_positives_single_look = 0

for _ in range(N_simS):
    control_outcomes = []
    treated_outcomes = []
    stopped_early = False

    for day in range(N_days):
        control_outcomes.extend(np.random.binomial(1, NULL_RATE, USERS_PER_ARM_PER_DAY))
        treated_outcomes.extend(np.random.binomial(1, NULL RATE, USERS_PER_arm_PER_DAY))

        # 如果每天都进行检测，就会出现假阳性结果：
        if len(control_outcomes) >= 10:
            _, p = stats.ttest_ind(treated_outcomes, control_outcomes)
            if p < 0.05 and not stopped_early:
                false_positives_peeking += 1
                stopped_early = True

    # 如果只在最后一次进行检测，就会得到更准确的结果：
    _, p_final = stats.ttest_ind(treated_outcomes, control_outcomes)
    if p_final < 0.05:
        false_positives_single_look += 1

print(f"每天进行检测时的假阳性率：{false_positives_peeking / N_SIMS:.1%}")
print(f"仅在最后一次进行检测时的假阳性率：{false_positives_single_look / N_simS:.1%}")

      预期输出：

      每天进行检测时的假阳性率：30.2%
仅在最后一次进行检测时的假阳性率：4.2%

      具体原理如下：

      每次模拟都会生成随机数据，其中两个组的数据都来自相同的60%成功率分布，因此任何观察到的“差异”实际上都是随机噪声。内层循环每天为每个组添加60个观测值，并对这些数据进行t检验。

      当p值首次低于0.05时，模拟就会认为出现了假阳性结果，并停止运行（这类似于在实际研究中，一旦发现显著性就立即终止实验的情况）。

      在第30天进行的检测属于“仅一次性检测”的情况，其结果更为准确。在这种情况下，假阳性率为4.2%，与理论预期值接近；而如果每天都进行检测，假阳性率会达到30.2%，这意味着有超过四分之一的“显著”实验实际上只是检测到了随机噪声。

      步骤2：实现mSPRT e值算法

      mSPRT算法会在每个时间点计算贝叶斯因子，即观察到的数据在各种假设下出现的可能性差异。对于每次试验只有两种结果的情景，并且已知每个组的成功率服从Beta(1,1)分布，那么可以使用对数贝塔函数来计算这个贝叶斯因子。

      from scipy.special import betaln

def compute_evalue_running(outcomes_treated, outcomes_control,
                           alpha_prior=1.0, beta_prior=1.0):
    """
    计算两个组在Bernoulli分布下的mSPRT e值。

    参数
    ----------
    outcomes_treated : 0/1类型的数组
    outcomes_control : 0/1类型的数组
    alpha_prior, beta_prior : Beta分布的先验参数（默认值为均匀分布）

    返回
    -------
    e_values : 数组，其中每个元素代表一次观测得到的e值
    """
    outcomes_treated = np.asarray(outcomes_treated, dtype=float)
    outcomes_control = npассив(outcomes_control, dtype(float)
    n = min(len(outcomes_treated), len(outcomes_control))

    cum_t = np.cumsum(outcomes_treated[:n])
    cum_c = np.cumsum(outcomes_control[:n])
    t_arr = np.arange(1, n + 1, dtype=float)

    # 备择假设：每个组都有独立的Beta分布先验
    log_ml_t = (betaln(alpha_prior + cum_t, beta_prior + t_arr - cum_t)
                - betaln(alpha_prior, beta_prior))
    log_ml_c = (betaln(alpha_prior + cum_c, beta_prior + t_arr - cum_c)
                - betaln(alpha_prior, beta_prior))

    # 假设零：两个组共享同一个Beta分布先验
    pooled_successes = cum_t + cum_c
    pooled_n = 2 * t_arr
    log_ml_h0 = (betaln(alpha_prior + pooled_successes,
                        beta_prior + pooled_n - pooled_successes)
                 - betaln(alphaprior, beta_prior))

    # 贝叶斯因子即为对数边际似然值的差
    log_bf = log_ml_t + log_ml_c - log_ml_h0

    return np.exp(log_bf)

      具体过程如下：该函数接收两个数组，其中每个数组包含按时间顺序排列的0或1结果。对于每一个时间点t，该函数会计算每个选项在 해당时间点的成功次数与总尝试次数的累积值。

      betaln函数返回beta函数的対数值，而beta函数正是Beta-Binomial边缘概率分布中的归一化常数。H1模型假设每个选项的失败率遵循独立的Beta分布；而H0模型则假设所有选项的失败率遵循同一个共享的Beta分布。

      对数贝叶斯因子就是这两个假设下计算结果的差异。将其取指数后得到的值即为e值。当某种干预措施确实有效时，e值会随时间逐渐增加；而在没有效果的情况下，e值会始终接近1，并且在H0假设下它满足非负超鞅的性质。

      通过对空数据进行的简单验证，可以确认上述理论预测的行为是正确的：

      np.random.seed(0)
null_t = np.random.binomial(1, 0.60, 500)
null_c = np.random.binomial(1, 0.60, 500)
ev_null = compute_evalue_running(null_t, null_c)
print(f"在零假设下，最终得到的e值应为接近1的值：{ev_null[-1]:.3f}")
print(f"在零假设下，最大的e值为：{ev_null.max():.3f}")

      预期结果：

      在零假设下，最终得到的e值应为接近1的值：0.078
最大的e值为：2.188

      具体过程如下：在零假设的情况下，最终的e值会接近1（在这个例子中为0.078，这是由于抽样误差导致的），而且在500次观测中得到的最大e值也远低于20这个停止阈值。根据Ville不等式，一个有效的零假设下的e值过程达到20的概率最多仅为5%，这一结果与5%的第一类错误率是一致的。在这次包含500次观测的实验中，得到的最大e值为2.188，这也是符合预期结果的。

      步骤3：将mSPRT方法应用于真实数据集

      现在，我们将这种检测方法应用到那些确实存在实际干预效果的人工生成数据上。你会逐日计算e值，并找出第一个超过停止阈值的日期。

      import matplotlib.pyplot as plt

np.random.seed(42)
treated_shuffled = treated.copy()
control_shuffled = control.copy()
np.random.shuffle(treated_shuffled)
np.random.shuffle(control_shuffled)

USERS_PER_arm_PER_DAY = 60
N_days_RUN = 30
n_per_arm = USERS_PER_ARM_PER_DAY * N_days_RUN  # 即1,800

treated_seq = treated_shuffled[:n_per_arm]
control_seq = control_shuffled[:n_per_arm]

e_values = compute_evalue_running(treated_seq, control_seq)

ALPHA = 0.05
THRESHOLD = 1 / ALPHA  # 即20

days = np.arange(1, len(e_values) + 1) / USERS_PER_ARM_PER_DAY
cross_indices = np.where(e_values >= THRESHOLD)[0]
if len(cross_indices) > 0:
    stopping_day = days[cross_indices[0]]
    print(f"mSPRT方法确定的停止日期为：{stopping_day:.1f}")
    print(f"在停止日期时得到的e值为：{e_values[cross_indices[0]]:.1f}")
else:
    stopping_day = None
    print("在这段时间内，mSPRT方法并未超过停止阈值")

print(f"在第N_days_RUN天得到的最终e值为：{e_values[-1]:.2f}")

      预期结果：

      mSPRT停止日期：25.9
停止时的E值：20.9
第30天的最终E值：75.64

      具体过程如下：》你会随机打乱实验组与对照组的顺序，以模拟用户每日到达的随机性（在实际实验中，用户的到来顺序通常是随机的），然后将每组前1,800个观测数据逐个输入到compute_evalue_running函数中进行分析。在第25.9天时，E值超过了20这一阈值，这意味着你可以在提前4天的时候就终止实验，而且所得到的分析结果仍然是有效的。到了第30天，E值已经上升到了75.64，远高于这个阈值。

      82ae7b10-c598-4597-80e6-375fa76b209d

      图2：在真实的50,000用户合成数据集上，mSPRT的E值变化轨迹（第1组为实验组，第2组为对照组）。蓝线表示对数尺度下的E值变化情况；红色虚线表示阈值1/α=20。

      绿色的虚线标记了第25.9这一天，也就是E值首次超过阈值的时刻。下方的图表显示，随着数据的积累，两组任务的完成率逐渐趋于一致。与图1中的模拟结果不同，这些数据来源于真实的共享数据集，因此实验组确实取得了4.85个百分点的效果提升。

      步骤4：与固定样本检验的效能进行比较

      mSPRT方法确实会带来一定的实际成本。当该方法能够产生显著效果时，它可以让实验提前终止；然而，如果实际效果低于预期，或者样本量较小，那么这种方法的效能就会大打折扣。这个模拟实验准确地量化了这种权衡关系。

      from scipy.stats import ttest_ind

np.random.seed(42)

N_SIMS = 1000
TRUE_EFFECT = 0.05
BASE_RATE = 0.60
N_PER_arm = 1800          # 每组30天，每天60名用户
DAILY_BATCH = 60
THRESHOLD = 20

msprt_stopping_days = []
ms prt Detected = 0
ttest_detected = 0

for sim in range(N_SIMS):
    t_obs = np.random.binomial(1, BASE_RATE + TRUE_EFFECT, N_PER_arm)
    c_obs = np.random.binomial(1, BASE_rate, N_PER.arm)

    e_vals = compute_evalue_running(t_obs, c_obs)
    days = np.arange(1, N_PER_ARM + 1) / DAILY_BATCH
    crosses = npWHERE(e_vals >= THRESHOLD)[0]
    if len(crosses) > 0:
        msprt Detected += 1
        ms prt_stopping_days.append(days[crosses[0]])
    else:
        msprt_stopping_days.append(30.0)

    _, p = ttest_ind(t_obs, c_obs)
    if p < 0.05:
        ttest_detected += 1

msprt_power = msprt Detected / N_SIMS
ttest_power = ttest Detected / N_simS
median_stop = np.median(ms prt_stopping_days)
pct_stopped_early = np.mean(np.array(msprt_stopping_days) < 30.0)

print(f"mSPRT方法的效能：{msprt_power:.1%}")
print(f"固定样本t检验的效能：{ttest_power:.1%}")
print(f>mSPRT方法的平均停止日期：{median_stop:.1f} / 30")
print(f>提前终止实验的比例：{pct_stopped_early:.1%}")

      预期输出：

      mSPRT的检测功效：34.7%  
固定样本t检验的检测功效：88.7%  
mSPRT停止试验的中位天数：第30天/第30天  
提前终止试验的比例：49.3%

      具体过程如下：你进行了1,000次模拟实验，其中实际存在的效应为5个百分点。对于mSPRT而言，会计算其连续检测值的e值，并记录下首次超过20的值出现的那一天；而对于固定样本检验来说，则在第30天结束时进行一次检测。实验结果显示，mSPRT在49.3%的实验中能够检测出这一效应，而固定样本检验的这一比例为88.7%。当效应大小为5个百分点，且每组数据有1,800个观测值时，mSPRT需要大约两倍多的观测次数才能达到与固定样本检验相同的检测功效。
      这就是“始终有效”的保证所需要付出的代价。当你每天进行检测时，mSPRT能够有效地控制I型错误的发生率；而如果使用固定样本检验，每天进行一次检测的话，其误报率会高达30.2%；但无论你在何时停止检测，mSPRT的误报率始终保持在5%左右。
      对于你的团队来说，究竟哪种方法代价更高，取决于哪种方法更不利于你完成研究任务：是让实验持续更长时间，还是让错误的结论被错误地传播出去。大多数团队在亲自进行这种模拟实验之前，都会低估检测功效带来的成本。

      根据真实情况验证效果

      np.random.seed(0)
t_full = treated_shuffled
c_full = control_shuffled[:len(t_full)]

e_full = compute_evalue_running(t_full, c_full)
days_full = np.arange(1, len(e_full) + 1) / USERS_PER_ARM_PER_DAY

cross_full = np.where(e_full >= THRESHOLD)[0]
if len(cross_full) > 0:
    print(f"mSPRT正确检测出了这一效应。)")
    print(f"本可以在第{days_full[cross_full[0]]:.1f}天停止实验。)")
    print(f"数据中的真实效应值为：{treated.mean() - control.mean():.4f}")
    print(f"在停止试验时的e值为：{e_full[cross_full[0]]:.1f}")
else:
    print("使用这些数据，mSPRT未能达到检测阈值。)")

      预期输出：

      mSPRT正确检测出了这一效应。
本可以在第27.1天停止实验。
数据中的真实效应值为：0.0485
在停止试验时的e值为：22.2

      具体过程如下：在对全部随机排列的数据进行检测时，mSPRT在第27.1天时其连续检测值的e值达到了检测阈值。数据中真实的效应值为4.85个百分点，这一数值与预设的5个百分点非常接近，因此mSPRT正确地检测出了这一效应；而设计为仅持续30天的固定样本检验，则会强制你在第30天才停止实验。如果每天每组数据有60名用户参与实验，那么mSPRT会允许你在第27.1天就结束实验，这样就能节省近3天的时间——毕竟这个特征最终还是会被检测出来的。

      步骤5：使用自助法计算置信区间

      停试日可以告诉你何时应该终止实验，但它无法告诉你这种效应的大小，也无法说明对该效应的估计值究竟有多精确。而自助法置信区间则能够同时提供这两方面的信息。rng = np.random.default_rng(7)
      point_est = treated.mean() - control.mean()

      boot_diffs = np.array([
      rng.choice(treated, size=len(treated), replace=True).mean() -
      rng.choice(control, size=len(control), replace=True).mean()
      for _ in range(500)
      ])

      lower = float(np.percentile(boot_diffs, 2.5))
      upper = float(nppercentile(boot_diffs, 97.5))

      print(f"点估计值（处理组 - 对照组）：{point_est:.4f} ({point_est*100:.2f}个百分点)"
      print(f"95% 自助法置信区间：[{lower:.4f}, {upper:.4f}] "
      f"([{lower*100:.2f}个百分点, {upper*100:.2f}个百分点])")
      print(f"根据真实情况，如果 lower ≤ 0.05 ≤ upper，则95%置信区间包含真实效应；否则不包含。"

      预期输出：
      点估计值（处理组 - 对照组）：0.0485 (4.85个百分点)
      95% 自助法置信区间：[0.0407, 0.0581] ([4.07百分点, 5.81个百分点])
      根据真实情况，95%置信区间包含真实效应。

      具体原理如下：你分别对处理组和对照组的样本进行500次有放回的重采样，每次计算两组均值之间的差异。这500个差异值的第2.5百分位数和第97.5百分位数构成了置信区间。该置信区间的范围为4.07百分点至5.81个百分点，因此它既能够确认这种效应确实存在，同时也表明这种效应的大小介于4.07百分点与5.81个百分点之间。由于每组有25,000名用户参与实验，因此这个置信区间具有较高的可靠性，它既能让你知道“这种干预措施是否有效”，也能告诉你“效果的具体幅度”。

      当mSPRT方法失效时

      即便采用了序贯检验的方法，实验设计仍然需要严谨性。在以下四种情况下，该方法的可靠性会受到破坏，或者使其在实际应用中变得毫无意义。

      先验分布设定错误

      mSPRT方法假设每组实验的完成率遵循Beta(1,1)先验分布，这一建模假设具有实际影响。但如果真实的效应值远远超出了该先验分布所预期的范围，那么这一假设就会失效。

      当真实效应值在3至10个百分点之间，且基础概率约为60%时，均匀分布的Beta(1,1)先验分布能够较为准确地估计效应大小。但如果真实效应仅为0.3个百分点，对于一个较为微小的AI功能改进来说，这种先验分布会导致计算得到的e值增长极其缓慢，从而使得在达到统计显著性阈值之前就耗尽所有样本资源。

      为了确保先验分布的合理性，你应该根据产品历史上进行的A/B测试数据来对其进行调整：利用最大似然法拟合历史数据中的效应大小分布，从而确定合适的Beta超参数；同时要验证所选先验分布是否能够将你的研究目标——即最小可检测效应值——置于其重点关注的范围之内。

      非平稳结果

      mSPRT方法的有效性要求在零假设条件下，e值过程必须是一个非负的超级鞅过程，而这意味着数据生成过程必须是平稳的。如果你的AI模型在实验进行过程中发生了变化，或者用户群体构成发生了改变（例如由于营销活动导致第12天的用户群体与之前不同），又或者任务难度存在星期几差异，这些因素都会导致e值受到环境噪声的影响，从而使你无法准确区分这种环境噪声与处理效应之间的真实差异。

      通过在保留样本上运行e值计算方法来诊断数据是否非平稳：如果原本应该接近1的零假设e值呈上升趋势，那就说明你的实验环境不够稳定，导致该方法无法产生可靠的结果。

      未进行多重性校正的多指标分析

      mSPRT方法能够控制单次比较中的I类错误。当你检测20个指标时，该方法本身并不会出现故障，因此每个单独计算的e值仍然有效。但问题在于整体错误率会升高：如果同时对20个指标应用mSPRT，并且一旦有某个指标的e值超过20就停止计算，那么出现至少1个假阳性的概率将会远远超过5%。

      为了解决多重性问题，可以在α=0.05、m=20的情况下，将阈值上调至1/(α/m) = 400；或者在实验结束后，对最终的e值应用Benjamini-Hochberg校正方法。

      多重性问题是与固定样本量检验中会遇到的问题完全相同的。mSPRT方法既不会使这个问题更加严重，也无法彻底解决它。这是一个值得明确指出的常见误解。

      最低检测时间要求仍然是有效的

      由于无论何时进行检测，这一保证都依然有效，因此人们可能会想要立即开始监测过程。但实际上并不应该这样做。虽然这个保证在任何时候都成立，但由于检验功效较低，即使真实存在效应，该方法也很少能够将其检测出来。

      第4步中的模拟实验清楚地说明了这一点：当每个组有1,800个观测值，效应大小为5个百分点时，mSPRT方法的检验功效仅为49.3%。在开始使用mSPRT进行监测之前，应先使用标准的功效计算工具，根据预期的效应大小计算出达到80%检测功效所需的最低样本量，并将这个数值作为起始检测的基准。只有当样本量达到这一要求后，才能开始计算e值。

      下一步该怎么做

      首先将mSPRT方法应用于你的主要评估指标，并确保使用的样本量至少能达到预期效应大小对应的80%检测功效。

      在开始正式实验之前，先使用历史保留数据进行A/A检验。这种校准检查不会耗费任何资源，而且能够在数据非平稳性影响实际实验之前及时发现问题。那些跳过这一步骤的团队，往往会在实际操作中发现校准错误，而这时再发现问题就为时已晚了。

      有关包含自助法置信区间在内的完整实现方法，请参阅配套代码库中的07_sequential_msprt/文件。

Tags:
Comments are closed.