在蒙特卡洛模拟中(Monte CarloMethod)中经常要用到服从各种分布的随机数发生器来估计模型参数,可是在一般的程序设计语言中只给出了服从均匀分布的随机数发生器,如C语言中的rand()函数,产生[0, RAND_MAX]中的均匀分布。在概率统计中我们可以由均匀分布构造服从任意分布的函数。
我们先看它的理论基础,下面是大学课本中概率统计中常见的一道证明题:
已知随机变量X的分布函数为F(X),证明Y = F(X)服从[0,1]的均匀分布。
由此可知,对于任意的概率密度 ,有 为X的分布函数,设X服从[0,1]的均匀分布,易得由(1)式我们可以得到y与x间的关系式,从而已知任意分布的概率密度函数,都可以由上式的均匀分布得到其相应的分布函数。
下面来看一个密度函数为指数分布的例子:
已知:
,X服从[0,1]上的均匀分布,由(1)式得: ,得到当然对于有些分布y是得不到关于x的显函数的,但是对于计算机来说这可以利用数值分析的方法来解。
好了,说了这么多,下面来看由C语言的rand()函数构造符合
的指数分布随机数发生器例子:#include <iostream>
#include <ctime>
#include <cstdlib>
#include <cmath>
using namespace std;
int main()
{
double x, y;
int i;
srand(time(NULL));
for(i = 0; i < 5000; ++i)
{
x = (double)rand() / RAND_MAX;//C语言的伪随机数发生器产生[0,1]的均匀分布;
y = log(1 - x) * -1;//由均匀分布构造生成服丛指数函数的y;
cout << y << endl;
}
return 0;
}