自定义似然函数变量提前引用报错的解决方案咨询

首先，你遇到的UnboundLocalError是直接的代码顺序问题：计算s1_s时用到的p1_s在后面才定义，但核心问题是s1_s、a1_s、p1_s形成了循环依赖——s1_s依赖p1_s，a1_s依赖s1_s，p1_s又依赖a1_s，这是一个隐式不动点问题，不能按线性顺序计算，必须通过迭代求解满足所有等式的稳定值。

不需要手动推导解析解，用迭代法就能实现，以下是修改后的代码：

import numpy as np
import pandas as pd
from scipy.optimize import minimize

def solve_fixed_point(p_hat, p1_bl, p2_bl, p3_bl, lamda_s, df, tol=1e-8, max_iter=1000):
    # 初始化p1_s/p2_s/p3_s为基准模型的概率，作为迭代起点
    p1_s, p2_s, p3_s = p1_bl.copy(), p2_bl.copy(), p3_bl.copy()
    
    for _ in range(max_iter):
        # 计算当前s1_s/s2_s/s3_s
        s1_s = p1_s * p_hat + p1_bl * (1 - p_hat)
        s2_s = p2_s * p_hat + p2_bl * (1 - p_hat)
        s3_s = p3_s * p_hat + p3_bl * (1 - p_hat)
        
        # 计算a1_s/a2_s/a3_s
        a1_s = 4 * s1_s + 1 * s2_s + 6 * s3_s
        a2_s = 0 * s1_s + 2 * s2_s + 8 * s3_s
        a3_s = 10 * s1_s + 4 * s2_s + 0 * s3_s
        
        # 计算新的p1_s/p2_s/p3_s
        new_p1_s = np.exp(lamda_s*(a1_s-a2_s)) / (np.exp(lamda_s*(a1_s-a2_s)) + 1 + np.exp(lamda_s*(a3_s-a2_s)))
        new_p2_s = 1 / (np.exp(lamda_s*(a1_s-a2_s)) + 1 + np.exp(lamda_s*(a3_s-a2_s)))
        new_p3_s = np.exp(lamda_s*(a3_s-a2_s)) / (np.exp(lamda_s*(a1_s-a2_s)) + 1 + np.exp(lamda_s*(a3_s-a2_s)))
        
        # 检查收敛：最大绝对误差小于阈值
        max_diff = max(np.max(np.abs(new_p1_s - p1_s)), 
                       np.max(np.abs(new_p2_s - p2_s)), 
                       np.max(np.abs(new_p3_s - p3_s)))
        if max_diff < tol:
            break
        
        # 更新p值，进入下一轮迭代
        p1_s, p2_s, p3_s = new_p1_s, new_p2_s, new_p3_s
    else:
        # 迭代到最大次数仍未收敛，抛出警告
        print("警告：迭代未达到收敛阈值")
    
    return s1_s, s2_s, s3_s, p1_s, p2_s, p3_s

def lik(params):
    p_s = params[0]
    lamda_s = params[1]
    lamda_bl = params[2]
    # phi_bl在当前代码中未使用，注意检查是否遗漏逻辑
    p_hat = params[4]

    # 计算基准模型的a和p
    a1_bl = 4 * df["s1"] + 1 * df["s2"] + 6 * df["s3"]
    a2_bl = 0 * df["s1"] + 2 * df["s2"] + 8 * df["s3"]
    a3_bl = 10 * df["s1"] + 4 * df["s2"] + 0 * df["s3"]

    denom_bl = np.exp(lamda_bl*(a1_bl-a2_bl)) + 1 + np.exp(lamda_bl*(a3_bl-a2_bl))
    p1_bl = np.exp(lamda_bl*(a1_bl-a2_bl)) / denom_bl
    p2_bl = 1 / denom_bl
    p3_bl = np.exp(lamda_bl*(a3_bl-a2_bl)) / denom_bl

    # 计算基准模型的f_bl（避免修改全局df，用局部变量）
    f_bl = p1_bl * y1 + p2_bl * y2 + p3_bl * y3

    # 求解循环依赖的不动点
    s1_s, s2_s, s3_s, p1_s, p2_s, p3_s = solve_fixed_point(p_hat, p1_bl, p2_bl, p3_bl, lamda_s, df)

    # 计算s模型的f_s
    f_s = p1_s * y1 + p2_s * y2 + p3_s * y3

    # 按pcode计算乘积
    pcodes = df["pcode"].unique()
    ff_bl = np.array([np.prod(f_bl[df['pcode'] == i]) for i in pcodes])
    ff_s = np.array([np.prod(f_s[df['pcode'] == i]) for i in pcodes])

    pp = p_s * ff_s + (1 - p_s) * ff_bl
    li = np.log(pp)
    LL = np.sum(li)

    return -LL

# 注意：需要确保df、y1、y2、y3已提前定义
# 添加参数边界（比如概率类参数在0-1，lambda非负）
bounds = [(0, 1), (0, None), (0, None), (None, None), (0, 1)]
results = minimize(lik, np.array([0.6, 0.6, 0.1, 0.3, 0.2]), method='L-BFGS-B', bounds=bounds)
print(results)

关键修改说明： #

• 新增solve_fixed_point函数，用迭代法求解循环依赖的变量，直到各p值的变化小于收敛阈值
• 避免修改全局df，改用局部变量存储f_bl和f_s，防止迭代过程中产生副作用
• 添加了参数边界约束，符合变量的实际意义（比如p_s是概率，范围0-1；lamda_s/lamda_bl通常非负）
• 注意到原代码中phi_bl参数未被使用，需检查是否遗漏了相关逻辑

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