Pythonで〇×ゲームのAIを一から作成する その248 グラフによる多腕バンディット問題のアルゴリズムの性質の視覚化と考察

目次と前回の記事

Python のバージョンとこれまでに作成したモジュール

本記事のプログラムは Python のバージョン 3.13 で実行しています。また、numpy のバージョンは 2.3.5 です。

リンク	説明
marubatsu.py	Marubatsu、Marubatsu_GUI クラスの定義
ai.py	AI に関する関数
mbtest.py	テストに関する関数
util.py	ユーティリティ関数の定義
tree.py	ゲーム木に関する Node、Mbtree クラスなどの定義
gui.py	GUI に関する処理を行う基底クラスとなる GUI クラスの定義

AI の一覧とこれまでに作成したデータファイルについては、下記の記事を参照して下さい。

今回の記事の内容

前回の記事では、多腕バンディット問題における探索と活用のトレードオフについて説明しました。前回の記事で行ったように、その性質を数式や数値だけで説明されても、直感的には理解しづらいという問題があります。そこで様々なアルゴリズムの紹介は次回以降の記事で行うことにし、今回の記事では、探索と活用を行うアルゴリズムの性質をグラフで視覚化するプログラムを実装します。

具体的には「原始モンテカルロ法」で多腕バンディット問題を解いた場合をグラフ化します。そして、前回の記事で説明した「探索に極振りした原始モンテカルロ法」の性質が、実際にグラフ上でどのように現れるかを視覚的にわかりやすく示します。

視覚化の方針

多腕バンディット問題における「活用」とは、「最も当たりが出やすいと推測されるスロットマシンのレバーを引く」ことです。この 活用で稼げるコインの期待値 は、前回の記事で説明したように 活用を行う前に収集した「経験（データ）」が多い程高くなります。

この経験は、スロットマシンのレバーを引くたびに蓄積されていきます。そのため、多腕バンディット問題を解く一般的なアルゴリズムでは、試行回数（レバーを引いた総数）多くなればなるほど、活用によって得られるコインの期待値も高くなります¹。つまり、多腕バンディット問題の目的関数（得られるコインの総数の期待値）は、回数を増やせば増やすほど理論上の最大値である $f_{max}$²に近づいていきます。

前回の記事のような数式だけの説明では「目的関数が $f_{max}$ に向かってどのように、どこまで近づいていくのか」をイメージするのは困難です。しかし「レバーを引いた回数」と「目的関数」の関係をグラフ化することで、その関係を直感的に捉えられるようになります。

また、アルゴリズムの性質をより深く理解するためには、目的関数の推移だけでなく、最高のスロットマシンの推測に使われる「各スロットマシンの当たりの平均値（推定確率）」がどのように真の確率へ収束していくかという性質も重要です。そこで今回の記事では、各スロットマシンの内部状態の推移についてもグラフで視覚化することにします。

具体的な検証の設定と用語の説明

グラフ化を行うにあたり、シミュレーションの条件を具体的に設定する必要があります。今回の記事では、スロットマシンのレバーを引く総数 $m$ を 1 回 から 1 万回 まで変化させたそれぞれのケースについて実際に問題を解き、結果をグラフ化します。

また、本記事では強化学習の用語にならって、以下のように定義します。

• エピソード：ある $m$ 回にわたってスロットマシンのレバーを選択して引くという一連の処理
• ステップ：エピソードの中で、レバーを実際に 1 回引く単一の処理。強化学習の慣例に従い、ステップ数は 0 から数える

原始モンテカルロ法の近似を行うアルゴリズムの視覚化

今回の記事では、多腕バンディット問題を「探索を行ってから最後に活用を行う」原始モンテカルロ法で解くアルゴリズムの視覚化を行い、前回の記事で説明した「探索に全振りした場合の性質」が、実際にグラフでどのように表されるかを確かめることにします。

原始モンテカルロ法をそのままシミュレートする際の問題点

レバーを引く総数（$m$）が 1 から 1万 までのそれぞれの場合に対して原始モンテカルロ法をそのまま実行しようとすると、計算量が膨大になりすぎるという問題に直面します。具体的には、合計で約 5,000 万回もの探索処理を行う必要があります。理由は以下の通りです。

• 原始モンテカルロ法の 1 エピソードは「$m-1$ 回の探索と、1 回の活用」で構成されるため、$m$ の値が変わると探索と活用の割り当て方法が完全に変わってしまう
• そのため、1 から 1 万までの すべての $m$ に対して、それぞれ独立して原始モンテカルロ法を実行する 必要がある
• 1 エピソードあたりの探索回数は $m - 1$ 回なので、すべての $m$ を足し合わせた探索の総数は、下記の式から約 $\frac{10000^2}{2} =$ 約 5,000 万回となる

$$0 + 1 + \dots + 9999 = \sum_{i=0}^{9999}i = \frac{{(0 + 9999)} \times 10000}{2} = \frac{9999 \times 10000}{2} \approx \frac{10000^2}{2}$$

1 エピソードで行われる活用の回数は最後の 1 回だけなので、活用の総数は 1 万回になります。

目的関数である「累積報酬の期待値に関するグラフ」を描画するためには、複数回のエピソードのシミュレーションを行い「累積報酬の平均値」を計算する必要があります。そのため、1 エピソードあたりのシミュレーションの処理時間が長いと、全体のシミュレーションに時間がかかりすぎる³という深刻な問題が発生します。

原始モンテカルロ法のほぼ正確な近似アルゴリズム

前回の記事で説明したように、最後以外の行動を探索に全振りする原始モンテカルロ法では、目的関数（得られたコインの総数の期待値）の値は、「最後の 1 回で行われる活用」からほとんど影響を受けません。

例えば $m = 10,000$ の場合では、最後の 1 回の活用が目的関数に与える影響はわずか 1 万分の 1 です。そのため、常にスロットマシンをランダムに選択し続ける 「探索のみを行うアルゴリズム」は、原始モンテカルロ法のほぼ正確な近似 とみなすことができます。

さらに、探索のみを行うアルゴリズムには「$m$ の値に関係なく、常にランダムな探索しか行わない」という大きなメリットがあります。これにより、1 万回の「探索のみを行うシミュレーション」を 1 回実行し、その途中（例えば 1〜5,000 回目まで）のデータを取り出すだけで、$m = 5,000$ の場合のシミュレーション結果としてそのまま流用できます。

つまり、$m=10,000$ として探索のみのシミュレーションをたった 1 回行うだけで、$m$ が 1〜10,000 の場合のすべての結果（最後の活用を除いた原始モンテカルロ法の近似）を同時に計算できるということです。

この場合の探索回数はわずか 1 万回 で済み、原始モンテカルロ法をシミュレートした場合（約 5,000 万回）と比べて5,000 分の 1 に激減します。その結果、プログラムの実行速度も約 5000 倍高速になります。そこで本記事では、この「探索のみを行うアルゴリズム」を原始モンテカルロ法の代わり（近似）として用いて視覚化を進めます。

原始モンテカルロ法の場合は、1 万回のシミュレーションのデータから 1 ~ 5,000 回目のデータを取り出しても $m = 5,000$ のシミュレーションの結果にはなりません。その理由は、取り出したデータには「5,000 回目の最後に活用を行う」という原始モンテカルロ法特有の処理が含まれておらず、単なる探索のみのデータになってしまうからです。

1 回のエピソードの処理を記録するプログラム

グラフ化を行うためのデータを集めるために、前回の記事と同じ下記の設定の多腕バンディット問題に対して、シミュレーションデータを収集するプログラムを実装します。

• スロットマシンの台数： $3$ 台
• 各スロットマシンの当たりの確率：$0.3, 0.5, 0.7$
• スロットマシンのレバーを引く総数：1 万回

収集するデータは、各ステップごとの下記の 5 つの値です。これらのデータをわざわざ記録しておく具体的な理由については、この後でグラフ化を行う際に詳しく説明します。

• 累積報酬：そのステップまでに得られたコインの総数
• 1 ステップあたりの累積報酬：累積報酬 ÷ その時点でのステップ数
• 各スロットマシンで得られたコインの総数
• 各スロットマシンを引いた回数
• 各スロットマシンで得られたコインの平均：得られたコインの総数 ÷ スロットマシンを引いた回数

プログラムの定義と実行

プログラムは前回の記事で作成した「原始モンテカルロ法のデータを記録するプログラム」とよく似ていますが、以下の点が大きく異なります。

• 試行（探索）回数を「1 〜 4 回」から 1 万回 へと大幅に増やす
• 探索の後の「活用」は行わない（探索に全振りする）
• 毎ステップの途中経過をすべて list に記録する

プログラムで用いる主な変数の一覧は下記の通りです。多腕バンディット問題を扱う MAB クラスについて忘れた方は以前の記事を復習して下さい。

変数	意味
mab	スロットマシンの当たりの確率を [0.3, 0.5, 0.7] に設定した MAB クラスのインスタンス
stepnum	シミュレーションを行う総ステップ数（1 万 回）
step	繰り返し処理内での現在処理中のステップ数（ステップの番号）
total_coin	そのステップまでに得られた「コインの総数」（最新値）
mean_coin	そのステップまでに得られたコインの「1 回あたりの平均枚数」（最新値）
total_coins	各ステップごとの total_coin の推移を記録する list
mean_coins	各ステップごとの mean_coin の推移を記録する list
slot	各スロットマシンの詳細データを管理する辞書 dict。下記のキーに 1 次元または 2 次元の list を記録する "total_coin[arm]"：arm のスロットマシンが得たコインの総数（最新値） "mean_coin[arm]"：arm のスロットマシンが得たコインの平均（最新値） "selectnum[arm]"：arm のスロットマシンが選択された回数（最新値） "total_coins[arm][step]"：arm のスロットマシンが step のステップまでに得たコインの総数 "mean_coins[arm][step]"：arm のスロットマシンが step のステップまでに得たコインの平均 "selectnums[arm][step]"：arm のスロットマシンが step のステップまで選択された回数

下記は 1 万回の探索を行いながらデータを収集、記録するプログラムです。

• 4 ~ 17 行目：初期設定を行う。14 ~ 16 行目では、スロットマシンの台数だけ空の list を要素とする 2 次元の list を作成して代入する。以前の記事で説明したように、ここで [[]] * mab.num と記述してしまうとすべての要素が同じ空リストのオブジェクトを共有してしまう（連動して書き換わってしまう）という問題が発生する点に注意すること
• 19 ~ 32 行目：stepnum（1 万）回の探索を行う繰り返し処理
  • 20 行目：探索ではランダムにスロットマシンを選択するので、常に np.random.choice(mab.num) でスロットマシンを選択する
  • 26 ~ 28 行目：レバーを選択したスロットマシンのデータ（総数、回数、平均）を更新する
  • 29 ~ 32 行目：そのステップにおけるすべてのスロットマシンの状態をリストに記録する。選ばれなかったマシンの状態も含めて、毎ステップすべてのスロットマシンのデータを記録し続ける必要がある点に注意すること

 1  from ai import MAB
 2  import numpy as np
 3  
 4  mab = MAB(p=[0.3, 0.5, 0.7])
 5  stepnum = 10000
 6  total_coin = 0
 7  total_coins = []
 8  mean_coin = 0
 9  mean_coins = []
10  slot = {
11      "total_coin": [0] * mab.num,
12      "mean_coin": [0] * mab.num,
13      "selectnum": [0] * mab.num,
14      "total_coins": [[] for _ in range(mab.num)],
15      "mean_coins": [[] for _ in range(mab.num)],
16      "selectnums": [[] for _ in range(mab.num)],
17  }
18  
19  for step in range(stepnum):
20      selectarm = np.random.choice(mab.num)
21      coin = mab.play(selectarm)
22      total_coin += coin
23      total_coins.append(total_coin)
24      mean_coin = total_coin / (step + 1)
25      mean_coins.append(mean_coin)
26      slot["total_coin"][selectarm] += coin
27      slot["selectnum"][selectarm] += 1
28      slot["mean_coin"][selectarm] = slot["total_coin"][selectarm] / slot["selectnum"][selectarm]
29      for arm in range(mab.num):
30          slot["total_coins"][arm].append(slot["total_coin"][arm])
31          slot["mean_coins"][arm].append(slot["mean_coin"][arm])
32          slot["selectnums"][arm].append(slot["selectnum"][arm])
33  
34  print(f"得られたコインの総数　 = {total_coin:5d}")
35  print(f"得られたコインの平均　 = {mean_coin:.3f}")
36  for arm in range(mab.num):
37      print(f"スロットマシン {arm} (p={mab.p[arm]:.2f})")
38      print(f"  得られたコインの総数 = {slot["total_coin"][arm]:5d}")
39      print(f"  レバーを引いた回数　 = {slot["selectnum"][arm]:5d}")
40      print(f"  得られたコインの平均 = {slot["mean_coin"][arm]:.3f}")

行番号のないプログラム

from ai import MAB
import numpy as np

mab = MAB(p=[0.3, 0.5, 0.7])
stepnum = 10000
total_coin = 0
total_coins = []
mean_coin = 0
mean_coins = []
slot = {
    "total_coin": [0] * mab.num,
    "mean_coin": [0] * mab.num,
    "selectnum": [0] * mab.num,
    "total_coins": [[] for _ in range(mab.num)],
    "mean_coins": [[] for _ in range(mab.num)],
    "selectnums": [[] for _ in range(mab.num)],
}

for step in range(stepnum):
    selectarm = np.random.choice(mab.num)
    coin = mab.play(selectarm)
    total_coin += coin
    total_coins.append(total_coin)
    mean_coin = total_coin / (step + 1)
    mean_coins.append(mean_coin)
    slot["total_coin"][selectarm] += coin
    slot["selectnum"][selectarm] += 1
    slot["mean_coin"][selectarm] = slot["total_coin"][selectarm] / slot["selectnum"][selectarm]
    for arm in range(mab.num):
        slot["total_coins"][arm].append(slot["total_coin"][arm])
        slot["mean_coins"][arm].append(slot["mean_coin"][arm])
        slot["selectnums"][arm].append(slot["selectnum"][arm])
        
print(f"得られたコインの総数　 = {total_coin:5d}")
print(f"得られたコインの平均　 = {mean_coin:.3f}")
for arm in range(mab.num):
    print(f"スロットマシン {arm} (p={mab.p[arm]:.2f})")
    print(f"  得られたコインの総数 = {slot["total_coin"][arm]:5d}")
    print(f"  レバーを引いた回数　 = {slot["selectnum"][arm]:5d}")
    print(f"  得られたコインの平均 = {slot["mean_coin"][arm]:.3f}")

実行結果（乱数が使われているため、実行結果は毎回下記の内容と若干異なります）

得られたコインの総数　 =  5024
得られたコインの平均　 = 0.502
スロットマシン 0 (p=0.30)
  得られたコインの総数 =  1019
  レバーを引いた回数　 =  3318
  得られたコインの平均 = 0.307
スロットマシン 1 (p=0.50)
  得られたコインの総数 =  1672
  レバーを引いた回数　 =  3371
  得られたコインの平均 = 0.496
スロットマシン 2 (p=0.70)
  得られたコインの総数 =  2333
  レバーを引いた回数　 =  3311
  得られたコインの平均 = 0.705

実行結果から、大数の法則に従った、下記のような期待通りのデータが収集できていることが確認できます。

• 1 ステップ当たりに得られるコインの平均は 3 つのスロットマシンの当たりの確率の平均である $p_{mean} = \frac{0.3 + 0.5 + 0.7}{3} = 0.5$ とほぼ同じ値になる
• 総ステップ数は 1 万回なので、得られたコインの総数はその 1 万倍の約 5000 枚になる
• 3 つのスロットマシンがいずれもほぼ均等に約 3300 回選択され、得られたコインの平均はスロットマシンの当たりの確率（p）とほぼ同じ値になっている

累積報酬のグラフ化

最初に、ステップごとの累積報酬（得られたコインの総数）の推移を表すグラフを描画します。その際にコインの総数の上限（理論上の最大値）である $f_{max}$ と、累積報酬の内訳（各スロットマシンで得られたコインの総数）の推移も併せて描画しました。

• 6 行目：スロットマシンの当たりの確率の最大値（maxp）を元に、各ステップごとの累積報酬の上限（upper bound）を計算してグラフを描画する。目出つように、このグラフを c="k" を指定して黒色の直線で描画した
• 8 行目：累積報酬の推移をグラフ化する
• 9、10 行目：各スロットマシンで得られたコインの総数をグラフ化する。その際に、個々のスロットマシンで得られた枚数であることが明確に区別できるように、linestyle="--"⁴ を指定してして点線で描画した

 1  import matplotlib.pyplot as plt
 2  import japanize_matplotlib
 3  
 4  plt.title(f"累積報酬の推移（原始モンテカルロ法の近似）")
 5  maxp = max(mab.p)
 6  upperbounds = [step * maxp for step in range(stepnum)]
 7  plt.plot(upperbounds, label="上限", c="k")
 8  plt.plot(total_coins, label="コインの総数")
 9  for arm in range(mab.num):
10      plt.plot(slot["total_coins"][arm], linestyle="--", label=f"arm {arm} (p={mab.p[arm]:.2f})")
11  plt.xlabel("ステップ数")
12  plt.ylabel("コインの枚数")
13  plt.legend()
14  plt.show()

行番号のないプログラム

import matplotlib.pyplot as plt
import japanize_matplotlib

plt.title(f"累積報酬の推移（原始モンテカルロ法の近似）")
maxp = max(mab.p)
upperbounds = [step * maxp for step in range(stepnum)]
plt.plot(upperbounds, label="上限", c="k")
plt.plot(total_coins, label="コインの総数")
for arm in range(mab.num):
    plt.plot(slot["total_coins"][arm], linestyle="--", label=f"arm {arm} (p={mab.p[arm]:.2f})")
plt.xlabel("ステップ数")
plt.ylabel("コインの枚数")
plt.legend()
plt.show()

実行結果

グラフから読み取れる性質と考察

グラフから、探索に全振りしたアルゴリズム（原始モンテカルロ法の近似）が持つ以下の性質が一目で分かります。これは、前回の記事で説明した「探索に全振りした場合の考察」をそのまま裏付ける結果となっています。

• 累積報酬が理論上の上限 $f_{max}$ に大きく及ばない
  常に最善のスロットマシン（確率 0.7）を引き続けられる「神様の視点（黒い直線）」と比べて、「コインの総数（青い直線）」の傾きが大幅に低くなっています。これは、実力が低いハズレ台（確率 0.3 や 0.5）に対しても、馬鹿正直に均等な回数だけレバーを割り当て続けて大損しているためです。

• 累積報酬や内訳のグラフが綺麗な「直線」になる
  すべてのグラフが完全な直線になっていることから、獲得コインの枚数が試行回数（ステップ数）に完全に比例していることがわかります。

• 各スロットマシンの傾きが「確率」に対応している
  点線の傾きに注目すると、最も当たりやすい arm 2（確率 0.7）が一番多くコインを稼いでおり、次いで arm 1（0.5）、arm 0（0.3）の順に綺麗に並んでいます。これは大数の法則により、各スロットマシンの獲得報酬が本来の確率通りの期待値へと収束していることを視覚的に表しています。

ステップ当たりの累積報酬のグラフ化

先程のグラフは、累積報酬も上限もステップ数が増えるにつれて右肩上がりに伸びていく性質があるため、両者の距離が「近づいているのか」それとも「引き離されているのか」を直感的に見極めることが困難です。

累積報酬の上限 $f_{max}$ は、下記のように最も当たりやすいスロットマシンの確率 $p_{max}$ にステップ数を掛けた値です。

$$f_{max} = p_{max} \times ステップ数$$

この数式をステップ数で割ると、「1 ステップあたりに得られるコインの上限（期待値）」は、常に $p_{max}$ という一定の定数 になります。このように、合計値ではなくステップ数で割った「1 回あたりの平均値」に変換して比較したほうが、アルゴリズムの本当の実力を直感的に評価しやすくなります。

下記は、ステップ当たりの累積報酬と、その上限と平均の推移を描画するプログラムです。先ほどのプログラムとの違いは以下の通りです。

• upperbounds （上限の基準線）のすべての要素を、定数である maxp（0.7）にした
• 比較基準として、3 台すべてのマシンの確率の平均値である $0.5$ を表す「当たりの平均」の線を加えた
• 1 ステップ当たりの獲得コインの推移を mean_coins で描画した
• グラフの縦軸の範囲を 0 ~ 1 に固定して変化を捉えやすくするため、plt.ylim(0, 1) を実行した
• 各スロットマシンの個別の平均推移（slot["mean_coins"][arm]）は、分母が「総ステップ数」ではなく「各マシンの選択回数」で計算されており、グラフの性質が異なるためここでは描画しない

plt.title(f"ステップ当たりの累積報酬の推移（原始モンテカルロ法の近似）")
maxp = max(mab.p)
upperbounds = [maxp for step in range(stepnum)]
plt.plot(upperbounds, label="上限", c="k")
mean = np.mean(mab.p)
means = [mean for step in range(stepnum)]
plt.plot(means, label="当たりの平均")
plt.plot(mean_coins, label="ステップ当たりのコインの枚数")
plt.ylim(0, 1)
plt.xlabel("ステップ数")
plt.ylabel("ステップ当たりのコインの枚数")
plt.legend()
plt.show()

実行結果

グラフから読み取れる性質と考察

グラフを観察すると、探索に全振りしたアルゴリズム（原始モンテカルロ法の近似）が持つ決定的な弱点が見えてきます。

• どれだけ試行回数を増やしても、実力が「当たりの平均 $0.5$」から全く向上しない
  グラフの「ステップ当たりのコインの枚数（オレンジ色の線）」は、スタート直後こそ大きくブレるものの、ステップ数が多くなると最低のラインである「当たりの平均（水色の基準線）」へ完全にぴったりと張り付いていきます。

• 最善の上限 $0.7$ に近づかない理由
  これは、過去のステップで「どの台が当たりやすいか」というデータ（経験）がどれだけ蓄積されても、それをその後の選択に全く活かさず、最後まで馬鹿正直にすべてのスロットマシンを完全ランダムに探索し続けているためです。

なお、今後の記事で紹介していく賢い強化学習アルゴリズムでは、ステップ数が増えるにつれてオレンジ色の線が「当たりの平均（$0.5$）」を突破し、神様の視点である「上限 $0.7$」に向かって右肩上がりに近づいていくという性質を示します。今回のグラフはその性能を測るための「最悪のベースライン（基準線）」として重要な意味を持ちます。

グラフの左端（ステップ数が 1〜数百回 のエリア）に注目すると、線が激しく上下にブレています。これは、試行回数が少ないうちはデータのばらつき（分散）が大きくなり、回数を増やすほど真の値に収束していくという 大数の法則 によるものです。

そのため、今後どのような天才的なアルゴリズムを導入したとしても、データが不足している最初の数ステップ〜数百ステップの間は、必ずこのようにグラフのブレが大きくなります。

スロットマシンの選択回数のグラフ化

探索と活用のアルゴリズムの性質をより深く理解するためには、目的関数（コインの枚数）の推移だけでなく、「AI がそれぞれのスロットマシンをどれくらい選んだか」という詳細なデータを把握することも重要です。

各スロットマシンが選択された回数の推移をグラフ化することで、ステップ数が増えて学習が進むにつれて、AI が優秀なスロットマシンを どのように「えこひいき」していったか（行動の偏り） を視覚的に知ることができます。

下記は、そのグラフを描画するプログラムです。

plt.title(f"スロットマシンの選択回数（原始モンテカルロ法の近似）")
for arm in range(mab.num):
    plt.plot(slot["selectnums"][arm], label=f"arm {arm} (p={mab.p[arm]:.2f})")
plt.xlabel("ステップ数")
plt.ylabel("選択回数")
plt.legend()
plt.show()

実行結果

グラフから読み取れる性質と考察

グラフを確認すると、探索に全振りしたアルゴリズム（原始モンテカルロ法の近似）の振る舞いについて、以下の特徴があることがわかります。

• 3 本の線が完全に重なり合って、綺麗な1本の直線に見える
  arm 0 から arm 2 までの選択回数を表すグラフは、ほぼぴったりと重なり合って右肩上がりに伸びています。これは、AI が「えこひいき」を一切行わず、1 万ステップの最後まですべてのスロットマシンを完全に平等（確率 $\frac{1}{3}$）に選び続けたことを意味しています。

• 情報が全く活かされない
  たとえ途中で「arm 2（確率 0.7）が当たりやすい」というデータが集まっていたとしても、このアルゴリズムにはその情報を活用する仕組みがありません。結果として、最悪のハズレ台である arm 0（確率 0.3）も、最高の当たり台と同じ回数（約 3,333回）だけ無駄に引き続けられてしまいます。

今後の記事で紹介する賢い強化学習アルゴリズムでは、最初は 3 本の線が横並びですが、学習が進むにつれて「当たりやすいスロットマシンの線だけが急激に上に跳ね上がり、ハズレ台の線は底を這うように横ばいになる」 という、えこひいきのグラフが現れるようになります。今回の完全な横並びのグラフは、それらと比較するための重要な基準点となります。

スロットマシンの当たりの確率の推定値のグラフ化

前回の記事で説明したように、強化学習の「活用」では 「各スロットマシンが得たコインの枚数の平均値（標本平均）」をそのスロットマシンの本当の当たりの確率の推定値 とみなし、その値が最大となるスロットマシンを選択します。そのため、この推定値がステップごとにどう変化していくかは、アルゴリズムの性質や推測の精度を評価するための極めて重要な要素となります。

下記は、各スロットマシンの当たりの確率の推定値を描画するプログラムです。推定値の精度が直観的に伝わるように、各スロットマシンの「本当の当たりの確率（真の確率）」を表す水平な基準線を黒色（c="k"）の点線（linestyle=":"）で併せて描画しました。

plt.title(f"スロットマシンの当たりの確率推定値（原始モンテカルロ法の近似）")
for arm in range(mab.num):
    plt.plot([mab.p[arm] for _ in range(stepnum)], linestyle=":", c="k")
    plt.plot(slot["mean_coins"][arm], label=f"arm {arm} (p={mab.p[arm]:.2f})")
plt.ylim(0, 1)
plt.xlabel("ステップ数")
plt.ylabel("当たりの確率")
plt.legend()
plt.show()

実行結果

グラフから読み取れる性質と考察

グラフを確認すると、ステップ数が進むにつれて AI の「推測の精度」が劇的に進化していく様子がはっきりとわかります。

• 大数の法則による収束
  スタート直後（ステップ数が少ない段階）は、どのマシンの推定値も激しく上下にブレています。しかし、ステップ数が 1,000、2,000 と増えていくにつれてブレが急速に収まり、黒い点線（真の確率：0.3, 0.5, 0.7）へと綺麗に収束 していきます。これこそが、以前の記事で紹介した「大数の法則」がもたらす数学的な恩恵です。

• 原始モンテカルロ法が「無駄な探索」をしている証拠
  グラフの「中盤から終盤（ステップ数 3,000〜10,000）」に注目してください。3,000 ステップを過ぎたあたりで、3 台の推定値は真の確率にほぼ重なっており、「どの台が一番優秀か（arm 2 であること）」はほぼ完全に明らかです。それにも関わらず、今回の「探索のみに全振りしたアルゴリズム」は、残りの 7,000 ステップもの間、わかりきっている答えを確かめるためだけに完全ランダムな探索を馬鹿正直に続けてしまっています。

この「すでに実力差が十分に見えているのに、ハズレ台を平等に引き続けてしまう無駄」をなくし、「確率の低いハズレ台の探索は早めに切り上げ、当たりそうな優秀な台の『活用』へとえこひいきしていく」 賢いアルゴリズムが、今後の記事で説明する ε-greedy法（イプシロン・グリーディ法）や UCB1 アルゴリズム になります。

複数回のエピソードの結果のグラフ化

ここまでは1回のエピソード（1万ステップ）で得られたデータのグラフ化を行ってきましたが、多腕バンディット問題の目的関数は、1回きりの偶然の結果ではなく、無限回のエピソードを繰り返したときの平均値を表す「累積報酬の期待値」 です。そのため、たった1回のエピソードの結果をグラフにするだけでは、そのアルゴリズムの真の性能を正しく分析するには不十分です。

とはいえ、現実的に無限回のエピソードを実行することは不可能なため、本記事では 1,000回 のエピソードを繰り返し実行し、その平均値を視覚化して分析することにします。

1 回のエピソードを行う関数の定義

任意の回数のエピソードのシミュレーションを簡単に実行できるようにするために、まずは 1 回のエピソードを行ってその結果を返り値として返す play_episode という関数を定義することにします。

play_episode は今後さまざまなバンディット問題に対応できるよう、下記の仮引数を持つ関数として定義します。また、返り値はエピソード内で収集した各データを代入した変数を要素とする tuple 型とします。

仮引数	意味
mab	スロットマシンの当たりの設定を行った MAB クラスのインスタンス
stepnum	エピソードのステップ数

下記は play_episode を定義するプログラムです。先ほどの 1 エピソードを行うプログラムを関数にまとめただけなので説明は省略します。

def play_episode(mab, stepnum):
    total_coin = 0
    total_coins = []
    mean_coin = 0
    mean_coins = []
    slot = {
        "total_coin": [0] * mab.num,
        "mean_coin": [0] * mab.num,
        "selectnum": [0] * mab.num,
        "total_coins": [[] for _ in range(mab.num)],
        "mean_coins": [[] for _ in range(mab.num)],
        "selectnums": [[] for _ in range(mab.num)],
    }

    for step in range(stepnum):
        arm = np.random.choice(mab.num)
        coin = mab.play(arm)
        total_coin += coin
        total_coins.append(total_coin)
        mean_coin = total_coin / (step + 1)
        mean_coins.append(mean_coin)
        slot["total_coin"][arm] += coin
        slot["selectnum"][arm] += 1
        slot["mean_coin"][arm] = slot["total_coin"][arm] / slot["selectnum"][arm]
        for arm in range(mab.num):
            slot["total_coins"][arm].append(slot["total_coin"][arm])
            slot["mean_coins"][arm].append(slot["mean_coin"][arm])
            slot["selectnums"][arm].append(slot["selectnum"][arm])
            
    return total_coin, total_coins, mean_coin, mean_coins, slot

定義した関数を使って、先程と同じ設定の多腕バンディット問題に対する 1 万ステップのエピソードの処理を行うと、下記の実行結果から先程とほぼ同じ結果が得られることが確認できました。

total_coin, total_coins, mean_coin, mean_coins, slot = play_episode(mab, stepnum)
print(f"得られたコインの総数　 = {total_coin:5d}")
print(f"得られたコインの平均　 = {mean_coin:.3f}")
for arm in range(mab.num):
    print(f"スロットマシン {arm} (p={mab.p[arm]:.2f})")
    print(f"  得られたコインの総数 = {slot["total_coin"][arm]:5d}")
    print(f"  レバーを引いた回数　 = {slot["selectnum"][arm]:5d}")
    print(f"  得られたコインの平均 = {slot["mean_coin"][arm]:.3f}")

実行結果

得られたコインの総数　 =  4967
得られたコインの平均　 = 0.497
スロットマシン 0 (p=0.30)
  得られたコインの総数 =  1013
  レバーを引いた回数　 =  3347
  得られたコインの平均 = 0.303
スロットマシン 1 (p=0.50)
  得られたコインの総数 =  1687
  レバーを引いた回数　 =  3415
  得られたコインの平均 = 0.494
スロットマシン 2 (p=0.70)
  得られたコインの総数 =  2267
  レバーを引いた回数　 =  3238
  得られたコインの平均 = 0.700

10 回のエピソードのステップ当たりの累積報酬の比較

1,000 回のエピソードの平均（期待値）を計算する前に、まずは小手調べとして 10 回のエピソードを実行し、それぞれの「ステップ当たりの累積報酬」の推移を重ねてグラフ化してみることにします。そうすることで、1 回ずつのシミュレーションにどれくらいの「ブレ（ばらつき）」があるのかを視覚的に確かめることができます。

下記はその処理を行うプログラムです。

plt.title(f"10 エピソードの累積報酬の推移（原始モンテカルロ法の近似）")
for _ in range(10):  
    total_coin, total_coins, mean_coin, mean_coins, slot = play_episode(mab, stepnum)
    plt.plot(mean_coins)
plt.xlabel("ステップ数")
plt.ylabel("ステップ当たりのコインの枚数")
plt.show()

実行結果

グラフから読み取れる性質と考察

グラフを確認すると、10 本すべての線がステップ数の増加とともに下記のような性質を持つことがわかります。

• 後半に向かって 1 本の「一直線」に束ねられていく
  スタート直後（ステップ数が少ない段階）では、10 本のグラフが上下に激しく散らばっており、エピソードごとの運・不運のブレが大きく表れています。しかし、ステップ数が 1,000、20,000と増えていくにつれて、すべてのグラフが「当たりの平均（0.5）」のラインへと綺麗に収束し、ほぼ 1 本の一直線上に重なり合います。

• 原始モンテカルロ法が持つ低品質だが安定した（ブレの少ない）性能
  この結果は、探索に全振りした原始モンテカルロ法（の近似）を実行した場合は、何回勝負（エピソード）をやり直したとしても、最終的に得られる 1 回あたりの報酬はスロットマシンの当たりの平均である $p_{mean} = 0.5$ 付近で完全に固定されてしまう という、低品質ではあるが、安定したブレのない結果が得られる という性質を裏付けています。

この「10 本の線が 1 つの同じゴール（0.5）に吸い込まれてしまう」という性質を視覚的に捉えられたのは、複数回のエピソードを重ねて描画したからです。

1,000 回のエピソードの「ステップ当たりの累積報酬」の平均と標準偏差

先ほどの 10 本の重ね合わせのグラフから、ステップ数が増える（試行回数が多くなる）につれて、ステップ当たりの累積報酬の「ばらつき（ブレ）」がほとんどなくなるという傾向があることが視覚的に確認できました。しかし、具体的に数値としてどれくらいのブレに収まっているのかまではわかりません。

そこで、1,000 回のエピソードのシミュレーションを実際に実行し、各ステップにおける「累積報酬の平均値」と、ブレの具体的な度合いを表す「標準偏差」を厳密に計算してグラフ化しててみることにします。

1,000 回のエピソードで得られたデータの記録

下記は 1,000 回のエピソードを行い、得られた膨大なデータ（ステップあたりの累積報酬と各スロットマシンの当たりの確率の推定値）を効率よく記録するプログラムです。

これまでのプログラムでは、空のリスト（[]）を用意して、データを 1 つずつ append メソッドで末尾に追加していく処理を行ってきました。しかし、今回のように「複数のエピソード」「複数のステップ」「複数のマシン」が絡み合う多次元のデータを append で記録しようとすると、リストが 3 つの入れ子状になり、プログラムの見た目が非常に複雑でわかりづらくなるという問題が発生します。

そこで今回は、あらかじめすべてのデータを記録できる大きさ（形状：shape）の多次元配列（ndarray）を numpy で最初に用意しておき、インデックスを指定してデータを代入する処理 を採用しました。

プログラムで使用する多次元配列を代入する変数一覧は以下の通りです。

変数名	意味
mean_coins_by_step	shape が (1000, 10000) の 2 次元の ndarray。[episode][step] の形で、あるエピソードの特定のステップにおける「ステップ当たりの累積報酬」を記録する
slotratio_by_step	shape が (3, 1000, 10000) の 3 次元の ndarray。[arm][episode][step] の形で、あるエピソードの特定のステップにおける、各スロットマシンの「当たり確率の推定値（標本平均）」を記録する

• 1 ~ 3 行目：各変数の初期化
  • 2、3 行目：np.zeros を利用して、指定した形状（shape）のすべての要素が 0 の ndarray を生成して代入する
• 5 ~ 9 行目：1,000 回のエピソードを実行し、得られたデータを ndarray の適切な要素に記録する
  • 7 行目：numpy の ndarray は Python の list と異なり、非常に強力な「スライス表記」を利用でき、mean_coins_by_step[episode, :] のように「すべての要素」を意味する : を記述することで、1 万ステップ分のデータ（mean_coins）を繰り返し処理を使わずにまとめて代入することができる。 忘れた方は以前の記事を復習すること
  • 8, 9 行目：同様の処理を slotratio_by_step の各スロットマシンに対して行う

1  episodenum = 1000
2  mean_coins_by_step = np.zeros((episodenum, stepnum))
3  slotratio_by_step = np.zeros((mab.num, episodenum, stepnum))
4  
5  for episode in range(episodenum):  
6      total_coin, total_coins, mean_coin, mean_coins, slot = play_episode(mab, stepnum)
7      mean_coins_by_step[episode, :] = mean_coins
8      for arm in range(mab.num):
9          slotratio_by_step[arm, episode, :] = slot["mean_coins"][arm]

行番号のないプログラム

episodenum = 1000
mean_coins_by_step = np.zeros((episodenum, stepnum))
slotratio_by_step = np.zeros((mab.num, episodenum, stepnum))

for episode in range(episodenum):  
    total_coin, total_coins, mean_coin, mean_coins, slot = play_episode(mab, stepnum)
    mean_coins_by_step[episode, :] = mean_coins
    for arm in range(mab.num):
        slotratio_by_step[arm, episode, :] = slot["mean_coins"][arm]

上記のプログラムを実行すると、筆者のパソコンでは約 2.5 分ほどで処理が完了しました。

ここで、今回の記事の冒頭で説明した、原子モンテカルロ法そのまま実行して計算する場合の処理時間を思い出してください。原子モンテカルロ法で上記の処理を愚直に計算した場合は 2.5 分の 5000 倍の 12500 分 = 約 9 日もの時間 がかかります。原始モンテカルロ法を、探索のみを行うアルゴリズムで近似するというアプローチが如何に重要であるかがこの、ことから実感できたのではないでしょうか。

1,000 回のエピソードの「ステップ当たりの累積報酬」の平均の計算とグラフ化

mean_coins_by_step に記録された膨大なデータを元に、各ステップごとの 「1,000エピソードにおける累積報酬の平均値（期待値）」 を計算します。

この計算は、numpy の平均を計算する関数 np.mean を利用することで、すべてのステップの平均をまとめて計算することができます。mean_coins_by_step の形状（shape）は (1000, 10000) であり、エピソード（試行回数）は最初の「0 番目の軸（axis）」に並んでいます。従って、np.mean の引数に axis=0 を指定して呼び出すだけで、1,000 個のエピソードを対象とした平均の計算処理を、すべてのステップ（0 番の軸（axis）の方向）に対して一括して実行することができます ⁵。

下記は、以下は、1,000 回のエピソードにおける「ステップ当たりの累積報酬の平均」を計算してグラフ化するプログラムです。

少しややこしいですが、変数名は「ステップ当たりの累積報酬（mean_coins）」の平均（mean）であることを表すために、mean_mean_coins_by_step という名前を採用しました。また、データが少ないスタート直後のブレ幅の激しさが視覚的にハッキリと伝わるように、今回はあえて y 軸の範囲を固定しない（plt.ylim を実行しない）工夫を行っています。

plt.title(f"ステップ当たりの累積報酬の平均\n{episodenum} エピソード（原始モンテカルロ法の近似）")
mean_mean_coins_by_step = np.mean(mean_coins_by_step, axis=0)
plt.plot(mean_mean_coins_by_step)
plt.xlabel("ステップ数")
plt.ylabel("ステップ当たりの累積報酬の平均")
plt.show()

実行結果

グラフから読み取れる性質と考察

1,000 回ものエピソードの平均を計算したことで、原始モンテカルロ法を用いた多腕バンディット問題の目的関数（累積報酬の期待値）が持つ本来の性質が綺麗なデータとして浮かび上がってきました。

• ステップ数に関わらず、常に $p_{mean} = 0.5$ 付近の値をとり続ける
  スタート直後の大きなブレが収まった後は、グラフの線が $0.5$ のラインに完全に固定されます。どれだけ試行回数（ステップ数）を重ねてデータを集めても、実力が 0.5 から右肩上がりに成長することはありません。

• ステップ数が増えると、ブレが完全に消えていく
  最初は激しく上下していた線が、ステップ数の増加とともに一直線へと収束していきます。1 回のエピソードのグラフよりもさらに線が滑らかになっており、無限回（期待値）の性質へ近づいている様子が視覚的に示されました。

1,000 回のエピソードの「ステップ当たりの累積報酬」の標準偏差の計算とグラフ化

続いて、「ステップ当たりの累積報酬」のブレの度合いを表す「標準偏差（standard deviation）」を計算してグラフ化します。

標準偏差は下記のプログラムのように np.std 利用することで、平均と同様にすべてのステップを一括で簡単に計算できます。エピソード（試行回数）の軸に沿って集計を行うため、こちらも引数に axis=0 を指定します。

なお、グラフから終盤の具体的な数値を読み取るのは難しいため、最終ステップ（インデックスが -1）の標準偏差を print で表示しました。

plt.title(f"ステップ当たりの累積報酬の標準偏差\n{episodenum} エピソード（原始モンテカルロ法の近似）")
std_mean_coins_by_step = np.std(mean_coins_by_step, axis=0)
plt.plot(std_mean_coins_by_step)
plt.xlabel("ステップ数")
plt.ylabel("ステップ当たりの累積報酬の標準偏差")
print(f"最終ステップの標準偏差 = {std_mean_coins_by_step[-1]:.5f}")

実行結果

最終ステップの標準偏差 = 0.00501

グラフから読み取れる性質と考察

グラフを確認すると、ステップ数（試行回数）がブレ幅を強力に抑え込んでいく様子が明らかになります。

• 標準偏差が右肩下がりに「急激に減少」していく
  スタート直後（ステップ数が少ない段階）は高い値を示しますが、ステップ数が増えるにつれてグラフは急激な下り坂を描き、ゼロに向かって収束していきます。

• 最終ステップでのブレは、わずか「0.005」
  最終ステップでの標準偏差は約 0.005 という極めて小さな値まで減少していることが実証されました。これは、1,000 回のエピソード間での結果の数値にほとんど差（ブレ）がなくなり、どのシミュレーションでもほぼ完全に同じ $p_{mean} = 0.5$ という結果に落ち着くことを表しています。

この結果から、「試行回数を増やせば増やすほど、平均値のブレ（標準偏差）はゼロに向かって小さくなっていく」という大数の法則の性質が確認できました。

特定のステップの「ステップ当たりの累積報酬」のヒストグラムの描画

「ステップ当たりの累積報酬」のばらつきがどのように推移しているのか、その具体的なデータの広がり（分布）を確認するために、特定のステップにおける 1,000 回のエピソードのヒストグラム（度数分布図）を描画することにします。

今回は試行回数の違いによる変化を明快に比較するため、0、9、99、999、9,999ステップ目のそれぞれの分布を以下のプログラムで順番に描画します。

for step in [0, 9, 99, 999, 9999]:
    plt.title(f"{step} ステップでのステップ当たりの累積報酬\n{episodenum} エピソード（原始モンテカルロ法の近似）")
    plt.hist(mean_coins_by_step[:, step], bins=20)
    plt.xlabel("ステップ当たりの累積報酬")
    plt.ylabel("回数")
    plt.show()

実行結果

グラフから読み取れる性質と考察

5 枚のヒストグラムをステップ数が小さい順から見比べていくと、統計学における極めて重要でな 2 つの現象がはっきりと確認できます。

• 横軸の「データの広がり」がみるみる縮んでいく
  最初の 0 ステップ（レバーを 1 回引いた場合）の時点では、獲得コインが 0 か 1 しかないので両端にデータが極端に分かれています。しかし、ステップ数が増えるにつれてデータのばらつきが急速に小さくなり、下記の表のようにグラフの横幅がどんどん狭くなっていきます。最後の 9,999 ステップ目では、すべてのデータが目標の平均値である $p_{mean} = 0.5$ のすぐ周り（0.485 ~ 0.525 の範囲）だけに凝縮されます。

ステップ数	グラフの x 軸の範囲	x 軸の幅
0	0 ~ 1	1
9	0.1 ~ 0.9	0.8
99	0.35 ~ 0.65	0.3
999	0.46 ~ 0.54	0.08
9999	0.485 ~ 0.525	0.04

• 綺麗な「左右対称の釣り鐘型（正規分布）」になる
  ステップ数が増えていくにつれて、グラフの形状がなだらかな左右対称の山形、すなわち正規分布（ガウス分布）にそっくりの形状へと近づいていくことがわかります。これは、以前の記事で説明した中心極限定理による減少です。

これまでの記事でに学んだ「大数の法則」によって平均値が 0.5 に収束し、さらに「中心極限定理」によってそのブレ方が綺麗な正規分布を描いて細くなっていくという、確率論の核心となる 2 大定理の正しさが、シミュレーションデータをグラフ化することで実際に確認することができました。

平均からのブレ（標準偏差）のグラフ化

以前の記事で説明したように、正規分布に従う確率分布から抽出した標本は、約 99.7 % の確率で 平均から標準偏差の 3 倍（± 3σ）の範囲 に収まるという性質があります。先ほどのヒストグラムの分析によって、ステップ数が多くなるとデータの分布が綺麗な正規分布の形状に近づいていくことが実証されました。

そこで、1,000 回分のエピソードのデータを用いて、「ステップ当たりの累積報酬の平均値」を中心に、データのほぼすべて（99.7%）が入り込む 平均 ± 3σ の信頼区間の広がり をグラフで視覚化することにします。。

統計学では、母集団の真の値（母平均など）が含まれることが、かなり確信 (confident) できる数値範囲のことを信頼区間と呼びます。ただし、1,000 回のエピソードから計算した「ステップ当たりの累積報酬の平均値」は母平均の推定値なので、グラフで描画した範囲は厳密には母平均が 99.7 % の確率で含まれる平均 ± 3σ の信頼区間とは微妙に異なるものです。

ただし、その違いはグラフ化した場合に肉眼で確認ほどではないので、下記のグラフで描画された範囲を平均 ± 3σ の信頼区間と考えても問題はありません。

参考までに信頼区間の Wikipedia の項目のリンクを下記に示します。

下記はそのグラフを描画するプログラムです。

• 3, 4 行目：各ステップにおける「ステップ当たりの累積報酬の平均 ± 標準偏差の 3 倍」の上限線と下限線を計算する
• 5 行目：np.fill_between を利用して、計算した上限と下限の間の領域を塗りつぶす。その際に、alpha=0.5 を指定して領域を 半透明 にすることで、2 行目で描画した「平均値の推移線」が影に隠れて見えなくならないように工夫した

1  plt.title(f"ステップ当たりの累積報酬の平均 ± 3σ の範囲\n{episodenum} エピソード（原始モンテカルロ法の近似）")
2  plt.plot(mean_mean_coins_by_step)
3  mean_minus_3sigma =  mean_mean_coins_by_step - 3 * std_mean_coins_by_step
4  mean_plus_3sigma =  mean_mean_coins_by_step + 3 * std_mean_coins_by_step
5  plt.fill_between(np.arange(stepnum), mean_minus_3sigma, mean_plus_3sigma, alpha=0.5)
6  plt.ylim(0, 1)
7  plt.xlabel("ステップ数")
8  plt.ylabel("ステップ当たりの累積報酬の平均")
9  plt.show()

行番号のないプログラム

plt.title(f"ステップ当たりの累積報酬の平均 ± 3σ の範囲\n{episodenum} エピソード（原始モンテカルロ法の近似）")
plt.plot(mean_mean_coins_by_step)
mean_minus_3sigma =  mean_mean_coins_by_step - 3 * std_mean_coins_by_step
mean_plus_3sigma =  mean_mean_coins_by_step + 3 * std_mean_coins_by_step
plt.fill_between(np.arange(stepnum), mean_minus_3sigma, mean_plus_3sigma, alpha=0.5)
plt.ylim(0, 1)
plt.xlabel("ステップ数")
plt.ylabel("ステップ当たりの累積報酬の平均")
plt.show()

実行結果

グラフから読み取れる性質と考察

実行結果のグラフから、ブレが 0 に近づいていくという 不確実性の収束 を確認できます。

• 帯の収束
  ステップ数が進むにつれて ±3σ の半透明の帯が急速に細くなり、ブレが減少します。

• 最終ステップの精度
  最終ステップの標準偏差は 0.005 であり、99.7% の確率で $0.5 \pm (0.005 \times 3)$、すなわち 0.485 〜 0.515 の範囲に収まります。

このグラフにより、試行回数を増やすことで平均値が一定値に収束し、不確実性が減少するという大数の法則が視覚的に表現できました。

1,000 回のエピソードの「スロットマシンの当たりの確率の推測値」の平均と標準偏差

次に、各スロットマシンの当たりの推定値（コインの総数の平均）に関しても、同様に 1,000 回のエピソードのデータを重ね合わせて平均（期待値）と標準偏差（ブレ幅）の推移をグラフで可視化することにします。また、「各スロットマシンの当たりの確率の推定値」の標準偏差の推移についてもグラフで可視化することにします。

下記は、それらのグラフを描画するプログラムです。

• 1, 2 行目：3 次元の ndarray が代入された slotratio_by_step の形状（shape）は (mab.num, episodenum, stepnum) であり、1 番目のインデックス（2 番目の軸）にエピソードが格納される。そのため、axis=1 を指定して np.mean と np.std を計算すれば良い。計算結果は、元の 3 次元の ndarray から episodenum の軸が取り除かれ、(mab.num, stepnum) の 2次元の ndarray、すなわち mean[arm][step] というデータ構造になる

 1  mean = np.mean(slotratio_by_step , axis=1) 
 2  std = np.std(slotratio_by_step , axis=1)
 3  
 4  plt.title(f"スロットマシンの当たりの確率の推測値の平均 ± 3σ の範囲\n{episodenum} エピソード（原始モンテカルロ法の近似）")
 5  for arm in range(mab.num):
 6      plt.plot(mean[arm], label=f"arm {arm} (p={mab.p[arm]:.2f})")
 7      mean_minus_3sigma = mean[arm] - 3 * std[arm]
 8      mean_plus_3sigma = mean[arm] + 3 * std[arm]
 9      plt.fill_between(np.arange(stepnum), mean_minus_3sigma, mean_plus_3sigma, alpha=0.5)
10      print(f"スロットマシン {arm} (p={mab.p[arm]:.2f})")
11      print(f"  最終ステップの当たりの確率の推測値の平均　　 = {mean[arm][-1]:.3f}")
12      print(f"  最終ステップの当たりの確率の推測値の標準偏差 = {std[arm][-1]:.5f}")
13  plt.ylim(0, 1)
14  plt.xlabel("ステップ数")
15  plt.ylabel("ステップ当たりの累積報酬の平均")
16  plt.legend()
17  plt.show()
18  
19  plt.title(f"スロットマシンの当たりの確率の推測値の標準偏差\n{episodenum} エピソード（原始モンテカルロ法の近似）")
20  for arm in range(mab.num):
21      plt.plot(std[arm], label=f"arm {arm} (p={mab.p[arm]:.2f})")
22  plt.legend()
23  plt.show()

行番号のないプログラム

mean = np.mean(slotratio_by_step , axis=1) 
std = np.std(slotratio_by_step , axis=1)

plt.title(f"スロットマシンの当たりの確率の推測値の平均 ± 3σ の範囲\n{episodenum} エピソード（原始モンテカルロ法の近似）")
for arm in range(mab.num):
    plt.plot(mean[arm], label=f"arm {arm} (p={mab.p[arm]:.2f})")
    mean_minus_3sigma = mean[arm] - 3 * std[arm]
    mean_plus_3sigma = mean[arm] + 3 * std[arm]
    plt.fill_between(np.arange(stepnum), mean_minus_3sigma, mean_plus_3sigma, alpha=0.5)
    print(f"スロットマシン {arm} (p={mab.p[arm]:.2f})")
    print(f"  最終ステップの当たりの確率の推測値の平均　　 = {mean[arm][-1]:.3f}")
    print(f"  最終ステップの当たりの確率の推測値の標準偏差 = {std[arm][-1]:.5f}")
plt.ylim(0, 1)
plt.xlabel("ステップ数")
plt.ylabel("ステップ当たりの累積報酬の平均")
plt.legend()
plt.show()

plt.title(f"スロットマシンの当たりの確率の推測値の標準偏差\n{episodenum} エピソード（原始モンテカルロ法の近似）")
for arm in range(mab.num):
    plt.plot(std[arm], label=f"arm {arm} (p={mab.p[arm]:.2f})")
plt.legend()
plt.show()

実行結果

スロットマシン 0 (p=0.30)
  最終ステップの当たりの確率の推測値の平均　　 = 0.300
  最終ステップの当たりの確率の推測値の標準偏差 = 0.00796
スロットマシン 1 (p=0.50)
  最終ステップの当たりの確率の推測値の平均　　 = 0.500
  最終ステップの当たりの確率の推測値の標準偏差 = 0.00870
スロットマシン 2 (p=0.70)
  最終ステップの当たりの確率の推測値の平均　　 = 0.700
  最終ステップの当たりの確率の推測値の標準偏差 = 0.00807

グラフから読み取れる性質と考察

グラフから、すべてのスロットマシンの「当たりの確率の推測値」の平均が真の確率に収束し、試行回数とともにすべてのスロットマシンの標準偏差が同じように 0 に向かって減少することが確認できます。

今後の記事で紹介する賢い強化学習アルゴリズムでは、当たりの確率が高いと推測されたスロットマシンをえこひいきして選択するため、当たりの確率が高い程スロットマシンの標準偏差だけが急速に 0 に近づいていくという現象が発生します。今回のグラフは、それらと比較するための重要な基準点となります。

上記のグラフでは 3 つのスロットマシンの標準偏差が完全に同じように減少しているように見えるかもしれませんが、実際には 1 番のスロットマシンと 0, 2 番のスロットマシンの標準偏差の値は微妙に異なります。

例えば、最終ステップの標準偏差はスロットマシン 0, 2 では約 0.8, スロットマシン 1 では約 0.87 となっています。

このような違いが生じる理由は、以前の記事で説明したように、多腕バンディット問題のスロットマシンの標準偏差が、スロットマシンの当たりの確率を $p$ とすると $\sqrt{p - p^2}$ になるからです。

今回の例に当てはめると、それぞれのスロットマシンの標準偏差は下記のようになります。

• スロットマシン 0：$\sqrt{0.3 - 0.3^2} = 約 0.459$
• スロットマシン 1：$\sqrt{0.5 - 0.5^2} = 約 0.5$
• スロットマシン 2：$\sqrt{0.7 - 0.7^2} = 約 0.459$

標本数が $m$ の標本平均の標準偏差は、母集団の標準偏差を $σ$ とすると、$\frac{σ}{\sqrt{m}}$ で計算することができます。最終ステップでは 3 つのスロットマシンが 1 万回の 1/3 の約 3333 回選択されるので、標本数が 3333 の標本平均になります。従って、各スロットマシンの当たりの推定値（標本数が 3333 の標本平均）の標準偏差は下記のように、先程の実行結果の最終ステップの標準偏差にかなり近い値になります。

• スロットマシン 0：$\frac{0.459}{\sqrt{3333}} = 約 0.00795$
• スロットマシン 1：$\frac{0.5}{\sqrt{3333}} = 約 0.00867$
• スロットマシン 2：$\frac{0.459}{\sqrt{3333}} = 約 0.00795$

今回の記事のまとめ

今回の記事では、多腕バンディット問題を解くアルゴリズムの性質を直観的に理解できるように、グラフで描画して視覚化するプログラムを実装しました。

具体的には、原始モンテカルロ法のほぼ正確な近似である「探索のみを行うアルゴリズム」をグラフ化し、前回の記事で説明した探索に極振りした場合の性質がグラフによってわかりやすく視覚化できたことを実際に確認しました。

次回の記事では、視覚化を行うプログラムを拡張して、多腕バンディッド問題を解く任意のアルゴリズムの視覚化を簡単に実装できるようにします。また、そのプログラムを利用して前回の記事で説明した、活用に全振りした場合のアルゴリズムの性質を視覚化し、今回の記事と同様にグラフから読み取れる性質について考察します。

本記事で入力したプログラム

リンク	説明
marubatsu.ipynb	本記事で入力して実行した JupyterLab のファイル

次回の記事

近日公開予定です

1. 前回の記事で説明したように、探索または活用に全振りした場合は例外的に $f_{max}$ に近づくことはありません ↩

2. 前回の記事で説明した、すべての回で当たりの確率最も高いスロットマシンのレバーを引いた場合に得られるコインの総数のことです。$f_{max}$ は、最高の当たりの確率を $p_{max}$、スロットマシンのレバーを引ける回数を $m$ とすると、$f_{max} = m \times p_{max}$ で計算することができます ↩

3. 今回の記事の後半で取り扱う具体例の処理では、一週間以上もの時間がかかります ↩

4. グラフの線の形状（スタイル）の指定方法については、こちらのMatplotlib の公式ドキュメントを参照して下さい ↩

5. axis を用いて ndarray を特定の軸に沿ってまとめて計算する方法について忘れた方は、以前の記事を復習して下さい。 ↩