NumPyの乱数生成について

はじめに

NumPy の乱数生成について、理解があいまいなまま使用していたせいでハマったため、ここでちゃんとまとめて供養しておきます。

• Version

Python: 3.10.17
NumPy: 2.2.6

Generator による「内部状態」管理

同じ seed を default_rng()（コンストラクタ） に渡して Generator を作成した場合、以下の rng1 と rng2 は別々の Generator オブジェクトですが、同じ内部状態 (state) を持ちます。

from numpy.random import default_rng
seed = 0
rng1 = default_rng(seed)
rng2 = default_rng(seed)

ここでは、この内部状態について詳しく見ていきます。
まず、rng1 と rng2 が持つ共通の内部状態を例えば、

state = 1

とします。（極端に単純化しています。実際ははるかに大きな状態空間と複雑な変換を使います。）
イメージとしては、「次の乱数を作るための現在状態」を表しています。
この内部状態を使って rng.random() メソッドが [0, 1) の乱数を生成するまでの過程を見ていきます。
説明用として、以下のような単純な疑似乱数生成器を考えます。

state = (5 * state + 1) % 10     # mod(10) で 0~9 の範囲に収める
output = state / 10              # [0, 1) にスケーリング

このとき、1回目の更新は、

state = (5 * 1 + 1) % 10
      = 6

となり、state の値が 1 ⇒ 6 に更新されます。
rng.random() メソッドは半開区間 [0, 1) の値を返したいので、この整数 6 を 10 で割った 0.6 を返す、という感じです。
まとめると、

state = (5 * state + 1) % 10
output = state / 10

となります。

さて、次にもう一度 rng.random() を呼び出す場合を考えます。
このとき、内部状態 state は呼び出しごとに初期化されているのではなく、先ほど更新した state = 6 から始まることになります。つまり、

state = (5 * 6 + 1) % 10
      = 1

となり、 state の値は 6 ⇒ 1 に更新されます。
そして、先ほどと同様に 1 を 10 で割った値 0.1 を返します。

以上をまとめると、

シード seed = 0　⇒ 初期状態 state = 1

1回目 random()
    state: 1 ⇒ 6　（決まったアルゴリズムで次の整数・ビット列をつくる + 内部状態を更新する）
    output: 6 / 10 = 0.6　（[0, 1) の小数にスケーリング）

2回目 random()
    state: 6 ⇒ 1
    output: 1 / 10 = 0.1

というように、シードによって一意に決まった Generator の内部状態 state は、Generator のメソッドを呼び出すたびに更新され、その状態から次の出力値が計算されます。

つまり、rng.random() は本当にランダムに小数を 「選んでいる」 のではなく、決まったアルゴリズムに従って 「計算」 しているだけなのです。なので、同じ seed から始めれば同じ内部状態になり、同じ内部状態から同じ計算をすれば、当然同じ値になる、というわけです。

また、rng.choice(10, 2) も同様に、seed によって内部状態が定まった後、決められたアルゴリズムに従って内部状態の更新と計算が行われ、最終的に2つの整数を返す、という仕組みです。
つまり、Generator が持つ内部状態から、疑似乱数アルゴリズムに従って乱数をつくり、それを各メソッドのアウトプットに応じた形へ変換している ということです。

random()
  内部状態 ⇒ 小数

choice()
  内部状態 ⇒ 整数選択

注意
厳密には、Generator 自体は分布メソッド（「どんな確率分布に従う乱数が欲しいか」を指定して乱数を返すメソッド）を提供し、Generator が内部に保持している BitGenerator (PCG64) が状態管理と乱数ビットの生成を担当します。

BitGenerator
    ↓
ランダムなビット列・整数列を作る
    ↓
Generator
    ↓
そのビット列を、ほしい確率分布の乱数に変換する

本記事では説明を簡単にするため、以降も「Generator が内部状態を持つ」と表現します。

同じ seed でも、呼び出しの回数や順序によって結果は変わる

重要な注意点として、「たとえ同じ seed でも、メソッド呼び出しのたびに内部状態は更新されるので、呼び出しの回数や順序によって乱数生成結果が異なる」 という点があります。

• 呼び出しの回数が違う場合

seed = 0
rng1 = default_rng(seed)
rng2 = default_rng(seed)
rng1.random() # rng1だけを1回余分に呼び出す
print(rng1.random())
print(rng2.random())

# 0.2697867137638703
# 0.6369616873214543

• 呼び出しの順序が違う場合

seed = 0
rng = default_rng(seed)
print(rng.random())
print(rng.choice(10, 2))

# 0.6369616873214543
# [5 2]

seed = 0
rng = default_rng(seed)
print(rng.choice(10, 2))
print(rng.random())

# [8 6]
# 0.2697867137638703

そのため、複数の実験条件を比較する場合には、1つの Generator (rng) を使いまわすのではなく、「実験単位」や「用途単位」で Generator を分けて管理すると、実行順序などに依存しない再現性のある比較がしやすくなります。

SeedSequence を使った実験管理

numpyの公式ドキュメント では、複数の実験間で再現性のある疑似乱数を生成するために SeedSequence を使った運用が紹介されています。SeedSequence を使うことで、例えば「実験IDに対応する値を親 seed として扱い、そこから用途別の子 Generator を作る」といった実験管理が可能になります。

SeedSequence は Generator を作るための「seed の設計図」のようなものであり、default_rng() に渡して初めて Generator となります。
SeedSequence に対して entropy と呼ばれる seed のようなものを与えることでインスタンス化します。

例えば、ある実験を seed (= 実験ID) を変えて繰り返し行う場合は以下のようになります。

from numpy.random import SeedSequence, default_rng
for seed in range(1000):
  ss = SeedSequence(seed)
  prob1_ss, prob2_ss, prob3_ss = ss.spawn(3)
  prob1_rng = default_rng(prob1_ss)
  prob2_rng = default_rng(prob2_ss)
  prob3_rng = default_rng(prob3_ss)

ss.spawn(3) によって、親の SeedSequence から 子の SeedSequence を3つ作成しています。子の SeedSequence は、同じ親 seed から生まれていますが、それぞれ用途別の Generator を作るために使うことができます。

最初の seed = 0 では、以下のイメージのようになり、

seed = 0
   ↓
SeedSequence(seed = 0)
   ↓ spawn(3) （目的別に seed 設計図を作る）
   |
   ├── prob1_ss  → prob1_rng（seed 設計図を default_rng() に渡して Generator をつくる）
   ├── prob2_ss  → prob2_rng
   └── prob3_ss  → prob3_rng

次の seed = 1 では、seed = 0 の時とは別の Generator セットが作られます。

seed = 1
   ↓
SeedSequence(seed = 1)
   ↓ spawn(3)
   |
   ├── prob1_ss  → prob1_rng（seed = 0 の時の prob1_rng とは異なる）
   ├── prob2_ss  → prob2_rng
   └── prob3_ss  → prob3_rng

このとき、3つの Generator は同じ親 seed 由来ですが、それぞれ異なる内部状態を持ちます。つまり、 実験管理用の seed を固定しつつ、内部状態の独立性が高い Generator を用途別に複数作ることができる、ということです。

なお、利便性を考慮し、直接 SeedSequence を使わずに 親 Generator から直接 spawn() することも可能です。

parant_rng = default_rng(seed)
prob1_rng, prob2_rng, prob3_rng = parent_rng.spawn(3)

今回は以上になります。