随机数:从蒙特卡洛法计算圆周率说起

蒙特卡洛法很流氓,思路却非常简单,往往暴力又好用。通过大量的随机样本,去估算并逼近得到一个渐近值,虽然永远不是确切值,但工程上有一定精度的结果就已经可以接受。
最近连续击败李世石的Google围棋AI AlphaGo,在估算围棋的势的时候,正是使用蒙特卡罗方法。对于一个棋面,要判断黑白子谁的势强,随机轮流落子,最后谁赢就是谁的势强,一盘下来得到的结论虽然比较扯,但大量样本得到的结果却非常有说服力。

那么如何用蒙特卡洛法计算圆周率pi呢?
在1X1的正方形上,以任意一顶点为圆心,1为半径画圆,被正方形的边截出1/4圆弧,其围成的扇形面积是(pi*r^2)/4,因为r=1,所以面积是pi/4。
取其顶点为原点,两边为x,y轴正方向,建立正交直角坐标系。在正方形内取大量随机点为样本,随机点落在扇形内的期望是扇形和正方形的面积之比,即pi/4比1,也即pi/4。
记总样本数为n, 落在扇形内的样本数为m, 当n趋向正无穷时,有m/n=pi/4。那么当n趋向无穷时,pi=4*m/n。

编程实现的话,

var m = 0, n = 10000;
for (var i = 0; i < n; ++i) {
  var x = Math.random();
  var y = Math.random();
  if (x * x + y * y <= 1) {
    m++;
  }
}
var pi = 4 * m / n;
console.log(pi);

在n=10000时,得到的pi大体在3.13到3.15之间浮动。
然而,在加大样本量到5000000,得到的结果也并不理想,小数点第三位后仍然不稳定,5000000的样本量也只是得到小数点后两位精度。

问题出在哪里?蒙特卡罗法事实上对样本的平均性要求很高。也即随机样本要真的均匀分布在指定范围内,而不是看起来好像是平均分布了。
具体到实现上,就是随机数的生成质量(均匀性)要非常好。

而反过来来说,上面的代码可以测试一个随机数生成算法的质量。

JavaScript的Math.random()得到的随机数是伪随机数,质量差强人意。
相信熟悉Linux的人都知道Linux上最先实现了真随机数生成器 (https://zh.wikipedia.org/wiki//dev/random)。
这是采集自硬件的热噪声并汇聚到系统熵池,根据熵产生的真随机数。当熵耗尽,取随机数便会被阻塞。这意味着,如果不采用专门的真随机数硬件,真随机数是稀少和高时间成本的。(大规模的熵收集可以通过采集大气噪声,这里有一家很有趣的公司:random.org)

折中方案是,利用真随机数作为种子,高质量伪随机数算法生成随机数序列。这种随机数既带有一定程度的熵,又能高效率产生,适合弱安全性场合。新一代浏览器中的crypto模块产生的随机数正是这种思路。
https://developer.mozilla.org/en-US/docs/Web/API/RandomSource/getRandomValues

我们再来计算一下圆周率

var m = 0, n = 10000;
for (var i = 0; i < n; ++i) {
  var arr = crypto.getRandomValues(new Uint32Array(2));
  var x = arr[0] / 0xFFFFFFFF;
  var y = arr[1] / 0xFFFFFFFF;
  if (x * x + y * y <= 1) {
    m++;
  }
}
var pi = 4 * m / n;
console.log(pi);

同样是5000000的样本,crypto.getRandomValues产生的随机数得到的精度可以达到小数点后4位,质量有了一定提高。

毫无疑问,AlphaGo内部肯定有一个高质量的随机数生成器,而且很可能是硬件生成器。

2016年3月10日 | 归档于 技术, 程序
  1. 2016年5月9日 01:33 | #1

    简单,直白,粗暴。上次才看到如何测试一个洗牌程序的好坏。。。似乎就是验证各种结果的概率,如果接近理论值就是比较好的算法。其中也是对如何“随机”特别讲究。
    话说。。。薰大,站点用了神马google插件,导致加载好慢,一直打不开页面。大陆用户悲剧了。

发表评论

XHTML: 您可以使用这些标签: <a href="" title=""> <abbr title=""> <acronym title=""> <b> <blockquote cite=""> <cite> <code> <del datetime=""> <em> <i> <q cite=""> <s> <strike> <strong>
:wink: :-| :-x :twisted: :) 8-O :( :roll: :-P :oops: :-o :mrgreen: :lol: :idea: :-D :evil: :cry: 8) :arrow: :-? :?: :!: