伪随机数与准随机数的比较.pdf
5页·32· 计算机与信息技术 软件纵横 伪随机数与准随机数的比较 王水花 张煜东 吴乐南 (东南大学 信息科学与工程学院,江苏 南京 210096) 摘 要 传统的随机数生成法主要采用逆转法,首先生成[0, 1]区间上的均匀分布U,然后令X=F-1(U) ,则X满足分布F然而,这样生成的均匀分布具有明显的差异性,在小样本或高维空间的情况下尤其严重因此,引进一种新的随机数生成方法,即准随机数生成器实验验证了准随机数生成器得到的随机数的差异性优于传统方法最后,提出一种基于准随机数生成器的蒙特卡罗积分方法,结果优于传统的蒙特卡罗积分 关键词 伪随机数;准随机数;Kolmogorov- Smirnov 假设检验 1 引言 随机数生成算法[1]是一类重要的算法, 广泛应用于仿真技术等场合然而,目前的伪随机数生成器(Pseudo-random number generator, PRNG)[2]存在一个重要缺陷,即样本分布与真实分布不一致,这主要发生在以下两种情况:①抽样代价过高,样本数目较少;②空间维数较高[3]。
因此,有必要寻找一类新的随机数发生器准随机数发生器(Quasi-random number generator, QRNG)[4]能够生成稳定、低差异性的(low-discrepancy)样本,而与样本数目或空间维数无关[5]故针对蒙特卡罗积分结果不稳定的情况,提出一种基于 QRNG 的蒙特卡罗积分,发现比传统方法性能有所提升 2 伪随机数介绍 伪随机数是由确定的算法生成的,其分布函数与相关性均能通过统计测试与真实随机数的差别在于,它们是由算法产生的,而不是一个真实的随机过程一般地,伪随机数的生成方法主要有以下 3 种[6]: (1) 直接法(Direct Method) ,根据分布函数的物理意义生成缺点是仅适用于某些具有特殊分布的随机数,如二项式分布、泊松分布 (2) 逆转法(Inversion Method) ,假设U服从[0,1]区间上的均匀分布, 令X=F-1(U) , 则X的累计分布函数 (CDF)为F该方法原理简单、编程方便、适用性广 (3)接受拒绝法(Acceptance-Rejection Method) :假设希望生成的随机数的概率密度函数(PDF)为f, 则首先找到一个 PDF 为g的随机数发生器与常数c,使得f(x)≤cg(x) ,然后根据接收拒绝算法求解。
由于算法平均运算c次才能得到一个希望生成的随机数,因此c的取值必须尽可能小显然,该算法的缺点是较难确定g与c 因此,伪随机数生成器(PRNG)一般采用逆转法,其基础是均匀分布,均匀分布 PRNG 的优劣决定了整个随机数体系的优劣[7]下文研究均匀分布的 PRNG 3 伪随机数生成器的缺点 重复做N=10000 次试验,每次产生S=20 与S=100 个随机分布的样本,同时采用 Kolmogorov- Smirnov 假设检验(hypothesis test)来确定样本是否满足均匀分布规定: ① 0 假设(null hypothesis)为样本服从均匀分布;② 1假设(alternative hypothesis)为样本不服从均匀分布 采用 P 值(∈[0, 1])衡量,P 值越趋近于 0,表示越有理由拒绝 0 假设, 即样本不服从均匀分布; P 值越趋近于 1,表示越有理由接受 0 假设,即样本服从均匀分布 如图 1 与图 2 所示:随着 P 值下降,样本也越来越不服从均匀分布实践中希望 P 值越大越好然而统计学的结论显示,P 值一定服从均匀分布,与N、S大小无关,这表明由于随机性,总会出现某次抽样得到的样本不服从、甚至远离均匀分布。
另外,样本大小的不同,造成检验标准的不同,直观上看 S=100 对应的均匀分布普遍比 S=20 对应的更均匀因此,小样本情况下均匀分布 PRNG 的差异性尤为严重 软件纵横 计算机与信息技术 ·33· 00.20.40.60.81020406080100120140p-valuesNumber of Tests00.20.40.60.81051015样本值频数(a) KS 假设检验 P 值直方图 (b) P>0.95 对应结果 00.20.40.60.81051015样 本 值频数00.20.40.60.8105101520样本值频数(c) 0.475≤P≤0.525 对应结果 (d) P0.95 对应结果 00.20.40.60.810123样本值频数00.20.40.60.8012345样本值频数(c) 0.475≤P≤0.525 对应结果 (d) P<0.05 对应结果 图 2 10000 次 S=20 的均匀分布 PRNG ·34· 计算机与信息技术 软件纵横 -1-0.500.51-1-0.500.51xy(a) 二维(P=π/4) (b) 三维(P=π/6) 图 3 维数灾难示意图 维数灾难(curse of dimensionality)是另一个造成差异性的潜在原因。
通常假设一个半径r维数为d的超球,被一个边长为 2r维数为d的超立方体所包围,假设超立方体内存在一个均匀分布的点, 则由于超球的体积为 2rdπd/2/[dΓ (d/2) ],超立方体的体积为(2r)d,因此该点同时也落在超球内的概率P为 /212(/ 2)ddPddπ−=Γ(1) 令维数d由 2 逐步增长到 20,则对应的概率P如图 4 所示显然,当d=20 时,P仅为 2.46×10-8因此,若在 2 维空间中 1 个样本在半径r的意义下能逼近一个正方形,则在20 维空间内,则需要 1/2.46×10-8=4.06×107个样本才能在半径r的意义下逼近超立方体因此,在高维空间中的大样本甚至不如低维空间中的小样本 246810121416182010-610-410-2d维数P概率图 4 概率 P 与维数 d 的关系 4 准随机数发生器 上节讨论了造成差异性的两个情况: 小样本与高维空间本节讨论如何构建一类新的随机数发生器,使其具有较低的差异性 PRNG 缺陷的根源在于“随机性”与“均匀性”的矛盾因此,不要求新的发生器模拟真实的均匀分布,而力求任意大小的样本(尤其是小样本)都能满足低差异性。
换言之,以牺牲随机性为代价,换来均匀性的提高,称其为准随机数发生器(QRNG) 均匀分布 QRNG 的优势在于,其生成的样本更趋于均匀分布在其基础上构建的各类分布(包括高斯分布)的 QRNG,其生成的样本也更趋于服从对应的分布 目前有 3 种准随机序列(Quasi-random sequency)可用来辅助生成均匀分布随机数,分别是 Halton 序列、Sobol 序列、Latin 超立方体序列以 Halton 序列为例[8],简单介绍其原理假定存在互质基底b1, b2, …, bs,则整数n可以展开为: 0( , )i ij inaj n b==∑(2) 若记第j个基底的逆函数为 1( )( , ) ji bij inaj n b− −Φ=∑(3) 则定义序列 ()12( ),( ),...,( ) snbbbxnnn= ΦΦΦ(4) 为 Halton 序列这样生成的序列有强相关性,为了消除相关性,一般采用“忽略(skip) ” 、 “跳跃(leap) ” 、与“置乱(scramble) ”3 个步骤忽略与跳跃的物理意义见图 5,另外常采用 RR2 算法对 Halton 序列进行置乱。
图 5 忽略与跳跃的物理意义 5 实验 实验在主频 3.2GHz、内存 2GB 的 IBM P4 电脑上运行,实验平台采用 Matlab2009b,并借助 Statistical Toolbox 软件纵横 计算机与信息技术 ·35· 5.1 P 值比较 对QRNG重新进行第0节的试验 (S=20) , 令Skip为1000,Leap 为 100图 6 显示了基于 QRNG 的 P 值直方图,此时 P值几乎全部落在点 1 上,放大可见 P 值分布在[0.994, 1]区 间内 对比图 6 (a) 与图 2 (a) , 可见 QRNG 效果优于 PRNG,其差异性更小 任取其中某次抽样结果,绘成直方图,示于图 6(b) ,可见较图 2(b) (c) (d)更为均衡 0.9950.9960.9970.9980.9991020040060080010001200140016001800p-valuesNumber of Tests00.20.40.60.8100.511.522.53样 本 值频数(a) KS 假设检验 P 值直方图 (b) 某次抽样结果直方图 图 6 10000 次 S=20 的均匀分布 QRNG 5.2 抱团现象 图 7 (a) 给出了另一个直观例子, 在二维空间中绘出 500个随机数的散点图。
可见基于 PRNG 的均匀分布有一部分点互相重叠,反之另一部分区域为空,这种“抱团现象”正是由于随机性引起图 7(b)给出了基于 QRNG 的散点图,可见此时“抱团现象”完全消除在高维空间中抱团现象能导致“子超立方体”为空,使得样本分布与均匀分布相差极大 (a) PRNG (b) QRNG 图 7 500 个均匀分布的样本点 6 应用 QRNG 可以应用在对随机数的低差异性要求极高的场合,下面给出一个例子传统的积分计算通过将自变量划分为一个个小块,计算对应的函数值,再累加求得这样对一维函数或二维函数比较适用,但是对高维函数则无能为力例如,对 10 维的函数求积分,假设每维划分 100 个区间,则一共需要计算大约 10010=1020,在通常的计算机上根本无法承受 因此, 往往采用蒙特卡罗方法[9]求解, 与 “插针法” 类似以二维为例,考虑平面上一个边长为 1 的正方形及其内部的一个形状不规则图形,为求出该图形的面积,蒙特卡罗方法向该正方形随机地投掷N个点,统计落于图形内的M个点,则该图形的面积近似为M/N 然而, 传统蒙特卡罗方法中的随机数均采用 PRNG 生成,差异性较高。
因此提出一种基于 QRNG 的蒙特卡罗方法,即算法中的随机数采用 QRNG 生成假设求解 10 维空间中一个半径为 1 的超球的体积,理论结果是 2.5502分别运行两种算法, 规定每运行 10 万次记录当前结果, 结果如图 8 所示很明显,传统方法得到的结果偏小,且与理想结果偏差较大;而本文方法的结果更趋于分布在理想结果两侧,且偏差更小小,明显优于传统方法 0.20.40.60.811.21.41.61.82x 1062.482.492.52.512.522.532.542.552.56计 算 次 数积分结果蒙 特 卡 罗 方 法 理 想 结 果0.20.40.60.811.21.41.61.82x 1062.52.512.522.532.542.552.562.572.582.59计 算 次 数积分结果蒙 特 卡 罗 方 法 理 想 结 果(a) 传统方法 (b) 本文方法 图 8 两种蒙特卡罗方法比较 ·36· 。





