scipy多项式pmf返回值nan

Question

我正在尝试使用（python）中的multinominal.pmf函数scipy.stats。

当我在输入中所有概率均大于零的情况下使用此函数时，它可以正常工作。问题是当我要使用其中一个概率为零的函数时。

以下示例说明了我的意思：

In [18]: multinomial.pmf([3, 3, 0], 6, [1/3.0, 1/3.0, 1/3.0])
Out[18]: 0.027434842249657095

In [19]: multinomial.pmf([3, 3, 0], 6, [2/3.0, 1/3.0, 0])
Out[19]: nan

可以看出，在第一次所有概率均大于0的情况下，使用该函数没有问题。但是，当我将其中一个概率更改为零时，该函数返回nan，即使通过该函数也应返回0.21948。

当其中一个概率为零时，是否有一种方法（在python中）来计算pmf？可以通过其他方式处理该问题，或者解决该功能。

附加信息

我在示例中使用Matlab中的mnpdf函数计算了示例中的函数应返回的值。但是，由于我的其余代码在python中，所以我更喜欢找到一种在python中进行计算的方法。

Answer 1

好地方！这是Scipy中的错误。源代码可以在这里找到。

3031至3051行：

def pmf(self, x, n, p):
    return np.exp(self.logpmf(x, n, p))

2997行到3017行：

def logpmf(self, x, n, p):
    n, p, npcond = self._process_parameters(n, p)

2939至2958行：

def _process_parameters(self, n, p):

    p = np.array(p, dtype=np.float64, copy=True)
    p[...,-1] = 1. - p[...,:-1].sum(axis=-1)

    # true for bad p
    pcond = np.any(p <= 0, axis=-1)  # <- Here is why!!!
    pcond |= np.any(p > 1, axis=-1)

    n = np.array(n, dtype=np.int, copy=True)

    # true for bad n
    ncond = n <= 0

    return n, p, ncond | pcond

行pcond = np.any(p <= 0, axis=-1)结果pcond是true如果任何值p<= 0。

然后在logpmf3029行中：

return self._checkresult(result, npcond_, np.NAN)

结果logpmf并pmf返回nan！

请注意，正确计算了实际结果（3020行，2994-2995）：

result = self._logpmf(x, n, p)

def _logpmf(self, x, n, p):
    return gammaln(n+1) + np.sum(xlogy(x, p) - gammaln(x+1), axis=-1)

用您的价值观：

import numpy as np
from scipy.special import xlogy, gammaln

x = np.array([3, 3, 0])
n = 6
p = np.array([2/3.0, 1/3.0, 0])

result = np.exp(gammaln(n+1) + np.sum(xlogy(x, p) - gammaln(x+1), axis=-1))
print(result)

>>>0.219478737997