単回帰分析~コンビニエンスストア

概要

Pythonを使った統計計算と図示の練習のため、コンビニエンスストアで単回帰分析をやってみた。

コンビニの店舗数は「商業動態統計年報」の2016年データを使い、説明変数として2015年の国勢調査人口、2018年の国土地理院による都道府県面積、2017年道路統計年報の道路実延長データを使った。

計算コードは以下の通り。

人口との関係

コンビニ店舗数と人口の散布図と回帰式を以下に示す。

相関係数が極めて高いのは当然で、やはりコンビニの出店計画には人口ファクターが強く作用していることがわかる。

定数項bが店舗数のオーダーに比べてほぼゼロというのも興味深い。

係数aの値からは5万人弱で1店舗ということになるが、人口10万くらいの都市で2店舗しかないことになり、ちょっと少ないような気がする。もしかすると、一定規模以上の市町村単位や都市単位くらいで層化して出店計画を立てているのかもしれない。

面積との関係

店舗数と面積の関係

次に店舗数と都道府県面積の関係を見てみた。

この結果はかなりはずれで、データを見ても人口が少なく面積が群を抜いて大きい北海道の影響を大きく受けている。

ここで、面積が極端に大きい北海道(83423.83㎢~1位、2906店舗~5位)と面積が小さいが集積度が極端に高い東京都(2193.96㎢~45位、7003店舗~1位)の2つを除いて計算してみる。

これはさらにおかしな結果で、面積が小さいほど店舗数が多いことになる。

考えてみれば、集客数を期待するなら人口が集積している地域が有利だから、人口密度に比例する可能性を考えた方がいいのかもしれない。もし面積が小さい県の方が集積度が高いと想定すると、面積だけを取り出したときに逆の関係になるとも考えられるが、相関係数や決定係数が小さすぎるので考察は難しい。

人口密度

以下は人口密度との関係。

今度はかなりきれいに相関の高さが出ている。

直接的な計算式に入れているかどうかわからないが、GISなどで出店計画を立てるとしたら、人口密度の高いエリアを選んでいくだろうことが想定される。

ただ、店舗数は人口などといった売り上げに直結するデータから導かれるのが普通で、人口密度が高くても人口が少なければ出店インセンティブにはならない。

人口と人口密度の関係

試しに人口を説明変数、人口密度を被説明変数として両者の関係を見てみると、驚くことに「人口が多いほど人口密度が高くなる(あるいはその逆)」という関係になる。

ここから先は人口論や地域論になりそうなので置いておくが、少なくとも日本においては、「狭いところほど人が集まっている傾向がある」ということになりそうである。

もちろんこれは他の国でも一般に当てはまることかもしれないが、朝のラッシュ時に特定の車両に無理やり乗り込んでいる割に離れた車両がすいているとか、1本電車を遅らせたらガラガラだったとか、そのあたりの行動パターンを見ていると、何となく日本に特有のような気がする。

道路延長との関係

最後に道路延長との関係を見てみる。

東京のように稠密な都市は例外とすると、概ね関係はありそうである。ただし相関係数、決定係数は高くない。

コンビニ店舗が道路の利便性に依っていることは推測できるが、やはり人口という売り上げ直結のデータに比べると関係は弱い。

高速道路の延長についても見てみたが、こちらはほとんど関係は見られなかった。

ただ、高速道路の延伸に伴ってコンビニエンスストアの店舗数が伸びているようであり、マクロな延長というよりも物流上のインパクトが大きいことは予想される。

 

相関係数

相関係数の定義

相関係数は、多数のデータの組がどの程度線形に近い性質を持つかを表す値で、以下で定義される。

(1)    \begin{eqnarray*} r &=& \frac{{\rm Cov}(X, Y)}{\sqrt{V(X) \cdot V(Y)}} \\ &=& \frac{E(X - \overline{X})(Y - \overline{Y})} {\sqrt{ ( E\left[ (X - \overline{X})^2 \right] E\left[ (Y - \overline{Y})^2 \right] }} \\ &=& \frac{\displaystyle \sum_{i=1}^{n}(x_i - \overline{x})(y_i - \overline{y})} {\left[ \displaystyle \sum_{i=1}^{n} (x_i - \overline{x})^2 \displaystyle \sum_{i=1}^{n} (y_i - \overline{y})^2 \right]^{1/2}} \end{eqnarray*}

相関係数の線形変換

変数を線形変換した場合の相関係数は分散・共分散の性質から、以下のようになり、元のままか符号が反転する。

(2)    \begin{eqnarray*} r' &=& \frac{{\rm Cov}(aX+b, cY+d)}{\sqrt{V(aX+b) \cdot V(cY+d)}} \\ &=& \frac{ac{\rm Cov}(X, Y)}{\sqrt{a^2 V(X) \cdot c^2 V(Y)}} \\ &=& \frac{ac}{|ac|} \cdot r \end{eqnarray*}

完全線形関係の場合の相関係数

XとYが完全な線形関係にある場合、相関係数は1または-1になる。このとき、傾きの大きさや、平行移動は相関係数に影響しない。

(3)    \begin{eqnarray*} r &=& \frac{ {\rm Cov}(X, aX+b) } { \left[ V(X) \cdot V(aX+b) \right]^{1/2} } \\ &=& \frac{ a {\rm Cov}(X, X) } { \left[ V(X) \cdot a^2 V(X) \right]^{1/2} } \\ &=& \frac{a}{\sqrt{a^2}} \cdot \frac{ {\rm Cov}(X, X) } { \left[ V(X) \cdot V(X) \right]^{1/2} } \\ &=& \frac{a}{|a|} \cdot \frac{ V(X) } { V(X) } \\ &=& 1 \; {\rm or} \; -1 \end{eqnarray*}

いろいろな分布の相関係数

完全な線形関係

以下のコードで確認。

相関係数は傾きや平行移動に対して影響を受けず、増加関数なら1、減少関数なら-1になることがわかる。

線形性が強い関係

線形関数の値に対して、乱数でばらつきを与えた場合の相関係数の違いを示す。ばらつきが大きい方が相関係数は小さくなる。

負の線形性が強い場合は、相関係数がマイナスになる。

放物線(強い関係があるのに相関が低いケース)

以下のような放物線では、XとYにきちんとした数学的関係があるのに、相関係数がゼロに近くなる。

相関係数はXとYが単調増加/単調減少の度合いが強いほど、また線形関係に近いほど1に近くなるが、それ以外の関係が強い場合にはそれを補足できない場合がある。

反比例(負の線形性に見えてしまう場合)

以下は反比例関数の場合。

関数の形状や範囲に寄るが、この場合は相関係数の絶対値が0.8以上と1に近く、これだけ見ると負の線形性が強そうに見える。

対数関数(正の線形性に見えてしまう場合)

対数関数の場合。この場合は0.9以上とかなり強い線形性を示唆している。

相関係数に関する注意

本来の関係との乖離

先にみたように、線形関係ではないが数学的な関係を持つ場合に、相関係数からは全く関係がない、元の関係とは異なり線形関係を持つ、といった解釈になることがある。

相関係数が高い場合に、線形回帰式などで物事を予測する際には注意が必要。

変数間に解析的な関係が見いだせるならそれを重視すべきであり、よしんばそれがわからないにしても、定義域の範囲で「ある程度は当たる」程度に考えておくべきか。

因果関係

堂々と間違えられるケースが、「科学的な」記事やマスメディアなどでよくみられる。

相関係数は「変数間の単調な増加/減少傾向が強いかどうか」だけを示すもので、必ずしも因果関係を示唆しない。

気温が高いとビールはよく売れるがおでんは売れない。その二つに負の相関があるからといって、「ビール好きはおでんが嫌い」と言えないが、形を変えてこのような解釈がなされる恐れがある。

もともとのメカニズムで因果関係が示唆されていて、その上で相関係数の大きさを論じるなら意義もあるが、その場合でも、事象に対する寄与度などをよく考えておかないと「それだけが原因」と考えるような間違いを犯すことになる。

 

回帰分析~単回帰

概要

(X,\: Y) = (x_i,\: y_i) \; (i=1 \ldots n)のデータからなる複数の標本に対して、線形式で最も説明性の高いものを求める。

そのために、各データのx_iに線形式を適用して得られた値\hat{y_i} = a x_i + by_iの差(残差)の平方和が最小となるような係数a,\: bを求める(最小二乗法)。

定式化

(1)    \begin{equation*} \sum_{i=1}^n (\hat{y_i} - y_i)^2 \; \rightarrow \; \min \quad {\rm where} \; \hat{y_i} = a x_i + b \end{equation*}

ここで残差を最小化するa,\: bを求めるために、それぞれで偏微分する。

(2)    \begin{equation*} \left\{ \begin{array}{l} \displaystyle \frac{\partial}{\partial a} \sum_{i=1}^n (a x_i + b - y_i)^2 = 0 \\ \displaystyle \frac{\partial}{\partial b} \sum_{i=1}^n (a x_i + b - y_i)^2 = 0 \end{array} \right. \quad \Rightarrow \quad \left\{ \begin{array}{l} \displaystyle \sum_{i=1}^n 2(a x_i + b - y_i)x_i = 0 \\ \displaystyle \sum_{i=1}^n 2(a x_i + b - y_i) = 0 \end{array} \right. \end{equation*}

展開して

(3)    \begin{equation*} \left\{ \begin{array}{l} \displaystyle a \sum_{i=1}^n x_i^2 + b \sum_{i=1}^n x_i - \sum_{i=1}^n x_i y_i = 0 \\ \displaystyle a \sum_{i=1}^n x_i + bn - \sum_{i=1}^n y_i = 0 \end{array} \right. \end{equation*}

各和の記号を新たに定義し、行列表示で表すと、

(4)    \begin{equation*} \left[ \begin{array}{cc} S_{xx} & S_x \\ S_x & n \end{array} \right] \left[ \begin{array}{c} a \\ b \end{array} \right] = \left[ \begin{array}{c} S_{xy} \\ S_y \end{array} \right] \end{equation*}

回帰式の係数

式(4)を解くと以下を得る。

(5)    \begin{equation*} \left[ \begin{array}{c} a \\ b \end{array} \right] = \left[ \begin{array}{c} \dfrac{n S_{xy} - S_x S_y}{n S_{xx} - {S_x}^2} \\ \dfrac{S_{xx} S_y - S_x S_{xy}}{n S_{xx} - {S_x}^2} \end{array} \right] \end{equation*}

あるいは式(5)については以下のようにも表せる。ただし以下の式で、\sigma_{XY} = {\rm Cov}(X, Y){\sigma_X}^2 = {\rm Var}(X)

(6)    \begin{equation*} a = \frac{\displaystyle \sum_{i=1}^n (x_i - \overline{X})(y_i - \overline{Y})}{\displaystyle \sum_{i=1}^n (x_i - \overline{X})^2} =\frac{\sigma_{XY}}{{\sigma _X}^2} \end{equation*}

(7)    \begin{equation*} b = \overline{Y} - a \overline{X} \end{equation*}

決定係数

全変動・回帰変動・残差変動

x_iの回帰式による推定値を\hat{y}_i = a x_i + bで表す。

ここで、y_i,\: \hat{y},\: \overline{Y}によって、以下の変動が定義される。

y_i - \overline{Y}
全変動。平均からのデータのばらつき。
\hat{y}_i - \overline{Y}
回帰変動。平均からのばらつきのうち、回帰式によって説明される部分。
y_i - \hat{y}_i
残差変動。平均からのばらつきのうち、回帰式によって説明されない部分。

当然、[全変動]=[回帰変動]+[残差変動]となることはすぐに確認できる。

また、これら3つの変動の2乗和をTSS(Total Sum of Squares)、ESS(Explained Sum of Squares、RSS(Residual Sum of Squares)として、以下のように定義する。

(8)    \begin{eqnarray*} TSS &=& \sum_{i=1}^n (y_i - \overline{Y})^2 \\ ESS &=& \sum_{i=1}^n (\hat{y}_i - \overline{Y})^2 \\ RSS &=& \sum_{i=1}^n (y_i - \hat{y}_i)^2 \end{eqnarray*}

このとき、それぞれについて以下のように表せる。

(9)    \begin{eqnarray*} TSS &=& n {\sigma _Y}^2 \\ ESS &=& n \frac{{\sigma_{XY}}^2}{{\sigma_x}^2} \\ RSS &=& n {\sigma_Y}^2 - n \frac{{\sigma_{XY}}^2}{{\sigma_X}^2} \end{eqnarray*}

これらから、全変動・回帰変動・残差変動の2乗和についても、以下の関係が成り立つことがわかる。

(10)    \begin{equation*} TSS = ESS + RSS \end{equation*}

補足~ESS

式(9) のESSの確認。以下のように、線形変換された分散から導ける。

(11)    \begin{eqnarray*} ESS &=& \sum_{i=1}^n \left( a x_i +b - (a \overline{X} + b) \right)^2 \\ &=& n {\rm Var}(a X + b) \\ &=& na^2{\rm Var}(X) \\ &=& n a^2 {\sigma_x}^2 \\ &=& n \left( \frac{\sigma_{XY}}{{\sigma_x}^2} \right)^2 {\sigma_x}^2 \\ &=& n \frac{{\sigma_{XY}}^2}{{\sigma_x}^2} \end{eqnarray*}

あるいは、以下のように地道に計算しても確認できる。

(12)    \begin{eqnarray*} \sum_{i=1}^n (\hat{y}_i - \overline{Y}) &=& \sum_{i=1}^n ({\hat{y}_i}^2 - 2 \hat{y}_i \overlne{Y} + \overline{Y}^2) \\ &=& \sum_{i=1}^n (\hat{y}_i^2 - 2(a x_i + b) \overline{Y} +\overline{Y}^2) \\ &=& \sum_{i=1}^n (\hat{y}_i^2 - 2a x_i \overline{Y} - 2b\overline{Y} + \overline{Y}^2) \\ &=& \sum_{i=1}^n \hat{y}_i^2 - 2an \overline{X} \; \overline{Y} -2nb \overline{Y} + n \overline{Y}\end{eqnarray*} \\

ここで

(13)    \begin{eqnarray*} \sum_{i=1}^n \hat{y}_i^2 &=& \sum_{i=1}^n (ax_i + b)^2 \\ &=&\sum_{i=1}^n (a^2 x_i^2 + 2ab x_i + b^2) \\ &=& a^2 \sum_{i=1}^n x_i^2 + 2abn \overline{X} + nb^2 \end{eqnarray*}

これを代入して、

(14)    \begin{eqnarray*} ESS &=& a^2 \sum_{i=1}^n x_i^2 + 2nab \overline{X} + nb^2 - 2na \overline{X} \; \overline{Y} - 2nb \overline{Y} + n \onerline{Y} \end{eqnarray*}

さらに\sum_{i=1}^n x_i^2 = \n \sigma_X^2 + n \overline{X}^2を考慮して、

(15)    \begin{eqnarray*} ESS &=& n a^2 {\sigma_X}^2 + n a^2 \overline{X}^2 + 2nab \overline{X} + nb^2 -2na\overline{X}\;\overline{Y} - 2nb\overline{Y} + n\overline{Y}^2 \\ &=& n a^2 {\sigma_X}^2 + n(a \overline{X} - \overline{Y})^2 + 2nab\overline{X} + nb^2 - 2nb\overline{Y} \\ &=& n a^2 {\sigma_X}^2 + nb^2 + 2nab\overline{X} + nb^2 - 2nb\overline{Y} \\ &=& n a^2 {\sigma_X}^2 + 2nb(b + a\overline{X} - \overline{Y}) \\ &=& n a^2 {\sigma_X}^2 \\ &=& n \left( \frac{\sigma_{XY}}{{\sigma_X}^2} \right)^2 {\sigma_X}^2 \\ &=& n \frac{{\sigma_{XY}}^2}{{\sigma_X}^2} \end{eqnarray*}

補足~RSS

式(9) のRSSの確認。

(16)    \begin{eqnarray*} RSS &=& \sum_{i=1}^n (y_i - \hat{y}_i)^2 = \sum_{i=1}^n (y_i - a x_i - b)^2 \\ &=& \sum_{i=1}^n (y_i - a x_i - \overline{Y} + a \overline{X})^2 \\ &=& \sum_{i=1}^n \left(y _i - \overline{Y} - a (x_i - \overline{X}) \right) ^2 \\ &=& \sum_{i=1}^n \left(y _i - \overline{Y} - \frac{\sigma_{XY}}{{\sigma_X}^2} (x_i - \overline{X}) \right) ^2 \\ &=& n {\sigma_Y}^2 - 2n \frac{\sigma_{XY}}{{\sigma_X}^2} \sigma_{XY} + n \frac{{\sigma_{XY}}^2}{{\sigma_X}^4} {\sigma_X}^2 \\ &=& n {\sigma_Y}^2 - 2n \frac{{\sigma_{XY}}^2}{{\sigma_X}^2} + n \frac{{\sigma_{XY}}^2}{{\sigma_X}^2} \\ &=& n {\sigma_Y}^2 - n \frac{{\sigma_{XY}}^2}{{\sigma_X}^2} \end{eqnarray*}

決定係数

決定係数は、全変動のうち回帰変動でどれだけ説明できるか(残差変動が少ないか)を表したもので、通常R^2で表される。

(17)    \begin{align*} R^2 = \frac{ESS}{TSS} = \frac{TSS - RSS}{TSS} = 1 - \frac{RSS}{TSS} \end{align*}

この値は、残差変動が小さいほど1に近づき、RSS = TSSすなわち全変動が回帰変動で全く説明できないときに0となる。

 

ロシアンルーレットの確率

基本

ロシアンルーレットは、賭の順番問題や条件付き確率の問題としてよく見かける。まず、一番基本的なルールで考える。

1人の相手とロシアンルーレットで命をかけることになった。6発入りの拳銃に1発だけ弾を込め、最初にシリンダーを高速で回転させる。シリンダーが止まったところで、相手と交互に引き金を引いていくとして、まず自分から始めるのと相手に譲るのと、どちらが生き残る可能性が高いか。

 

n回目に実弾に当たる確率をP_n \; (n = 1~6)とすると、

1発目で当たる確率は、6つの弾倉のどこかに実弾が入っているので、

(1)    \begin{equation*} P_1 = \frac{1}{6} \end{equation*}

2発目で当たる確率は、1発目で実弾に当たらず、2発目は残り5つの弾倉のどこかに実弾が入っているので、

(2)    \begin{equation*} P_2 = \frac{5}{6} \times \frac{1}{5} = \frac{1}{6} \end{equation*}

以下、6発目まで同じように考えて、

(3)    \begin{eqnarray*} P_3 &=& \frac{5}{6} \times \frac{4}{5} \times \frac{1}{4} = \frac{1}{6} \\ P_4 &=& \frac{5}{6} \times \frac{4}{5} \times \frac{3}{4} \times \frac{1}{3} = \frac{1}{6} \\ P_5 &=& \frac{5}{6} \times \frac{4}{5} \times \frac{3}{4} \times \frac{2}{3} \times \frac{1}{2} = \frac{1}{6} \\ P_6 &=& \frac{5}{6} \times \frac{4}{5} \times \frac{3}{4} \times \frac{2}{3} \times \frac{1}{2} \times \frac{1}{1}= \frac{1}{6} \end{eqnarray*}

すなわち、通常のロシアンルーレットでは、1~6回目のどこで実弾に当たる確率も等しく1/6で、順番の後先に有利性はない。

毎回シリンダーを回転させる場合

1発だけ弾を込めた拳銃で、2人でロシアンルーレットで命をかける。ただし、それぞれが自分の番になる毎にシリンダーを回転させて、ランダムに止まったところで引き金を引く。自分が最初に打つのと、相手から打たせるのと、どちらを選ぶか。

 

毎回シリンダーを回転させるので、そこまで発砲されていなければ、実弾である確率は毎回1/6。n回目に実弾に当たる確率をp_nとすると、n-1回は当たらないことを考慮して、

(4)    \begin{equation*} p_n = \left( \frac{5}{6} \right) ^{n-1} \times \frac{1}{6} \end{equation*}

まず、この確率はnが大きくなるほどゼロに近づく。たとえば1回めで当たる確率は約16.7%、10回目で当たる確率は約3.2%、20回目で当たる確率は0.5%、40回目だと約0.01%程度となる。

次に、n回目までに当たる確率をP_nとすると、これは各回の確率を累積していくことで求められるので、等比数列の和より、

(5)    \begin{equation*} P_n = \sum_{i=1}^{n} \left( \frac{5}{6} \right) ^{n-1} \times \frac{1}{6} = \frac{1}{6} \frac{1-\left( \dfrac{5}{6} \right)^n}{1 - \dfrac{5}{6}} = 1-\left( \dfrac{5}{6} \right)^n \end{equation*}

この累積確率はn \rightarrow \inftyのとき1に近づくが、10回目で約83.8%、20回目で約97.4%、30回目で約99.6%、40回目で約99.9%となり、100回も打ってその時点で生き延びている確率はほとんどゼロに近い。

次に交互に打っていく場合を考える。

先攻で打つ方になった場合、1回目、3回目・・・と奇数階に打つことになるので、その累積確率をP_{1 \cdot n}とすると、

(6)    \begin{eqnarray*} P_{1\cdot 2n-1} &=& \sum_{k=1}^{n} \left( \frac{5}{6} \right)^{2(k-1)} \times \frac{1}{6} = \sum_{k=1}^{n} \frac{1}{6} \left( \frac{25}{36} \right)^{k-1} = \frac{1}{6} \cdot \frac{1 - \left( \dfrac{25}{36}\right) ^n }{1 - \dfrac{25}{36}} \\ &=& \frac{6}{11} \left( 1 - \left( \dfrac{25}{36}\right) ^n \right) \end{eqnarray*}

後攻で打つ方になった場合は、2回目、4回目・・・と偶数階に打つことになるので、その累積確率をP_{2 \cdot n}とすると、

(7)    \begin{eqnarray*} P_{2\cdot 2n} &=& \sum_{k=1}^{n} \left( \frac{5}{6} \right)^{2k-1} \times \frac{1}{6} = \sum_{k=1}^{n} \frac{5}{36} \left( \frac{25}{36} \right)^{k-1} = \frac{5}{36} \cdot \frac{1 - \left( \dfrac{25}{36}\right) ^n }{1 - \dfrac{25}{36}} \\ &=& \frac{5}{11} \left( 1 - \left( \dfrac{25}{36}\right) ^n \right) \end{eqnarray*}

同じnに対して常にP_{1 \cdot n} / P_{2 \cdot n} = 1.2となり、先に打つのを選んだほうが実弾に当たる確率は高く、この比率はn \rightarrow \inftyとしたときの結果と同じ。

したがって、「お先にどうぞ」と相手に最初を譲ったほうが、少しでも生き残る確率が高くなる。また、3人以上で順番に打っていく場合は、できるだけ後の順番で打つ方が実弾に当たる確率が低くなる。

条件付きの問題

2発の弾を並べて込める

ロシアンルーレットと少し違うが、条件付きの問題。

筋悪の金貸しから借りたかなりの借金を返せなくなり、命を取られようとしているが、「チャンスをやる」と言われた。彼は6発入りの拳銃の弾倉に2 発の弾を並べて込め、シリンダーを回す。シリンダーがランダムな位置で止まったところで、彼は言う。「今から空に向かって私が1回引き金を引く。その後にお前のこめかみに銃口を当てて打ってもらうが、そのまま打ってもいいし、もう一度シリンダーを回して止まった位置で打ってもいい」

 

1発目が空砲の場合

金貸しが銃を空に向けて引き金を引くと空砲だった。そのまま打つか、シリンダーを回転させるか、どちらを選択すべきか。

 

まず、シリンダーをもう一度回転させてから打つ場合は、6弾倉中2発残っているので、当たる確率は1/3。逆に助かる率は2/3となる。

一方シリンダーを回転させない場合は、手渡された拳銃の弾の状態を、円形に並んだ弾倉を一列で表し、下表のように整理する。ただし、左欄の●は実弾がある場所、○は実弾が入っていない場所で、1番目は金貸しが打ったときの位置、2番目は手渡された時の位置。右欄は2発目が実弾の場合は●、空砲の場合は○とする。

●●○○○○ 対象外
○●●○○○
○○●●○○
○○○●●○
○○○○●●
●○○○○● 対象外

金貸しが打ったときは空砲だったので、1番目が実弾の事象は対象外。このとき、実弾に当たらない確率は3/4となる。

すなわち、もう一度シリンダーを回転させるよりも、そのまま続けて打った方が助かる可能性が高い。

1発目が実弾の場合

もし、金貸しが最初に空へ向けて打ったときに実弾だった場合はどうか?

 

この場合に続けて打つと、1発目が実弾である2ケースに対して、空砲で助かる率は1/2。

シリンダーを回転させてから打つと、実弾は残り1発なので、当たらない確率は5/6で、回転させて打った方が生存確率は高くなる。

条件に応じた選択と確率

以上の結果を整理すると下表の通りとなり、1発目が空砲の場合と実弾の場合で、その後の選択による生存確率の高い方が逆転している。

空砲も実弾もその位置が連続しているため、1発目が空砲の場合は続けて空砲の確率が高く、1発目が実弾の場合は次も実弾の確率が高そうだが、1発目が実弾の場合は弾が1つ消費されるので、その影響で結果が逆転している。

1発目 選択 生存確率
空砲 シリンダー回転 1/3 (0.33)
空砲 続けて打つ 3/4 (0.75)
実弾 シリンダー回転 5/6 (0.83)
実弾 続けて打つ 1/2 (0.50)

実弾2発をランダムに込める

上記の金貸しの問題では実弾2発を並べて込めたが、これをランダムな位置に2発込めた場合はどうなるか。

この場合、2発の実弾の装填状況は下表のようになる。

1発目が→ 空砲 実弾
●●○○○○ 対象外  ●
●○●○○○ 対象外
○●●○○○ 対象外
●○○●○○ 対象外
○●○●○○ 対象外
○○●●○○ 対象外
●○○○●○ 対象外
○●○○●○ 対象外
○○●○●○ 対象外
○○○●●○ 対象外
●○○○○● 対象外
○●○○○● 対象外
○○●○○● 対象外
○○○●○● 対象外
○○○○●● 対象外

1発目が空砲だった場合

シリンダーを回転させた場合の生存確率は、4/6 = 2/3 (0.667)。

続けて打つ場合は、1発目が実弾のケースを除いて2発目が空砲の場合なので、6/10 = 3/5 (0.6)。

1発目が実弾だった場合

シリンダーを回転させた場合の生存確率は、5/6 (0.833)。

続けて打つ場合は、1発目が空砲のケースを除いて2発目が空砲の場合なので、4/5 (0.8)。

1発目が空砲でも実弾でも、シリンダーを回転させた方が生き残る確率が高くなる。

条件付確率~男女問題

概要

条件付き確率の代表的な問題として、2人の子の性別問題をよく見かけるが、前提条件によって答えが異なることを明確にしている例は少ない。

また、この問題がベン図を用いた分析や、ベイズの定理の親和性についても考えてみる。

兄妹の場合

まず準備として、以下の問題を考える。

2人兄弟の上の子が男であった場合、下の子が女である確率はいくらか。

 

普通に考えれば、男か女どちらかなのだから、確率は1/2と考えるが、まずは生真面目にすべての事象を並べてみる。

conditional-probability-brother-cases

上の子が兄か姉かに関係なく下の子が女である確率は1/2と直感的に考えられるが、それが確認できた。

一応、以下の様なベン図を考えてみる。図中、色付けしたところが条件下での事象で、濃い色が事後確率として求めたい部分。

conditional-probability-brother-venn-diagram

この場合、上の子が男(事象EM)のもとでの下の子が女(事象EMかつYF)ということになり、EMの要素が(男, 男)、(男, 女)の2つに対してEMかつYFの要素は(男,女)だけなので、確率は1/2となる。

なお、ベン図の代わりに場合分けの表で整理してみると、以下の様になる。
conditional-probability-brother-table

これを敢えてベイズの定理で書くと以下の様になるが、まあ当たり前という感じ。

(1)    \begin{equation*} P(YF|EM) = \frac{P(EM \cap YF)}{P(EM)} = \frac{1}{2} \end{equation*}

これをベイズの定理の別の表現で書くと以下の様になる。

(2)    \begin{equation*} P(YF|EM) = \frac{P(EM|YF)P(YF)}{P(EM)} = \frac{\dfrac{1}{2} \cdot \dfrac{1}{2}}{\dfrac{1}{2}}=\frac{1}{2} \end{equation*}

ただし P(EM|YF)=1/2としているが、この項の考え方は結局求めたい確率と同じ道筋で考えなければならない。

いずれにしてもこの問題は、第一子が男であろうが女であろうが、第二子が男か女かと問われれば、確率1/2でどちらか、という自然な結果となり、ベン図や場合分け表、ベイズの定理を使ったとしても、あまり恩恵は得られない。

区別なしの場合

次に以下の様な問題を考える。

2人の子のうち、1人が男の場合、もう1人が女である確率はいくらか。

 

これも普通に考えれば、男か女どちらかなのだから、確率は1/2と考えてしまいそうになるが、まじめに事象を考えてみる。

conditional-probability-2children-cases

この場合は1人が男とわかっているので、(姉, 妹)の事象は存在しなくなり、もう1人が女というだけで姉か妹かの区別はないことから、確率は2/3となる。

これをベン図で描いてみると以下のようになるが、男と女の組み合わせが2通りあることを確認する点では、上のような整理のほうがわかりやすい。

conditional-probability-2children-venn-diagram-1

場合分けの表は以下の通りになる。

conditional-probability-2children-table-1

これをベイズの定理で表現すると以下のようになる。

(3)    \begin{eqnarray*} P(F|M) &=& \frac{P(F \cap M)}{P(M)} = \frac{P(FM) + P(MF)}{P(MM) + P(MF) + P(FM)} \\ &=& \frac{\dfrac{1}{4} + \dfrac{1}{4}}{\dfrac{1}{4} + \dfrac{1}{4} + \dfrac{1}{4}} = \frac{2}{3} \end{eqnarray*}

条件付き確率から求める方法は、式(2)と同じ理由でメリットがない。

順番がある場合

次に以下の様な問題を考える。

2人の子がいる家庭を訪ねた。最初に挨拶に出てきたのが男の子だった場合、もう1人が女の子である確率はいくらか。

 

この場合は上の子と下の子を区別していないので、一つ前の問題と同じ答えになりそうだが、事象を列挙してみる。

conditional-probability-2children-by-order-cases

この問題でも、2人とも女のケースが除かれるが、兄弟構成に対して挨拶の順番がそれぞれ2通りずつあり、2人とも男の場合は兄、弟のいずれも最初に挨拶しうるので前提条件の事象に2つ含まれる。結果として、2人目が女の子である確率は1/2となる。

このように、2人兄弟の性別に関する問題でも、条件の付し方によって答えが異なってくるが、いずれの問題もベイズの定理を駆使するというものではなさそうである。

3人兄弟の場合

問題を以下の様に拡張する。

3人の子のうち2人が男の子の場合、残る1人が女の子である確率はいくらか。

 

これについても以下のように事象を列挙し、確率3/4を得る。

conditional-probability-3children-cases

 

ベイズの定理で表すと以下のようになるが、この場合もまじめに事象のパターンを列挙した考えるほうが素直。

(4)    \begin{eqnarray*} P(F|MM) &=& \frac{P(F \cap MM)}{P(MM)} \\ &=& \frac{P(MMF) + P(MFM) + P(FMM)}{P(MMM) + P(MMF) + P(MFM) + P(FMM)} \\ &=& \frac{\dfrac{1}{8} + \dfrac{1}{8} + \dfrac{1}{8}}{\dfrac{1}{8} + \dfrac{1}{8} + \dfrac{1}{8} + \dfrac{1}{8}} = \frac{3}{4} \end{eqnarray*}

 

ベイズの定理~難病検査

概要

ベイズの定理で典型例として説明される問題。ある病気の検査で陽性が出た場合に罹患している確率を求める。類似の問題として、鉱脈を探るのに試掘結果から鉱脈が存在する確率を求める問題がある。

ロジックとしては、「誤判定の可能性がある調査・検査の結果から、目的とする事実の確率を求める」という構図で共通している。

難病問題

病気の検査の問題を例示する。

ある難病の検査結果が陽性だったとき、その病気に罹患している確率を求めよ。ただし以下の情報が与えられている。

  • その難病は1000人に1人が罹患している
  • その難病にかかっている場合、90%の確率で検査結果は陽性になる
  • その難病にかかっていなくても、20%の確率で検査結果は陽性になる

 

ここで以下の記号を定義する。

  • 全体としての難病の罹患率→P(D)=0.001
  • 難病にかかっている場合に陽性になる確率→P(P|D)=0.9
  • 難病にかかっていない場合に陽性になる確率→P(P| \overline{D})=0.2

このとき、ベイズの定理により以下の確率を得る。

(1)    \begin{eqnarray*} P(D|P) &=& \frac{P(P|D) P(D)}{P(P|D) P(D) + P(P|\overline{D}) P(\overline{D})} \\ &=& \frac{0.9 \times 0.001}{0.8 \times 0.001 + 0.2 \times 0.999} \\ &\approx& 0.004484 \end{eqnarray*}

検査結果が陽性であっても、この難病にかかっている確率は0.5%未満ということになる。検査しなければ1000人に1人の罹患率なので0.1%なので、検査によってその確からしさが約4.5倍となるが、事前確率が小さい場合、検査の精度を上げても劇的な判定効果は望めない。

なお、この約4.5倍という値はP(P|D) / P(P|\overline{D}) = 4.5に相当し、事前確率P(D)の値が小さいほど4.5に近づく(ベイズの定理の解釈参照)。

 

ベイズの定理の解釈

ベイズの定理の変形

ベイズの定理は、一般に以下のように表現される。

(1)    \begin{equation*} P(T_k |C) = \frac{P(C|T_k)P(T_k)}{\displaystyle \sum_{i=1}^n P(C|T_i)P(T_i)} \end{equation*}

またはT\overline{T}の2事象の場合は、

(2)    \begin{equation*} P(T|C) = \frac{P(C|T)P(T)}{P(C|T)P(T) + P(C|\overline{T})P(\overline{T})} \end{equation*}

ここで次のように確率の記号を定義する。また、イメージを掴みやすくするために、難病検査の問題の場合を対比させる。

P_t = P(T)
ターゲットの事象が真である確率(true)あるいは事前確率。たとえば日本全体での難病の罹患率など。
P_c = P(T|C)
判定が真の条件下で事象が真である確率(under condition)あるいは事後確率。たとえば検査結果が陽性の場合に難病に罹患している確率。
P_h = P(C|T)
事象が真の場合に判定が真となる確率(hit)。たとえば難病に罹患している場合に検査結果が陽性となる確率。
P_f = P(C| \overline{T})
事象が偽なのに判定が正となる確率(fail)。たとえば難病に罹患していないのに検査結果が陽性となる確率。

P(\overline{T}) = 1 - P(T)であることに留意して、これらの記号でベイズの定理を表すと以下の通り。

(3)    \begin{equation*} P_c = \frac{P_h P_t}{P_h P_t + P_F (1 - P_t)} = \frac{1}{1 + \dfrac{P_f}{P_h} \left( \dfrac{1}{P_t} - 1 \right) } \end{equation*}

いくつかのケース

誤判定を絶対しない場合

(8)においてP_f = 0とおくと、P_t = 1となる。

難病の問題で言えば、「罹患していないときには必ず陰性判定となる(陽性とならない)」という検査に相当し、このような検査の結果が陽性判定なら必ず罹患していることになる。この時、罹患時に陽性となる確率P_hがいくらかは関係ない。

ただしこの逆は必ずしも保証されておらず、陰性判定でも罹患している場合がある。

完璧な検査

難病の問題で、上記に加えて「罹患しているときは必ず陽性判定となる(P_h = 1」という条件を加える。

まず、(2)をCの余事象について書くと以下の様になる。

(4)    \begin{eqnarray*} P(T|\overline{C}) &=& \frac{P(\overline{C}|T)P(T)}{P(\overline{C}|T)P(T) + P(\overline{C}|\overline{T})P(\overline{T})} \\ &=& \frac{(1 - P_h) P_t}{(1 - P_h) P_t + (1 - P_f)(1 - P_t)} \\ \end{eqnarray*}

ここで、P(\overline{T} | \overline{C}) = 1 - P(T| \overline{C})を考える。これは、陰性反応の時に難病にかかっていない確率に相当する。

(5)    \begin{eqnarray*} P(\overline{T} | \overline{C}) = 1 - P(T| \overline{C}) &=& 1 - \frac{(1 - P_h) P_t}{(1 - P_h) P_t + (1 - P_f)(1 - P_t)} \\ &=& \frac{(1 - P_f)(1 - P_t)}{(1 - P_h) P_t + (1 - P_f)(1 - P_t)} \end{eqnarray*}

(5)においてP_f = 0と置くと、

(6)    \begin{eqnarray*} P(\overline{T} | \overline{C}) = 1 - P(T| \overline{C}) = \frac{(1 - P_t)}{(1 - P_h) P_t + (1 - P_t)} \end{eqnarray*}

ここでP_h  = 1の場合はP(\overline{T} | \overline{C}) = 1となる。

この結果を難病の問題で言えば、検査結果が陰性の場合は必ず罹患していない、ということになる。

以上をまとめると、検査の性質が「罹患していなければ必ず陰性、罹患していれば必ず陽性となる」という場合、検査結果が陽性なら必ず罹患しており、陰性なら罹患していないという、至極当然の結果となる。

意味のない検査

(8)においてP_h = P_fと置くと、P_c = P_tとなる。

難病の例で言えば、罹患している時に陽性となる確率と、罹患していないときに陽性となる誤判定率が同じ場合、検査の意味がない(一般的な罹患率でしか罹患の程度がわからない)ということになる。

正しい判定の率が高くても、誤判定率が同じ場合は意味をなさない。

意味のある検査

(8)をP_f / P_hについて解く。

(7)    \begin{equation*} \frac{P_f}{P_h} = \frac{\dfrac{1}{P_c} - 1}{\dfrac{1}{P_t}- 1} \end{equation*}

この式で、P_f < P_hのとき、P_c > P_tとなり、検査が意味を持つことがわかる。

稀な事象

事前確率P_tがとても小さい場合を考える。たとえば難病の罹患率が1万人や10万人に1人といった場合。

P_t \approx 0とすると(8)は以下のように変形できる。

(8)    \begin{eqnarray*} P_c &=& \frac{1}{1 + \dfrac{P_f}{P_h} \left( \dfrac{1}{P_t} - 1 \right) } = \frac{P_t}{P_t + \dfrac{P_F}{P_h} (1 - P_t)} \\ &=& \frac{\dfrac{P_h}{P_F} P_t}{1 + \left( \dfrac{P_h}{P_f} - 1 \right) P_t} \approx \frac{P_h}{P_F} P_t} \end{eqnarray*}

すなわち、事前確率が十分に小さい場合は、検査の正誤判定確率の倍率分だけ事後確率が大きくなる。

たとえば難病検査の問題の場合、P_h = 0.9 , \; P_f = 0.2として、P_tが小さくなるほど、事後確率は正誤判定確率の倍率0.9/0.2 = 4.5倍に近づく。

(9)    \begin{eqnarray*} P_t = 0.1 & \rightarrow & P_c = \frac{0.9 \times 0.1}{0.9 \times 0.1 + 0.2 \times 0.9} \approx 0.333 \\ P_t = 0.01 & \rightarrow & P_c = \frac{0.9 \times 0.01}{0.9 \times 0.01 + 0.2 \times 0.99} \approx 0.0435 \\ P_t = 0.001 & \rightarrow & P_c = \frac{0.9 \times 0.001}{0.9 \times 0.01 + 0.2 \times 0.999} \approx 0.00448 \\ P_t = 0.0001 & \rightarrow & P_c = \frac{0.9 \times 0.0001}{0.9 \times 0.01 + 0.2 \times 0.9999} \approx 0.000449 \end{eqnarray*}

 

条件付き確率~ベイズの定理

基本

条件付き確率~ベイズの定理(Baye’s theorem)はとっつきにくいが、記号の使い方や図で理解することでわかりやすくなる。

事象Cが起こったときの事象Tの条件付き確率は、以下で計算される。

(1)    \begin{equation*} P(T|C) = \frac{P(T \cap C)}{P(C)} \end{equation*}

事象の記号ににABを使うのが一般的だが、私にはどちらが条件で、どちらを最終的に求めたいのかわかりにくいので、ここでT:Target、C:Conditionの記号を使った

P(C)は事前確率、P(T|C)は事後確率、P(T \cap C)は同時確率と呼ばれる。

ここでP(T \cap C) = P(T|C)P(C) = P(C|T)P(T)から、以下のようにも表現される。

(2)    \begin{equation*} P(T|C) = \frac{P(C|T)P(T)}{P(C)} \end{equation*}

一般表現

より一般的には、事象T_i (i=1,2,\ldots ,n)は互いに排反で、C \subset T_1 \cup T_2 \ldots \cup T_nとするとき、

(3)    \begin{equation*} P(T_k |C) = \frac{P(C|T_k)P(T_k)}{\displaystyle \sum_{i=1}^n P(C|T_i)P(T_i)} \end{equation*}

たとえばT_kが背反する2事象、すなわちT, \overline{T}の場合は以下のようになる。

(4)    \begin{equation*} P(T|C) = \frac{P(C|T)P(T)}{P(C|T)P(T) + P(C|\overline{T})P(\overline{T})} \end{equation*}

具体例として、癌などの難病の検査に関する問題が見られる。

また、ベイズの定理についてこちらでもう少し詳しい解釈をしている。

 

中心極限定理

概要

中心極限定理(central limit theorem: CLT)は、一言で言えば次のようになる。「母集団がどのような確率分布に従うとしても、標本の数を十分大きくしたときには、その合計値あるいは標本平均は、正規分布に従う」

具体的には、母集団の平均を\mu、標準偏差を\sigmaとし、nが十分に大きいとき、

  • 標本の合計S_n = \sum X_{i}は正規分布N(n \mu,n\sigma^2)に従う
  • 標本平均\overline{X}_n = \frac{1}{n} \sum X_{i}は正規分布N(\mu, \frac{\sigma^2}{n})に従う

 

表現

中心極限定理は、一般には以下のように表される。

(1)    \begin{equation*} \lim_{n \rightarrow \infty} \Pr \left( \frac{S_n - n \mu}{\sqrt{n} \sigma} \leq \alpha \right) = \int_{-\infty}^{\alpha} \frac{1}{\sqrt{2} \pi} e^{- \frac{x^2}{2}} dx \end{equation*}

これを少し変形すると、

(2)    \begin{equation*} \lim_{n \rightarrow \infty} \Pr \left( \frac{\overline{X}_n - \mu}{\sigma / \sqrt{n}} \leq \alpha \right) = \int_{-\infty}^{\alpha} \frac{1}{\sqrt{2} \pi} e^{- \frac{x^2}{2}} dx \end{equation*}

実用

たとえば、サイコロをn回振った目の合計を考える。全て1(合計がn)や全て6(合計が6n)というケースは稀なので、その間の値になりそうだと予想される。

中心極限定理を用いると、n個のサイコロの目の平均と分散より、n個のサイコロの目の合計は、N( \frac{7}{2} , \frac{35}{12n})に従うことになる。

これをRの下記コードで試してみた。一回の試行でサイコロを投げる回数をn.dicesに設定して、その平均を求める試行を1000回繰り返す。

n.dicesの回数を変化させた実行結果は以下の通りで、このケースの場合は、n=10程度でもかなり平均の周りに尖った分布となる。

CLT_dice_n=01CLT_dice_n=02n=5n=10

 

 

単純な事象の平均と分散

コイン

コイントスで表→1、裏→0としたときの平均、分散。分布は{表, 裏]の一様分布。

平均

(1)    \begin{equation*} E(X) = 0\times \frac{1}{2} + 1 \times \frac{1}{2} = \frac{1}{2} \end{equation*}

分散

(2)    \begin{equation*} V(X) = 0^2 \cdot \frac{1}{2} + 1^2 \cdot \frac{1}{2} - \left( \frac{1}{2} \right)^2 = \frac{1}{4} \end{equation*}

サイコロ

サイコロ1つを投げたときの目の数の平均、分散。分布は{1,2, 3, 4, 5, 6}の一様分布。

平均

(3)    \begin{equation*} E(X) = \sum_{k=1}^6 k \times \frac{1}{6} = \frac{1+2+3+4+5+6}{6} = \frac{21}{6} = \frac{7}{2} = 3.5 \end{equation*}

分散

(4)    \begin{eqnarray*} V(X) &=& \sum_{k=1}^6 k^2 \times \frac{1}{6} - \left( \frac{7}{2} \right)^2 = \frac{1 + 4 + 9 +16 + 25 + 36}{6} - \frac{49}{4} \\ &=& \frac{91}{6} - \frac{49}{4} = \frac{364-294}{24} = \frac{70}{24} \\ &=& \frac{35}{12} \simeq 2.92 \end{eqnarray*}

トランプ

トランプを一枚引いたときの数の平均、分散。ここでAは1とする。分布は{1, …, 13}の一様分布。

平均

(5)    \begin{equation*} E(X) = \sum_{k=1}^{13} k\times \frac{1}{13} = \frac{13 \times 14}{2} \cdot \frac{1}{13} = 7 \end{equation*}

分散

(6)    \begin{eqnarray*} V(X) &=& \sum_{k=1}^{13} k^2 \times \frac{1}{13} - 7^2 \\ &=& \frac{13 (13+1) (2 \cdot 13 +1)}{6} \frac{1}{13} -49 = \frac{14 \cdot 27}{6} -49 = 7 \cdot 9 - 49 \\ &=& 14 \end{eqnarray*}