SciPy:从 PMF 生成自定义随机变量 [英] SciPy: generating custom random variable from PMF

问题描述

我试图在 Python 中根据某个丑陋的分布生成随机变量.我有一个 PMF 的显式表达式，但它涉及一些产品，这使得获取和反转 CDF 令人不愉快(请参阅下面的代码以了解 PMF 的显式形式).

I'm trying to generate random variables according to a certain ugly distribution, in Python. I have an explicit expression for the PMF, but it involves some products which makes it unpleasant to obtain and invert the CDF (see below code for explicit form of PMF).

本质上，我试图通过其 PMF 在 Python 中定义一个随机变量，然后让内置代码完成从分布中采样的艰苦工作.如果 RV 的支持是有限的，我知道如何做到这一点，但这里的支持是可数无限的.

In essence, I'm trying to define a random variable in Python by its PMF and then have built-in code do the hard work of sampling from the distribution. I know how to do this if the support of the RV is finite, but here the support is countably infinite.

我目前尝试按照以下@askewchan 的建议运行的代码是:

The code I am currently trying to run as per @askewchan's advice below is:

import scipy as sp
import numpy as np

class x_gen(sp.stats.rv_discrete):
    def _pmf(self,k,param):
        num = np.arange(1+param, k+param, 1)
        denom = np.arange(3+2*param, k+3+2*param, 1)

        p = (2+param)*(np.prod(num)/np.prod(denom))

        return p

pa_limit = limitrv_gen()
print pa_limit.rvs(alpha,n=1)

但是，这在运行时返回错误:

However, this returns the error while running:

File "limiting_sim.py", line 42, in _pmf
    num = np.arange(1+param, k+param, 1)
TypeError: only length-1 arrays can be converted to Python scalars

基本上，np.arange() 列表似乎在 def _pmf() 函数内部不起作用.我不知道为什么.任何人都可以在这里启发我和/或指出修复方法吗?

Basically, it seems that the np.arange() list isn't working somehow inside the def _pmf() function. I'm at a loss to see why. Can anyone enlighten me here and/or point out a fix?

编辑 1: 清除了 askewchan 的一些问题，上面反映了编辑.

EDIT 1: cleared up some questions by askewchan, edits reflected above.

编辑 2: askewchan 使用阶乘函数提出了一个有趣的近似值，但我正在寻找更多的精确解决方案，例如我正在尝试使用 np.arange 的解决方案.

EDIT 2: askewchan suggested an interesting approximation using the factorial function, but I'm looking more for an exact solution such as the one that I'm trying to get work with np.arange.

推荐答案

你应该能够像这样子类 rv_discrete :

You should be able to subclass rv_discrete like so:

class mydist_gen(rv_discrete):
    def _pmf(self, n, param):
        return yourpmf(n, param)

然后您可以使用以下命令创建分发实例:

Then you can create a distribution instance with:

mydist = mydist_gen()

并使用以下方法生成样本:

And generate samples with:

mydist.rvs(param, size=1000)

或者您可以使用以下命令创建一个冻结的分发对象:

Or you can then create a frozen distribution object with:

mydistp = mydist(param)

最后生成样本:

mydistp.rvs(1000)

<小时>

在您的示例中，这应该可以工作，因为 factorial 会自动广播.但是，对于足够大的 alpha，它可能会失败:

---

With your example, this should work, since factorial automatically broadcasts. But, it might fail for large enough alpha:

import scipy as sp
import numpy as np
from scipy.misc import factorial

class limitrv_gen(sp.stats.rv_discrete):
    def _pmf(self, k, alpha):
        #num = np.prod(np.arange(1+alpha, k+alpha))
        num = factorial(k+alpha-1) / factorial(alpha)
        #denom = np.prod(np.arange(3+2*alpha, k+3+2*alpha))
        denom = factorial(k + 2 + 2*alpha) / factorial(2 + 2*alpha)

        return (2+alpha) * num / denom

pa_limit = limitrv_gen()
alpha = 100
pa_limit.rvs(alpha, size=10)

这篇关于SciPy:从 PMF 生成自定义随机变量的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！