使用 scipy.optimize 和对数似然法查找 beta 二项分布的 alpha 和 beta答案

【问题标题】：Finding alpha and beta of beta-binomial distribution with scipy.optimize and loglikelihood使用 scipy.optimize 和对数似然法查找 beta 二项分布的 alpha 和 beta
【发布时间】：2019-02-03 16:46:40
【问题描述】：

如果 p（成功概率）在二项分布中具有形状参数 α > 0 和 β 的 beta 分布，则该分布是 beta-binomial > 0。形状参数定义了成功的概率。我想从 beta 二项分布的角度找到最能描述我的数据的 α 和 β 值。我的数据集players 包含有关命中数 (H)、击球数 (AB) 和转换 (H / AB) 的数据）很多棒球运动员。我在Beta Binomial Function in Python Beta Binomial Function in PythonJulienD 的回答的帮助下估计了 PDF

from scipy.special import beta
from scipy.misc import comb

pdf = comb(n, k) * beta(k + a, n - k + b) / beta(a, b)

接下来，我编写了一个我们将最小化的对数似然函数。

def loglike_betabinom(params, *args):
   """
   Negative log likelihood function for betabinomial distribution
   :param params: list for parameters to be fitted.
   :param args:  2-element array containing the sample data.
   :return: negative log-likelihood to be minimized.
   """

   a, b = params[0], params[1]
   k = args[0] # the conversion rate
   n = args[1] # the number of at-bats (AE)

   pdf = comb(n, k) * beta(k + a, n - k + b) / beta(a, b)

   return -1 * np.log(pdf).sum()

现在，我想写一个函数来最小化loglike_betabinom

 from scipy.optimize import minimize
 init_params = [1, 10]
 res = minimize(loglike_betabinom, x0=init_params,
                args=(players['H'] / players['AB'], players['AB']),
                bounds=bounds,
                method='L-BFGS-B',
                options={'disp': True, 'maxiter': 250})
 print(res.x)

结果是 [-6.04544138 2.03984464]，这意味着 α 是负数，这是不可能的。我的脚本基于以下 R-sn-p。他们得到 [101.359, 287.318]..

 ll <- function(alpha, beta) { 
    x <- career_filtered$H
    total <- career_filtered$AB
    -sum(VGAM::dbetabinom.ab(x, total, alpha, beta, log=True))
 }

 m <- mle(ll, start = list(alpha = 1, beta = 10), 
 method = "L-BFGS-B", lower = c(0.0001, 0.1))

 ab <- coef(m)

谁能告诉我我做错了什么？非常感谢您的帮助！

【问题讨论】：

你现在如何最小化你的损失函数？您是否编写了自己的方法，还是使用了包中的某些东西？无论哪种方式，细节是什么？
我用from scipy.optimize import minimize
@HJA24 什么是players？你能分享一下这些数据吗？这是我测试答案所缺少的部分。
您在哪里指定要传递给minimize 的bounds？
我删除了边界部分

标签： python scipy distribution beta binomial-theorem

【解决方案1】：

需要注意的一点是，对于数据集中的 n 和 k 的值，对数似然中的 comb(n, k) 在数值上可能表现不佳。您可以通过将comb 应用于您的数据并查看是否出现infs 来验证这一点。

修改事物的一种方法可能是按照https://stackoverflow.com/a/32355701/4240413 中的建议重写负对数似然，即作为伽玛函数的对数函数，如

from scipy.special import gammaln
import numpy as np

def loglike_betabinom(params, *args):

    a, b = params[0], params[1]
    k = args[0] # the OVERALL conversions
    n = args[1] # the number of at-bats (AE)

    logpdf = gammaln(n+1) + gammaln(k+a) + gammaln(n-k+b) + gammaln(a+b) - \
     (gammaln(k+1) + gammaln(n-k+1) + gammaln(a) + gammaln(b) + gammaln(n+a+b))

    return -np.sum(logpdf)

然后您可以使用

最小化对数似然

from scipy.optimize import minimize

init_params = [1, 10]
# note that I am putting 'H' in the args
res = minimize(loglike_betabinom, x0=init_params,
            args=(players['H'], players['AB']),
            method='L-BFGS-B', options={'disp': True, 'maxiter': 250})
print(res)

这应该会给出合理的结果。

如果您想进一步修改代码，可以查看How to properly fit a beta distribution in python? 以获取灵感。

【讨论】：