拟合离散数据：负二项式、泊松、几何分布

Question

在 scipy 中，不支持使用数据拟合离散分布。我知道这方面有很多话题。

例如，如果我有一个如下所示的数组：

x = [2,3,4,5,6,7,0,1,1,0,1,8,10,9,1,1,1,0,0]

我无法申请这个数组：

from scipy.stats import nbinom
param = nbinom.fit(x)

但是我想问你最新的情况，有没有什么方法可以拟合这三个离散分布，然后选择最适合离散数据集的？

Answer 1

您可以使用Method of Moments 来适应任何特定的分布。

基本思想：获取经验的第一、第二等矩，然后从这些矩中推导出分布参数。

所以，在所有这些情况下，我们只需要两个时刻。让我们得到它们：

import pandas as pd
# for other distributions, you'll need to implement PMF
from scipy.stats import nbinom, poisson, geom

x = pd.Series(x)
mean = x.mean()
var = x.var()
likelihoods = {}  # we'll use it later

注意：我使用 pandas 而不是 numpy。这是因为 numpy 的 var() 和 std() 不适用 Bessel's correction，而 pandas 则适用。如果您有 100 多个样本，则应该不会有太大差异，但对于较小的样本，这可能很重要。

现在，让我们获取这些分布的参数。 Negative binomial 有两个参数：p、r。让我们估计它们并计算数据集的可能性：

# From the wikipedia page, we have:
# mean = pr / (1-p)
# var = pr / (1-p)**2
# without wiki, you could use MGF to get moments; too long to explain here
# Solving for p and r, we get:

p = 1 - mean / var  # TODO: check for zero variance and limit p by [0, 1]
r = (1-p) * mean / p

UPD： Wikipedia 和 scipy 使用 p 的不同定义，一种将其视为成功的概率，另一种将其视为失败的概率。因此，为了与 scipy 概念保持一致，请使用：

p = mean / var
r = p * mean / (1-p)

UPD 结束

计算可能性：

likelihoods['nbinom'] = x.map(lambda val: nbinom.pmf(val, r, p)).prod()

Poisson同理，只有一个参数：

# from Wikipedia,
# mean = variance = lambda. Nothing to solve here
lambda_ = mean
likelihoods['poisson'] = x.map(lambda val: poisson.pmf(val, lambda_)).prod()

Geometric distribution 也一样：

# mean = 1 / p  # this form fits the scipy definition
p = 1 / mean

likelihoods['geometric'] = x.map(lambda val: geom.pmf(val, p)).prod()

最后，让我们找到最合适的：

best_fit = max(likelihoods, key=lambda x: likelihoods[x])
print("Best fit:", best_fit)
print("Likelihood:", likelihoods[best_fit])

如果您有任何问题，请告诉我