关于随机数的笔记
随机数
在实际生活中,随机数的应用范围非常广,包括物理仿真、统计采样、密码学、博彩等。获 得随机数的方法一般有2种,一种是基于物理现象由硬件产生,得到真随机数;另一种是由数值算法产生,得到伪随机数。对于很多应用,如仿真来说,数值方法是 最好的,因为数值方法得到的随机数并不是随机的,但是这些数看起来足够随机,并且随机数之间无关联,符合中心极限定理。更重要的是,这种方法很快并且占用 内存很少。在很多场合,产生不随机的数也是一种优势,因为可以产生相同的随机数序列,这样便于debug和确认程序,也可用来比较不同的系统。
随机数生成质量的评价
- 精确性:所需要产生的随机数的分布函数是确定的, 那么生成的方法也是确定的。
- 数学上的有效性:是否满足某些条件。
- 速度:生成随机数的速度必须是可以接受的。
- 空间:生成随机数所需的计算机内存。尤其是需要应用一些外部的表格时。
- 简单:算法和实现都必须满足此条件。
- 参数稳定性:对所有的参数而言都是速度均匀的。
真随机数的产生
真随机数只能用某些随机物理过程来产生。例如:放射性衰变、电子设备的热噪音、宇宙射线的触发时间等等。如果采用随机物理过程来产生蒙特卡洛计算用的随机数,理论上不存在问题,但是实际应用中,要做出速度很快而又准确的随机物理过程产生器是很困难的。Intel810RNG的原理大概是:利用热噪声(是由导体中电子的热震动引起的)放大后,影响一个由电压控制的振荡器,通过另一个高频振荡器来收集数据。
伪随机数的产生
实际应用的随机数通常都是通过某些数学公式计算而产生的伪随机数。这样的伪随机数从数学意义上讲已经一点不是随机的了。但是,只要伪随机数能够通过随机数的一系列的统计检验,我们就可以把它当作真随机数而放心地使用。这样我们就可以很经济地、重复地产生出随机数。理论上要求伪随机数产生器要具备以下特征:良好的统计分布特性、高效率的伪随机数产生、伪随机数产生的循环周期长,产生程序可移植性好和伪随机数可以重复产生。其中满足良好的统计特性是最重要的。
最早的伪随机数产生器可能是冯●诺依曼平方取中法。现在比较流行且用的最多的是线性同余法。
平方取中法代码如下:
public class CustomRandom {
static final int FIGURES = 10000;
static long mRandom;
public static void main(String[] args) {
long seed = System.currentTimeMillis();
mRandom = seed % FIGURES;
for (int i = 0; i < 10; i++)
System.out.println(getRandom(seed));
}
private static long getRandom(long seed) {
return mRandom = (mRandom * mRandom / (long) Math.pow(10, 5/2)) % FIGURES;
}
}
线性同余方法如下:
x0称为种子,改变它的值就得到基本序列的不同区段。a是乘子,c是增量,m是模。要仔细挑选m和a,使得产生的伪随机数的循环周期尽可能地长。
对于伪随机数的统计检验包括:均匀性检验、独立性检验、组合规律检验、无连贯性检验、参数检验等等。最基本的是均匀性和独立性。
生成特定分布随机数
得到【0,1】之间的均匀分布的随机数是很容易的,但是如何通过这样的随机数得到任意其他分布的随机数呢?
其中一种方法是Inverse Transform Method,基本思想是将【0,1】均匀分布的随机数u当做概率,返回随机数x,使得CDF函数F(x)=u。
其证明方法如下:
算法步骤如下:
以指数分布函数为例子:
对于离散的分布函数而言,计算会相对简单一些,算法如下:
证明如下:
这种方法需要可逆的F(x),并且速度可能会比较慢,因为需要一些比较。当然也有一些优点,如保留了单调性和相关性,这对于其他的应哟公式很有帮助的。
对于正太分布而言,这种方法是无效的,但是可以采用另一种方法:Box–Muller transform。其基本形式是:
正如以上所述,为x的cdf找到逆函数并不总是可能的,并且即使这是可能的,这也并不是效率很高,存在效率比较高的其他方法,比如接下来将要介绍的acceptance-rejection method。算法如下:
注意事项:
证明:
对其工作原理的一般性解释:
生成几种分布的 C++代码如下:
double generate()
{
static bool seeded = false;
if(!seeded)
{
srand(time(0));
seeded = true;
}
return static_cast<double>(rand()) / static_cast<double>(RAND_MAX);
}
//
// Uniform random numbers
//
class UniformRandomGenerator
{
int lowerBound;
int upperBound;
public:
UniformRandomGenerator(int lowerBound, int upperBound) : lowerBound(lowerBound), upperBound(upperBound) {}
int draw();
};
int UniformRandomGenerator::draw()
{
return generate() * (upperBound-lowerBound) + lowerBound;
}
//
// Exponentially distributed random numbers
//
class ExponentialRandomGenerator
{
int offset;
int lambda;
public:
ExponentialRandomGenerator(int offset, int lambda) : offset(offset), lambda(lambda) {}
int draw();
};
int ExponentialRandomGenerator::draw()
{
// found this somewhere on the web
double r = -log(generate());
return r*lambda + offset;
}
class NormalRandomGenerator
{
int mean;
int sigma;
public:
NormalRandomGenerator(int mean, int sigma) : mean(mean), sigma(sigma) {}
int draw();
};
int NormalRandomGenerator::draw()
{
// again, found somewhere on the web
double u1, u2, z;
u1 = generate();
u2 = generate();
z = sqrt(-2.0 * log(u1)) * cos(2.0 * PI * u2);
return mean + sigma*z;
}