Ch.8 「質的データの統計的分析」の章末問題の解答例 〜自然科学の統計学(東京大学出版会)〜

当記事は「基礎統計学Ⅲ 自然科学の統計学(東京大学出版会)」の読解サポートにあたってChapter.8の「質的データの統計的分析」の章末問題の解説について行います。
※ 基本的には書籍の購入者向けの解答例・解説なので、まだ入手されていない方は下記より入手をご検討いただけたらと思います。また、解説はあくまでサイト運営者が独自に作成したものであり、書籍の公式ページではないことにご注意ください。(そのため著者の意図とは異なる解説となる可能性はあります)

下記などの演習も並行で抑えておくと良いと思います。
https://www.hello-statisticians.com/practice/stat_practice7.html

章末の演習問題について

問題8.1の解答例

i)
下記を実行することで$\alpha_0, \alpha_1$の推定値の$\hat{\alpha}_0, \hat{\alpha}_1$を得ることができる。

import numpy as np
from scipy import linalg

X = np.ones([5,2])
X[:,1] = np.arange(0.,5.,1)
y = np.array([0.067, 0.267, 0.500, 0.767, 0.900]).reshape([5,1])
X_ = linalg.inv(np.dot(X.T,X))
alpha = np.dot(X_,np.dot(X.T,y))

print(alpha)

・実行結果

> print(alpha)
[[ 0.067 ]
 [ 0.2166]]

上記より、$\hat{\alpha}_0=0.067, \hat{\alpha}_1=0.2166$が得られる。

ⅱ)
下記を実行することで結果を得ることができる。

import numpy as np
from scipy import linalg

X = np.ones([5,2])
X[:,1] = np.arange(0.,5.,1)
y = np.array([0.067, 0.267, 0.500, 0.767, 0.900]).reshape([5,1])
X_ = linalg.inv(np.dot(X.T,X))
alpha = np.dot(X_,np.dot(X.T,y))

print(np.dot(X,alpha))

・実行結果

> print(np.dot(X,alpha))
[[ 0.067 ]
 [ 0.2836]
 [ 0.5002]
 [ 0.7168]
 [ 0.9334]]

ⅲ)
省略

iv)
$0.2166X+0.067>1, 0.2166X+0.067<0$をそれぞれ解けばよい。
・$0.2166X+0.067>1$
$X>4.3074…$

・$0.2166X+0.067<0$
$X<-0.3093…$

問題8.2の解答例

i)
$\displaystyle P = F(\beta_{0}+\beta_{1}X) = \frac{1}{2}$は$\beta_{0}+\beta_{1}X=0$の時に成立する。よってこのときの$X$の値は下記のように導出できる。
$$
\large
\begin{align}
\beta_{0} + \beta_{1}X &= 0 \\
\beta_{1}X &= -\beta_{0} \\
X &= – \frac{\beta_{0}}{\beta_{1}}
\end{align}
$$

ⅱ)
それぞれ下記のように考えることができる。
・プロビット
$$
\large
\begin{align}
P &= F(\beta_{0}+\beta_{1}X) \\
F(Z) &= \int_{-\infty}^{Z} \frac{1}{\sqrt{2 \pi}} e^{-\frac{x^2}{2}} dx
\end{align}
$$

上記を元に$\displaystyle \frac{dP}{dX}$を計算する。
$$
\large
\begin{align}
\frac{dP}{dX} &= \frac{dP}{dZ} \cdot \frac{dZ}{dX} \\
&= \frac{1}{\sqrt{2 \pi}} e^{-\frac{Z^2}{2}} \cdot \beta_{1}
\end{align}
$$

上記を$Z$の関数と見ると、$Z=0$のとき最大値を取ることがわかる。よって、$\displaystyle \left| \frac{dP}{dX} \right|$は$\displaystyle Z=0, X = – \frac{\beta_{0}}{\beta_{1}}$のとき下記のような最大値を取る。
$$
\large
\begin{align}
\frac{1}{\sqrt{2 \pi}} e^{-\frac{0^2}{2}} \cdot \beta_{1} &= \frac{1}{\sqrt{2 \pi}} \cdot \beta_{1} \\
&= 0.3989 \beta_{1}
\end{align}
$$

・ロジット
$$
\large
\begin{align}
P &= F(\beta_{0}+\beta_{1}X) \\
F(Z) &= \frac{e^{Z}}{1+e^{Z}} \\
&= \frac{1}{1+e^{-Z}}
\end{align}
$$

$Z=\beta_{0}+\beta_{1}X$とおき、上記を元に$\displaystyle \frac{dP}{dX}$を計算する。
$$
\large
\begin{align}
\frac{dP}{dX} &= \frac{dP}{dZ} \cdot \frac{dZ}{dX} \\
&= -\frac{1}{(1+e^{-Z})^2} \cdot -e^{-Z} \cdot \beta_{1} \\
&= \frac{e^{-Z}}{(1+e^{-Z})^2} \cdot \beta_{1}
\end{align}
$$

ここで上記を最大にする$Z$がわかれば$Z=\beta_{0}+\beta_{1}X$より対応する$X$の値も計算できる。よって、以下では$\displaystyle g(Z) = \frac{e^{-Z}}{(1+e^{-Z})^2}$を最大にする$Z$を考える。
$$
\large
\begin{align}
\frac{dg(Z)}{dZ} &= \frac{d}{dZ} ((1+e^{-Z})^{-2}e^{-Z}) \\
&= -2(1+e^{-Z})^{-3} \times -e^{-Z} \times e^{-Z} + (1+e^{-Z})^{-2} \times -e^{-Z} \\
&= \frac{e^{-Z}(e^{-Z}-1)}{2(1+e^{-Z})^{3}}
\end{align}
$$

上記は単調減少かつ$Z=0$で$\displaystyle \frac{dg(Z)}{dZ}=0$なので$Z=0$のとき$g(Z)$は最大値を取る。このことは$\displaystyle \left| \frac{dP}{dX} \right|$は$\displaystyle Z=0, X = – \frac{\beta_{0}}{\beta_{1}}$のとき下記のような最大値を取ることを意味する。
$$
\large
\begin{align}
\frac{e^{-0}}{(1+e^{-0})^2} \cdot \beta_{1} &= \frac{1}{(1+1)^2} \cdot \beta_{1} \\
&= 0.25 \beta_{1}
\end{align}
$$

問題8.3の解答例

例8.2と例8.5に関してそれぞれ確認を行う。
・例8.2

print(2.324/1.086)
print(1.162/0.543)

・実行結果

> print(2.324/1.086)
2.13996316759
> print(1.162/0.543)
2.13996316759

・例8.5

print(50.26/27.56)
print(2.526/1.358)
print(0.612/0.337)

・実行結果

> print(50.26/27.56)
1.8236574746
> print(2.526/1.358)
1.86008836524
> print(0.612/0.337)
1.81602373887

$\displaystyle \frac{\pi}{\sqrt{3}} \simeq 1.81$より、例8.5で概ね成立していることが確認できる。

問題8.4の解答例

下記を実行することで結果を得ることができる。

import numpy as np

X = np.array([[0.,89.,0.], [0.,85.,0.], [1.,84.,1.], [0.,83.,0.], [1.,81.,1.], [0.,88.,1.], [1.,87.,0.], [1.,86.,0.], [0.,86.,1.], [1.,84.,10.], [0.,87.,0.], [1.,83.,1.], [1.,82.,1.], [0.,83.,0.], [0.,82.,0.], [0.,86.,0.], [0.,78.,0.], [0.,77.,1.], [0.,85.,1.], [1.,75.,1.], [1.,80.,0.], [1.,82.,0.], [1.,79.,1.], [0.,83.,1.], [0.,82.,1.], [0.,89.,0.], [1.,86.,1.], [1.,80.,1.], [1.,81.,0.], [0.,78.,0.]])
y = np.array([1., 1., 0., 1., 0., 1., 0., 0., 1., 0., 1., 0., 0., 1., 0., 1., 0., 0., 1., 0., 0., 0., 0., 0., 0., 1., 1., 0., 1., 0.])

y_probit = -27.56 - 1.358*X[:,0] + 0.337*X[:,1] - 0.447*X[:,2]
y_probit[y_probit>0] = 1
y_probit[y_probit<=0] = 0
y_logi = -50.26 - 2.526*X[:,0] + 0.612*X[:,1] - 0.610*X[:,2]
y_logi[y_logi>0] = 1
y_logi[y_logi<=0] = 0

print(np.float(np.sum(y==y_probit))/y.shape[0])
print(np.float(np.sum(y==y_logi))/y.shape[0])

・実行結果

> print(np.float(np.sum(y==y_probit))/y.shape[0])
0.833333333333
> print(np.float(np.sum(y==y_logi))/y.shape[0])
0.9

上記より正答率はプロビットに対して$83.3$%、ロジットに対して$90.0$%であることがわかる。

問題8.5の解答例

まとめ

https://www.amazon.co.jp/dp/4130420674