～幾何分布ってなに？～統計検定2級を勉強してみよう！【経営工学を専門にしている大学生の日記】

はじめに

こんにちは！しゅんです！

今回は幾何分布について解説していきます！

幾何分布は確率・統計で登場する数学用語です。この記事では幾何分布について図を使って分かりやすく解説していきます。

また今回はpythonを使って幾何分布を扱っていこうと思います。

それでは解説していきましょう！

統計検定2級に関する他の記事はこちらから見れます！
ぜひ他の記事も読んでみてください！

もっとみる

このブログの簡単な紹介はこちらに書いてあります。
興味があったら見てみてください。

このブログについて

このブログでは経営工学を勉強している現役理系大学生が、経営工学に関することを色々話していきます！

ぼくが経営工学を勉強している中で感じたことや、興味深かったことを皆さんと共有出来たら良いなと思っています。

そもそも経営工学とは何なのでしょうか。Wikipediaによると

経営工学（けいえいこうがく、英: engineering management）は、人・材料・装置・情報・エネルギーを総合したシステムの設計・改善・確立に関する活動である。そのシステムから得られる結果を明示し、予測し、評価するために、工学的な分析・設計の原理・方法とともに、数学、物理および社会科学の専門知識と経験を利用する。

引用元 : 経営工学 – Wikipedia

長々と書いてありますが、要は経営、経済の課題を理系的な観点から解決する学問です。

幾何分布ってなに？

幾何分布を理解するためにはベルヌーイ分布を知っておく必要があります。ベルヌーイ分布とはコイントスのように一定の確率で2種類の結果のいずれかが必ず起こる試行におけるその結果の確率分布のことです。

ベルヌーイ分布については前回の記事で解説しているのでより詳しいことが知りたい場合はそちらもぜひ見てみてください！

2023年5月7日～ベルヌーイ分布ってなに？～統計検定2級を勉強してみよう！【経営工学を専門にしている大学生の日記】

幾何分布はこのベルヌーイ試行を初めて成功するまで繰り返したときの確率分布です。例えば下の図を見てみましょう。

おもてが出るまでコインを投げる試行を考えてみましょう。上の図のように5回目で初めておもてが出る確率はいくつになるでしょうか。下が計算式となります。

\((\frac{1}{2}) × (\frac{1}{2})^4 =\frac{1}{32}=0.03125\)

この式を簡単に説明すると、\(\frac{1}{2}\)の確率でおもてが出て\(\frac{1}{2}\)の確率でうらが出ます。そのためおもてが1回、うらが4回でる確率が\((\frac{1}{2}) × (\frac{1}{2})^4\)となります。

計算式が二項分布と似ていますが、計算式は違います。二項分布は何回目におもてが出るかも考えていましたが、幾何分布では最後しかおもてが出ないことが確定しているのでコンビネーションの計算は必要ありません。

以上より上記の式が出てきます。より詳しい説明は高校数学で確率を勉強してみてください！

この計算式からコインを繰り返し投げるとき5回目にはじめておもてが出る確率は0.03125ということが分かりました。例えば3回目で初めておもてが出る確率は下のようになります。

\((\frac{1}{2}) × (\frac{1}{2})^2 =\frac{1}{8}=0.125\)

より一般化してみましょう。コインを繰り返し投げるとき、\(x\)回目に初めておもてが出る確率は以下のように計算できます。

\((\frac{1}{2}) × (\frac{1}{2})^{x-1}\)

この式は成功確率\(\frac{1}{2}\)のベルヌーイ試行を初めて成功するまで繰り返す幾何分布の確率関数の式となります。\(x\)に初めて成功した回数を代入することによってそれぞれの成功回数となる確率を求めることができます。

さらに一般化しましょう。成功確率が\(p\)のベルヌーイ試行を初めて成功するまで繰り返す幾何分布の確率関数は以下のようになります。

\(p × (1-p)^{x-1}\)

この式が一般の幾何分布の確率関数となります。

幾何分布の確率関数

\(P(X=x) \equiv p × (1-p)^{x-1}\)　\((x=1,2,…)\)

幾何分布の期待値と分散

幾何分布の期待値と分散は下のようになります。

幾何分布の期待値と分散

期待値 : \(\mu = \frac{1}{p}\)
　分散 : \(\sigma^2 = \frac{1-p}{p^2}\)

期待値はその確率分布が取る値の平均を表します。例えば成功確率が\(\frac{1}{5}\)の試行において5回目でおもてになるだろうということを表しています。分散はその確率分布が取る値にどれだけばらつきがあるかを表します。

日常で見る幾何分布の例

二項分布はベルヌーイ試行を複数回行う場合を考えているので、ベルヌーイ分布の記事と同じように日常でも結構見かけます。日常での具体例があった方がイメージがつかみやすいと思うのでこの章ではいくつか日常の例を紹介していきます。

コイントス

これは先ほども説明しましたがコイントスはベルヌーイ試行です。コイントスは成功確率\(p\)が\(\frac{1}{2}\)となるようなベルヌーイ試行となります。したがってコイントスをして初めて表が出るまでに投げた回数は\(p=\frac{1}{2}\)の幾何分布に従います。

コイントスをこれまでに説明した式で表すと以下のようになります。

確率関数 :
\(P(X=x) \equiv \frac{1}{2} × (\frac{1}{2})^{x-1}\)　\((x=0,1,2,…)\)

期待値 : \(\mu = 2\)
　分散 : \(\sigma^2= 2\)

サイコロ

サイコロを1回投げて3,6が出る場合を成功、それ以外が出る場合を失敗としたときにこれは成功確率\(p\)が\(\frac{1}{3}\)のベルヌーイ試行となります。したがって3,6が初めて出るまでにサイコロを投げた回数は\(p=\frac{1}{3}\)の幾何分布に従います。

これまでに説明した式で表すと以下のようになります。

確率関数 :
\(P(X=x) \equiv \frac{1}{3} × (\frac{2}{3})^{x-1}\)　\((x=1,2,…)\)

期待値 : \(\mu = 3\)
　分散 : \(\sigma^2= 6\)

サッカーのPK

例えばあるサッカー選手がPKを決める確率が\(\frac{1}{10}\)だとします。このときこのサッカー選手が初めてPKを決めるまでにボールを蹴った回数は\(p=\frac{1}{10}\)の幾何分布に従います。

これまでに説明した式で表すと以下のようになります。

確率関数 :
\(P(X=x) \equiv \frac{1}{10} × (\frac{9}{10})^{x-1}\)　\((x=1,2,…)\)

期待値 : \(\mu = 10\)
　分散 : \(\sigma^2= 90\)

幾何分布をpythonで扱う

最後にpythonを使って幾何分布を扱いたいと思います。プログラミングの知識が必要なので、分からなければプログラムの所は飛ばして結果だけ確認してみてください。

事前準備

pythonで確率・統計を扱うときはscipyというライブラリをよく使います。その中でも今回は幾何分布を扱いたいのでgeomをインポートします。geomは幾何分布(geometric distribution)のことです。また数値計算でよく使うnumpy、グラフ作成でよく使うmatplotlibというライブラリもインポートしておきます。

import scipy.stats as stats
from scipy.stats import binom
import matplotlib.pyplot as plt
import numpy as np

幾何分布の確率を計算する

ということでまずは成功確率\(p\)を0.5として幾何分布の確率関数を計算します。例えば4回目で初めて成功する確率を求めてみましょう。

p = 0.5 # 成功確率
k = 4 # 4回目で初めて成功する確率を求める

geo_dist = geom(p) # 幾何分布のインスタンスを生成

prob = geo_dist.pmf(k) # 4回目で初めて成功する確率を計算
print(prob) # 結果を出力

ということで成功確率が0.5の幾何分布において、4回目で初めて成功する確率は0.0625となることが分かりました。

幾何分布の確率関数をグラフ化する

次にこの幾何分布をグラフ化してみましょう。つまり初めて成功する回数が1,2,…回のときそれぞれに対して確率がいくつになるかを棒グラフで表します。

p = 0.5 # 成功確率
k = np.arange(1, 11) # 初めて成功する回数
geo_dist = geom(p) # 幾何分布のインスタンスを生成
pmf = geo_dist.pmf(k) # 確率質量関数を計算

# 棒グラフを描画
plt.bar(k, pmf, color = "blue", alpha = 0.5)
plt.xticks(k)
plt.xlabel("Number of trials until first success")
plt.ylabel("Probability")
plt.title("Geometric distribution with p = 0.5")
plt.show()

ということでこの図をみると回数が増えるにつれ確率が小さくなっていますね。幾何分布は初めて成功する回数が増えるにつれてその確率は小さくなっていきます。

このグラフは成功確率0.5の試行で1回目に成功する確率は0.5、2回目に成功する確率は0.25、3回目に成功する確率は0.125、…ということを表しています。当然初めて成功する回数が大きくなるにつれてその確率はどんどん小さくなっていくのでこのような分布になっています。

グラフでは10回までしか表示していません。11回目以降も存在しますが、ほぼ0で違いが分かりづらいので10回目まで表示しました。

幾何分布の累積分布関数をグラフ化する

累積分布関数とは、それまでの確率を全部足し算した関数です。例えば\(x=3\)のとき、幾何分布の累積分布関数の値は1回目に初めて成功する確率、2回目に初めて成功する確率、3回目に初めて成功する確率を足し算した値になります。

それではpythonで累積分布関数をグラフ化してみましょう。\(p=0.5\)の幾何分布の累積分布関数を作成します。

p = 0.5
x = np.arange(1, 11)

# 幾何分布の累積分布関数を計算する
cumulative_dist = geom.cdf(x, p)

# グラフを描画する
plt.plot(k, cumulative_dist, color = "blue", alpha = 0.5, marker="o")
plt.yticks(np.arange(0,1.1,0.1))
plt.xticks(np.arange(0,11,1))
plt.xlabel('Number of trials until success')
plt.ylabel('Cumulative distribution function')
plt.title('Geometric distribution with p=0.5')
plt.show()

このグラフを見ると4回までに初めて成功する確率が90%を超えているのが分かりますね。つまりだいたい4回くらいやったら90%の確率で成功するということです。

シミュレーションしてみる

それでは実際にコイントスをシミュレーションしてみましょう。成功確率\(p=0.5\)として、初めて成功する（おもてが出る）まで投げ続ける試行を1回として、この試行を10回、100回、1000回、10000回おこなったときに結果がどのような分布になるかをpythonでシミュレーションしていきたいと思います。

プログラムは試行回数が10回のときのものを表しています。5行目の数字を変えることで好きな試行回数にすることができます。

# 幾何分布のパラメータ
p = 0.5

# 試行回数（ここを好きな数値に変更する）
num_trials = 10

# 初めて成功した回数を格納するリスト
successes = []

# 幾何分布に従って試行を行う
for i in range(num_trials):
    # 成功までの試行回数
    num_failures = np.random.geometric(p)
    
    # 初めて成功した回数をリストに追加する
    successes.append(num_failures)

# ヒストグラムを描画する
plt.hist(successes, bins=range(1, max(successes) + 2), align="left", rwidth=0.8, density=True, color = "blue", alpha = 0.5)
plt.xlabel("Number of trials until success")
plt.ylabel("Frequency")
plt.title("Geometric distribution with p=0.5 and {} trials".format(num_trials))
plt.xticks(np.arange(0,max(successes)+1,1))
plt.show()

このグラフは横軸が初めて成功するまでにコイントスをした回数、縦軸がその回数が結果として現れた確率を表しています。

左上から順番に10,100,1000,10000回となっています。10回の場合は思っている形とは全然違いますが、回数を増やすにつれてグラフの形が確率分布と同じようになっているのが分かると思います。

つまり確率分布とは、とある試行をたくさんやっていったら各結果の出現確率はこのような分布になっていくんだよってことを表しているんですね。

期待値をグラフで確認する

最後に期待値についてもグラフを使って確認してみましょう。先ほどのコイントスをやる回数を10,100,1000,10000回と増やしていったときに、各回数で得られた結果の平均値がどのように推移していくかをグラフを使って確認してみましょう。

p = 0.5 # 成功確率
n_sim = [10,100,1000,10000] # 試行回数を設定
result=[] # 平均値を格納するための空リストを作成

for n in n_sim:
  simulated_data = np.random.geometric(p, n) # 幾何分布に従う乱数を生成
  mean = simulated_data.mean() # 乱数の平均値を計算
  result.append(mean) # 結果に平均値を格納

# グラフを描画する
plt.plot(n_sim, result, marker ="o", color = "blue", alpha = 0.5, label = "Simulated")
plt.hlines(2, 10,10000, color = "orange", label = "Theoretical")
plt.xscale('log')
plt.xlabel('Number of trials')
plt.ylabel('Average number of trials until success')
plt.title('Geometric distribution with p=0.5')
plt.legend()
plt.show()

橙線は期待値で2のラインを表しています。つまり橙線は理論的にはコイントスしたときに初めておもてが出る回数は2回になるはずということを言っています。

一方青線は実際にシミュレーションした結果平均何回の試行回数で初めておもてが出たかを表しています。試行回数が少ないときは期待値との誤差がありますが、試行回数が多くなるにつれてどんどん期待値に近づいているのが分かりますね。

ランダムに試行されるのでプログラムを実行するたびにグラフは変わりますが、基本的には試行回数が大きくなるにつれて誤差は小さくなっていくはずです。