Python3中寻找与R VGAM::pbetabinom()结果一致的等价函数

你得到差异结果的根本原因是两个接口的参数化规则完全不同，直接按位置把VGAM的prob、rho参数传给scipy的betabinom.cdf是错误的：

• VGAM的pbetabinom采用均值-离散参数化：输入prob是二项分布概率的均值，rho是衡量过离散程度的相关系数参数
• scipy的betabinom采用标准Beta分布形状参数化：第三、第四个入参分别是Beta分布的形状参数a（即α）、b（即β），和prob、rho没有直接对应关系，不能直接传值。

VGAM官方文档给出的两类参数转换关系为，当给定size = n、prob = μ、rho = ρ时，对应的Beta形状参数为：

α = μ * (1 - ρ) / ρ
β = (1 - μ) * (1 - ρ) / ρ

基于上述转换规则，封装调用scipy的接口即可得到和VGAMpbetabinom()完全一致的计算结果：

from scipy.stats import betabinom

def pbetabinom_vgam(q, size, prob, rho):
    """
    等价于R中VGAM::pbetabinom(q, size, prob, rho)的累积分布计算
    """
    alpha = prob * (1 - rho) / rho
    beta = (1 - prob) * (1 - rho) / rho
    return betabinom.cdf(q, n=size, a=alpha, b=beta)

# 测试样例
if __name__ == "__main__":
    result = pbetabinom_vgam(
        q=884.0,
        size=2425,
        prob=0.374491636284026,
        rho=9.192015e-05
    )
    print(result)  # 输出0.1849791，和R侧结果完全匹配

你之前的写法betabinom.cdf(884.0, 2425, 0.374491636284026, 9.192015e-05)等价于传入a=0.3745、b=9.19e-5，对应的Beta分布极度右偏，和目标参数下的分布形态差异极大，因此得到了0.000189的错误结果。

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