チョコボール統計スマートフォン用の表示で見るチョコボール統計チョコボールの秘密を統計解析で明らかにしていく。おもちゃのカンヅメ欲しい。終わりに参考文献(お世話になった書籍)広告はじめに当ブログでは、日々のチョコボール開封結果に基づいて金のエンゼルと銀のエンゼルの出現確率を推定しています。このモデルは、現実を正しく表現していません。チョコボールのエンゼルは金のエンゼルと銀のエンゼルが同時に現れることは無いですよね。「金のエンゼル」、「銀のエンゼル」、または、「ハズレ」の3状態のいずれかが実現します。これも簡単のために、この期間のデータを除外して推定してきました。単純にデータが減るだけなので、この期間のデータも積極的に利用していきたいと思います。現実に即したモデルでエンゼルの出現確率を推定していこうと思います。【トップに戻る】アプローチ具体的に3つの状態(金のエンゼル、銀のエンゼル、ハズレ)を合わせて推定するモデルを考えていきます。モデルを考える前に、前提、仮定を明示しておきます。(金のエンゼル2倍キャンペーンは除く)。フレーバーによって出現確率は変わらないと仮定します。(そのように公表されているので信じることにします)。モデル次に、チョコボールのエンゼルの発生過程をモデル化します。1章で述べたとおり実現値は3値{金、銀、ハズレ}です。そのため、多項分布を利用します。全試行回数Mのうち、{金、銀、ハズレ}のそれぞれの個数を表すベクトルです。{金、銀、ハズレ}の出現確率を表す3次元のベクトルで、このパラメータが推定対象です。今回は、最尤推定ではなく、ベイズ推論を行いたいので、事前分布を設定する必要があります。そこで事前分布には、多項分布の共役事前分布であるディリクレ分布を利用します。ディリクレ分布のハイパーパラメータです。ここまでを合わせて、Kruschkeのダイアグラム*4にすると以下のようなモデルとなります。図1.単純な多項分布のモデル.多項分布のパラメータの事前分布にはディリクレ分布を利用する.次に、金のエンゼル2倍キャンペーンのデータを利用する方法を考えます。実はこのアイディアはネタが丸かぶりしている以下の記事を参考にさせていただきました(参考にというかほぼそのままです...)。bayesmax.sblo.jp上記の記事では、このキャンペーン期間のデータには確率に重みが付くというモデルにされています。以下のベクトルとなります。2.0データ今回利用するデータは、2017年11月~2019年7月までに当ブログ内で開封した566箱が対象です。グアムで購入したチョコボールについては、金のエンゼルも銀のエンゼルも共に存在しないため、対象外としています*5。データをまとめると以下の通りです。以下のように2つの期間に分けたデータとしました。1実験パラメータ推論2.2節に示した多項分布モデルのパラメータを推論します。2.2節で述べたとおり今回の実験では、事前分布には共役事前分布であるディリクレ分布を利用します。そのため、ベイズの定理に従って事後分布を計算すると以下の通りディリクレ分布になります*6。ディリクレ事前分布のハイパーパラメータで、今回は適当な値を設定します。今回はあえてPyMC3*7を利用し、サンプルによる近似事後分布を求めます(MCMC)*8。今回は近似事後分布を求めました*9。具体的なコードは、以下を参照ください。実験結果3章で示したデータを利用して、金のエンゼルと銀のエンゼルの出現確率を推定した結果を示します。2章で述べたとおり、金のエンゼル2倍キャンペーンを含めないモデルと含めるモデルをそれぞれ推定しました。金のエンゼル2倍キャンペーンを除いた場合まず、問題を単純にするために金のエンゼル2倍キャンペーンを除いた場合の結果です。x1は、多項分布のパラメータの事後分布をプロットしたものです。3つの図はそれぞれ左から、ハズレ,銀のエンゼル,金のエンゼルの出現確率を表しています。それぞれの図には、期待値(mean)と95%HPD*10が記載されています。図x1から、95%の確率でハズレを引き、銀のエンゼルは約4.6%(2.8%~6.6%の間)で出現するであろうと推定できていることがわかります。金のエンゼルは0.2%(0%~0.7%の間)であると推定しています*11。ここで利用したデータでは金のエンゼルは一つも出現していません。事後分布では0.7%以下になるだろうと柔軟な予測ができています。次に、2章で述べたように、エンゼルの出現確率に重みを載せて、全てのデータを利用して推定した結果を示します。x2は全てのデータを利用して、多項分布のパラメータの事後分布を推定した結果です。図のそれぞれの位置関係は、先の図x1と同じです。図x1と図x2を比較すると、ほぼ同じような結果なのですが、金のエンゼルの事後分布を見ていただくと、分布の形状が微妙に異なっていることがわかると思います。今回のデータでは、金のエンゼル2倍キャンペーンの期間で金のエンゼルが一つ出ています*12。確率0%では無いことははっきりしているため、0の付近の山が少し下がり、分布が右にシフトしている様子が確認できます。推定結果には影響を与えていないということがわかります。買ったらエンゼルが当たるのか?エンゼルの出現確率が推定できたのなら、次は、何個買えばエンゼルが当たるのかの見積もりが気になりますよね。この分布を利用することで、エンゼルを得るまでに必要なチョコボールの個数を見積もることができそうです。4.2節の推定で得られた結果は確率分布(事後分布)でした。この確率分布にはデータがまだ十分でないための曖昧性が表現されているため、代表点で推定することは避けたいです。そのため、以下の積分を計算することで、事後分布を利用した予測結果を得ることができます。4.2節で推定した事後分布です。期待値を計算するということですね。サンプル集合として得られていることを思い出します。サンプル集合で事後分布を予測できているため、サンプルごとの平均で積分を計算することができます。金のエンゼル1枚と銀のエンゼルを5枚出すために必要なチョコボールの購入数を見積もって見ます。図x3.銀のエンゼルを5つ手に入れるまでに必要なチョコボールの個数の分布.図x3は銀のエンゼルを5つ得るまでに必要な個数の分布(累積確率)です。事後分布を使って推定した結果(青線)と事後分布の期待値を使って推定した結果(赤線)を載せています。この図から、100個程度のチョコボールを買うことで、銀のエンゼルが5個得られる確率が50%を超えそうだということがわかります。銀のエンゼルの予測は、期待値を使った場合も事後分布を使った場合も概ね同じ程度であることがわかります。図x4.金のエンゼルを1つ手に入れるまでに必要なチョコボールの個数の分布.次に、図x4は金のエンゼルを1枚得るまでに必要な個数の分布(累積確率)です。こちらの図でも事後分布を使って推定した結果(青線)と事後分布の期待値を使って推定した結果(赤線)を載せています。この図から、金のエンゼルを得るためには、250個ほど買うことで50%を超えるということがわかります。1,000個も買えば80%の確率で金のエンゼルが当たるという予想になっています。期待値を使って予測した方がポジティブな予測になっているのがわかります。期待値が少し高めなのだろうということがわかります。終わりに以上本記事は、金のエンゼルと銀のエンゼルを合わせて推定してみました。金のエンゼルは0.3%、銀のエンゼルは4.6%という推定値となりました。事後分布を活用したメリットに触れることができました。正確に推定するにはもっと多くのデータが必要です。そのため、今後も計測を続けていこうと思います。見るベイズ的な手法を軸にした機械学習の参考書。単行本この商品を含むブログを見るpymc3を使った入門書。理論面も多少書かれており、こちらで手を動かしながら学習するのが良いと思う。広告Amazonの欲しいものリスト作りました。チョコボールのカンパ募集中です。Simerstメディア:この商品を含むブログを見る【トップに戻る】*1:独立に計算するということは、金のエンゼルと銀のエンゼルが共に出現する状態を許容してしまいます。いうか、金のエンゼルと銀のエンゼルは互いに関係ない別のものと考えることになっちゃいます。*2:たまたま、金のエンゼルの出現確率が極端に低いためにこれまでに問題にはならなかっただけです*3:ただただ単純にサボっていただけという。知りたいパラメータが95%の確率で含まれる区間*11:正確には、95%HPDの下限はは0%ではありません。チョコボール統計*13:MAP値もよく使われます*14:共役事前分布を使うことで、事後分布をディリクレ分布として推定することができるため、解析的に計算することは可能です。チョコボール計測(ピーナツ)こんばんは、チョコボール統計研究所です。チョコボール計測(ピーナツ)こんばんは、チョコボール統計研究所です。チョコボール計測(ピーナツ)こんばんは、チョコボール統計研究所です。チョコボール計測(ミルクココア)こんばんは、チョコボール統計研究所です。チョコボール計測(ピーナツ)こんばんは、チョコボール統計研究所です。週間)おもちゃのカンヅメ欲しい内容に間違いや気づいたことがありましたら、お気軽にツッコミお願いします。(1)月別アーカイブカンパ募集中カンパ募集中(Amazon欲しいものリスト)リンクはてなブログブログをはじめる週刊はてなブログはてなブログPro広告はてなブログをはじめよう!hippy-hikkyさんは、はてなブログを使っています。ブログを報告するスターをつけました引用をストックしましたストック一覧を見る閉じる引用するにはまずログインしてくださいログイン閉じる引用をストックできませんでした。再度お試しください閉じる限定公開記事のため引用できません。