如何在PyMC3中实现带条件自定义似然函数的线性回归

首先修正初始代码中的变量名笔误，再按需求实现自定义似然逻辑，完整实现如下：

实现步骤 #

• 首先将指示变量delta和观测目标y也注册为模型数据
• 修正线性预测项mu的变量名错误
• 自定义对数似然函数，通过pm.math.switch实现分条件的似然规则：
  • delta=0时直接使用标准正态分布的对数似然（对应最小二乘规则）
  • delta=1时仅当观测值大于预测值时计入正态似然贡献，否则该样本不贡献似然
• 使用pm.Potential将自定义似然加入模型

完整代码 #

import pymc3 as pm
import numpy as np
import pandas as pd

# 假设你的数据已存入df变量
df = pd.DataFrame({
    "x1": [10, 20, 5, 10],
    "x2": [1, 2, -1, -2],
    "delta": [0, 0, 1, 1],
    "observed_target": [100, 50, 200, 100]
})

with pm.Model() as ttf_model:
    # 先验定义
    intercept = pm.Normal('param_intercept', mu=0, sd=5)
    beta_0 = pm.Normal('param_x1', mu=0, sd=5)
    beta_1 = pm.Normal('param_x2', mu=0, sd=5)
    std = pm.HalfNormal('param_std', beta=0.5)

    # 注册模型输入数据
    x_1 = pm.Data('var_x1', df['x1'])
    x_2 = pm.Data('var_x2', df['x2'])
    delta = pm.Data('var_delta', df['delta'])
    y_obs = pm.Data('var_y', df['observed_target'])

    # 线性预测项（修正原代码变量名笔误）
    mu = intercept + beta_0 * x_1 + beta_1 * x_2

    # 自定义对数似然
    def custom_logp(y, mu, std, delta):
        # 计算所有样本的标准正态对数似然
        normal_ll = pm.Normal.dist(mu=mu, sigma=std).logp(y)
        # 按规则筛选似然贡献
        ll = pm.math.switch(
            delta == 0,
            normal_ll,
            pm.math.switch(
                y > mu,
                normal_ll,
                0.0  # 该样本不计入似然贡献，对数似然为0
            )
        )
        return ll.sum()
    
    # 将自定义似然加入模型
    pm.Potential('custom_likelihood', custom_logp(y_obs, mu, std, delta))

    # 后续可正常执行采样、推断等操作
    trace = pm.sample(2000, tune=1000, cores=2)

注意事项 #

如果实际运行中出现参数收敛异常，可将delta=1且y<=mu时的对数似然值从0.0调整为极小的负数（如-1e3），避免模型无限制增大预测值mu来规避似然惩罚。

内容的提问来源于stack exchange，提问作者Mike Tauber