有没有方法用SciPy的随机变量进行算术运算？

Question

SciPy的stats模块里有一种叫“随机变量”的对象（他们称之为rv_frozen）。这个对象可以让我们很方便地绘制某种分布的随机变量的累积分布函数（cdf）。下面是一个非常简单的例子：

import scipy.stats as stats
n = stats.norm()
x = linspace(-3, 3)
y = n.cdf(x)
plot(x, y)

我在想，是否可以对这些随机变量进行一些基本的算术运算。下面这个例子其实是个美好的幻想（它并不能正常工作）。

du_list = [stats.randint(2, 5) for _ in xrange(100)]
du_avg = sum(du_list) / len(du_list)
x = linspace(0, 10)
y = du_avg.cdf(x)
plot(x, y)

这个美好的幻想例子应该能生成一个图，显示的是一个随机变量的累积分布函数，这个随机变量是100个独立同分布的随机变量的平均值，每个随机变量的值均匀分布在{2,3,4}这个集合里。

Answer 1

你可以手动计算这个问题。

这里的X是一个随机变量，它是多个随机变量Xi的总和，这些Xi是均匀分布在2到5之间的。你可以从这个X生成的分布中抽样，得到概率密度函数（pdf），然后再进行积分，得到累积分布函数（cdf）。

或者你也可以尝试找到这个问题的解析解。

可以查看一下Irwin-Hall分布和Math Stackexchange上相关的讨论。

Answer 2

完全符合你描述的方法并不存在。不同分布的累积分布函数（CDF）都是在 **/scipy/stats/distributions.py` 这个源文件中定义的。例如：

玻尔兹曼分布的累积分布函数（第7675行）：

def _cdf(self, x, lambda_, N):
    k = floor(x)
    return (1-exp(-lambda_*(k+1)))/(1-exp(-lambda_*N))

你可以先估计最大似然估计（MLE），然后再调用累积分布函数的方法，看看这个示例：

import scipy.stats as ss
unknown=np.random.normal(loc=1.1, scale=2.0, size=100)
Loc, Scale=ss.norm.fit_loc_scale(unknown) #making a MLE fit
unknown_cdf=lambda x: ss.norm.cdf(x, loc=Loc, scale=Scale) #the cdf of the MLE to the data
plt.plot(np.linspace(-10, 10), unknown_cdf(np.linspace(-10, 10)), '-')

Answer 3

我知道这个回答有点晚，但我想说一下，以防将来还有其他人需要这个信息。最近我也需要类似的功能，甚至考虑过扩展scipy的rv_discrete来实现这个功能，但后来我发现了PaCAL。

PaCAL是一个用于处理随机变量的Python软件包，它支持很多种分布，包括连续分布。它甚至还支持一些二元联合分布。这个软件包可以在PyPI上找到。不过，它只适用于Python 2.x。

补充一下：现在PaCAL的Github仓库也支持Python 3.x了。