利用倾向得分进行产品实验：基于Python的大语言模型特征的因果推断

如果其中任何一个假设不成立，即使另外两个假设成立，估计结果也会出现偏差。本教程的最后会详细分析这些可能导致的错误情况。

先决条件

你需要使用Python 3.11或更高版本，并且需要对pandas和scikit-learn有基本的了解，同时还需要对逻辑回归有一定的熟悉程度。

请安装本教程所需的软件包：

pip install numpy pandas scikit-learn matplotlib

具体工作流程如下：这四个软件包共同完成了整个数据分析流程。pandas用于加载数据，NumPy负责处理权重计算及数组运算，scikit-learn用于构建倾向得分模型并执行匹配操作，而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

具体操作步骤如下：克隆代码仓库后，运行generate_data.py命令即可生成本教程中会用到的合成数据集。使用参数“seed 42”可以确保数据集的可重复性，而50,000名用户的样本量可以为本教程中的各种分析方法提供可靠的数据支持。最终生成的CSV文件会被保存在data/synthetic_llm_logs.csv路径下。

设置测试示例

这个合成数据集模拟了一个SaaS产品，在该产品中，用户可以选择使用成本更高的模型。在拥有五万名用户的情况下，不同用户群体的选择比例存在显著差异：重度使用者的选择比例为65%，中度使用者为35%，而轻度使用者的选择比例仅为12%。

数据生成器中设定的真实因果效应表现为：选择使用高级模型的用户的任务完成率会提高8个百分点。但如果进行简单的比较，这一数值会被夸大到约21个百分点，因为选择偏差使得被纳入分析的群体中包含了那些参与度最高的用户。

只有了解真实的因果关系，才能验证你所使用的倾向得分方法是否能够准确反映这一关系。

现在来加载数据，并观察其中存在的选择偏差问题：

import pandas as pd

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

print(df.groupby("engagement_tier").opt_in_agent_mode.mean().round(3))

naive_effect = (
    df[df.opt_in_agent_mode == 1].task_completed.mean()
    - df[df/opt_in_agent_mode == 0].task_completed.mean()
)
print(f"\n朴素选择效应：{naive_effect:+.4f}")

预期输出结果：

engagement_tier
heavy     0.647
light     0.120
medium    0.353
Name: opt_in_agent_mode, dtype: float64

Naive opt-in effect: +0.2106

实际分析结果如下：你加载了50,000条数据，按照用户参与程度对它们进行分组，然后计算了每组中选择使用高级模型的用户的比例。结果显示，重度使用者的选择比例远高于轻度使用者，这正是数据中存在的选择偏差现象。朴素方法计算出的效应值为+0.2106，几乎是真实值+0.08的三倍。正是这种差距，需要通过倾向得分方法来消除。

步骤1：估算倾向得分

倾向得分是一种模型输出的结果，该模型能够根据可观察的特征预测用户是否会选择使用某种功能。逻辑回归是一个合适的起点，因为它既易于解释，计算速度也很快；但在步骤4中要关注平衡诊断结果——如果任何加权SMD值超过0.1，那就说明逻辑回归模型忽略了某些交互作用因素，此时应该改用梯度提升算法。

对于这个数据集来说，相关的可观察特征包括用户参与程度和查询的置信度。在现实产品中，你还可以纳入其他可能影响用户选择行为的变量，比如设备类型、使用时长、套餐等级以及历史使用习惯等。

from sklearn.linear_model import LogisticRegression
from sklearn.metrics import roc_auc_score

X = pd.get_dummies(
    df[["engagement_tier", "query_confidence"],
    drop_first=True
).astype(float)
y_treat = df.opt_in_agent_mode

ps_model = LogisticRegression(max_iter=1000).fit(X, y_treat)
df["propensity"] = ps_model.predict_proba(X)[:, 1]

# 基本验证步骤
print(df.groupby("engagement_tier").propensity.mean().round(3))
print(
    f"\n处理组倾向得分范围：  "
    f"{df[df.opt_in_agent_mode == 1].propensity.min():.3f} - "
    f"{df[df/opt_in_agent_mode == 1].propensity.max():.3f}"
)
print(
    f"对照组倾向得分范围：  "
    f"{df[df(opt_in_agent_mode == 0).propensity.min():.3f} - "
    f"{df[df.opt_in_agent_mode == 0].propensity.max():.3f}"
)
print(f"倾向得分模型的AUC值：{roc_auc_score(y_treat, df.propensity):.3f}")

预期输出：

engagement_tier
heavy     0.646
light     0.120
medium    0.353
名称：倾向得分，数据类型：float64

倾向得分范围（实验组）：0.114 – 0.675
倾向得分范围（对照组）：0.114 – 0.673
倾向得分模型的AUC值：0.744

具体处理过程如下：你将用户参与程度划分为不同的类别，并将这些类别转换为虚拟变量；同时保持查询结果的置信度为连续值，然后使用逻辑回归模型进行拟合。该模型预测出的概率值即为每个用户的倾向得分。

Scikit-learn中的`LogisticRegression`默认会应用L2正则化（`C=1.0`），这种正则化方式会使倾向得分的值略微向0.5靠近。如果在实际生产环境中使用，你可以将`penalty=None`设置为无正则化处理。

在每个参与程度类别中，计算出的平均倾向得分几乎能够准确反映该类别中用户的真实选择比例，因此该模型已经经过了校准。0.744的AUC值说明该模型在区分用户是否选择参与时，其表现远优于随机猜测的结果（0.5）。

此外，实验组与对照组的倾向得分范围存在重叠（两者的范围都大致在0.11到0.67之间），这正是我们所要求的视觉上的重叠条件。

图2：在同一包含50,000个用户的合成数据集上进行的两次正向性检验结果。

在上图中，顶部面板展示了各组用户倾向得分的平滑核密度曲线。三个峰值分别与三个参与程度类别相对应（轻度参与的用户的倾向得分约为0.12，中度参与的约为0.35，重度参与的约为0.65），这一结果符合预期，因为用户的参与程度确实是由其所处类别决定的。底部面板将同样的分布转换为各类别中的用户数量：每个类别中都包含数千名选择参与和未选择参与的用户，这正是进行正向性检验所需要的数据结构。

虽然图1只是以示意图的形式说明了这一概念，但该图证实这一概念在真实数据中也同样成立，因此后续进行的加权处理和匹配操作将能够基于实际的数据结果来进行。

步骤2：逆概率加权

逆概率加权方法会根据用户的倾向得分为其分配一个权重，这个权重的大小与用户的倾向得分成反比。例如，一个倾向得分为0.12的用户属于罕见情况（即尽管参与程度较低但仍选择参与的用户），因此他的权重相当于1 / 0.12 ≈ 8，这意味着他的数据能够代表群体中大约8个类似用户的情况；而一个倾向得分为0.12的对照组用户，则属于典型情况，因此他的权重为1 / (1 – 0.12) ≈ 1.14。

import numpy as np

# ATE加权：实验组的权重为1/P(treat)，对照组的权重为1/P(no treat)
df["ipw"] = np.where(
    df.opt_in_agent_mode == 1,
    1 / df.propensity,
    1 / (1 - df.propensity)
)

t = df[df(opt_in_agent_mode == 1]
c = df[df/opt_in_agent_mode == 0]
ate_ipw = (
    (t.task_completed * t.ipw).sum() / t.ipw.sum()
    - (c.task_completed * c.ipw).sum() / c.ipw.sum()
)
print(f"IPW平均处理效应（ATE）：{ate_ipw:+.4f}")

# ATT：选择参与的用户实际获得的收益
df["ipw_att"] = np.where(
    df.opt_in_agent_mode == 1,
    1,
    df.propensity / (1 - df.propensity)
)
t = df[df(opt_in_agent_mode == 1]   # 在计算ATT时重新划分数据集
c = df[df/opt_in_agent_mode == 0]
treated_mean = t.task_completed.mean()
control_w_mean = (c.task_completed * c.ipw_att).sum() / c.ipw_att.sum()
att_ipw = treated_mean - control_w_mean
print(f"IPW对实验组的平均处理效应（ATT）：{att_ipw:+.4f}")

预期输出：

IPW平均处理效应：+0.0851
接受处理的用户的IPW平均处理效应：+0.0770

具体计算过程如下：首先，为每位用户计算ATE权重，然后比较选择参与实验的用户与未选择参与实验的用户在任务完成度上的差异。接着计算ATT权重，这些权重仅用于重新调整对照组的分布，使其与接受处理用户的协变量分布相匹配，从而计算出接受处理用户的平均处理效应。
ATE用于回答“对于那些可能选择参与实验、也可能没有选择的随机用户来说，这一措施会产生什么影响？”这个问题；而ATT则用于回答“选择参与实验的用户实际上获得了哪些收益？”这个问题。在这个数据集上，ATE的值为+0.0851，ATT的值为+0.0770，这两个数值都接近真实值+0.08，并且远高于简单的+0.2106这一估计结果。
在实际应用中，这两种方法的区别非常重要。决定是否向未选择参与实验的用户推出某项功能需要使用ATE来进行评估；而要了解选择参与实验的用户所获得的实际收益，则需要使用ATT。

步骤3：最近邻匹配

最近邻匹配采用另一种计算方法：将每位选择参与实验的用户与倾向得分最接近的未选择参与实验的用户进行配对，然后计算这些配对对之间的结果差异平均值，从而得到ATT的估计值。

from sklearn.neighbors import NearestNeighbors

treated_ps = df[df.opt_in_agent_mode == 1][["propensity"]].values
control_ps = df[df(opt_in_agent_mode == 0][["propensity"]].values)

nn = NearestNeighbors(n_neighbors=1).fit(control_PS)
_, idx = nn.kneighbors(treated_ps)

treated_outcomes = df[df.opt_in_agent_mode == 1].task_completed.values
matched_control_outcomes = (
    df[df/opt_in_agent_mode == 0][["task_completed"]].values[idx.flatten()]
)

att_match = (treated_outcomes - matched_control_outcomes).mean()
print(f"1-NN匹配得到的ATT值：{att_match:+.4f}")

预期输出：

1-NN匹配得到的ATT值：+0.0752

具体计算过程如下：首先提取每组用户的倾向得分，然后对对照组使用最近邻匹配算法进行建模，从而为每位接受处理的用户找到倾向得分最接近的对照用户。由于NearestNeighbors算法允许同一个对照用户被多个接受处理的用户所选中，因此这是一种带有替换机制的匹配方法。
接下来，获取每位接受处理用户的实际结果以及与其配对的对照用户的结果，计算每对数据之间的差异值，并对这些差异值进行平均处理，从而得到ATT的估计值。这个数值反映了选择参与实验的用户相比那些未选择参与实验但情况相似的用户所获得的实际收益。
最终得到的+0.0752这一结果接近真实值+0.08，但略低于使用IPW方法计算得出的ATT值。这是1-NN匹配方法的典型特征，因为使用最近邻匹配得到的估计值通常具有较高的方差。

值得注意的是，这种匹配方法还存在两种变体。其中一种变体允许同一个对照用户被多个接受处理的用户所选中，这样可以在缺乏合适的对照对象时减少偏差，但也会导致方差的增加。

不进行替换的匹配方式会确保每个对照组用户至多只与一个实验组用户配对，这样虽然能降低方差，但当实验组的人数远超过可用的对照组人数时，就会导致配对质量下降。

对于大多数实际分析来说，使用k=3-5的k最近邻匹配方法，并允许进行替换，是一个合理的默认选择。

步骤4：检查协变量平衡性

只有当倾向得分方法真正实现了各组之间协变量的平衡时，它们才能发挥应有的作用。你必须验证这一点，因为如果平衡性没有得到保障，你的估计结果就会出错。

判断协变量平衡性的标准方法是计算每个协变量的标准化均值差（SMD）。SMD是通过将实验组的平均值与对照组的平均值相减，然后除以合并后的标准差来得出的。

在加权之前，SMD可以反映原始分组之间的不平衡程度；而在加权之后，这些数值应该会变小（通常认为|SMD|小于0.1表示平衡性良好）。

def smd(treated_vals, control_vals, treated_w=None, control_w=None):
    """标准化均值差，可选项包括使用权重。"""
    if treated_w is None:
        treated_w = np.ones(len(treated_vals))
    if control_w is None:
        control_w = np.ones(len(control_vals))
    t_mean = np.average(treated_vals, weights=treated_w)
    c_mean = np_average(control_vals, weights=control_w)
    pooled_std = np.sqrt((treated_vals.var() + control_vals.var()) / 2)
    return (t_mean - c_mean) / pooled_std

engagement_heavy = (df.engagement_tier == "heavy").astype(float).values
qc = df.query_confidence.values
tr = (df.opt_in_agent_mode == 1).values

covariates = {
    "engagement_tier_heavy": engagement_heavy,
    "query_confidence": qc,
}

print(f"{'Covariate':<30} {'Raw SMD':>10} {'Weighted SMD':>15}")
for name, vals in covariates.items():
    smd_raw = smd(vals[tr], vals[~tr])
    smd_weighted = smd(
        vals[tr], vals[~tr],
        treated_w=df[tr].ipw.values,
        control_w=df[~tr].ipw.values,
    )
    print(f"{name:<30} {smd_raw:>+10.3f} {smd_weighted:>+15.3f}")

预期输出结果：

Covariate                         Raw SMD    Weighted SMD
engagement_tier_heavy              +0.742          +0.002
query_confidence                   -0.032          -0.003

具体计算过程如下：该辅助函数会为任意协变量计算标准化均值差，用户也可以选择使用权重进行计算。

随后，你会看到每个协变量的原始SMD值和加权后的SMD值。对于engagement_tier_heavy这一协变量来说，原始SMD值为+0.742，这意味着重度使用用户的参与比例远高于其他用户；而加权后的SMD值降到了+0.002，说明这种不平衡在加权处理后得到了消除。对于query_confidence这一协变量而言，其在原始数据上就已经接近平衡状态，加权处理后依然保持了这种平衡性。如果任何加权后的SMD值的绝对值超过了0.1，那就说明你的倾向得分模型存在问题；通常解决方法是在逻辑回归中添加更多的特征或交互项。

从视觉上看，上图已经证实了那些通过数值分析得出的结论：重叠条件确实成立，因此实现平衡是完全可能的。

步骤5：自助法置信区间

仅凭点估计值是无法完全说明问题的。任何你需要向产品团队报告的估算结果，都应当附带一个置信区间，以便他们能够判断+0.08这个数值与+0.03或+0.12这两个数值是否有显著差异。由于需要使用估算得到的倾向得分，因此对于IPW和匹配分析来说，计算标准误差会显得相当复杂；所以，最简单且最可靠的方法就是使用非参数自助法。

def estimate_all(sample):
    """在自助抽样样本的基础上返回(ATE_IPW, ATT_IPW, ATT_match)这三个数值。"""
    s = sample.copy()
    X_s = pd.get_dummies(
        s[["engagement_tier", "query_confidence"]], drop_first=True
    ).astype(float)
    ps = LogisticRegression(max_iter=1000).fit(X_s, s.opt_in_agent_mode)
    s["p"] = ps.predict_proba(X_s)[:, 1]

    s["w_ate"] = np.where(
        s.opt_in_agent_mode == 1, 1 / s.p, 1 / (1 - s.p)
    )
    s["w_att"] = npwhere(
        s.opt_in_agent_mode == 1, 1, s.p / (1 - s.p)
    )
    t, c = s[s(opt_in_agent_mode == 1], s[s.opt_in_agent_mode == 0])

    ate = (
        (t.task_completed * t.w_ate).sum() / t.w_ate.sum()
        - (c.task_completed * c.w_ate).sum() / c.w_ate.sum()
    )
    att = t.task_completed.mean() - (
        (c.task_completed * c.w_att).sum() / c.w_att.sum()
    )
    nn_b = NearestNeighbors(n_neighbors=1).fit(c[["p"]].values)
    _, idx_b = nn_b.kneighbors(t[["p"]].values)
    match = (
        t.task_completed.values
        - c.task_completed.values[idx_b.flatten()]
    ).mean()
    return ate, att, match

rng = np.random.default_rng(7)
n_reps = 500
results = np.zeros((n_reps, 3))
for i in range(n_reps):
    boot = df.iloc[rng.integers(0, len(df), size=len(df))]
    results[i] = estimate_all(boot)

for name, col in zip(["IPW ATE", "IPW ATT", "1-NN ATT"], range(3)):
    lo, hi = np.percentile(results[:, col], [2.5, 97.5])
    print(f"{name:<10} 95%置信区间: [{lo:+.4f}, {hi:+.4f}]")

预期输出结果：

IPW ATE    95%置信区间: [+0.0745, +0.0954]
IPW ATT    95%置信区间: [+0.0687, +0.0865]
1-NN ATT   95%置信区间: [+0.0659, +0.0940]

具体计算过程如下：首先对数据集进行500次有放回的重采样，然后重新拟合倾向得分模型，并针对每次重采样结果重新计算各项估算值；最后将自助法分布中第2.5百分位和第97.5百分位的数值作为95%置信区间。这三个置信区间都包含了真实值+0.08这个范围，而明显排除了+0.21这个不合理的数值。其中，IPW ATT的置信区间最为狭窄，因为这种估算方法仅对对照组进行了加权处理；而1-NN匹配方法的置信区间最为宽泛，因为这种匹配方式会排除那些不在匹配组中的对照用户。在笔记本电脑上运行这段代码大约需要90秒的时间。在撰写面向利益相关者的报告时，应该以点估计值作为主要内容，并同时列出置信区间，这样团队就能清楚地看到这个数值所包含的不确定性。

当倾向得分方法失效时

在假设成立的情况下，倾向得分能使基于用户自愿参与的比较变得严谨。然而，一旦这些假设不成立，它们就会产生有偏的估计结果。

四种常见的失效情况与前面提到的三个识别假设有关。

1. 未测量的混杂因素（违反无混杂性假设）

如果某些因素同时影响了用户是否自愿参与以及最终的结果，但这些因素并未被纳入倾向得分模型中，那么使用IPW或匹配方法所得到的估计结果就会产生偏差。这在实际应用中是最常见的失效原因。

举个例子：那些选择使用“代理模式”的用户，往往也是会关注你的技术博客并阅读发布说明的用户。如果阅读博客这一行为本身就能独立影响任务完成率，而这一因素被忽略的话，那么由此产生的效果就会被错误地归因于“代理模式”，从而导致估计结果失真。

应对这种问题的有效方法包括：深入了解哪些因素会促使用户选择自愿参与、在倾向得分模型中加入更多相关的特征变量，以及使用正式的敏感性分析工具（如Rosenbaum界限、E值）来评估未测量混杂因素的影响程度。

2. 正性重叠问题（违反正向重叠假设）

如果某些用户的自愿参与概率接近于零或接近于1，那么我们就无法为这些用户找到合适的对照组来进行比较。

在这种情况下，倾向得分方法会生成极端的权重值（例如1 / 0.001 = 1,000），从而导致某个极端值主导最终的估计结果。因此，匹配过程也会产生质量较差的配对结果。

如果存在极端值，建议先检查倾向得分的分布情况，并剔除那些位于[0.05, 0.95]区间之外的数据，然后再进行加权处理。

3. 指定错误的倾向得分模型（在实际应用中降低无混杂性）

线性逻辑回归模型无法捕捉非线性的关系。例如，如果用户是否自愿参与取决于其参与程度与查询复杂度之间的交互作用，那么仅使用主效应模型就无法准确反映这一关系，从而导致估计结果失真。

应使用更灵活的模型（比如基于倾向得分的梯度提升算法，或在加权之后再进行回归调整），并且务必在加权后检查最终的平衡性。如果加权后的平衡性不佳，那就说明模型设定有误。

4>用户之间的影响传递（违反独立处理单元假设）

倾向得分方法假定所有用户是相互独立的。但如果某个用户选择使用“代理模式”会影响到其他用户的任务完成情况（比如团队成员在共享工作空间中同时采用这一功能），那么所估计的效果就会包含这种影响传递效应。

这种情况违反了“独立处理单元假设”。要妥善处理这类问题，需要采用不同的方法：对于在工作空间层面采用的共同功能，可以使用集群随机化方法；而对于用户层面的影响传递现象，则需要设计基于网络结构的实验设计。

这些错误机制在回归系数中是无法被察觉的。它们表现为那些在理论上看起来很合理的估计值，但当这些功能真正推广给更广泛的用户群体时，这些问题就会暴露出来。

进行平衡性诊断分析，查看数据重叠图，并记录下你可能忽略的任何问题——这些才是你真正的防护措施。

下一步该怎么做

当你的功能是通过用户主动选择来使用的，并且你有丰富的协变量可以用于模型选择时，倾向得分方法就是最适合的工具。

如果用户的参与是按照明确的规则来决定的（比如根据查询复杂度设置门槛，或者需要支付费用才能使用），那么使用回归分析来处理这种差异会更为合适。如果你怀疑存在未被观察到的干扰因素，并且拥有外部随机化手段（例如随机部署实验或基于速率限制的路由机制），那么工具变量方法会更加有效。

为了防止由于倾向得分模型设定错误而导致的估计偏差，双重稳健估计方法结合了倾向得分加权与回归调整，只要这两个组成部分中的任何一个模型设定正确，最终的估计结果就会保持一致性。

本教程配套的笔记本文件可以在这里找到。你可以克隆这个代码仓库，生成合成数据集，然后运行psm_demo.ipynb（或psm_demo.py），从而重现本教程中的所有代码示例、计算结果和图表。

当某个人工智能功能是通过用户主动选择来使用的机制进行推广时，单纯使用“选择使用该功能的用户”与“未选择的用户”进行对比通常会得到错误的结果。倾向得分方法能为你提供“与那些选择了该功能的用户相似的用户群体”作为对照组，而自助法则能帮你得出一个区间范围，这样当利益相关者询问你的结论有多可靠时，你就可以用这个区间来回应他们。