摘要

生成随机变量的逆变换方法

几个缩写

• pmf(probability mass function)：概率质量函数。离散随机变量在各特定取值上的概率。只有离散型随机变量才有概率质量函数。

• PDF/pdf(probability density function)：概率密度函数，简称密度函数。描述随机变量的输出值，在某个确定的取值点附近的可能性的函数

• CDF/cdf(cumulatative distributionfunction)：累积分布函数，简称分布函数。是概率密度函数的积分，能完整描述一个实验随机变量X的概率分布。

连续型分布情况下

存在一个定理：

如果X是一个连续型随机变量，它的分布函数为
$$F_X(X)\sim U(0,1)$$。
逆变换方法生成随机数要求了概率积分的变化，定义逆变换为：
F X − 1 ( u ) = i n f { x : F X ( x ) = u } , 0 < u < 1 F_X^{-1}(u)=inf\{x:F_X(x)=u\},0<u<1 FX−1​(u)=inf{ x:FX​(x)=u},0<u<1
如果U~Uniform(0,1)，则对所有的x∈R
P ( F X − 1 ( U ) ≤ x ) = P ( i n f { t : F X ( t ) = U } ≤ x ) = P ( U ≤ F X ( x ) ) = F U ( F X ( x ) ) = F X ( x ) P(F_X^{-1}(U)\le x)=P(inf\{t:F_X(t)=U\}\le x) \\ =P(U\le F_X(x))\\ =F_U(F_X(x))\\ =F_X(x) P(FX−1​(U)≤x)=P(inf{ t:FX​(t)=U}≤x)=P(U≤FX​(x))=FU​(FX​(x))=FX​(x)
并且$F_X^{-1}(u)$和X有相同的分布。并且和X有相同的分布。因此，要想生成X的随机观测值，首先生成一个服从Uniform(0,1)的随机变量，并且传参给$F_X^{-1}(u)$得到相应值，即为X的观测值。只要提供的$F_X^{-1}$易于计算，那么这个方法就简单。同时，这个方法能够应用在连续型和离散型随机变量的生成。这个方法简单概括如下：

1. 推导出逆函数$F_X^{-1}(u)$
2. 写一个命令或函数来计算$F_X^{-1}(u)$
3. 对每一随机变量都要求：
   （1）生成的随机数u来自于Uniform(0, 1)
   （2）得到的x = $F_X^{-1}(u)$

例1. -------- 服从密度分布为$F_X(X)=3x^2(0<x<1)$的随即样本。

n <- 1000  #样本大小设置1000
U <- runif(n)   #生成均匀分布随机数
X <- U^(1/3)   #逆变换
hist(X, prob = TRUE)  #样本密度直方图， main = expression(f(x)==3*x^2)加上此参数可出现表达式标题
y <- seq(0, 1, .01)
lines(y, 3*y^2)    #理论密度曲线f(x)

结论：生成的随机数直方图与理论密度分布曲线保持一致。那就符合预期。
练习3.2-------

#练习 3.2   方法一  推导出逆函数的定义域后，即u的区间，将逆函数分段分别生成随机数
U1 = runif(5000,0.5,1)
X1 = -log(2-2*U1);   #u>1/2
U2 = runif(5000,0,0.5)
X2 = log(2*U2);    #u<1/2
X = append(X1, X2)    # x的样本
hist(X, prob = TRUE, ylim = c(0,1), breaks = 100)
x = seq(-60,60,0.1)
lines(x,0.5*exp(-abs(x)))
plot(x,0.5*exp(-abs(x)))

#练习 3.2   方法二
n = 10000  #样本大小
u = runif(n)
x3 = rep(0,n)
ind = u>.5   #判断不同定义域时对应的逆函数
x3[ind] = -log(2*(1-u[ind]))
x3[!ind]=log(2*u[!ind])
hist(x3, breaks=100, prob=TRUE, ylim = c(0,1))
y = seq(-60,60,.05)
lines(y,.5*exp(-abs(y)))

离散型分布情况下

逆变换法同样可以应用在离散型分布上。如果X是一个离散型随机变量并且
$$...<x_{i-1}<x_i<x_{i+1}<...$$
是$F_X(x)$的不连续点，则逆变换是当$F_X(x_{i-1})$<u<$F_X(x_{i+1})$时，使$F_X^{-1}(u)=x_i$。
对每个随机变量要求：

1. 生成的随机样本u服从U(0,1)
2. 当$F_X(x_{i-1})$<u<$F_X(x_{i+1})$，输出$x_i$

例 3.4 （两点分布）
这个例子利用逆变换来生成p=0.4的伯努利变量的一个随机样本。尽管在R中有更简单的办法生成两点分布的随机数，但这个简单的例子仅仅是用来阐述离散随机变量分布函数的逆变换计算方法。
在这个例子中，$F_X(0)=f_X(0)=1-p并且F_X(1)=1$ .因此，如果 u>0.6 ，$F_X^{-1}(u)=1$ .如果u≤0.6，$F_X^{-1}(u)=0$ 。这个生成器就把逻辑表达式u>0.6的值进行传递给X。

n <- 1000 
p <- 0.4 
u <- runif(n) 
x <- as.integer(u > 0.6)    #(u > 0.6) is a logical vector 
> mean(x) 
[1] 0.41 
> var(x) 
[1] 0.2421421

还有两种简单的二项分布随机数生成办法是利用R内置的函数直接的进行样本抽取：

n=100   #样本容量
rbinom(n, size = 1, prob = p) 
sample(c(0,1), size = n, replace = TRUE, prob = c(.6,.4))

练习3.5

#生成离散型分布函数随机数
u = runif(1000)
y1=integer(1000)   #1.初始化一个y1，用来放置离散概率函数的随机数，大小与均匀分布随机数保持一致
x = seq(0,4,1)     #2. 离散型分布随机变量的可能取值
p = c(0.1,0.2,0.2,0.2,0.3)   #3. 写出离散型概率函数分布列
f = c(0,cumsum(p))          #4. 写出离散型概率分布函数
m=length(f)         #判断对于均匀分布随机数进行的区间判断个数，准备找到对应离散型随机数的取值
#--------------------------------------------
y1 [u>0 &u<=0.1] = 0      #5. 区间判断
y1 [u>0.1 &u<=0.3] = 1
y1 [u>0.3 & u<=0.5] = 2
y1 [u>0.5 & u<= 0.7] = 3
y1 [u>0.7 & u<1] = 4
# ----------上下表达意思相同，选择一段即可-------
for(i in 1:m-1){
                      #5. 区间判断
  ind= u<=f[i+1] & u>f[i]
  y1[ind]=x[i]
}
#----------------------------------------------
b = table(y1)        #6. 随机数频数统计，加上"/n"统计频率

#sample函数生成离散型分布函数随机数,非逆变换法
k = sample(seq(0,4,1), size = 1000, replace = TRUE, prob = c(.1,.2,.2,.2,.3))
c = table(k)/n
b = table(y1)/n
rbind(p,b,c)
>
      0     1     2     3     4
p 0.100 0.200 0.200 0.200 0.300
b 0.111 0.205 0.211 0.198 0.275
c 0.097 0.225 0.199 0.192 0.287

可以看出，两种方法生成的样本还是服从题中分布的。