標準正規分布から（いい感じの）乱数をとってみた

目次

1. はじめに
2. 一様分布と標準正規分布
3. ボックス＝ミュラー法
4. 実装してみた
5. まとめ

1. はじめに

(擬似)乱数といえば、みなさん一度は使ったことがあると思います。
各言語で準備されている関数などを使ってランダムな値を取得できます。
特によくあるのは0以上1未満の値を返す関数です。

しかしたまにあるのが「こういう乱数じゃないんだよな」という場面です。
僕の場合は、具体的にはニューラルネットワークの重みを初期化する時です。
ニューラルネットワークの重みの初期値は、正規分布からとるといいと言われています。

そこで今回は、各言語で準備されている乱数から、
正規分布に従う乱数を抽出する方法を紹介・実装しました。

github.com

まずは乱数にまつわる用語や分布を紹介します。

2. 一様分布と標準正規分布

各言語でよくある乱数のメソッドから考えてみましょう。
0以上1未満の値を返す関数はほとんどどの言語でありますよね。
JavaScriptでいう Marh.random() です。

このメソッドはどんな大きさの値もほぼ同じ確率で抽出できます。
このように一様に値がでてくる分布を一様分布といいます。

図にすると次のようなものです。

次に標準正規分布を紹介します。
標準正規分布も図で確認してみましょう。

この正規分布は自然界でもたまに目にすることもできるようです。
また、中心極限定理でも登場するので目にすることも多い分布です。

このように有用な標準正規分布の乱数を取得する方法を紹介します。

3. ボックス＝ミュラー法

一様分布の乱数から標準正規分布に従う乱数を取得する方法があります。
それがボックス＝ミュラー法という方法です。
こちらはWikipediaに十分な説明があります。

ja.wikipedia.org

数式だけこちらで紹介します。
が互いに独立に 上の一様分布に従うとします。

$$ Z_1 = \sqrt{-2 log X} \cdot cos(2 \pi Y) \\ Z_2 = \sqrt{-2 log X} \cdot sin(2 \pi Y) $$

こんなに簡単な数式で標準正規分布に従う乱数 を取得できます。

4. 実装してみた

では上記の数式に従って、標準正規分布に従う乱数を計算する関数を実装します。
言語は TypeScript を採用します。

export default class RandomUtil {

  static extractFromStandardNormalDistribution(): number {
    const x: number = Math.random()
    const y: number = Math.random()
    const coefficient: number = Math.sqrt(-2 * Math.log(x))
    const radian: number = 2 * y * Math.PI
    const z1: number = coefficient * Math.cos(radian)
    // const z2: number = coefficient * Math.sin(radian)
    return z1
  }
}

最後に、本当に正規分布になっているのか確認したいと思います。
実装した乱数を10万個生成して、その分布を可視化します。

import RandomUtil from './utils/RandomUtil'
import Numbers from './utils/Numbers'

const numbers: number[] = Numbers.range(0, 1000 * 1000).map(() => {
  return RandomUtil.extractFromStandardNormalDistribution()
})

const min: number = -5
const max: number = 5
const interval: number = 0.1
const maxCharacters: number = 100
const countList: number[] = []
const legends: string[] = []

for (let x = min; x <= max; x += interval) {
  const count: number = numbers.filter((value) => {
    return value > x - interval / 2 && value < x + interval / 2
  }).length
  countList.push(count)
  legends.push((x>0 ? ' ' : '') + x.toFixed(1).toString())
}

const maxValue: number = Math.max(...countList)

countList.forEach((count: number, index: number) => {
  let bar: string = ''
  Numbers.range(0, count / maxValue * maxCharacters).forEach((numCharacters: number) => {
    bar += '='
  })
  console.log(legends[index], bar)
})

上記のプログラムを実行した結果は下記の通りです。
しっかり正規分布の形になっていますね！

-5.0 
-4.8 =
-4.6 =
-4.4 =
-4.2 =
-4.0 =
-3.8 =
-3.6 =
-3.4 =
-3.2 =
-3.0 =
-2.8 ==
-2.6 ==
-2.4 ===
-2.2 =====
-2.0 =======
-1.8 ===========
-1.6 ==============
-1.4 ====================
-1.2 =========================
-1.0 ===============================
-0.8 =====================================
-0.6 ==========================================
-0.4 ===============================================
-0.2 ==================================================
 0.0 ==================================================
 0.2 =================================================
 0.4 ===============================================
 0.6 ==========================================
 0.8 =====================================
 1.0 ===============================
 1.2 =========================
 1.4 ===================
 1.6 ==============
 1.8 ===========
 2.0 =======
 2.2 =====
 2.4 ===
 2.6 ==
 2.8 =
 3.0 =
 3.2 =
 3.4 =
 3.6 =
 3.8 =
 4.0 =
 4.2 =
 4.4 =
 4.6 =
 4.8 

5. まとめ

今回はボックス＝ミュラー法で標準正規分布の乱数を計算しました。
ニューラルネットワークの重みの初期値などに使う予定です。

また、他の分布の乱数を取得するには、逆関数法が有名ですね。
いつかそちらの方法で指数分布などの乱数も取得できるようにしたいです。

以上です。