概率统计Python计算：离散型2-维随机向量的联合分布律及边缘分布

此处，我们假定2-维离散型随机向量(X,Y)(X, Y)(X,Y)的联合分布律为

即随机变量XXX取mmm个值，YYY取nnn个值，将(X,Y)(X, Y)(X,Y)的联合分布中的概率值构成一个m×nm\times nm×n的矩阵，记为PXYP_{XY}PXY，即
PXY=(p11p12⋯p1np21p22⋯p2n⋮⋮⋯⋮pm1pm2⋯pmn)P_{XY}=\begin{pmatrix} p_{11}&p_{12}&\cdots&p_{1n}\\ p_{21}&p_{22}&\cdots&p_{2n}\\ \vdots&\vdots&\cdots&\vdots\\ p_{m1}&p_{m2}&\cdots&p_{mn} \end{pmatrix}PXY=⎝⎜⎜⎜⎛p11p21⋮pm1p12p22⋮pm2⋯⋯⋯⋯p1np2n⋮pmn⎠⎟⎟⎟⎞

1. 联合分布律的表示

Python的scipy.stats包并未提供2-维分布，但numpy包的array数组类对象却能很好地表示这样的2-维离散型随机向量的联合分布律。
例1 从含有3个正品，2个次品的5个产品中依次无放回地抽取两个。设XXX表示第1次取到的次品个数，YYY表示第2次取到的次品个数。求(X,Y)(X, Y)(X,Y)的联合分布律。
解：显然，XXX和YYY的所有可能取值均为{0,1}\{0, 1\}{0,1}。由于是无放回抽取，
P(X=i,Y=j)=P(X=i)P(Y=j∣X=i),i,j=0,1P(X=i, Y=j)=P(X=i)P(Y=j|X=i), i, j=0, 1P(X=i,Y=j)=P(X=i)P(Y=j∣X=i),i,j=0,1
P(X=0,Y=0)=(3/5)(2/4)=3/10P(X=0, Y=0)=(3/5)(2/4)=3/10P(X=0,Y=0)=(3/5)(2/4)=3/10，P(X=0,Y=1)=(3/5)(2/4)=3/10P(X=0, Y=1)=(3/5)(2/4)=3/10P(X=0,Y=1)=(3/5)(2/4)=3/10，
P(X=1,Y=0)=(2/5)(3/4)=3/10P(X=1, Y=0)=(2/5)(3/4)=3/10P(X=1,Y=0)=(2/5)(3/4)=3/10，P(X=1,Y=1)=(2/5)(1/4)=1/10P(X=1, Y=1)=(2/5)(1/4)=1/10P(X=1,Y=1)=(2/5)(1/4)=1/10。
于是，(X,Y)(X, Y)(X,Y)的联合分布律中由概率构成的矩阵PXY=(310310310110)P_{XY}=\begin{pmatrix}\frac{3}{10}&\frac{3}{10}\\\frac{3}{10}&\frac{1}{10}\end{pmatrix}PXY=(103103103101)。下列代码在Python中表示PXYP_{XY}PXY。

import numpy as np                      #导入numpy
from sympy import Rational as R         #导入Rational
Pxy=np.array([[R(3,10), R(3,10)],       #创建2-维数组Pxy[R(3,10), R(1,10)]])
print(Pxy)                              #输出2-维数组

第3~4行创建一个名为Pxy的array类对象，将其设置为两个等长的数组的数组，从而构成一个矩阵。Pxy中的每一个元素设置为表示有理数的Rational对象（第2行导入，别名为R）。运行该程序，输出

[[3/10 3/10][3/10 1/10]]

2. 边缘分布的计算

由上例可见，numpy的array类对象可将矩阵表示为2-维数组——数组的数组。2-维数组有两个“轴”：纵向记为axis=0，横向记为axis=1，如下图所示。

为计算变量XXX及YYY的边缘分布律，可调用array类对象Pxy的sum函数，指定按行对列标j（axis=1）相加得到XXX的边缘分布律，这是一个具有mmm个元素的数组，记为PXP_XPX；按列对行标i（axis=0）相加得YYY的边缘分布律，是一个具有nnn个元素的数组，记为PYP_YPY。下列程序定义了按此方法根据联合分布律PXYP_{XY}PXY计算边缘分布律PXP_XPX和PYP_YPY的Python函数。

import numpy as np                  #导入numpy
def margDist(Pxy):                  #定义计算边缘分布的函数Px=Pxy.sum(axis=1)              #按行相加得X分布律Py=Pxy.sum(axis=0)              #按列相加得Y分布律return (Px.reshape(Px.size, 1), #返回Px，PyPy.reshape(1, Py.size))

程序中第2~6行定义了用于根据(X,Y)(X,Y)(X,Y)的联合分布律计算XXX和YYY的边缘分布的函数margDist。参数Pxy是组织为2-维数组的(X,Y)(X,Y)(X,Y)的联合分布律。第3行、第4行分别对Pxy按行相加和按列相加得到XXX、YYY的边缘分布律存于Px和Py。第5~6行将Px，Py作为返回值返回。需要提及的是，为将array类对象表示的1-维数组，设置为一个列向量或行向量，以便与2-维数组表示的矩阵进行统一的运算，要调用该数组的reshape函数。因此，第5行返回的Px为列向量(p1⋅p2⋅⋮pm⋅)\begin{pmatrix} p_{1\cdot}\\ p_{2\cdot}\\\vdots\\p_{m\cdot}\end{pmatrix}⎝⎜⎜⎜⎛p1⋅p2⋅⋮pm⋅⎠⎟⎟⎟⎞，Py为行向量(p⋅1,p⋅2,⋯,p⋅n)(p_{\cdot1},p_{\cdot2},\cdots,p_{\cdot n})(p⋅1,p⋅2,⋯,p⋅n)。
例2 计算例1中2-维随机向量(X,Y)(X, Y)(X,Y)的边缘分布律。
解：由例1知，(X,Y)(X, Y)(X,Y)的联合分布律中概率矩阵为PXY=(310310310110)P_{XY}=\begin{pmatrix}\frac{3}{10}&\frac{3}{10}\\\frac{3}{10}&\frac{1}{10}\end{pmatrix}PXY=(103103103101)。下列代码在Python中表示PXYP_{XY}PXY，按行对列标求和得XXX得边缘分布概率序列PX=(3525)P_X=\begin{pmatrix}\frac{3}{5}\\\frac{2}{5}\end{pmatrix}PX=(5352)，按列对行标求和得YYY的边缘分布概率序列PY=(35,25)P_Y=(\frac{3}{5},\frac{2}{5})PY=(53,52)。
下列代码计算例1的随机向量(X,Y)(X, Y)(X,Y)中XXX和YYY的边缘分布。

import numpy as np                  #导入numpy
from sympy import Rational as R     #导入Rational
Pxy=np.array([[R(3,10), R(3,10)],   #创建联合分布律Pxy[R(3,10), R(1,10)]])
Px, Py=margDist(Pxy)                #计算边缘分布律
print('Px:%s'%Px)
print('Py:%s'%Py)

第3~4行设置(X,Y)(X, Y)(X,Y)的联合分布律的概率矩阵Pxy。第5行调用上列程序定义的函数margDist，计算结果赋予Px，Py。运行程序，输出

Px:[[3/5][2/5]]
Py:[[3/5 2/5]]

写博不易，敬请支持：
如果阅读本文于您有所获，敬请点赞、评论、收藏，谢谢大家的支持！
返回《导引》