モンティ・ホール問題を一般化して解いてみる

####モンティ・ホール問題
モンティ・ホール問題とは確率論に関する有名な問題で、ご存じの方も多いでしょうが、設定を簡単に説明すると次のようになります。
1) 3つの選択肢があり、その中に当たりが1つ存在する。
2) まず、挑戦者が3つの中から1つを選ぶ。
3) 次に、当たりの選択肢を知っている司会者が、残りの2つの中からハズレのものを1つ選んで挑戦者に知らせる。
4) このとき、挑戦者は最初に選んだ選択を維持するべきか、それとも残るもう1つの選択肢に変更すべきか？
この状況において、選択肢を変更した方が有利であるというのがモンティ・ホール問題の答えなのですが、この問題設定をより一般化してみようというのが、この記事の目的です。
####モンティ・ホール問題を一般化する
一般化と言っても、そのやり方はいろいろと考えられますが、今回は次のような設定にしてみたいと思います。
1) 選択肢の数は$\,N(\ge3)$個、この中に当たりが1つ存在する。
2) 挑戦者はこの中から1つを無作為に選ぶ。
3) 当たりを知っている司会者は残りの$\,N - 1\,$個の中から、ハズレの選択肢$\,M(\ge1)\,$個を無作為に選んで、挑戦者に開示する。
4) このとき、挑戦者は最初の選択肢のままにすべきか、それとも残りの$\,N - M - 1\,$個の中から1つを無作為に選んで、その選択肢に変更すべきか？
ここで、残される選択肢は1つ以上でなければならないので、$N - M - 1 \ge1$。つまり、$1\le M \le N - 2$ となります。
この設定において、オリジナルのモンティ・ホール問題は$\,N = 3$、$M = 1\,$の状況に対応します。
さて、挑戦者は選択肢を維持した方がよいのか、変更した方がよいのか、実際に計算して確かめてみたいと思います。
####確率を計算する
挑戦者が最初に選んだ選択肢を$\,F$、変更した場合の選択肢を$\,S\,$と呼ぶことにします。
まず求めたい確率は、
$D_M$ : 特定の$\,M\,$個のハズレが開示された状況の下で、$W_F$ : $\,F\,$が当たりである条件付き確率$\,P(W_F\mid D_M)\,$です。
ベイズの定理より、
$$
P(W_F|D_M) = \frac{P(D_M \mid W_F)\,P(W_F)}{P(D_M)}
$$
右辺のそれぞれの確率を求めていくと、
まず、$P(W_F)$ : $F\,$が当たりである事前確率は、$N\,$個の中に当たりが1つなので、
$$
P(W_F) = \frac {1}{N}
$$
次に、$P(D_M \mid W_F)$ : $F\,$が当たりのときに、特定の$\,M\,$個のハズレが開示される条件付き確率は、当たりである$\,F\,$以外の$\,N -1\,$個の中から、特定の$\,M\,$個の組が選ばれる確率なので、
$$
P(D_M \mid W_F) = \frac {1}{{_{N-1}}{C}{_M}}
$$
最後に、$P(D_M)$ : 特定の$\,M\,$個のハズレが開示される事前確率は、$W_F$ : $F\,$が当たりである場合と、$\overline{W_F}$ : $F\,$がハズレである場合に分けることで
$$
P(D_M) = P(D_M \mid W_F)\,P(W_F)\,+\,P(D_M \mid \overline{W_F})\,P(\overline{W_F})  
$$
となります。右辺第一項の各確率はすでに分かっているので、第二項の各確率を求めると、
$P(\overline{W_F})$ : $F\,$がハズレである事前確率は、$W_F\,$の余事象の確率ですから、

$$\begin{align*}
P(\overline{W_F}) &= 1 - P(W_F) \newline
                  &= 1 - \frac{1}{N} \newline
                  &= \frac{N- 1}{N}
\end{align*}$$

次に、$P(D_M \mid \overline{W_F})$ : $F\,$がハズレであるときに、特定の$M\,$個のハズレが開示される条件付き確率は、
1) $F\,$以外の$\, N- 1\,$個の中で、その$M\,$個を除いた$N - M -1\,$個の中に当たりが存在し、
2) $F\,$と、その当たり以外の$N - 2\,$個の中から、特定の$M\,$個の組が選ばれる確率なので、

$$
P(D_M \mid \overline{W_F}) = \frac{N - M -1}{N - 1} \cdot \frac{1}{{_{N-2}}{C}{_M}} 
$$

となります。以上で求めた各確率と、
$$\begin{align*}
_{N-1}{C}{_M} &= \frac{(N-1)!}{(N-M-1)!\,M!} \newline
                &= \frac{N-1}{N-M-1} \cdot \frac{(N-2)!}{(N-M-2)!\,M!} \newline
                &= \frac{N-1}{N-M-1} \cdot\, _{N-2}{C}{_M}  
\end{align*}$$

の関係を利用すれば、

$$\begin{align*}
P(W_F \mid D_M) &= \frac{P(D_M \mid W_F)\,P(W_F)}{P(D_M)} \newline \newline
                &= \frac{P(D_M \mid W_F)\,P(W_F)}{P(D_M \mid W_F)\,P(W_F) + P(D_M \mid \overline{W_F})\,P(\overline{W_F})} \newline \newline
                &= \frac{\frac{1}{\_{N-1}{C}{\_M}} \cdot \frac{1}{N}}{\frac{1}{\_{N-1}{C}{\_M}} \cdot \frac{1}{N} \, + \frac{N-M-1}{N-1} \cdot \frac{1}{\_{N-2}{C}{\_M}} \cdot \frac{N-1}{N}} \newline \newline
                &= \frac{1}{1 \, + (N-M-1)\cdot \frac{\_{N-1}{C}{\_M}}{\_{N-2}{C}{\_M}}} \newline \newline
               &= \frac{1}{1 + (N-M-1) \cdot \frac{N-1}{N-M-1}} \newline \newline
               &= \frac{1}{1 + (N-1)} \newline
               &= \frac{1}{N}
\end{align*}$$

結局、ハズレが$M\,$個開示された後でも、$F\,$が当たりである確率は$\,1/N\,$のままであることが分かります。
次に、挑戦者が選択肢を$F\,$から$S\,$に変更するとどうなるでしょうか。
ここで求めたい確率は、$D_M$ : 特定の$M$個のハズレが開示された状況の下で、$W_S$ : $S\,$が当たりである条件付き確率$\,P(W_S \mid D_M)\,$です。
ベイズの定理より、
$$
P(W_S \mid D_M) = \frac{P(D_M \mid W_S)\,P(W_S)}{P(D_M)}
$$

右辺の確率ですが、
まず、$P(W_S)$ : $S\,$が当たりである事前確率は、$P(W_F)\,$と等しく$\,1/N\,$です。

$$
P(W_S) = P(W_F)
$$

次に、$P(D_M \mid W_S)$ : $S\,$が当たりのときに、特定の$\,M\,$個のハズレが開示される条件付き確率は、当たりである$\,S\,$と、挑戦者が最初に選んだ$\,F\,$以外の$\,N -2\,$個の中から、特定の$\,M\,$個の組が選ばれる確率なので、

$$
P(D_M \mid W_S) = \frac{1}{_{N-2}{C}{_M}}
$$

となります。ここで再び

$$\begin{align*}
_{N-1}{C}{_M} &= \frac{N-1}{N-M-1} \cdot _{N-2}{C}{_M}  
\end{align*}$$

の関係を使うことによって、$S\,$が当たりである確率は

$$\begin{align*}
P(W\_S \mid D\_M) &= \frac{P(D\_M \mid W\_S)\,P(W\_S)}{P(D\_M)} \newline \newline
                &= \frac{\frac{1}{\_{N-2}{C}{\_M}} \cdot P(W\_F)}{P(D\_M)} \newline \newline 
                &=  \frac{N-1}{N-M-1} \cdot \frac{\frac{1}{\_{N-1}{C}{\_M}} \cdot P(W_F)}{P(D\_M)} \newline \newline 
               &= \frac{N-1}{N-M-1} \cdot \frac{P(D_M \mid W_F)\,P(W_F)}{P(D_M)} \newline \newline
               &= \frac{N-1}{N-M-1} \cdot P(W_F \mid D_M) \newline \newline
               &= \frac{N-1}{N-M-1} \cdot \frac{1}{N}
\end{align*}$$

となりますが、$M\,\ge1\,$なので、$(N-1)/(N-M-1) > 1\,$。
よって、

$$
P(W_S \mid D_M) > P(W_F \mid D_M)
$$

つまり、挑戦者は選択肢を変更した方がよいということが分かります。
ちなみに、もともとのモンティ・ホール問題である$\,N=3$、$M = 1\,$を代入すると、

$$\begin{align*}
P(W_F \mid D_M) &= \frac{1}{3} \newline \newline
P(W_S \mid D_M) &= \frac{2}{3}
\end{align*}$$

というよく知られた結果が得られます。

####シミュレーションで確認
さて、以上の計算結果が正しいのかをシミュレーションで確認してみたいと思います。
python で以下のように実装しました。

```python
import random
import matplotlib.pyplot as plt

class GeneralizedMontyHallProblem:
    def __init__(self, num_choices, num_disclosed_choices):
        assert num_choices > 2, 'The number of choices must be greater than 2.'
        assert num_disclosed_choices < num_choices - 1, f'The number of choices to be disclosed must be less than the number of choices - 1, less than {num_choices - 1} in this case.'
        self.N = num_choices
        self.M = num_disclosed_choices
        self.p_f_theory = 1 / self.N
        self.p_s_theory = self.p_f_theory * (self.N - 1) / (self.N - self.M - 1)
        self.list_pf, self.list_ps = [], []

    def plot_result(self):
        num_h = len(self.list_pf)
        h_axis = [i + 1 for i in range(num_h)]
        pf_theory_line = [self.p_f_theory for i in range(num_h)]
        ps_theory_line = [self.p_s_theory for i in range(num_h)]
        plt.plot(h_axis, self.list_pf, 'b', h_axis, self.list_ps, 'r')
        plt.plot(h_axis, pf_theory_line, 'b--', h_axis, ps_theory_line, 'r--')
        plt.legend(['P_F simulation', 'P_S simulation', 'P_F theory', 'P_S theory'])
        plt.xlabel('trials')
        plt.ylabel('probability')
        plt.show()        

    def simulate(self, num_trials):
        assert num_trials > 0, 'The number of trials must be greater than 0.'
        win_f, win_s = 0, 0
        p_f, p_s = 0, 0
        for trial in range(num_trials):
            # the winning choice is selected randomly
            w = random.randint(0, self.N - 1)
            # the player selects a choice randomly
            f = random.randint(0, self.N - 1)
            # make a list of choices other than w and f 
            remainig_choices = [c for c in range(self.N) if c != w]
            remainig_choices = [c for c in remainig_choices if c != f]
            # M loosing choices are selected randomly and disclosed to the player
            disclosed_choices = random.sample(remainig_choices, self.M)
            # make a list of choices for the second choice
            second_choices = [c for c in range(self.N) if c != f]
            for o in disclosed_choices:
                second_choices = [c for c in second_choices if c != o]
            # the player selects the second choice randomly
            s = random.sample(second_choices, 1)[0]
            # check if f or s is a winning one
            if f == w:
                win_f += 1
            if s == w:
                win_s += 1
            # update each probability
            p_f = win_f / (trial + 1)
            p_s = win_s / (trial + 1)
            self.list_pf.append(p_f)
            self.list_ps.append(p_s)
        # show results
        print(f'P_F\ntheory : {self.p_f_theory}\nsimulation : {self.list_pf[-1]}\n')
        print(f'P_S\ntheory : {self.p_s_theory}\nsimulation : {self.list_ps[-1]}')
        self.plot_result()
```

試しに、$N = 50$、$M = 17\,$という設定にしてみましょう。

```python
N = 50
M = 17
mh = GeneralizedMontyHallProblem(N, M)
```

先ほど求めた結果によれば、挑戦者が最初に選んだ選択肢が当たりである確率$\,P_F\,$と、選択肢を変更したときに当たりである確率$\,P_S\,$は、それぞれ、
$$\begin{align*}
P_F &= \frac{1}{50} \newline
    &= 0.02 \newline \newline
P_S &= \frac{49}{1600} \newline
    &= 0.030625
\end{align*}$$

に近い値になるはずです。
というわけで、試行回数10000回でシミュレーションしてみます。

```python
trials = 10000
mh.simulate(trials)
```

結果の値と、試行ごとのプロットは以下のようになりました。

```python
P_F
theory : 0.02
simulation : 0.021

P_S
theory : 0.030625
simulation : 0.0309
```

![Figure_1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1911288/fef949a1-64be-80ab-377c-4753a53bf0fa.png)


確かに正しいようです。
ご覧いただきありがとうございました。