随机变量的函数

概率论08 随机变量的函数

作者：Vamei 出处：http://www.cnblogs.com/vamei 欢迎转载，也请保留这段声明。谢谢！

 

随机变量的函数

在前面的文章中，我先将概率值分配给各个事件，得到事件的概率分布。

通过事件与随机变量的映射，让事件“数值化”，事件的概率值转移到随机变量上，获得随机变量的概率分布。

我们使用随机变量的函数，来定制新的随机变量。随机变量的函数是从旧有的随机变量到一个新随机变量的映射。通过函数的映射功能，原有随机变量对应新的随机变量。通过原有随机变量的概率分布，我们可以获知新随机变量的概率分布。事件，随机变量，随机变量函数的关系如下:

 

 

一个简单的例子是掷硬币。出现正面的话，我赢1个筹码，负面的话，我输1个筹码。那么，投掷一次，赢的筹码数是一个随机变量X，X可能取值为1和-1。因此X的分布为:

P(1)=0.5

P(−1)=0.5

换一个角度来思考，我们将正负面“换算”成输赢的钱。如果一个筹码需要10元钱买，那么投掷一次硬币，赢的钱是一个随机变量Y，且 Y=10X 。Y的分布为:

P(10)=0.5

P(−10)=0.5

Y实际上是随机变量X的一个函数。X的1对应Y的10，X的-1对应Y的-10。即 Y=10X

小总结，在上面的实验中，硬币为正面为一个事件。赢得的筹码数为一个随机变量X。赢得的钱是X的函数Y，它也是一个随机变量。 

 

随机变量的函数还可以是多变量函数， Y=g(X1,X2,...,Xn) 。Y的值y对应的是多维空间的点 (x1,x2,...,xn) 。比如掷硬币，第一次赢的筹码为 X1 ，第二次赢的筹码为 X2 。我们可以构成一个新的随机变量 Y=X1+X2 ，即两次赢得的筹码的总和。

 

获得新概率分布的基本方法

一个核心问题是，如何通过X的概率分布，来获得 Y=g(X) 的概率分布。基本的思路是，如果我们想知道Y取某个值y的概率，可以找到对应的X值x的概率。这两个概率相等。

 

因此，我们使用如下方法来获得Y的概率。如果有函数关系 Y=g(X1,X2,...,Xn) ，获得Y分布的基本方法是:

1. 通过 Y=g(X1,X2,...,Xn) ，找到对应 {Y≤y} 的 (x1,x2,...,xn) 区间I。

2. 在区间I上，积分 f(x1,x2,...,xn) ，获得 P(Y≤y)

3. 通过微分，获得密度函数。

 

如果有函数关系 Y=X2 ， 而X满足下面的分布:

f(x)=12π−−√e−x2/2

对于任意 y≥0 来说，

F(y)=P(Y≤y)=P(X2≤y)=P(−y√≤X≤y√)

F(y)=∫y√−y√12π−−√e−x2/2dx=2∫y√012π−−√e−x2/2dx

对上面的F(y)微分，即获得密度函数

f(y)=12π−−√y−1/2e−y/2,0≤y≤∞

 

绘制密度函数 

 

import numpy as np
import matplotlib.pyplot as plt

pi = np.pi

x = np.linspace(-10, 10, 200)
y = np.linspace(0.1, 10, 100)

fx = 1/np.sqrt(2*pi)*np.exp(-x**2/2)
fy = 1/np.sqrt(2*pi)*(y**(-1/2))*np.exp(-y/2)

plt.plot(x, fx, color = "red", label="X distribution")
plt.plot(y, fy, label="Y distribution")

plt.title("Y = X*X")
plt.xlabel("RV")
plt.ylabel("pdf")

plt.legend()

plt.show()

 

上面的例子展示的是单变量函数，我们看一个多变量函数的例子。即 Y=g(X1,X2,...,Xn) ，且已知 X1,X2,...,Xn 的联合分布为 f(x1,x2,...,xn) 。我们需要找到满足 g(x1,x2,...,xn)≤y 的区间。

比如， Y=X1+X2 ，且 X1,X2 满足如下分布:

f(x1,x2)=12πexp(−12(x21+x22))

为了让 x1+x2≤y ，我们可以让 x1 任意取值，而让 x2≤y−x1

FY(y)=∫∞−∞∫y−x1−∞f(x1,x2) dx2dx1

让x_2 = v - x_1，有

FY(y)=∫∞−∞∫y−∞f(x1,v−x1)dvdx1= ∫y−∞∫∞−∞f(x1,v−x1)dvdx1

微分，可得y的分布为:

fY(y)=∫∞−∞f(x1,y−x1)dx1 =∫∞−∞ 12πexp(−12(x21+(y−x1)2))dx1

上述方程也可以使用数值方法求解:

代码如下:

# By Vamei

import numpy as np
import scipy.integrate
import matplotlib.pyplot as plt
pi = np.pi

'''
core of the integral
'''
def int_core(y):
    f = lambda x: 1.0/(2*pi)*np.exp(-0.5*(x**2 + (y-x)**2))
    return f

'''
calculate f(y)
'''
def density(y):
    rlt = scipy.integrate.quad(int_core(y), -np.inf, np.inf)
    return rlt[0]

# get distribution
y  = np.linspace(-10, 10, 100)
fy = map(density, y)

plt.plot(y, fy)
plt.title("PDF of X1+X2")
plt.ylabel("f(y)")
plt.xlabel("y")
plt.show()

上面的int_core()函数是一个闭包，它表示积分核部分。density()函数用于求某个y值下的积分结果。

 

(我们也可以利用解析的方法，推导出f(y)满足分布 N(0,2√) 。如果有微积分基础，可以将此作为练习。)

 

单变量函数的通用公式

上面求新的随机变量分布的步骤较为繁琐。在一些特殊情况下，我们可以直接代入通用公式，来获得新的分布。

(通用公式实际上是从基本方法推导出的数学表达式)

 

对于单变量函数来说，如果 Y=g(X) ，g是一个可微并且单调变化的函数 (在该条件，存在反函数 g−1 ，使得 X=g−1(Y)) 。那么我们可以使用下面的通用公式，来获得Y的分布:

fY(y)=fX(g−1(y))⋅ddyg−1(y)

 

假设X为标准分布，即 N(0,1) ，且 Y=5X+1 ，那么 g−1(y)=(y−1)/5 ，因此:

fY(y)=fX((y−1)/5)⋅(1/5)=152π−−√e−(y−1)2/(2×25)

可以看到，新的分布是一个 μ=1,σ=5 的正态分布，即 N(1,5)

 

并不是所有的函数都有反变换，所以这里的“通用”公式并不能适用于所有的情况。

 

多变量函数的通用公式

在一些特殊情况下，我们可以使用多变量函数的通用公式。

 

如果 U=g1(X,Y),V=g2(X,Y) ，且存在反变换，使得

X=h1(U,V)

Y=h2(U,V)

那么，我们可以通过如下公式，从X,Y的分布获得U,V的联合分布:

fUV(u,v)=fXY(h1(u,v),h2(u,v))|J|

J表示雅可比变换(Jacobian tranformation)，表示如下

J=∣∣∣∣∣∂x∂u∂y∂u∂x∂v∂y∂v∣∣∣∣∣=∂x∂u∂y∂v−∂x∂v∂y∂u

 

如果X和Y是独立的随机变量，且有相同的分布

f(x)=e−x,x≥0

。如果 U=X+Y,V=Y ，求U和V的联合分布。

由于X和Y独立，所以

fXY(x,y)=f(x)f(y)=e−xe−y

根据 U=X+Y ， V=Y ，可以得到 u≥0,v≥0 , 且有:

X=U−V

Y=V

因此

f(u,v)=e−(u−v)e−v=e−u,u≥0,v≥0

总结

通过随机变量的函数，我们可以利用已知随机变量，创建新的随机变量，并获得其分布。

欢迎继续阅读“数据科学”系列文章