Box-Muller 算法, 每次生成一对正态随机数
设 `X, Y overset"iid"~ N(0, 1)`, 其联合密度为
`p(x, y) = 1/(2pi) "e"^(-(x^2+y^2)/2)`.
作变量替换 `X = r cos theta`, `Y = r sin theta`, 则
`1 = int_(-oo)^oo int_(-oo)^oo 1/(2pi) "e"^(-(x^2+y^2)/2) dx dy`
`= int_0^(2pi) int_0^oo 1/(2pi) "e"^(-r^2/2) r "d"r "d"theta`.
从而得到 `r, theta` 的分布函数
`P(r le r_0) = int_0^(2pi) int_0^(r_0) 1/(2pi) "e"^(-r^2/2) r "d"r "d"theta`
`= 1 - "e"^(-r_0^2/2)` `quad (r_0 ge 0)`,
`P(theta le theta_0) = int_0^(theta_0) int_0^oo 1/(2pi) "e"^(-r^2/2) r "d"r "d"theta`
`= theta_0/(2pi)` `quad (0 le theta_0 le 2pi)`.
如何生成满足上面这两个分布的随机变量 `r, theta`?
注意到 `theta ~ U(0, 2pi)`; 至于 `r`, 可以利用反函数:
令 `z = 1 - "e"^(-r^2/2)`, 反解得 `r = sqrt(-2 ln(1-z))`.
于是当 `z ~ U(0, 1)` 时, `r` 的分布函数为 `1 - "e"^(-r^2/2)`.
总之, 当 `U_1, U_2 ~ U(0, 1)` 时, 令
`X = r cos theta = sqrt(-2 ln U_1) cos 2pi U_2`,
`Y = r sin theta = sqrt(-2 ln U_1) sin 2pi U_2`,
则 `x, y` 服从标准正态分布.