Python科学计算（第2版）
上QQ阅读APP看书，第一时间看更新
上一章目录下一章

3.4.1 连续概率分布

可以使用下面的语句获得stats模块中所有的连续随机变量：

连续随机变量对象都有如下方法：

●rvs：对随机变量进行随机取值，可以通过size参数指定输出的数组的大小。

●pdf：随机变量的概率密度函数。

●cdf：随机变量的累积分布函数，它是概率密度函数的积分。

●sf：随机变量的生存函数，它的值是1-cdf(t)。

●ppf：累积分布函数的反函数。

●stat：计算随机变量的期望值和方差。

●fit：对一组随机取样进行拟合，找出最适合取样数据的概率密度函数的系数。

下面以正规分布为例，简单地介绍随机变量的用法。下面的语句获得默认正规分布的随机变量的期望值和方差，我们看到默认情况下它是一个均值为0、方差为1的随机变量：

    stats.norm.stats( )
    (array(0.0), array(1.0))

norm可以像函数一样调用，通过loc和scale参数可以指定随机变量的偏移和缩放参数。对于正态分布的随机变量来说，这两个参数相当于指定其期望值和标准差，标准差是方差的算术平方根，因此标准差为2.0时，方差为4.0：

    X = stats.norm(loc=1.0, scale=2.0)
    X.stats( )
    (array(1.0), array(4.0))

下面调用随机变量X的rvs( )方法，得到包含一万次随机取样值的数组x，然后调用NumPy的mean( )和var( )计算此数组的均值和方差，其结果符合随机变量X的特性：

    x = X.rvs(size=10000) # 对随机变量取10000个值
    np.mean(x), np.var(x) # 期望值和方差
    (1.0043406567303883, 3.8899572813426553)

也可以使用fit( )方法对随机取样序列x进行拟合，它返回的是与随机取样值最吻合的随机变量的参数：

    stats.norm.fit(x) # 得到随机序列的期望值和标准差
    (1.0043406567303883, 1.9722974626923433)

在下面的例子中，计算取样值x的直方图统计以及累积分布，并与随机变量的概率密度函数和累积分布函数进行比较。

❶其中histogram( )对数组x进行直方图统计，它将数组x的取值范围分为100个区间，并统计x中的每个值落入各个区间的次数。histogram( )返回两个数组pdf和t，其中pdf表示各个区间的取样值出现的频数，由于normed参数为True，因此pdf的值是正规化之后的结果，其结果应与随机变量的概率密度函数一致。

❷t表示区间，由于其中包括区间的起点和终点，因此t的长度为101。计算每个区间的中间值，然后调用X.pdf(t)和X.cdf(t)计算随机变量的概率密度函数和累积分布函数的值，并与统计值比较。❸计算样本的累积分布时，需要与区间的大小相乘，这样才能保证其结果与累积分布函数相同。

    pdf, t = np.histogram(x, bins=100, normed=True)  ❶
    t = (t[:-1] + t[1:]) *
 0.5 ❷
    cdf = np.cumsum(pdf) *
 (t[1] - t[0]) ❸
    p_error = pdf - X.pdf(t)
    c_error = cdf - X.cdf(t)
    print "max pdf error: {}, max cdf error: {}".format(
        np.abs(p_error).max( ), np.abs(c_error).max( ))
    max pdf error: 0.0217211429624, max cdf error: 0.0209887986472

图3-11（左）显示了概率密度函数和直方图统计的结果，可以看出二者是一致的。下右图显示了随机变量X的累积分布函数和数组pdf的累加结果。

图3-11 正态分布的概率密度函数（左）和累积分布函数（右）

有些随机分布除了loc和scale参数之外，还需要额外的形状参数。例如伽玛分布可用于描述等待k个独立随机事件发生所需的时间，k就是伽玛分布的形状参数。下面计算形状参数k为1和2时的伽玛分布的期望值和方差：

    print stats.gamma.stats(1.0)
    print stats.gamma.stats(2.0)
    (array(1.0), array(1.0))
    (array(2.0), array(2.0))

伽玛分布的尺度参数θ和随机事件发生的频率相关，由scale参数指定：

    stats.gamma.stats(2.0, scale=2)
    (array(4.0), array(8.0))

根据伽玛分布的数学定义可知其期望值为kθ，方差为kθ2。上面的程序验证了这两个公式。

当随机分布有额外的形状参数时，它所对应的rvs( )、pdf( )等方法都会增加额外的参数来接收形状参数。例如下面的程序调用rvs( )对k=2、θ=2的伽玛分布取4个随机值：

    x = stats.gamma.rvs(2, scale=2, size=4)
    x
    array([ 2.47613445,  1.93667652,  0.85723572,  9.49088092])

接下来调用pdf( )查看与上面4个随机值对应的概率密度：

    stats.gamma.pdf(x, 2, scale=2)
    array([ 0.17948513,  0.18384555,  0.13960273,  0.02062186])

也可以先创建将形状参数和尺度参数固定的随机变量，然后再调用其pdf( )计算概率密度：

    X = stats.gamma(2, scale=2) 
    X.pdf(x)
    array([ 0.17948513,  0.18384555,  0.13960273,  0.02062186])
上一章目录下一章