ブログ

解答

[1] 解答

$\boxed{ \ \mathsf{13}\ }$ ： ⑤

$\mu_0=0, \sigma_0=1, \sigma=2$のとき観測値$X=2$が得られた場合、選択肢の事後分布の平均値$\tilde\mu$の値を求めてみると、
・①,②,③$$\tilde\mu=\frac{\sigma^2X+\sigma_0^2\mu_0}{\sigma^2+\sigma_0^2}=\frac{2^2\times2+1^2\times0}{2^2+1^2}=\frac{8}{5}=1.6$$
・④,⑤$$\tilde\mu=\frac{\sigma_0^2X+\sigma^2\mu_0}{\sigma^2+\sigma_0^2}=\frac{1^2\times2+2^2\times0}{2^2+1^2}=\frac{2}{5}=0.4$$
グラフを見ると事後分布が最大になるところが概ね$0.4$付近でなので④か⑤が正解候補。
また、グラフから事前分布の最大値より事後分布の最大値が大きいことから、事後分布の分散$\tilde\sigma$は事前分布の分散$\sigma_0$より小さくなることから、正解は⑤となる。

[2] 解答

$\boxed{ \ \mathsf{14}\ }$ ： ②

事前分布が正規分布なので、事後分布も正規分布になることから、③と④は除外される。事後分散は
$$\left(\frac{n}{\sigma^2}+\frac{1}{\sigma_0^2}\right)^{-1}=\left(\frac{10}{3.07^2}+\frac{1}{2.7^2}\right)^{-1}\fallingdotseq0.83$$
となり、事前分布の分散より小さくなるので②が正解である。

[3] 解答

$\boxed{ \ \mathsf{15}\ }$ ： ③

$\mathrm{(A)}$と$\mathrm{(C)}$は正しい。$\mathrm{(B)}$は分散が未知のときの分散に対する共役事前分布としては、逆ガンマ分布が知られている。

解説

事前確率分布と事後確率分布

$n$個の同時に成立しない互いに排反な原因事象 $H_1, H_2, \ldots, H_n$ について、この原因の発生する確率$P(H_1), P(H_2), \ldots, P(H_n)$は事前にわかっているものとする（事前確率）。この原因によって結果となる事象$A$が発生するが、このとき各原因$H_i$によって結果$A$が生じる条件付き確率を$P(A|H_i)$とすると、結果$A$によって各原因の発生する確率は更新され、その確率はベイズの定理より、
$$P(H_i|A)=\frac{P(A|H_i)P(H_i)}{\sum_{k=1}^nP(A|H_k)P(H_k)}$$
となり、この$P(H_i|A)$を事後確率という。

これを連続型確率分布に適用する。データ$x$（結果）を観測する以前に持っている母数$\mathbf\theta$（原因）の確率分布を事前(確率)分布といい、これを$\pi(\mathbf\theta)$とする。この母数$\pi(\mathbf\theta)$のもとで観測されたデータ$x$の条件付き確率分布を$f(x|\mathbf\theta)$とすると、ベイズの定理より、
$$\pi(\mathbf\theta|x)=\frac{f(x|\mathbf\theta)\pi(\mathbf\theta)}{\displaystyle\int_\Theta f(x|\mathbf\theta)\pi(\mathbf\theta)}$$
となる。ここに、$\theta$は母数$\theta$の取りうる値の集合（母数空間）で、$\pi(\mathbf\theta|x)$はデータの観測によって更新された母数$\mathbf\theta$の確率分布で事後(確率)分布という。上式の分母$$\displaystyle\int_\Theta f(x|\mathbf\theta)\pi(\mathbf\theta)$$は周辺尤度といい、$\mathbf\theta$によらない値である。よって、事後分布は$$\pi(\mathbf\theta|x)\propto f(x|\mathbf\theta)\pi(\mathbf\theta)$$と表わすことがある。

正規分布の事後分布

母平均$\mu$、母分散$\sigma^2$（母分散は既知）に従う正規母集団から大きさ$n$の標本を抽出し標本平均$\bar x$を得たとする。母平均$\mu$の事前分布として平均$\mu_0$、分散$\sigma_0^2$の正規分布としたとき、平均$\mu$の事後分布を求める。

事前分布及び母集団の確率密度関数は、
$$
\begin{eqnarray}
\mu\sim N(\mu_o,\sigma_0^2)\Rightarrow\pi(\mu)=\frac1{\sqrt{2\pi\sigma_o^2}}\exp\left[-\frac{(\mu-\mu_0)^2}{2\sigma_0^2}\right]\\
X\sim N(\mu,\sigma^2)\Rightarrow f(x)=\frac1{\sqrt{2\pi\sigma^2}}\exp\left[-\frac{(x-\mu)^2}{2\sigma^2}\right]
\end{eqnarray}
$$
となる。抽出された標本データ$x_1, x_2, \cdots, x_n$ は正規分布に独立同一に従う$(i.i.d.)$ので、その同時確率分布は、
$$
f(x|\mu)=\prod_{i=1}^n\frac1{\sqrt{2\pi\sigma^2}}\exp\left[-\frac{(x_i-\mu)^2}{2\sigma^2}\right]=\left(\frac1{\sqrt{2\pi\sigma^2}}\right)^n\exp\left[-\sum_{i=1}^n\frac{(x_i-\mu)^2}{2\sigma^2}\right]
$$
となる。この式のうち、指数部分についてみると、標本平均$\bar x=\sum_{i=1}^nx_i/n$から、
$$\begin{align}
-\sum_{i=1}^n\frac{(x_i-\mu)^2}{2\sigma^2}&=-\frac1{2\sigma^2}{\sum_{i=1}^n(x_i-\mu)^2}\\
&=-\frac1{2\sigma^2}\left(n\mu^2-2\sum_{i=1}^nx_i\mu+\sum_{i=1}^nx_i^2\right)\\
&=-\frac1{2\sigma^2}\left(n\mu^2-2n\bar x\mu+\sum_{i=1}^nx_i^2\right)\\
&=-\frac1{2\sigma^2}\left\{n(\mu-\bar x)^2+\sum_{i=1}^n(x_i-\bar x)^2\right\}\\
\end{align}$$
したがって、
$$\begin{align}
f(x|\mu)&=\left(\frac1{\sqrt{2\pi\sigma^2}}\right)^n\exp\left[-\frac1{2\sigma^2}\left\{n(\mu-\bar x)^2+\sum_{i=1}^n(x_i-\bar x)^2\right\}\right]\\
&=\exp\left[-\frac{n(\mu-\bar x)^2}{2\sigma^2}\right]\cdot\left(\frac1{\sqrt{2\pi\sigma^2}}\right)^n\exp\left[-\frac{\sum_{i=1}^n(x_i-\bar x)^2}{2\sigma^2}\right]
\end{align}$$
となり、$\mu$の関数は第一因数だけとなる。
$$\begin{align}
\therefore\ \pi(\mu|x)\propto f(x|\mu)\pi(\mu)&\propto\exp\left[-\frac{n(\mu-\bar x)^2}{2\sigma^2}\right]\exp\left[-\frac{(\mu-\mu_0)^2}{2\sigma_0^2}\right]\\
&=\exp\left[-\frac{n(\mu-\bar x)^2}{2\sigma^2}-\frac{(\mu-\mu_0)^2}{2\sigma_0^2}\right]\\
&=\exp\left[-\frac{n\sigma_0^2(\mu-\bar x)^2+\sigma^2(\mu-\mu_0)^2}{2\sigma^2\sigma_0^2}\right]\\
&=\exp\left[-\frac{n\sigma_0^2+\sigma^2}{2\sigma^2\sigma_0^2}\left(\mu-\frac{n\sigma_0^2\bar x+\sigma^2\mu_0}{n\sigma_0^2+\sigma^2}\right)^2+\frac{n(\mu_0-\bar x)^2}{2(n\sigma_0^2+\sigma^2)}\right]\\
&=\exp\left[-\frac1{2\left(\frac{n}{\sigma^2}+\frac{1}{\sigma_0^2}\right)^{-1}}\cdot\left(\mu-\frac{n\sigma_0^2\bar x+\sigma^2\mu_0}{n\sigma_0^2+\sigma^2}\right)^2\right]\cdot\exp\left[\frac{n(\mu_0-\bar x)^2}{2(n\sigma_0^2+\sigma^2)}\right]\\
\end{align}$$
以上より、事後分布$\pi(\mu|x)$は正規分布
$$N\left(\frac{n\sigma_0^2\bar x+\sigma^2\mu_0}{n\sigma_0^2+\sigma^2},\ \left(\frac{n}{\sigma^2}+\frac{1}{\sigma_0^2}\right)^{-1}\right)$$
となる。

ここで、平均の分子分母を$\sigma^2\sigma_0^2$で除すると、
$$\frac{n\sigma_0^2\bar x+\sigma^2\mu_0}{n\sigma_0^2+\sigma^2}=\frac{(n/\sigma^2)\bar x+(1/\sigma_0^2)\mu_0}{(n/\sigma^2)+(1/\sigma_0^2)}$$
となり、事後分布の平均は、標本平均$\bar x$と事前分布の平均$\mu_0$の重み付き平均となっている。また、事後平均の分散は
（問題の[1]のケースでは、上式で$n=1,\bar x=X$を代入すると、解答と同じ結果が得られる。）

「統計学実践ワークブック」演習問題中心第4章「変数変換」

投稿日: 2021-07-032022-10-05 投稿者: yoichi_t

統計検定準$1$級対応の公式テキストである「統計学実践ワークブック」を1章から順に演習問題を中心に解説していきます。
今回は第$4$章「変数変換」です。

日本統計学会公式認定統計検定準1級対応統計学実践ワークブック

3,080円(11/10 13:44時点)

重要ポイント

確率変数にある関数をかけて変換したいことがあります。例えば、対数正規分布を考えた場合など。

Xの確率密度関数$f(X)$とした場合に、新たな確率変数$Y=g(X)$とした場合のYの確率密度関数は、単純に$f(g(X))$ではありません。以下の式となります。

$$
\begin{align}
h(y) &= f(g^{-1}(y))\frac{d(g^{-1}(y))}{dy}
\end{align}
$$

この辺りの詳説は別途記事を作成予定です。また、参考文献に挙げている書籍がわかりやすいと思います。

変数変換と置換積分

「問題演習で理解する統計学【$15$】」の「変数変換と置換積分」で取り扱いました。

ヤコビアンの図形的解釈

「問題演習で理解する統計学【$15$】」の「行列式$\det \mathbf{A}$の図形的解釈と平行四辺形の面積」と「ヤコビ行列$\mathbf{J}$とヤコビアン$\det \mathbf{J}$」で取り扱いました。

演習問題解説

問4.1

対数正規分布に従う確率変数Yの期待値、分散、確率密度関数を導出せよという問題です。

正規分布に従う確率変数$X$に指数関数をかけて変換したものを$Y$という確率変数にすることで、対数正規分布が得られます。ということで、変数変換を使って導出していこうというのが問題の趣旨になります。

(1) $Y$の期待値を求めよ

正規分布に従う確率変数$X$を利用して、$\exp(X)$の期待値を考えます。定義に従って愚直に計算していけば導出できますが、次の問では分散を導出することになっているため、モーメント母関数を利用します。

モーメント母関数$m(\theta)$は、$m(\theta) = E[\exp \left\{ \theta X \right\}]$です。つまり、$\theta=1$とすれば$Y$の期待値が導出できることがわかります。

$$
\begin{align}
E[\exp \left\{ \theta X \right\}] &= \int \exp\left\{ \theta X \right\} \frac{1}{\sqrt{2 \pi \sigma^2}}\exp \left\{ – \frac{(x-\mu)^2}{2\sigma^2} \right\} dx \\
&= \frac{1}{ \sqrt{2 \pi \sigma^2} } \int \exp\left\{ \theta x – \frac{(x-\mu)^2}{2\sigma^2} \right\} dx \\
&= \frac{1}{ \sqrt{2 \pi \sigma^2} } \int \exp \left\{ – \frac{x^2 – 2\mu x + \mu^2 – 2 \sigma^2 \theta x}{2 \sigma^2} \right\} dx \\
&= \frac{1}{ \sqrt{2 \pi \sigma^2} } \int \exp \left\{ – \frac{1}{2\sigma^2} (x – \mu – \sigma^2 \theta)^2 + \left( \theta \mu + \frac{\theta^2}{2}\sigma^2 \right) \right\} dx \\
&= \frac{1}{ \sqrt{2 \pi \sigma^2} } \exp \left\{ \theta \mu + \frac{\theta^2}{2}\sigma^2 \right\} \int \exp \left\{ – \frac{1}{2\sigma^2} (x – \mu – \sigma^2 \theta)^2 \right\} dx
\end{align}
$$

上記の3行目から4行目への変換は平方完成をしただけです。これで、$x$に関する項と$x$が関係しない項に分解できるので、積分の中の一部を外に出すことができシンプルになります(5行目)。
最後に残った積分ですが、$(x – \mu – \sigma^2 \theta) = t$とおくことで、ガウス積分が利用できることがわかります。

$$
\begin{align}
\int \exp \left\{ – \frac{1}{2\sigma^2} (x – \mu – \sigma^2 \theta)^2 \right\} dx &= \int \exp\left\{ – \frac{1}{2\sigma^2} t^2 \right\} dt \\
&= \sqrt{2 \pi \sigma^2}
\end{align}
$$

よって、正規分布のモーメント母関数$m(\theta)$は次のようになることがわかります。

$$
\begin{align}
m(\theta) &= E[\exp \left\{ \theta X \right\}] \\
&= \exp \left\{ \theta\mu + \frac{\theta^2}{2} \sigma^2 \right\}
\end{align}
$$

ここから、$\theta=1$として、$Y = \exp \{X\}$の期待値が導出できます。

$$
E[\exp \left\{ X \right\}] = \exp \left\{ \mu + \frac{1}{2} \sigma^2 \right\}
$$

(2) $Y$の分散を求めよ

分散は以下のように変形できます(参考)。

$$
\begin{align}
V[Y] = E[Y^2] – (E[Y])^2
\end{align}
$$

$E[Y]$については(1)で導出しているので、$E[Y^2]$を考えます。ここで、$Y^2 = \exp \left\{ 2X \right\}$となるため、$E[Y^2]$は$E[Y^2] = E[\exp \left\{ 2X \right\}]$となります。これは(1)で求めた正規分布のモーメント母関数で$\theta=2$とおいたものになります。

$$
\begin{align}
E[Y^2] &= E[\exp \left\{ 2X \right\}] \\
&= \exp\left\{ 2\mu + 2\sigma^2 \right\}
\end{align}
$$

よって分散$V[Y]$は以下の通り導出できます。

$$
\begin{align}
V[Y] &= E[Y^2] – (E[Y])^2 \\
&= \exp\left\{ 2\mu + 2\sigma^2 \right\} – \left( \exp \left\{ \mu + \frac{1}{2} \sigma^2 \right\} \right) \\
&= \exp\left\{ 2\mu + 2\sigma^2 \right\} – \exp \left\{ 2\mu + \sigma^2 \right\}
\end{align}
$$

(3) $Y$の確率密度関数を求めよ

密度関数は変数変換による密度関数の変換式(下式)を利用します。

$$
\begin{align}
h(y) &= f(g^{-1}(y))\frac{d(g^{-1}(y))}{dy}
\end{align}
$$

ここで、変数を変換する関数$g$は$\exp$なので、以下のようになります。

$$
\begin{align}
g^{-1}(y) &= \ln y \\
\frac{d(g^{-1}(y))}{dy} &= \frac{1}{y}
\end{align}
$$

よって、確率密度関数は以下の通り導出できます。

$$
\begin{align}
h(y) &= f(g^{-1}(y))\frac{d(g^{-1}(y))}{dy} \\
&= \frac{1}{\sqrt{2\pi \sigma^2}} \exp \left\{ – \frac{(\ln y – \mu)^2}{2\sigma^2} \right\} \frac{1}{y} \\
&= \frac{1}{\sqrt{2\pi \sigma^2}y} \exp \left\{ – \frac{(\ln y – \mu)^2}{2\sigma^2} \right\}
\end{align}
$$

問4.2

指数分布に従う独立な確率変数$X, Y$について、$X+Y$の確率密度関数を求めよとの問題です。

$X+Y$の確率密度関数を求めよ

$X+Y = z$とおいて、$z$の密度関数を考えます。($P(X+Y) = k(z)$)

ここで、$X+Y=z$となるのは、$X=x$で$Y=z-x$の形式における$x$の全ての場合です。$X,Y$は独立なので、確率密度の積でこれを表し、$x$の全ての場合についての考慮として、$x$について積分をします。

また、$z$が与えられれば、$x$の範囲は$0$〜$z$の範囲に限定されるため、積分範囲がこの範囲に限定されます。

$$
\begin{align}
k(z) &= \int^{\infty}_{-\infty} g(x)h(z-x) dx \\
&= \int^{\infty}_{-\infty} \lambda \exp(-\lambda x) \lambda \exp(-\lambda (z-x)) dx \\
&= \lambda^2 \int^{z}_{0} \exp \left\{ -\lambda z \right\} dx \\
&= \lambda^2 \exp\left\{ – \lambda z \right\} \left[ x \right]^z_0 \\
&= \lambda^2 z \exp \left\{ – \lambda z \right\}
\end{align}
$$

こちらの回答は参考文献に挙げた書籍から引用しました。ワークブックの記述よりもわかりやすいかなと思います。

確率変数の変換による同時確率密度関数の導出の詳細

$Z=X+Y, W=Y$のように確率変換の変換を行うことを考える。このとき$x=z-w, y=w$のように対応するので、変数変換におけるヤコビ行列式$\det{\mathbf{J}}$は下記のように計算できる。
$$
\large
\begin{align}
\det{\mathbf{J}} &= \left| \begin{array}{cc} \frac{\partial x}{\partial z} & \frac{\partial x}{\partial w} \\ \frac{\partial y}{\partial z} & \frac{\partial y}{\partial w} \end{array} \right| \\
&= \left| \begin{array}{cc} 1 & -1 \\ 0 & 1 \end{array} \right| \\
&= 1
\end{align}
$$

よって、変数変換後の同時確率密度関数を$g(z,w)=g(z)g(w)$とおくとき、下記のように表すことができる。
$$
\large
\begin{align}
g(z,w) &= g(z)g(w) = f(x)f(y) \left| \begin{array}{cc} 1 & -1 \\ 0 & 1 \end{array} \right| \\
&= f(x)f(y) = f(z-w)f(w)
\end{align}
$$

以下、ワークブックの内容に沿って$w$の積分を行うことで確率密度関数を計算できる。

参考文献

ワークブック以外の参考資料として以下のものがおすすめです。

松原ら, 統計学入門, 1991, 東京大学出版会

統計学入門 (基礎統計学Ⅰ)

確率分布の再生性とその導出（二項分布・ポアソン分布・正規分布）

投稿日: 2021-07-022022-11-18 投稿者: lib-arts

確率分布の再生性は二項分布、ポアソン分布、正規分布において成立しますが、導出にあたっては大きく分けると「①モーメント母関数を用いる手法」と「②畳み込み(convolution)を用いる手法」の二つが存在します。
当記事では二項分布、ポアソン分布、正規分布のそれぞれの分布に対し、二通りの導出を確認します。畳み込みを用いる手法ではかなりトリッキーな変形などを行うので、モーメント母関数を考える利点についても実感できるようなトピックであるのではないかと思います。
モーメント母関数については下記の記載も参照してみてください。
https://www.hello-statisticians.com/explain-terms-cat/prob_generating.html

二項分布の再生性の導出

モーメント母関数を用いた導出

二項分布のモーメント母関数を考えるにあたっては、ベルヌーイ分布の母関数を考えてそこから導出するのがシンプルでわかりやすい。二項分布の確率変数を$X=X_1+X_2+…+X_n$のように定義し、$X_i$はベルヌーイ分布の確率変数とする。このとき、$X_i$は確率$p$で$X_i=1$、確率$1-p$で$X_i=0$となるので、ベルヌーイ分布のモーメント母関数$m_{X_i}(t)$は下記のようになる。
$$
\large
\begin{align}
m_{X_i}(t) &= E\left[ e^{tX_i} \right] \\
&= e^{1 \times t}P(X_i=1|p) + e^{0 \times t}P(X_i=0|p) \\
&= e^t \times + e^0 \times (1-p) \\
&= pe^t +1 – p
\end{align}
$$

このとき二項分布のモーメント母関数$m_{X_i}(t)$は下記のように計算できる。（二項分布はベルヌーイ試行の和と対応するが、ベルヌーイ試行はそれぞれの試行がそれぞれ独立であることを仮定する）
$$
\large
\begin{align}
m_{X}(t) &= E\left[ e^{tX} \right] \\
&= E\left[ e^{t(X_1+X_2+…+X_n)} \right] \\
&= E\left[ e^{tX_1}e^{tX_2}…e^{tX_n} \right] \\
&= E\left[ e^{tX_1} \right] E\left[ e^{tX_2} \right]…E\left[ e^{tX_n} \right] \\
&= E\left[ e^{tX_i} \right]^n \\
&= (pe^t +1 – p)^n
\end{align}
$$

$n$回の試行に対応する二項分布のモーメント関数は上記のように$m_{X}(t)=(pe^t +1 – p)^n$となる。このとき確率$p$が等しい際に、$m$回と$n$回に対応する二項分布のモーメント母関数をそれぞれ$m_{X(m)}(t)$、$m_{X(n)}(t)$とすると、確率変数$X(m,n)=X(m)+X(n)$に対応するモーメント母関数は下記のようになる。
$$
\large
\begin{align}
m_{X(m,n)}(t) &= E\left[ e^{tX(m,n)} \right] \\
&= E\left[ e^{t(X(m)+X(n))} \right] \\
&= E\left[ e^{tX(m)} \right] E\left[ e^{tX(n)} \right] \\
&= (pe^t +1 – p)^m \times (pe^t +1 – p)^n \\
&= (pe^t +1 – p)^{m+n}
\end{align}
$$

よって、$X(m)$が二項分布$\mathrm{Binom}(m,p)$、$X(n)$が二項分布$\mathrm{Binom}(n,p)$に従うとき、$X(m,n)=X(m)+X(n)$は二項分布$\mathrm{Binom}(m+n,p)$に従う。（確率分布とモーメント母関数は$1$対$1$対応で、かつ$X(m,n)$のモーメント母関数は$\mathrm{Binom}(m+n,p)$のモーメント母関数に一致したことによる）

畳み込みを用いた導出

https://www.hello-statisticians.com/explain-books-cat/toukeigakunyuumon-akahon/ch7_practice-html.html
上記の問題$7.9$の解答と一致する。

ポアソン分布の再生性の導出

モーメント母関数を用いた導出

ポアソン分布のモーメント母関数を$m(t)$とすると下記のように求めることができる。
$$
\large
\begin{align}
m(t) &= E\left[ e^{tX} \right] \\
&= \sum_{k=0}^{\infty} e^{tk} \times \frac{\lambda^k e^{-\lambda}}{k!} \\
&= e^{-\lambda}\sum_{k=0}^{\infty} e^{tk} \times \frac{\lambda^k}{k!} \\
&= e^{-\lambda}\sum_{k=0}^{\infty} \frac{(\lambda e^{t})^k}{k!}
\end{align}
$$

ここで$e^x$のマクローリン展開の$\displaystyle e^{x}=1+\frac{x}{1!}+\frac{x^2}{2!}+…$と対応づけて考えると、下記のように変形することができる。
$$
\large
\begin{align}
m(t) &= e^{-\lambda}\sum_{k=0}^{\infty} \frac{(\lambda e^{t})^k}{k!} \\
&= e^{-\lambda} e^{\lambda e^{t}} \\
&= e^{(\lambda e^{t}-\lambda)} \\
&= e^{\lambda(e^{t}-1)}
\end{align}
$$

上記がポアソン分布のモーメント母関数となる。ここで確率変数$X_1$はパラメータ$\lambda_1$のポアソン分布に従い、確率変数$X_2$はパラメータ$\lambda_2$のポアソン分布に従うとする。また、このときそれぞれのモーメント母関数を$m_{X_1}(t)$、$m_{X_2}(t)$でおくとする。このときの確率変数$X=X_1+X_2$に関するモーメント母関数を$m_{X}(t)$とすると、$m_{X}(t)$は下記のように求めることができる。
$$
\large
\begin{align}
m_{X}(t) &= E\left[ e^{tX} \right] \\
&= E\left[ e^{t(X_1+X_2)} \right] \\
&= E\left[ e^{t(X_1+X_2)} \right] \\
&= E\left[ e^{tX_1} \right] E\left[ e^{tX_2} \right] \\
&= e^{\lambda_1(e^{t}-1)}e^{\lambda_2(e^{t}-1)} \\
&= e^{(\lambda_1+\lambda_2)(e^{t}-1)}
\end{align}
$$

よって、$X_1$がポアソン分布$\mathrm{Poisson}(\lambda_1)$、$X_2$がポアソン分布$\mathrm{Poisson}(\lambda_2)$に従うとき、$X=X_1+X_2$はポアソン分布$\mathrm{Poisson}(\lambda_1+\lambda_2)$に従う。（確率分布とモーメント母関数は$1$対$1$対応で、かつ$X$のモーメント母関数が$\mathrm{Poisson}(\lambda_1+\lambda_2)$のモーメント母関数に一致したことによる）

畳み込みを用いた導出

https://www.hello-statisticians.com/explain-books-cat/toukeigakunyuumon-akahon/ch7_practice-html.html
上記の問題$7.9$の解答と一致する。

正規分布の再生性の導出

モーメント母関数を用いた導出

平均$\mu$、分散$\sigma^2$の正規分布$\mathcal{N}(\mu,\sigma^2)$のモーメント母関数を$m(t)$とし、簡易的に導出を行う。
$$
\large
\begin{align}
m(t) &= E[e^{tX}] \\
&= \int_{-\infty}^{\infty} e^{tx} \times \frac{1}{\sqrt{2 \pi \sigma^2}}e^{-\frac{(x-\mu)^2}{2 \sigma^2}} dx \\
&= \int_{-\infty}^{\infty} \frac{1}{\sqrt{2 \pi \sigma^2}} \exp \left[ -\frac{(x-\mu)^2}{2 \sigma^2}+tx \right] dx \\
&= \int_{-\infty}^{\infty} \frac{1}{\sqrt{2 \pi \sigma^2}} \exp \left[ -\frac{x^2 – 2\mu x + \mu^2 – 2 \sigma^2 tx}{2 \sigma^2} \right] dx \\
&= \int_{-\infty}^{\infty} \frac{1}{\sqrt{2 \pi \sigma^2}} \exp \left[ -\frac{x^2 – 2(\mu + \sigma^2 t)x + (\mu + \sigma^2 t)^2 – (\mu + \sigma^2 t)^2 + \mu^2}{2 \sigma^2} \right] dx \\
&= \exp \left[ \frac{(\mu + \sigma^2 t)^2 – \mu^2}{2 \sigma^2} \right] \int_{-\infty}^{\infty} \frac{1}{\sqrt{2 \pi \sigma^2}}e^{-\frac{(x – (\mu + \sigma^2 t))^2}{2 \sigma^2}} dx \\
&= \exp \left[ \frac{\cancel{\mu^2} + 2 \mu \sigma^2 t + \sigma^4 t^2 – \cancel{\mu^2}}{2 \sigma^2} \right] \int_{-\infty}^{\infty} \frac{1}{\sqrt{2 \pi \sigma^2}}e^{-\frac{(x – (\mu + \sigma^2 t))^2}{2 \sigma^2}} dx \\
&= e^{\mu t + \frac{\sigma^2 t^2}{2}} \int_{-\infty}^{\infty} \frac{1}{\sqrt{2 \pi \sigma^2}}e^{-\frac{(x – (\mu + \sigma^2 t))^2}{2 \sigma^2}} dx \\
&= e^{\mu t + \frac{\sigma^2 t^2}{2}}
\end{align}
$$

上記のように、$e^{f(x)}$を$x$に関して平方完成し、$x$に関係しない定数項を積分の外に出し、確率分布の定義から積分の中を$1$とすることでモーメント母関数を求めることができる。

ここで確率変数$X_1$はパラメータ$\mu_1$、$\sigma_1^2$の正規分布$\mathcal{N}(\mu_1,\sigma_1^2)$に従い、確率変数$X_2$はパラメータ$\mu_2$、$\sigma_2^2$の正規分布$\mathcal{N}(\mu_2,\sigma_2^2)$に従うとする。また、このときそれぞれのモーメント母関数を$m_{X_1}(t)$、$m_{X_2}(t)$でおくとする。このときの確率変数$X=X_1+X_2$に関するモーメント母関数を$m_{X}(t)$とすると、$m_{X}(t)$は下記のように求めることができる。
$$
\large
\begin{align}
m_{X}(t) &= E[e^{tX}] \\
&= E[e^{t(X_1+X_2)}] \\
&= E[e^{t(X_1+X_2)}] \\
&= E[e^{tX_1}]E[e^{tX_2}] \\
&= e^{\mu_1 t + \frac{\sigma_1^2 t^2}{2}}e^{\mu_2 t + \frac{\sigma_2^2 t^2}{2}} \\
&= e^{(\mu_1+\mu_2)t + \frac{(\sigma_1^2+\sigma_2^2)t^2}{2}}
\end{align}
$$

よって、$X_1$が正規分布$\mathcal{N}(\mu_1,\sigma_1^2)$、$X_2$が正規分布$N(\mu_2,\sigma_2^2)$に従うとき、$X=X_1+X_2$は正規分布$\mathcal{N}(\mu_1+\mu_2,\sigma_1^2+\sigma_2^2)$に従う。（確率分布とモーメント母関数は$1$対$1$対応で、かつ$X$のモーメント母関数が$\mathcal{N}(\mu_1+\mu_2,\sigma_1^2+\sigma_2^2)$のモーメント母関数に一致したことによる）

畳み込みを用いた導出

https://www.hello-statisticians.com/explain-books-cat/toukeigakunyuumon-akahon/ch7_practice-html.html
上記の問題$7.9$の解答と一致する。

参考
・基礎統計学Ⅰ 統計学入門(東京大学出版会)

統計学入門 (基礎統計学Ⅰ)

-> 7.4章にモーメント母関数を用いた簡単な導出と、畳み込みを用いたポアソン分布の導出が記載されている。

Ch.5 「確率変数」の章末問題の解答例　〜基礎統計学Ⅰ 統計学入門(東京大学出版会)〜

投稿日: 2021-07-012022-12-31 投稿者: lib-arts

当記事は基礎統計学Ⅰ 統計学入門(東京大学出版会)」の読解サポートにあたってChapter.5の「確率変数」の章末問題の解説について行います。
※ 基本的には書籍の購入者向けの解説なので、まだ入手されていない方は下記より入手をご検討いただけたらと思います。また、解説はあくまでサイト運営者が独自に作成したものであり、書籍の公式ページではないことにご注意ください。（そのため著者の意図とは異なる解説となる可能性はあります）
https://www.amazon.co.jp/dp/4130420658

統計学入門 (基礎統計学Ⅰ)

章末の演習問題について

問題5.1の解答例

i)
$[0,6]$の一様分布なので、密度関数は$0 \leq x \leq 6$で$\displaystyle f(x) = \frac{1}{6}$、それ以外の区間で$f(x)=0$となる。また、期待値$E[X]$と分散$V[X]$は下記のように計算できる。
$$
\begin{align}
E[X] &= \int_{0}^{6} xf(x) dx \\
&= \int_{0}^{6} \frac{1}{6}x dx \\
&= \left[ \frac{1}{12}x^2 \right]_{0}^{6} \\
&= \left( \frac{6^2}{12} – \frac{0^2}{12} \right) \\
&= 3 \\
V[X] &= E[X^2] – (E[X])^2 \\
&= \int_{0}^{6} x^2f(x) dx – 3^2 \\
&= \left[ \frac{1}{18}x^3 \right]_{0}^{6} – 9 \\
&= 12-9 \\
&= 3
\end{align}
$$

ⅱ)
$$
\begin{align}
P\left(|X-E[X]| \geq k\sqrt{V[X]}\right) &\leq \frac{1}{k^2} \\
P\left(\frac{|X-E[X]|}{\sqrt{V[X]}} \geq k\right) &\leq \frac{1}{k^2} \\
P\left(\frac{|X-3|}{\sqrt{3}} \geq k\right) &\leq \frac{1}{k^2}
\end{align}
$$
チェビシェフの不等式が成立することを示すには、上記が$k>0$すべての$k$に対して成立することを示せば良い。まず、$\sqrt{3} \leq k$の際はここでの確率密度関数の定義より$\displaystyle P(|X-3| \geq \sqrt{3}k) = 0 \leq \frac{1}{k^2}$となるためチェビシェフの不等式は成立する。
また、$0<k<\sqrt{3}$のときは$\displaystyle P(|X-3| \geq \sqrt{3}k) = 1-\frac{k}{\sqrt{3}}$となり、これの$\displaystyle \frac{1}{k^2}$との大小関係を比較し、チェビシェフの不等式が成立することがわかる。（区間$0<k<\sqrt{3}$における$k$に関する3次関数の最小値問題を考えることで示すことができる。）

ⅲ)
$[0,1]$上の一様分布なので密度関数は$0 \leq x \leq 1$で$\displaystyle f(x) = 1$、それ以外の区間で$f(x)=0$となる。
このとき、期待値$E[X]$と分散$V[X]$は下記のように計算できる。
$$
\begin{align}
E[X] &= \int_{0}^{1} xf(x) dx \\
&= \int_{0}^{1} x dx \\
&= \left[ \frac{1}{2}x^2 \right]_{0}^{1} \\
&= \frac{1}{2} \\
V[X] &= \int_{0}^{1} (x-E[X])^2f(x) dx \\
&= \int_{0}^{1} \left( x-\frac{1}{2} \right)^2 dx \\
&= \left[ \frac{1}{3}\left( x-\frac{1}{2} \right)^3 \right]_{0}^{1} \\
&= \frac{1}{3} \left( \left( \frac{1}{2} \right)^3 – \left( -\frac{1}{2} \right)^3 \right) \\
&= \frac{1}{12}
\end{align}
$$
よって、歪度$\alpha_3$、尖度$\beta_4=\alpha_4-3$は下記のように計算できる。
$$
\begin{align}
\alpha_3 &= \frac{1}{V[X]^{3/2}} \int_{0}^{1} (x-E[X])^3f(x) dx \\
&= \frac{1}{V[X]^{3/2}} \int_{0}^{1} \left( x-\frac{1}{2} \right)^3 dx \\
&= \frac{1}{V[X]^{3/2}} \left[ \frac{1}{4}\left( x-\frac{1}{2} \right)^4 \right]_{0}^{1} \\
&= \frac{1}{V[X]^{3/2}} \frac{1}{4} \left( \left(\frac{1}{2}\right)^4 – \left(-\frac{1}{2}\right)^4 \right) \\
&= 0 \\
\alpha_4 &= \frac{1}{V[X]^{2}} \int_{0}^{1} (x-E[X])^4f(x) dx \\
&= \frac{1}{V[X]^{2}} \int_{0}^{1} \left( x-\frac{1}{2} \right)^4 dx \\
&= \frac{1}{V[X]^{2}} \left[ \frac{1}{5}\left( x-\frac{1}{2} \right)^5 \right]_{0}^{1} \\
&= \frac{1}{V[X]^{2}} \frac{1}{5} \left( \left(\frac{1}{2}\right)^5 – \left(-\frac{1}{2}\right)^5 \right) \\
&= \frac{12^2}{5} \frac{2}{32} \\
&= \frac{9}{5} \\
\beta_4 &= \alpha_4-3 \\
&= \frac{9}{5}-3 \\
&= -\frac{6}{5}
\end{align}
$$

問題5.2の解答例

求める期待値を$E[X]$とすると、$E[X]$は下記のように求めることができる。
$$
\begin{align}
E[X] &= \frac{1}{13000000}(40000000 \times 7 + 10000000 \times 14 + 200000 \times 903 \\
&+ 10000000 \times 5 + 100000 \times 645 + 1000000 \times 130 + 140000 \times 130 + 10000 \times 1300 \\
&+ 1000 \times 26000 + 200 \times 1300000) \\
&= 89.40769…
\end{align}
$$

問題5.3の解答例

i)
$$
\begin{align}
P(X=2^k) = \frac{1}{2^k} (k=1, 2, 3, …)
\end{align}
$$
得られる額$X$の確率分布は上記のようになる。

ⅱ)
i)で求めた全ての$k$に対して、それぞれの期待値は$1$である。これにより、$E[X]=1+1+….1+…$となり、$\infty$に発散する。

問題5.4の解答例

$E[X-a]^2$は下記のように変形できる。
$$
\begin{align}
E[X-a]^2 = (E[X]-a)^2
\end{align}
$$
上記の値が最小になるのは$a=E[X]$のときであり、最小値は$E[X-E[X]]^2$となる。

問題5.5の解答例

正$n$面体から発生させた乱数の期待値$E[X]$と分散$V[X]=E[X^2]-E[X]^2$を求めるにあたっては、$E[X]$と$E[X^2]$がわかれば良い。以下それぞれについて計算する。
$$
\begin{align}
E[X] &= \frac{1}{n}(1+2+3+…) \\
&= \sum_{x=1} \frac{1}{n} x \\
&= \frac{n+1}{2} \\
E\left[X^2\right] &= \frac{1}{n} \left( 1^2 + 2^2 + 3^2 + … \right) \\
&= \frac{1}{n} \sum_{x=1} \frac{1}{n} x^2 \\
&= \frac{1}{6}(n+1)(2n+1) \\
V[X] &= E\left[X^2\right] – E[X]^2 \\
&= \frac{1}{6}(n+1)(2n+1) – \frac{1}{2^2}(n+1)^2 \\
&= \frac{1}{6}(2n^2+3n+1) – \frac{1}{4}(n^2+2n+1) \\
&= \frac{1}{12}(4n^2+6n+2) – \frac{1}{12}(3n^2+6n+3) \\
&= \frac{1}{12}(n^2-1)
\end{align}
$$
上記より、$\displaystyle E[X] = \frac{n+1}{2}$、$\displaystyle V[X] = \frac{1}{12}(n^2-1)$となる。

問題5.6の解答例

確率変数を$Y=X^2$のようにおくと、$Y$の累積分布関数$F_{Y}(y)$は下記のように表すことができる。
$$
\begin{align}
F_{Y}(y) &= P(Y \leq y) \\
&= P(X^2 \leq y) \\
&= P(0 \leq X \leq \sqrt{y}) \\
&= \sqrt{y}
\end{align}
$$
確率密度関数$f_{Y}(y)$は$f_{Y}(y)=F’_{Y}(y)$のように表せるので、下記のように計算できる。
$$
\begin{align}
f_{Y}(y) &= F’_{Y}(y) \\
&= \left( \sqrt{y} \right)’ \\
&= \frac{1}{2 \sqrt{y}}
\end{align}
$$
また、期待値$E[X^2]=E[Y]$と分散$V[X^2]=V[Y]$は下記のように計算できる。
$$
\begin{align}
E[X^2] &= \int_{0}^{1} x^2 dx \\
&= \left[ \frac{1}{3}x^3 \right]_{0}^{1} \\
&= \frac{1}{3} \\
V[X^2] &= E\left[ (X^2)^2 \right] – (E[X^2])^2 \\
&= \int_{0}^{1} x^4 dx – (E[X^2])^2 \\
&= \frac{1}{5} – \frac{1}{3^2} \\
&= \frac{4}{45}
\end{align}
$$

問題5.7の解答例

確率変数を$Y=X^2$のようにおくと、$Y$の累積分布関数$F_{Y}(y)$は下記のように表すことができる。
$$
\begin{align}
F_{Y}(y) &= P(Y \leq y) \\
&= P(X^2 \leq y) \\
&= P(-\sqrt{y} \leq X \leq \sqrt{y}) \\
&= 2P(0 \leq X \leq \sqrt{y}) \\
&= 2\left( \Phi(\sqrt{y})-\frac{1}{2} \right)
\end{align}
$$
上記では正規分布が左右対称であることから、$P(-\sqrt{y} \leq X \leq \sqrt{y})=2P(0 \leq X \leq \sqrt{y})$が成立することを用いている。また、$\Phi(z)$は標準正規分布の累積分布関数を表すとすると、$\displaystyle \Phi(z) = \int_{-\infty}^{z} \frac{1}{\sqrt{2 \pi}}e^{-x^2/2} dx$を満たす（121ページの(6.22）式の周辺を参照）。
また、$\displaystyle \Phi(0)=\frac{1}{2}$より$\displaystyle P(0 \leq X \leq \sqrt{y})=\Phi(\sqrt{y})-\frac{1}{2}$を導出することができる。

確率密度関数は$f_{Y}(y)=F’_{Y}(y)$のため下記のように計算できる。
$$
\begin{align}
f_{Y}(y) &= F’_{Y}(y) \\
&= 2\left( \Phi(\sqrt{y})-\frac{1}{2} \right)’ \\
&= 2\Phi'(\sqrt{y})(\sqrt{y})’ \\
&= 2 \times \frac{1}{\sqrt{2 \pi}}e^{-y/2} \frac{1}{2 \sqrt{y}} \\
&= \frac{1}{\sqrt{2 \pi y}}e^{-y/2}
\end{align}
$$

また、期待値$E[X^2]=E[Y]$と分散$V[X^2]=V[Y]$は下記のように計算できる。
$$
\begin{align}
E[X^2] &= \int_{-\infty}^{\infty} x^2 \times \frac{1}{\sqrt{2 \pi}}e^(-x^2/2) dx \\
&= \frac{1}{\sqrt{2 \pi}}\int_{-\infty}^{\infty} x^2 \times e^(-x^2/2) dx \\
&= 1 \\
V[X^2] &= E\left[ (X^2)^2 \right] – (E[X^2])^2 \\
&= E[X^4] – (E[X^2])^2 \\
&= 3 – 1 \\
&= 2
\end{align}
$$

問題5.8の解答例

二項分布を含む離散型の確率分布では累積分布関数が94ページの図5.9のような階段状となるが、$<$を用いるか$\leq$を用いるかで左連続か右連続かが変わる。$\leq$を用いる場合は右から連続(右連続)であり、$<$を用いる場合は左から連続(左連続)である。

まとめ

Ch.5で取り扱った確率変数は推測統計などを考えるにあたっての重要なトピックである一方で、なかなか具体的なイメージがつかみにくい概念だと思います。基本的なトピックが少々抽象的であるのに加え、モーメント母関数、チェビシェフの不等式、確率変数の変換など高度なトピックも含まれており、他の章と比較しても難易度が高い章であるという印象でした。

統計学入門 (基礎統計学Ⅰ)

「統計学実践ワークブック」演習問題中心第3章分布の特性値

投稿日: 2021-06-282022-08-03 投稿者: yoichi_t

統計検定準1級対応の公式テキストである「統計学実践ワークブック」を1章から順に演習問題を中心に解説していきます。
今回は第3章「分布の特性値」です。

日本統計学会公式認定統計検定準1級対応統計学実践ワークブック

3,080円(11/10 13:44時点)

重要ポイント

第2章では確率関数を扱いました。確率変数は確率関数に従います。確率関数は確率分布とも呼ばれますが、この確率分布を代表する値について本章では扱っています。

期待値と分散

確率分布を代表する値として、期待値、分散が代表的です。期待値と分散については、こちらを参照ください。

期待値（平均）と分散は、単峰で対称な確率分布でよく使われますが、偏った分布の場合には別の指標を使ったほうが分布の位置や形状を把握するのに有効な場合があります。期待値のような分布の位置を表す値として、モード(mode、最頻値)や中央値(median)があります。単峰で対称な場合、期待値とこれらの指標は一致します。一方、分布の広がりを表現する値として四分位範囲があります。

共分散

確率変数が一つではない場合、分布の広がりを表す値として共分散があります。共分散についてもこちらを参照ください。

その他、テキストには細かい指標があります。ご確認ください。

演習問題解説

演習問題の全文は掲載しません。テキストは各自で用意するか、以下の抜粋から想像してください。

問3.1

ある動物群の平均体重が60kg、標準偏差が12kgだった。これらの動物の1ヶ月後の平均体重は65kgだったけど変動係数は変化なかった。

(1) 変動係数を求めよ

「変動係数」とは、標準偏差を平均値で割ったものです。標準偏差は「ばらつきの大きさ」の意味があり、値が大きいデータに対しては標準偏差は大きめになりがちです。そのため、平均値で割って正規化したという意味があります。
定義通りに以下のように導出できます。

$$
\begin{eqnarray}
\frac{\sqrt(V[X])}{E[X]} = \frac{12}{60} = 0.2
\end{eqnarray}
$$

(2) 1ヶ月後の体重の標準偏差を求めよ

問題には「変動係数は変化がなかった」と書かれています。そこで、以下の方程式が成り立ちます。

$$
\begin{eqnarray}
\frac{S}{65} &=& \frac{1}{5} \\
S &=& \frac{65}{5} = 13
\end{eqnarray}
$$

これも「変動係数」という量の定義を覚えていれば容易に解けます。

問3.2

以下の(1)~(3)の問いを「加重平均」、「幾何平均」、「調和平均」を使って算出せよ。

(1) 片道48kmの道のりを行きは8km/h、帰りは12km/hで往復した。往復の平均速度を求めよ。

速度は単位時間あたりに進む距離なので、割合と見ることができます。速度の平均なので、調和平均を適用すると平均速度を算出できます。

なお$x_1, \cdots, x_n$の調和平均$m_h$とは、以下の通り割合の平均の逆数で定義されています。

$$
m_h^{-1} = \frac{1}{n} \sum^{n}_{i=1} \left( \frac{1}{x_i} \right)
$$

なお、$m^{-1}$はmの逆数を表しています。

では、問題の回答に移ります。

$$
\begin{eqnarray}
m_h = \left\{ \frac{1}{2} \left( \frac{1}{8} + \frac{1}{12}\right) \right\}^{-1} = \frac{48}{5} = 9.6 [km/h]
\end{eqnarray}
$$

(2) 3種類の定食があり。それぞれの売り上げ数から、定食１食あたりの平均金額を求める。

テキストでは、A定食は600円、B定食は700円、C定食は500円であり、それぞれ、500食、700食、800食の売り上げがあったとされています。
愚直に計算しても良さそうですが、題意から「加重平均」を使って計算します。

重みをどう考えるかですが、それぞれの定食の売り上げ割合を重みとします。全体で2,000食売れていますので、A, B, Cの定食のそれぞれの売上数の割合は以下のようになります。

A定食: $\frac{500}{2000} = \frac{1}{4}$
B定食: $\frac{700}{2000} = \frac{7}{20}$
C定食: $\frac{800}{2000} = \frac{2}{5}$

ということで以下のように平均金額を算出できます。

$$
m_w = 600\frac{1}{4} + 700\frac{7}{20} + 500 \frac{2}{5} = 150+245+200 = 595 [円]
$$

(3) 消費者物価指数の平均伸び率を求めよ。

消費者物価指数の4年間の伸び率が1.15, 0.98, 1.03, 0.99であったとされています。

「伸び率」は今年の値と翌年の値を割ったものです（前年の値$x_0$に伸び率$\alpha_0$をかければ翌年の値$x_1$になる, $x_1 = \alpha_0 x_0$、下図も参照）。そのため、4年後の値は初年度に各年の伸び率の総積をかけることで算出されます。

$$
\begin{eqnarray}
x_4 &=& \alpha_3 x_3 \\
x_3 &=& \alpha_2 x_2 \\
&\cdots&
\end{eqnarray}
$$

なので、初年度と4年後の関係は$x_4 = (\alpha_0 \cdot \alpha_1 \cdot \alpha_2 \cdot \alpha_3)x_0$となります。

「平均伸び率」を利用すると、上記の総積の部分が平均伸び率$\alpha$の4乗になります。ということで、幾何平均を用いれば良いということがわかります。

$$
\begin{eqnarray}
m_g = \sqrt[4]{\prod \alpha_i} = \sqrt[4]{1.15 \cdot 0.98 \cdot 1.03 \cdot 0.99} \simeq 1.04
\end{eqnarray}
$$

問3.3

以下の厚さのパンとハムを使ってサンドイッチを作る。

パンの厚さX[mm]: E[X]=10, V[X] = $0.5^2$
ハムの厚さY[mm]: E[Y]=3, V[Y] = $0.4^2$

サンドイッチを作る際に上記の材料を2通りの方法で作ることが考えられる。

[方法1]
- 予めパンを二つに切断して貯めておく
- 溜め込んだパン片をランダムに抽出して、ランダムに選んだハムを挟んでサンドイッチを作る
[方法2]
- パンを選びそれを二つに切断する
- 二つに切断したパンを使ってサンドイッチを作成する

二つの方法の違いですが、方法1は、サンドイッチを構成する2枚のパン片の厚さが異なります。一方、方法2では同じパンを使ってサンドイッチを作るので必ず2枚のパン片は同じ厚さになります。（下図参照）

それぞれの方法で作成するサンドイッチの厚さの分散を求めよ

問題に書かれている二つのサンドイッチ作成方法から、方法1ではパン片が平均$E[X]$で分散$V[X]$にそれぞれ独立に従うことがわかります。厚さ$(E[X], V[X])$に従うパンを分けて溜め込んでおき、そこからランダムに選択しますので。一方、方法2ではパンを二つに分けてサンドイッチを作成するので、2つのパン片は同じ厚さになります。

方法1でのサンドイッチの厚さの分散は、次のように計算できます。

$$
\begin{eqnarray}
V[x_1 + x_2 + y] &=& V[x] + V[x] + V[y] \\
&=& 0.5^2 + 0.5^2 + 0.4^2 = 0.66
\end{eqnarray}
$$

一つ目の式変形は、$x_1, x_2$がそれぞれ独立であることを利用しています。

次に方法2での分散です。

$$
\begin{eqnarray}
V[2x + y] &=& 2^2 \cdot V[x] + V[y] \\
&=& 4 \cdot 0.5^2 + 0.4^2 = 1.16
\end{eqnarray}
$$

式変形は分散の性質を利用しています。

参考文献

ワークブック以外の参考資料として以下のものがおすすめです。

松原ら, 統計学入門, 1991, 東京大学出版会

Ch.7 「多次元の確率分布」の章末問題の解答例　〜基礎統計学Ⅰ 統計学入門(東京大学出版会)〜

投稿日: 2021-06-252022-11-09 投稿者: lib-arts

当記事は基礎統計学Ⅰ 統計学入門(東京大学出版会)」の読解サポートにあたってChapter.$7$の「多次元の確率分布」の章末問題の解説について行います。
基本的には書籍の購入者向けの解説なので、まだ入手されていない方は購入の上ご確認ください。また、解説はあくまでサイト運営者が独自に作成したものであり、書籍の公式ページではないことにご注意ください。

統計学入門 (基礎統計学Ⅰ)

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

章末の演習問題について

問題7.1の解答例

i)
$V(X+Y)=E[(X+Y)^2]-E[X+Y]^2$を用いて導出する。
$$
\begin{align}
V(X+Y) &= E[(X+Y)^2] – E[X+Y]^2 \\
&= E[X^2+Y^2+2XY] – (E[X]+E[Y])^2 \\
&= E[X^2]+E[Y^2]+2E[XY] – (E[X]^2+E[Y]^2+2E[X]E[Y]) \\
&= (E[X^2]-E[X]^2) + (E[Y^2]-E[Y]^2) + 2(E[XY]-E[X]E[Y]) \\
&= V(X) + V(Y) + 2\mathrm{Cov}(X,Y)
\end{align}
$$

上記より、$V(X+Y) = V(X) + V(Y) + 2\mathrm{Cov}(X,Y)$が成立する。

ⅱ)
i)と同様に$V(aX+bY)=E[(aX+bY)^2]-E[aX+bY]^2$を用いて導出する。
$$
\begin{align}
V(aX+bY) &= E[(aX+bY)^2] – E[aX+bY]^2 \\
&= E[a^2X^2+b^2Y^2+2abXY] – (E[aX]+E[bY])^2 \\
&= a^2E[X^2]+b^2E[Y^2]+2abE[XY] – (a^2E[X]^2+b^2E[Y]^2+2abE[X]E[Y]) \\
&= a^2(E[X^2]-E[X]^2) + b^2(E[Y^2]-E[Y]^2) + 2ab(E[XY]-E[X]E[Y]) \\
&= a^2V(X) + b^2V(Y) + 2ab\mathrm{Cov}(X,Y)
\end{align}
$$

上記より、$V(aX+bY) = a^2V(X) + b^2V(Y) + 2ab\mathrm{Cov}(X,Y)$が成立する。

問題7.2の解答例

i)
$R_p = xR_1+(1-x)R_2$の期待値E(R_p)と分散V(R_p)は下記のようになる。
$$
\begin{align}
E(R_p) &= E(xR_1+(1-x)R_2) \\
&= xE(R_1) + (1-x)E(R_2) \\
&= xe_1 + (1-x)e_2 \\
V(R_p) &= V(xR_1+(1-x)R_2) \\
&= x^2V(R_1) + (1-x)^2V(R_2) + 2x(1-x)\mathrm{Cov}(R_1,R_2) \\
&= x^2 \sigma_1^2 + (1-x)^2 \sigma_2^2 + 2x(1-x)\rho\sigma_1\sigma_2 \\
&= (\sigma_1^2-2\rho\sigma_1\sigma_2+\sigma_2^2)x^2 -2(\sigma_2^2-\rho\sigma_1\sigma_2)+\sigma_2^2
\end{align}
$$

問題7.3の解答例

確率変数$X$の取りうる値は$X=0,1,2,3$、確率変数$Y$の取りうる値は$Y=1,2$のため、$2 \times 4 = 8$パターンについて計算することで同時確率分布を求めることができる。
$$
\begin{align}
P(X=0,Y=1) &= \frac{1}{2^3} = \frac{1}{8} \\
P(X=1,Y=1) &= 0 \\
P(X=2,Y=1) &= 0 \\
P(X=3,Y=1) &= \frac{1}{2^3} = \frac{1}{8} \\
P(X=0,Y=2) &= 0 \\
P(X=1,Y=2) &= \frac{3}{2^3} = \frac{3}{8} \\
P(X=2,Y=2) &= \frac{3}{2^3} = \frac{3}{8} \\
P(X=3,Y=2) &= 0
\end{align}
$$

ここで$\displaystyle P(X=0)=\frac{1}{8}$、$\displaystyle P(Y=1)=\frac{1}{2}$であるが、$\displaystyle P(X=0, Y=1) = \frac{1}{8} \neq \frac{1}{32} = P(X=0)P(Y=1)$なので、$X$と$Y$は独立ではない。
また、$X$、$Y$、$XY$の確率分布を作成すると下記のようになる。
$$
\begin{align}
P(X=0) &= \frac{1}{8} \\
P(X=1) &= \frac{3}{8} \\
P(X=2) &= \frac{3}{8} \\
P(X=3) &= \frac{1}{8} \\
P(Y=1) &= \frac{1}{4} \\
P(Y=2) &= \frac{3}{4} \\
P(XY=0) &= P(X=0,Y=1)+P(X=0,Y=2) = \frac{1}{8} \\
P(XY=1) &= P(X=1,Y=1) = 0 \\
P(XY=2) &= P(X=2,Y=1)+P(X=1,Y=2) = \frac{3}{8} \\
P(XY=3) &= P(X=3,Y=1) = \frac{1}{8} \\
P(XY=4) &= P(X=2,Y=2) = \frac{3}{8} \\
P(XY=6) &= P(X=3,Y=2) = 0
\end{align}
$$

上記より、$E[X]$、$E[Y]$、$E[XY]$は下記のように計算できる。
$$
\begin{align}
E[X] &= 0 \cdot \frac{1}{8} + 1 \cdot \frac{3}{8} + 2 \cdot \frac{3}{8} + 3 \cdot \frac{1}{8} \\
&= \frac{3}{8} + \frac{6}{8} + \frac{3}{8} \\
&= \frac{3}{2} \\
E[Y] &= 1 \cdot \frac{1}{4} + 2 \cdot \frac{3}{4} \\
&= \frac{7}{4} \\
E[XY] &= 0 \cdot \frac{1}{8} + 1 \cdot 0 + 2 \cdot \frac{3}{8} + 3 \cdot \frac{1}{8} + 4 \cdot \frac{3}{8} + 6 \cdot 0 \\
&= \frac{21}{8}
\end{align}
$$

上記より$\displaystyle E[XY] = \frac{21}{8} = \frac{3}{2} \cdot \frac{7}{4} = E[X]E[Y]$が成立するので$\mathrm{Cov}(X,Y)=0$となり無相関である。

問題7.4の解答例

(Ⅰ)の方法で$m_A$と$m_b$を測った際の計測値をそれぞれ$X_A$、$X_B$とする。このとき、天秤の測定誤差が$\sigma^2$であることから、$V(X_A)=V(X_B)=\sigma^2$となる。
また、(Ⅱ)の方法で計測した重さの和を$Y$、差を$Z$とすると、$m_A$と$m_B$の推定値はそれぞれ$\displaystyle \frac{Y+Z}{2}$、$\displaystyle \frac{Y-Z}{2}$となる。このとき、$\displaystyle V\left( \frac{Y+Z}{2} \right) = V\left( \frac{Y-Z}{2} \right) = \frac{1}{2}\sigma^2$となり、(Ⅱ)の方が優れた方法であるといえる。

問題7.5の解答例

$$
\large
\begin{align}
V[U] &= V[aX+b] \\
&= a^2V[X] \\
V[V] &= V[cY+d] \\
&= c^2V[Y] \\
\mathrm{Cov}(U,V) &= \mathrm{Cov}(aX+b, cY+d) \\
&= ac\mathrm{Cov}(X,Y)
\end{align}
$$

上記を利用して、下記のように導出することができる。
$$
\large
\begin{align}
\rho_{UV} &= \frac{\mathrm{Cov}(U,V)}{\sqrt{V[U]V[V]}} \\
&= \frac{\mathrm{Cov}(aX+b, cY+d)}{\sqrt{V[aX+b]V[cY+d]}} \\
&= \frac{\cancel{ac}\mathrm{Cov}(X,Y)}{\cancel{ac}\sqrt{V[X]V[Y]}} \\
&= \frac{\mathrm{Cov}(X,Y)}{\sqrt{V[X]V[Y]}} \\
&= \rho_{XY}
\end{align}
$$

問題7.6の解答例

i)
$$
\large
\begin{align}
X_{1}, X_{2} & \sim N(0,1) \quad i.i.d., \\
Y_{1} &= aX_{1} + bX_{2} \\
Y_{2} &= cX_{1} + dX_{2}
\end{align}
$$

$7.3$節の議論により、上記のように定義した$Y_{1}, Y_{2}$の相関係数$\rho$は下記のように表すことができる。
$$
\large
\begin{align}
\rho &= \frac{\mathrm{Cov}(Y_1,Y_2)}{\sqrt{V[Y_1]}\sqrt{V[Y_2]}} \\
&= \frac{ac+bd}{\sqrt{a^2+b^2}\sqrt{c^2+d^2}}
\end{align}
$$

上記に$a=1, b=0, c=c, d=1$を代入し、$\rho=0.5$を$c$に関して解く。
$$
\large
\begin{align}
\rho &= 0 \\
\frac{1 \times c + 0 \times 1}{\sqrt{1^2+0^2}\sqrt{c^2+1^2}} &= 0.5 \\
c &= 0.5 \sqrt{c^2+1^2} \\
c^2 &= 0.25(c^2+1) \\
c^2 &= \frac{1}{3} \\
c &= \frac{1}{\sqrt{3}}
\end{align}
$$

ⅱ)
i)と同様に考え、$c^2 = \rho^2(c^2+1)$を$c$に関して解く。
$$
\large
\begin{align}
c^2 &= \rho^2(c^2+1) \\
(1-\rho^2)c^2 &= \rho^2 \\
c &= \frac{\rho}{\sqrt{1-\rho^2}}
\end{align}
$$

ⅲ)

問題7.7の解答例

i)
並列のシステムは全てが故障するまで継続すると考えられるため、$Y=\max(X_1,X_2)$が成立する。このとき確率変数$Y$の確率密度関数を$f_Y(y)$、累積分布関数を$F_Y(y)$とするとそれぞれ下記のように求めることができる。
$$
\large
\begin{align}
F_Y(y) &= F_Y(Y \leq y) \\
&= F_Y(X_1 \leq y, X_2 \leq y) \\
&= F_Y(X_1 \leq y)F_Y(X_2 \leq y) \\
&= (1-e^{-\lambda y})^2
\end{align}
$$
$$
\large
\begin{align}
f_Y(y) &= F_Y'(Y \leq y) \\
&= ( (1-e^{-\lambda y})^2 )’ \\
&= 2\lambda e^{-\lambda y}(1-e^{-\lambda y})
\end{align}
$$

ⅱ)
直列のシステムは一つが故障するまで継続すると考えられるため、$Y=\min(X_1,X_2)$が成立する。このとき確率変数$Y$の確率密度関数を$f_Y(y)$、累積分布関数を$F_Y(y)$とするとそれぞれ下記のように求めることができる。
$$
\large
\begin{align}
F_Y(y) &= F_Y(Y \leq y) \\
&= 1 – F_Y(Y > y) \\
&= 1 – F_Y(X_1 > y)F_Y(X_2 > y) \\
&= 1 – (e^{-\lambda y})^2 \\
&= 1 – (e^{-2\lambda y})
\end{align}
$$
$$
\large
\begin{align}
f_Y(y) &= F_Y'(Y \leq y) \\
&= ( 1 – (e^{-2\lambda y}) )’ \\
&= 2\lambda(e^{-2\lambda y})
\end{align}
$$

問題7.8の解答例

問題7.9の解答例

i) 二項分布

$X_1+X_2=x$を考えた際に$X_1=k$とすると、$X_2=x-k$となる。このときの確率分布を$P(X_1+X_2=x|m,n,p)$を考えると、下記のように計算することができる。（$X_1$に対応する試行回数を$m$、$X_2$に対応する試行回数を$n$とする）
$$
\begin{align}
P(X_1+X_2=x|m,n,p) &= \sum_{k=0}^{x} P(X_1=k|m,p)P(X_2=x-k|n,p) \\
&= \sum_{k=0}^{x} {}_m C_k p^{k}(1-p)^{m-k} {}_n C_{x-k} p^{x-k}(1-p)^{n-(x-k)} \\
&= p^{k}(1-p)^{m-k}p^{x-k}(1-p)^{n-(x-k)} \sum_{k=0}^{x} {}_m C_k {}_n C_{x-k} \\
&= p^{k+x-k}(1-p)^{m-k+n-(x-k)} \sum_{k=0}^{x} {}_m C_k {}_n C_{x-k} \\
&= p^{x}(1-p)^{n+m-x} \sum_{k=0}^{x} {}_n C_k {}_m C_{x-k}
\end{align}
$$
このときヴァンデルモンドの畳み込みより、$\displaystyle {}_{m+n} C_x = \sum_{k=0}^{x} {}_m C_k {}_n C_{x-k}$が成立するので、$P(X_1+X_2=x|m,n,p)$は下記のように変形できる。
$$
\begin{align}
P(X_1+X_2=x|m,n,p) &= p^{x}(1-p)^{n+m-x} \sum_{k=0}^{x} {}_n C_k {}_m C_{x-k} \\
&= p^{x}(1-p)^{n+m-x} {}_{m+n} C_x \\
&= {}_{m+n} C_x p^{x}(1-p)^{n+m-x}
\end{align}
$$
上記により二項分布の再生性を示すことができた。（ヴァンデルモンドの畳み込みの式は、$m+n$個の中から$x$個取り出す際に、$m$個から$k$個、$n$個から$x-k$個取り出す方法を考えることで示すことができる）

ⅱ) ポアソン分布

$X_1+X_2=x$を考えた際に$X_1=k$とすると、$X_2=x-k$となる。このときの確率分布を$P(X_1+X_2=x|\lambda_1,\lambda_2)$を考えると、下記のように計算することができる。（$X_1$に対応するポアソン分布のパラメータを$\lambda_1$、$X_2$に対応するポアソン分布のパラメータを$\lambda_2$とする）
$$
\begin{align}
P(X_1+X_2=x|\lambda_1,\lambda_2) &= \sum_{k=0}^{x} P(X_1=k|\lambda_1)P(X_2=x-k|\lambda_2) \\
&= \sum_{k=0}^{x} \frac{\lambda_1^k e^{-\lambda_1}}{k!} \frac{\lambda_2^{x-k} e^{-\lambda_2}}{(x-k)!} \\
&= e^{-\lambda_1}e^{-\lambda_2} \sum_{k=0}^{x} \frac{1}{x!} \frac{x!}{k!(x-k)!} \lambda_1^k \lambda_2^{x-k} \\
&= \frac{e^{-(\lambda_1+\lambda_2})}{x!} \sum_{k=0}^{x} {}_x C_k \lambda_1^k \lambda_2^{x-k} \\
&= \frac{e^{-(\lambda_1+\lambda_2})}{x!} (\lambda_1+\lambda_2)^x \\
&= \frac{(\lambda_1+\lambda_2)^x e^{-(\lambda_1+\lambda_2)}}{x!}
\end{align}
$$
上記のたたみこみ(convolution)演算によってポアソン分布の再生性が示された。

ⅲ) 正規分布

$X_1+X_2=x$を考えた際に$X_1=k$とすると、$X_2=x-k$となる。このときの確率分布を$P(X_1+X_2=x|\mu_1,\mu_2,\sigma_1^2,\sigma_2^2)$を考えると、下記のように計算することができる。（$X_1$に対応する正規分布のパラメータを$\mu_1$と$\sigma_1^2$、$X_2$に対応する正規分布のパラメータを$\mu_2$と$\sigma_2^2$とする）
$$
\begin{align}
P(X_1 &+ X_2=x|\mu_1,\mu_2,\sigma_1^2,\sigma_2^2) \\
&= \sum_{k=0}^{x} P\left( X_1=k|\mu_1,\sigma_1^2\right) P\left(X_2=x-k|\mu_2,\sigma_2^2\right) \\
&= \sum_{k=0}^{x} \frac{1}{\sqrt{2 \pi \sigma_1^2}}\exp \left( -\frac{(k-\mu_1)^2}{2\sigma_1^2} \right) \times \frac{1} {\sqrt{2 \pi \sigma_2^2}}\exp \left( -\frac{((x-k)-\mu_2)^2}{2\sigma_2^2} \right)
\end{align}
$$
上記に対してガウス積分などを用いて導出を行う。

まとめ

Chapter.$7$「多次元の確率分布」について確認する内容でした。他のChapterに比較しても高度な話題が多いため、徐々に理解すれば十分な内容だと思います。

統計学入門 (基礎統計学Ⅰ)

Ch.4 「最尤法」の章末問題の解答例　〜基礎統計学Ⅲ 自然科学の統計学(東京大学出版会)〜

投稿日: 2021-06-212022-08-07 投稿者: lib-arts

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

自然科学の統計学 (基礎統計学)

章末の演習問題について

問題4.1の解答例

下記を実行すればよい。

import numpy as np
import matplotlib.pyplot as plt

x = np.array([3.03, 3.09, 3.15, 3.20, 3.25, 3.30, 3.34, 3.37, 3.41, 3.45, 3.48])
p = np.array([0., 0., 0.11, 0.28, 0.48, 0.71, 0.83, 0.91, 0.95, 0.98, 1.])

plt.plot(np.log(x),p)
plt.show()

・実行結果

問題4.2の解答例

i)
$i$回目の実験における表が出るまでの回数を$x_i$とする。この確率を$P(X=x_i|p)$とすると、$P(X=x_i|p)=f(x_i)$より下記のように表せる。
$$
\begin{align}
P(X=x_i|p) &= f(x_i) \\
&= p(1-p)^{x_i}
\end{align}
$$
このとき、$n$回の実験に関する同時確率(joint probability)は$P(X=x_1,X=x_2,…,X=x_n|p)$と表すことができ、$X_1$〜$X_n$が独立であるので下記のように式変形ができる。
$$
\begin{align}
P(X=x_1,X=x_2,…,X=x_n|p) &= P(X=x_1|p)P(X=x_2|p)…P(X=x_n|p) \\
&= f(x_1)f(x_2)…f(x_n) \\
&= p^n(1-p)^{\sum_{i=1}^{n} x_i-1} \\
&= p^n(1-p)^{T-n}
\end{align}
$$
上記では$\displaystyle T = \sum_{i=1}^{n} x_i$を利用した。また、このときの尤度を$L(p)$とすると、$L(p)=P(X=x_1,X=x_2,…,X=x_n|p)$なので、$L(p)=p^n(1-p)^{T-n}$となる。以下、対数尤度の微分$\displaystyle \frac{\delta \log{L(p)}}{\delta p} = 0$を満たす$p$を求める。
$$
\begin{align}
\log{L(p)} &= \log{p^n(1-p)^{T-n}} \\
&= n\log{p}+(T-n)\log{(1-p)}
\end{align}
$$
$$
\begin{align}
\frac{\delta \log{L(p)}}{\delta p} &= \frac{n}{p}-\frac{T-n}{1-p} = 0 \\
\frac{n}{p} &= \frac{T-n}{1-p} \\
n(1-p) &= (T-n)p \\
n-np &= Tp-np \\
n &= Tp \\
p &= \frac{n}{T}
\end{align}
$$
上記より、$\displaystyle \tilde{p} = \frac{n}{T}$を示すことができる。
また、ここでは$n=100$、$T=192$より、$\displaystyle \tilde{p} = \frac{100}{192} = 0.52083…$となる。

問題4.3の解答例

i)
平均$\mu=0$、分散$\sigma^2$の正規分布は下記のように書ける。
$$
\begin{align}
P(x|\mu=0,\sigma^2) = \frac{1}{\sqrt{2\pi\sigma^2}} \exp\left( -\frac{x^2}{2\sigma^2} \right)
\end{align}
$$
このとき、$P(x_1,x_2,…,x_n|0,\sigma^2)$の同時確率は下記のように計算できる。
$$
\begin{align}
P(x_1,x_2,…,x_n|0,\sigma^2) &= P(x_1|0,\sigma^2)P(x_2|0,\sigma^2)…P(x_n|0,\sigma^2) \\
&= \prod_{i=1}^{n} \frac{1}{\sqrt{2\pi\sigma^2}}\exp\left( -\frac{x_i^2}{2\sigma^2} \right)
\end{align}
$$
上記の同時確率は$\sigma^2$に着目すると尤度と見なせるので、$L(\sigma^2)=P(x_1,x_2,…,x_n|0,\sigma^2)$が成立する。このとき積の形を和の形に直すにあたって、$L(\sigma^2)$の対数の$\log{L(\sigma^2)}$を下記で計算する。
$$
\begin{align}
\log{L(\sigma^2)} &= \log\left( \prod_{i=1}^{n} \frac{1}{\sqrt{2\pi\sigma^2}}\exp\left( -\frac{x_i^2}{2\sigma^2} \right) \right) \\
&= \log\left(\prod_{i=1}^{n} \frac{1}{\sqrt{2\pi\sigma^2}} \right) + \log\left(\prod_{i=1}^{n} \exp\left( -\frac{x_i^2}{2\sigma^2} \right) \right) \\
&= \log\left( \frac{1}{\sqrt{2\pi\sigma^2}} \right)^n + \log\left(\exp\left( \sum_{i=1}^{n} -\frac{x_i^2}{2\sigma^2} \right)\right) \\
&= \log\left( 2\pi\sigma^2 \right)^{-n/2} – \sum_{i=1}^{n} \frac{x_i^2}{2\sigma^2} \\
&= – \frac{n}{2}\log\left( 2\pi \right) – \frac{n}{2}\log\left( \sigma^2 \right) – \sum_{i=1}^{n} \frac{x_i^2}{2\sigma^2}
\end{align}
$$
次に上記で求めた$\log{L(\sigma^2)}$を$\sigma^2$で微分する。
$$
\begin{align}
\frac{\delta \log{L(\sigma^2)}}{\delta \sigma^2} = – \frac{n}{2\sigma^2} + \sum_{i=1}^{n} \frac{x_i^2}{2(\sigma^2)^2}
\end{align}
$$
ここで$\displaystyle \frac{\delta \log{L(\sigma^2)}}{\delta \sigma^2}=0$を考える。
$$
\begin{align}
& \frac{\delta \log{L(\sigma^2)}}{\delta \sigma^2} = 0 \\
& – \frac{n}{2\sigma^2} + \sum_{i=1}^{n} \frac{x_i^2}{2(\sigma^2)^2} = 0 \\
& \sum_{i=1}^{n} \frac{x_i^2}{2(\sigma^2)^2} = \frac{n}{2\sigma^2} \\
& \sum_{i=1}^{n} x_i^2 = \frac{n \times 2(\sigma^2)^2}{2\sigma^2} \\
& \sum_{i=1}^{n} x_i^2 = n\sigma^2 \\
& \sigma^2 = \frac{\sum_{i=1}^{n} x_i^2}{n}
\end{align}
$$
上記より、$\displaystyle \tilde{\sigma} = \sqrt{\frac{\sum_{i=1}^{n} x_i^2}{n}}$を求めることができる。

問題4.4の解答例

i)
指数分布の確率密度関数を$f(t)=P(t|\lambda)$としたとき、下記のように表すことができる。
$$
\begin{align}
P(t|\lambda) = \lambda e^{-\lambda t} \\
(t \geq 0)
\end{align}
$$
この際に観測されたサンプルを$t_1$、$t_2$、…、$t_{10}$、尤度を$L(\lambda)$とすると、尤度はサンプルが観測された同時確率をパラメータの関数と読み替えたものなので下記が成立する。
$$
\begin{align}
L(\lambda) &= P(t_1|\lambda)P(t_2|\lambda)…P(t_{10}|\lambda) \\
&= \prod_{i=1}^{10} \lambda e^{-\lambda t_i}
\end{align}
$$
ここで上記の対数をとった対数尤度の$\log{L(\lambda)}$について計算すると下記のようになる。
$$
\begin{align}
\log{L(\lambda)} &= \log\left( \prod_{i=1}^{10} \lambda e^{-\lambda t_i} \right) \\
&= \log\left( \lambda^{10} \prod_{i=1}^{10}e^{-\lambda t_i} \right) \\
&= \log\left( \lambda^{10} e^{\sum_{i=1}^{10} -\lambda t_i} \right) \\
&= 10\log{\lambda} – \sum_{i=1}^{10} \lambda t_i \\
&= 10\log{\lambda} – \lambda \sum_{i=1}^{10} t_i
\end{align}
$$
上記を$\lambda$について微分し、$0$になる$\lambda$が$\tilde{\lambda}$となる。

$$
\begin{align}
& \frac{\delta \log{L(\lambda)}}{\delta \lambda} = 0 \\
& \frac{10}{\lambda} – \sum_{i=1}^{10} t_i = 0 \\
& \frac{10}{\lambda} = \sum_{i=1}^{10} t_i \\
& \lambda = \frac{10}{\sum_{i=1}^{10} t_i}
\end{align}
$$

上記より、$\displaystyle \tilde{\lambda} = \frac{10}{\sum_{i=1}^{10} t_i}$が成立する。これを元に下記のように$\tilde{\lambda}$が推定できる。
$$
\begin{align}
\tilde{\lambda} &= \frac{10}{523+971+878+130+696+201+1643+491+665+431} \\
&= \frac{10}{6629} \\
&= 0.0015085…
\end{align}
$$

ⅱ)
指数分布の累積分布関数を$F(t)=F(X \leq t)$とすると、下記のように表すことができる。
$$
\begin{align}
F(X \leq t) &= f(t) \\
&= 1-e^{-\lambda t} \\
F(X \geq t) &= 1-F(X \leq t) \\
&= 1-(1-e^{-\lambda t}) \\
&= e^{-\lambda t}
\end{align}
$$
よって、$F(X \geq 1000)$は下記のように表すことができる。
$$
\begin{align}
F(X \geq 1000) &= e^{-(10 \times 1000)/6629} \\
&= e^{-(10000)/6629} \\
&= 0.221236…
\end{align}
$$

問題4.5の解答例

$\tilde{\theta}=\max{x_i}$が最尤推定量となる。

まとめ

Chapter.$4$で取り扱った最尤法は得られたサンプルの背後に仮定した確率分布のパラメータの推定を行う手法です。同時確率を尤度と読み替えるところがなかなかイメージがつかないかもしれませんが、最重要項目なのでこの辺の理解を中心に理解を進めていくと良いと思います。

自然科学の統計学 (基礎統計学)

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

Ch.1 「確率の基礎」の章末問題の解答例〜基礎統計学Ⅲ 自然科学の統計学(東京大学出版会)〜

投稿日: 2021-06-212022-08-07 投稿者: lib-arts

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

自然科学の統計学 (基礎統計学)

章末の演習問題について

問題1.1の解答例

二項分布$\mathrm{Binom}(n,p)$の確率変数を$X$としたとき、$X=k$となる確率$P(X=k|n,p)$は下記を用いて表すことができる。
$$
\begin{align}
P(X=k|n,p) = {}_n C_k p^k (1-p)^k
\end{align}
$$
上記に基づいて、下記のように期待値$E[X]$を計算することができる。
$$
\begin{align}
E[X] &= \sum_{k=0}^{n} k {}_n C_k p^k (1-p)^k \\
&= \sum_{k=1}^{n} k {}_n C_k p^k (1-p)^k \\
&= \sum_{k=1}^{n} k \frac{n(n-1)…(n-k+1)}{k!} p^k (1-p)^{n-k} \\
&= \sum_{k=1}^{n} np \frac{(n-1)…(n-k+1)}{(k-1)!} p^{k-1} (1-p)^{n-k} \\
&= np \sum_{k=1}^{n} {}_{n-1} C_{k-1} p^{k-1} (1-p)^{n-k} \\
&= np \sum_{k=0}^{n-1} {}_{n-1} C_{k} p^{k} (1-p)^{n-(k+1)} \\
&= np \sum_{k=0}^{n-1} {}_{n-1} C_{k} p^{k} (1-p)^{(n-1)-k)} \\
&= np (1+(1-p))^{n-1} \\
&= np
\end{align}
$$
また、$V[X]$は下記のように考えることができる。（確率母関数を用いた計算と同様な計算式となることは押さえておく方が良い）
$$
\begin{align}
V[X] &= E[X^2] – (E[X])^2 \\
&= E[X(X-1)+X] – (E[X])^2 \\
&= E[X(X-1)] + E[X] – (E[X])^2
\end{align}
$$
$E[X]$はすでに求めたので、以下では$E[X(X-1)]$について求めることを考える。
$$
\begin{align}
E[X] &= \sum_{k=0}^{n} k(k-1) {}_n C_k p^k (1-p)^k \\
&= \sum_{k=2}^{n} k(k-1) {}_n C_k p^k (1-p)^k \\
&= \sum_{k=2}^{n} k(k-1) \frac{n(n-1)…(n-k+1)}{k!} p^k (1-p)^{n-k} \\
&= \sum_{k=2}^{n} n(n-1)p^2 \frac{(n-2)…(n-k+1)}{(k-2)!} p^{k-2} (1-p)^{n-k} \\
&= n(n-1)p^2 \sum_{k=2}^{n} {}_{n-2} C_{k-2} p^{k-2} (1-p)^{n-k} \\
&= n(n-1)p^2 \sum_{k=0}^{n-2} {}_{n-2} C_{k} p^{k} (1-p)^{n-(k+2)} \\
&= n(n-1)p^2 \sum_{k=0}^{n-2} {}_{n-2} C_{k} p^{k} (1-p)^{(n-2)-k)} \\
&= n(n-1)p^2 (1+(1-p))^{n-2} \\
&= n(n-1)p^2
\end{align}
$$
よって、$V[X]$は下記のように求めることができる。
$$
\begin{align}
V[X] &= E[X(X-1)] + E[X] – (E[X])^2 \\
&= n(n-1)p^2 + np – n^2p^2 \\
&= np((n-1)p+1-np) \\
&= np(np-p+1-np) \\
&= np(1-p)
\end{align}
$$
ここまでの議論より、期待値$E[X]=np$、分散$V[X]=np(1-p)$が導出できる。

問題1.2の解答例

確率$p$で事象$A$が観測されるベルヌーイ試行を考えた際に、事象$A$が最初に観測されるまでに要した試行回数を確率変数$X$で定めると、確率$P(X=r|p)$は下記のように表すことができる。
$$
\begin{align}
P(X=r|p) = p(1-p)^{r-1}
\end{align}
$$
このとき上記は幾何分布(geometric distribution)と呼ばれる。ここで幾何分布の期待値を$E[X]$とおくと、$E[X]$は下記のように求めることができる。
$$
\begin{align}
E[X] &= \sum_{r=1}^{\infty} rP(X=r|p) \\
&= \sum_{r=1}^{\infty} rp(1-p)^{r-1} \\
&= p\sum_{r=1}^{\infty} r(1-p)^{r-1}
\end{align}
$$
このとき$\displaystyle \sum_{r=1}^{\infty} r(1-p)^{r-1}$を考えるにあたって、以下のようにマクローリン展開を考える。
$$
\begin{align}
\frac{1}{1-x} &= \sum_{n=0}^{\infty} x^n
\end{align}
$$
上記の両辺を微分すると下記が得られる。
$$
\begin{align}
\frac{1}{(1-x)^2} &= \sum_{n=1}^{\infty} nx^{n-1}…(1)
\end{align}
$$
(1)式の右辺を、$n=r$、$x=(1-p)$で置き換えると$\displaystyle \sum_{r=1}^{\infty} r(1-p)^{r-1}$が得られる。よって、下記が成立する。
$$
\begin{align}
\sum_{r=1}^{\infty} r(1-p)^{r-1} &= \frac{1}{(1-(1-p))^2} \\
&= \frac{1}{p^2}
\end{align}
$$
したがって、$E[X]$は下記のように求められる。
$$
\begin{align}
E[X] &= p\sum_{r=1}^{\infty} r(1-p)^{r-1} \\
&= p \frac{1}{p^2} \\
&= \frac{1}{p}
\end{align}
$$

問題1.3の解答例

事象が$k$回起こるまでの回数を確率変数$X$としたとき、$X$は負の二項分布に従う。ここで$P(X=n|p)$は下記のように表すことができる。
$$
\begin{align}
P(X=n|p) &= {}_{n-1} C_{k-1} p^{k-1} (1-p)^{n-k} p \\
&= {}_{n-1} C_{k-1} p^k (1-p)^{n-k}
\end{align}
$$
これを期待値の定義に基づいて計算するのは大変なので、このときの確率変数$X$は同じくパラメータ$p$の幾何分布の確率変数を$Y_i$とした際に、$X=Y_1+Y_2+…++Y_k$と表すことができる。
ここで問題1.2の解答より、$\displaystyle E[Y_i] = \frac{1}{p}$なので、負の二項分布の期待値$E[X]$は下記のように求めることができる。
$$
\begin{align}
E[X] &= E[Y_1+Y_2+…++Y_k] \\
&= \frac{1}{p} + \frac{1}{p} + … + \frac{1}{p} \\
&= \frac{k}{p}
\end{align}
$$
よって求める期待値は$\displaystyle E[X] = \frac{k}{p}$となる。

問題1.4の解答例

i)
$(a+b)^N=(a+b)^{R}(a+b)^{N-R}$の$a^{n}b^{N-n}$の係数を左辺と右辺に関して考える。
・左辺の$a^{n}b^{N-n}$の係数
$$
\begin{align}
{}_N C_{n}
\end{align}
$$

・右辺の$a^{n}b^{N-n}$の係数
$$
\begin{align}
\sum_{r=0}^{\min(n,R)} {}_{R} C_{r} \times {}_{N-R} C_{n-r}
\end{align}
$$

ここで「左辺の$a^{n}b^{N-n}$の係数=右辺の$a^{n}b^{N-n}$の係数」より、下記が成立する。
$$
\begin{align}
{}_N C_{n} &= \sum_{r=0}^{\min(n,R)} {}_{R} C_{r} \times {}_{N-R} C_{n-r} \\
1 &= \frac{\sum_{r=0}^{\min(n,R)} {}_{R} C_{r} \times {}_{N-R} C_{n-r}}{{}_N C_{n}} = p_r \\
p_r &= 1
\end{align}
$$
ⅱ)
$p_r$は下記のように変形できる。
$$
\large
\begin{align}
p_r &= \frac{{}_R C_r {}_{N-R} C_{n-r}}{{}_N C_n} \\
&= \frac{R(R-1)…(R-r+1)}{r!} \times \frac{(N-R)(N-R-1)…(N-R-(n-r)+1)}{(n-r)!} \times \frac{n!}{N(N-1)…(N-n+1)} \\
&= \frac{n!}{r!(n-r)!} \times \frac{R(R-1)…(R-r+1) \times (N-R)(N-R-1)…(N-R-(n-r)+1)}{N(N-1)…(N-n+1)} \\
&= {}_n C_r \times \frac{R/N(R/N-1/N)…(R/N-(r-1)/N) \times (1-R/N)(1-R/N-1/N)…(1-R/N-(n-r)/N+1/N)}{(1-1/N)…(1-n/N+1/N)} \\
&= {}_n C_r \times \frac{p(p-1/N)…(p-(r-1)/N) \times (1-p)(1-p-1/N)…(1-p-(n-r)/N+1/N)}{(1-1/N)…(1-n/N+1/N)}
\end{align}
$$

上記に対して$N \to \infty$を考えると、下記のように二項分布$Bin(n,p)$の確率関数の式が導出できる。
$$
\large
\begin{align}
\lim_{N \to \infty} p_r &= {}_n C_r \times \lim_{N \to \infty} \frac{p(p-1/N)…(p-(r-1)/N) \times (1-p)(1-p-1/N)…(1-p-(n-r)/N+1/N)}{(1-1/N)…(1-n/N+1/N)} \\
&= {}_n C_r p^{r} (1-p)^{n-r}
\end{align}
$$

問題1.5の解答例

確率変数$X$がポアソン分布の$Po(\lambda)$に従うとき$X=k$となる確率を$P(X=k|\lambda)$とすると、$P(X=k|\lambda)$は下記のように表すことができる。
$$
\large
\begin{align}
P(X=k|\lambda) = \frac{\lambda^k e^{-\lambda}}{k!}
\end{align}
$$
よって、ポアソン分布の期待値を$E[X]$とすると、$E[X]$は下記のように計算することができる。
$$
\large
\begin{align}
E[X] &= \sum_{k=0}^{\infty} kP(X=k|\lambda) \\
&= \sum_{k=1}^{\infty} k \times \frac{\lambda^k e^{-\lambda}}{k!} \\
&= e^{-\lambda} \sum_{k=1}^{\infty} \frac{\lambda^k}{(k-1)!} \\
&= e^{-\lambda} \sum_{k=0}^{\infty} \frac{\lambda^{k+1}}{k!} \\
&= \lambda e^{-\lambda} \sum_{k=0}^{\infty} \frac{\lambda^k}{k!}
\end{align}
$$
上記の計算にあたっては$f(x)=e^x$のマクローリン展開の式を考える。
$$
\large
\begin{align}
e^x = \sum_{n=0}^{\infty} \frac{x^n}{n!}
\end{align}
$$
上記の式を$n=k$、$x=\lambda$で置き換えると、$\displaystyle e^{\lambda} = \sum_{k=0}{\infty} \frac{\lambda^k}{k!}$とすることができる。よって、ポアソン分布の期待値$E[X]$は下記のように計算できる。
$$
\large
\begin{align}
E[X] &= \lambda e^{-\lambda} \sum_{k=0}^{\infty} \frac{\lambda^k}{k!} \\
&= \lambda e^{-\lambda} e^{\lambda} \\
&= \lambda
\end{align}
$$
次にポアソン分布の分散を$V[X]$とする。定義より下記のように計算できる。（ここでの計算は確率母関数を用いた計算と同様な流れになる。）
$$
\large
\begin{align}
V[X] &= E[X^2]-E[X]^2 \\
&= E[X(X-1)]+E[X]-E[X]^2
\end{align}
$$
$E[X]$はすでに求めているので、ここでは$E[X(X-1)]$を求めることを考える。
$$
\large
\begin{align}
E[X(X-1)] &= \sum_{k=0}^{\infty} k(k-1)P(X=k|\lambda) \\
&= \sum_{k=2}^{\infty} k(k-1) \times \frac{\lambda^k e^{-\lambda}}{k!} \\
&= e^{-\lambda} \sum_{k=2}^{\infty} \frac{\lambda^k}{(k-2)!} \\
&= e^{-\lambda} \sum_{k=0}^{\infty} \frac{\lambda^{k+2}}{k!} \\
&= \lambda^2 e^{-\lambda} \sum_{k=0}^{\infty} \frac{\lambda^k}{k!} \\
&= \lambda^2 e^{-\lambda} e^{\lambda} \\
&= \lambda^2
\end{align}
$$
よって、$V[X]=E[X(X-1)]+E[X]-E[X]^2$は下記のように計算できる。
$$
\large
\begin{align}
V[X] &= E[X(X-1)]+E[X]-E[X]^2 \\
&= \lambda^2 + \lambda – (\lambda)^2 \\
&= \lambda
\end{align}
$$
ここまでの導出によって、ポアソン分布の期待値は$E[X]=\lambda$、分散は$V[X]=\lambda$となる。

問題1.6の解答例

i)
下記を実行することで当てはまり具合を調べる。

import numpy as np

x = np.array([0., 1., 2., 3., 4.])
observed = np.array([109., 65., 22., 3., 1.])
lamb = np.sum(x*observed)/np.sum(observed)

fact_x = np.array([1., 1., 2., 2.*3., 2.*3.*4.])
poisson_p = lamb**x * np.e**(-lamb)/fact_x

print(lamb)
print(poisson_p*np.sum(observed))
print(poisson_p*np.sum(observed)-observed)

・実行結果

> print(lamb)
0.61
> print(poisson_p*np.sum(observed))
[ 108.67017381   66.28880603   20.21808584    4.11101079    0.62692915]
> print(poisson_p*np.sum(observed)-observed)
[-0.32982619  1.28880603 -1.78191416  1.11101079 -0.37307085]

実行結果を確認すると、誤差が$2$以内であることが確認できる。

ⅱ)
下記を実行することで当てはまり具合を調べる。

import numpy as np

x = np.array([0., 1., 2., 3., 4., 5., 6., 7.])
observed = np.array([74., 98., 90., 61., 21., 13., 7., 1.])
lamb = np.sum(x*observed)/np.sum(observed)

fact_x = np.array([1., 1., 2., 2.*3., 2.*3.*4., 2.*3.*4.*5., 2.*3.*4.*5.*6., 2.*3.*4.*5.*6.*7.])
poisson_p = lamb**x * np.e**(-lamb)/fact_x

print(lamb)
print(poisson_p*np.sum(observed))
print(poisson_p*np.sum(observed)-observed)

・実行結果

> print(lamb)
1.8054794520547945
> print(poisson_p*np.sum(observed))
[60.00440052  108.33671217   97.79985386   58.85854219   26.56697212
    9.59322445    2.88672827    0.74456123]
> print(poisson_p*np.sum(observed)-observed)
[-13.99559948  10.33671217   7.79985386  -2.14145781   5.56697212
  -3.40677555  -4.11327173  -0.25543877]

実行結果を確認すると、i)に比べて当てはまりが良くないことがわかる。これは交通事故が天候などによることからポアソン分布の前提のベルヌーイ試行が成立しないことに起因すると考えることができる。

問題1.7の解答例

期待値を$E[X]$、分散を$V[X]$とおく。それぞれ下記のように計算できる。
$$
\begin{align}
E[X] &= \int_{-\infty}^{\infty} xf(x) dx \\
&= \int_{a}^{b} \frac{1}{b-a}x dx \\
&= \left[ \frac{1}{2(b-a)} x^2 \right]_{a}^{b} \\
&= \left( \frac{1}{2(b-a)}b^2 – \frac{1}{2(b-a)}a^2 \right) \\
&= \frac{b^2-a^2}{2(b-a)} \\
&= \frac{(b-a)(b+a)}{2(b-a)} \\
&= \frac{a+b}{2} \\
V[X] &= \int_{-\infty}^{\infty} (x-E[X])^2f(x) dx \\
&= \int_{a}^{b} \frac{1}{b-a}(x-E[X])^2 dx \\
&= \int_{a}^{b} \frac{1}{b-a} \left( x-\frac{a+b}{2} \right)^2 dx \\
&= \left[ \frac{1}{3(b-a)} \left( x-\frac{a+b}{2} \right)^3 \right]_{a}^{b} \\
&= \frac{1}{3(b-a)} \left( \left( b-\frac{a+b}{2} \right)^3 – \left( b-\frac{a+b}{2} \right)^3 \right) \\
&= \frac{1}{3(b-a)} \left( \left(b-\frac{a+b}{2}\right)-\left(a-\frac{a+b}{2}\right) \right) \left( \left(b-\frac{a+b}{2}\right)^2 + \left(b-\frac{a+b}{2}\right)\left(a-\frac{a+b}{2}\right) + \left(a-\frac{a+b}{2}\right)^2 \right) \\
&= \frac{1}{3(b-a)}(b-a)\left( \frac{(a-b)^2}{4} + \frac{(a-b)(b-a)}{4} + \frac{(b-a)^2}{4} \right) \\
&= \frac{1}{3}\frac{(b-a)^2}{4} \\
&= \frac{(b-a)^2}{12}
\end{align}
$$
上記では$V[X]$の計算にあたって、$a^3-b^3=(a-b)(a^2+ab+b^2)$が成立することを用いた。

問題1.8の解答例

$X \sim Po(\lambda)$のモーメント母関数を$m_{X}(t)$とおくと、確率密度関数が$\displaystyle f(x|\lambda) = \frac{\lambda^{x}e^{-\lambda}}{x!}$であることより、$m_{X}(t)$は下記のように求めることができる。
$$
\large
\begin{align}
m_{X}(t) &= E[e^{tX}] = \sum_{k=0}^{\infty} e^{tk} f(x|\lambda) \\
&= \sum_{k=0}^{\infty} e^{tk} \times \frac{\lambda^{k}e^{-\lambda}}{k!} \\
&= e^{-\lambda} \sum_{k=0}^{\infty} \frac{(\lambda e^{t})^{k}}{k!} \\
&= e^{-\lambda} \times e^{\lambda e^{t}} \\
&= e^{\lambda(e^{t}-1)}
\end{align}
$$

上記の結果を活用することで、$X_1+X_2$のモーメント母関数$m_{X_1+X_2}(t)$は下記のように導出できる。
$$
\large
\begin{align}
m_{X_1+X_2}(t) &= E[e^{t(X_1+X_2)}] = E[e^{tX_1}]E[e^{tX_2}] \\
&= e^{\lambda_1(e^{t}-1)} \times e^{\lambda_2(e^{t}-1)} \\
&= e^{(\lambda_1+\lambda_2)(e^{t}-1)}
\end{align}
$$

上記が$Po(\lambda_1+\lambda_2)$に一致するので、モーメント母関数と確率分布の1対1対応より、$X_1+X_2 \sim Po(\lambda_1+\lambda_2)$であることを示すことができる。

問題1.9の解答例

i)
確率変数$X$が$X \sim Ex(\lambda)$のときのモーメント母関数を$m_{X}(t)$とおくと、$Ex(\lambda)$の確率密度関数が$f(x|\lambda) = \lambda e^{-\lambda x}, x \geq 0$であることから$m_{X}(t)$は下記のように計算できる。
$$
\large
\begin{align}
m_{X}(t) &= E[e^{tX}] = \int_{0}^{\infty} e^{tx} f(x|\lambda) dx \\
&= \int_{0}^{\infty} e^{tx} \times \lambda e^{-\lambda x} dx \\
&= \lambda \int_{0}^{\infty} e^{-(\lambda-t) x} dx \\
&= \frac{-\lambda}{\lambda-t} \left[ e^{-(\lambda-t) x} \right]_{0}^{\infty} \\
&= \frac{\lambda}{\lambda-t}
\end{align}
$$

上記の計算は$\lambda-t>0$の前提の元に行なったが、前提が成立しない時はモーメント母関数は存在しない。よって、$\lambda>t$の場合のみを以下では考える。

ここまでの結果を活用することで、$S_n$のモーメント母関数$m_{S_n}(t)$は下記のように導出できる。
$$
\large
\begin{align}
m_{S_n}(t) &= E[e^{t(X_1+X_2+…+X_n)}] = E[e^{tX_1}]…E[e^{tX_n}] \\
&= \prod_{i=1}^{n} \frac{\lambda}{\lambda-t} \\
&= \left( \frac{\lambda}{\lambda-t} \right)^n
\end{align}
$$

ⅱ)
「自然科学の統計学」の表記ではガンマ分布$Ga(\alpha,\lambda)$の確率密度関数$f(x|\alpha,\lambda)$は下記のように表される。
$$
\large
\begin{align}
f(x|\alpha,\lambda) = \frac{\lambda^{\alpha}}{\Gamma(\alpha)} x^{\alpha-1} e^{-\lambda x} \quad x \geq 0
\end{align}
$$
「現代数理統計学」では二つ目のパラメータが逆数で表されるなど、表記が異なる場合もあるので注意が必要である。

ここで確率変数$X \sim Ga(\alpha,\lambda)$を考えた時のモーメント母関数を$m_{X}(t)$とおくと、$m_{X}(t)$は下記のように計算できる。
$$
\large
\begin{align}
m_{X}(t) &= E[e^{tX}] = e^{tx} f(x|\alpha,\lambda) dx \\
&= \int_{0}^{\infty} e^{tx} \times \frac{\lambda^{\alpha}}{\Gamma(\alpha)} x^{\alpha-1} e^{-\lambda x} dx \\
&= \frac{\lambda^{\alpha}}{\Gamma(\alpha)} \int_{0}^{\infty} x^{\alpha-1} e^{-(\lambda-t) x} dx \\
&= \frac{\lambda^{\alpha}}{\Gamma(\alpha)} \times \left( \left[ \frac{-1}{\lambda-t} x^{\alpha-1} e^{-(\lambda-t) x} \right]_{0}^{\infty} + \frac{\alpha-1}{\lambda-t} \int_{0}^{\infty} x^{\alpha-2} e^{-(\lambda-t) x} dx \right) \\
&= \frac{\lambda^{\alpha}}{\Gamma(\alpha)} \times \frac{\alpha-1}{\lambda-t} \times \int_{0}^{\infty} x^{\alpha-2} e^{-(\lambda-t) x} dx \\
&= \frac{\lambda^{\alpha}}{\Gamma(\alpha)} \times \frac{\alpha-1}{\lambda-t} \times \frac{\alpha-2}{\lambda-t} \times \int_{0}^{\infty} x^{\alpha-3} e^{-(\lambda-t) x} dx \\
&= … \\
&= \frac{\lambda^{\alpha}}{\Gamma(\alpha)} \times \frac{\alpha-1}{\lambda-t} \times … \times \frac{1}{\lambda-t} \times \int_{0}^{\infty} x^{0} e^{-(\lambda-t) x} dx \\
&= \frac{\lambda^{\alpha}}{(\lambda-t)^{\alpha-1}} \times \int_{0}^{\infty} e^{-(\lambda-t) x} dx \\
&= \frac{\lambda^{\alpha}}{(\lambda-t)^{\alpha-1}} \times \left[ \frac{-1}{\lambda-t} e^{-(\lambda-t) x} \right]_{0}^{\infty} \\
&= \frac{\lambda^{\alpha}}{(\lambda-t)^{\alpha}} = \left( \frac{\lambda}{\lambda-t} \right)^{\alpha}
\end{align}
$$

上記の計算にあたっては$\Gamma(\alpha) = (\alpha-1)!$が成立することを用いた。

ⅲ)
ⅱ)で$\alpha=n$のようにおくと、i)の結果に一致する。モーメント母関数が一致することより、$S_n$の分布はガンマ分布$Ga(n,\lambda)$に等しいことが確認できる。

問題1.10の解答例

$X_1, X_2, …, X_n \sim Po(\lambda), i.i.d., \quad \lambda=1$であり、ポアソン分布の再生性より$S_n = X_1+X_2+…+X_n \sim Po(n)$が成立する。よって、下記のように考えることができる。
$$
\large
\begin{align}
e^{-n} \sum_{k=0}^{n} \frac{n^{k}}{k!} &= \sum_{k=0}^{n} \frac{n^{k} e^{-n}}{k!} \\
&= \sum_{k=0}^{n} P(S_n=k) \\
&= P(S_n \leq n) \\
&= P(S_n-n \leq 0) \\
&= P \left( \frac{S_n-n}{\sqrt{n \lambda}} \leq 0 \right) \\
\lambda &= 1
\end{align}
$$

上記に対して中心極限定理を適用すると、下記のように導出できる。
$$
\large
\begin{align}
\lim_{n \to \infty} P \left( \frac{S_n-n}{\sqrt{n \lambda}} \leq 0 \right) &= \Phi(0) = \frac{1}{2} \\
\Phi(x) &= \int_{-\infty}^{x} \frac{1}{\sqrt{2 \pi}} e^{-\frac{z^2}{2}} dz
\end{align}
$$

まとめ

Chapter.$1$では確率分布を中心に確率の基本について確認する内容でした。基本的に入門書でも記載のある内容ですが、重要事項のためなるべく問題演習を行なって慣れておくと良い内容だと思われました。

自然科学の統計学 (基礎統計学)

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

「統計学実践ワークブック」演習問題中心第2章確率分布と母関数

投稿日: 2021-06-202022-10-05 投稿者: yoichi_t

統計検定準$1$級対応の公式テキストである「統計学実践ワークブック」を$1$章から順に演習問題を中心に解説していきます。
今回は第$2$章「確率分布と母関数」です。

日本統計学会公式認定統計検定準1級対応統計学実践ワークブック

3,080円(11/10 13:44時点)

重要ポイント

本章では、第$1$章でも扱った確率関数(probability function)について深めていきます。確率関数はデータの素性、対象のモデルを仮定するために必要な概念で統計的な分析や推論を行う上で非常に重要です。確率関数の具体例として$5,6$章で代表的な確率分布が紹介されます。また、それらを利用した推論などは$8$章以降となります。

累積分布関数(cumulative distribution)

確率密度関数は、「密度」なので範囲を指定することで「確率」となります。累積分布関数$F(x)$は$-\inf$から$x$までの確率を表します。

$$
F(t) = \int^{t}_{-\infty} f(x) dx
$$

離散確率変数の場合も同様です。

同時確率密度関数(joint probability density function)

二つの確率変数$x, y$の確率密度関数を同時確率密度関数$p(x, y)$と呼びます。

条件付き確率密度関数(conditional probability density function)

同時確率において、ある変数の値が特定の値となると設定するときの確率密度関数を条件付き確率密度関数$p(x|y)$と呼びます。$p(x|y)$の時は確率変数yの値が特定の値と決められている時のxについての確率密度関数という意味です。

条件付き確率密度関数と同時確率密度関数は以下の性質を持っています。

$$
p(x | y) = \frac{p(x, y)}{p(y)} = \frac{p(x, y)}{\int p(x, y) dx}
$$

母関数(generating function)

母関数は、統計の分野では、確率密度関数の期待値や分散など確率関数の性質を知るために用いられることが多いです。離散確率変数の場合には「確率母関数」、連続確率変数の場合には「モーメント母関数」が用いられることが多いです。

なお、母関数については、こちらも参照してください。

演習問題解説

演習問題の全文は掲載しません。テキストは各自で用意するか、以下の抜粋から想像してください。

問$2.1$

$xy$平面の単位正方形（$0 \le x \le 1, 0 \le y \le 1$）上の以下の確率密度関数$f(x, y)$について。

$$
f(x, y) = c(x+y)
$$

$(1)$ 規格化定数$c$は？

確率密度関数は、積分して$1$になる必要があります。積分するということは、面積になるので、x, yの範囲から$1$であることがすぐにわかります。

また、定義に従って積分しても結局$1.0$となることがわかります。
$$
\begin{eqnarray}
\int_0^1 \int^1_0 (x+y) dxdy &=& \int^1_0 x dx + \int^1_0 y dy \\
&=& \left[ \frac{1}{2}x^2 \right]^0_1 + \left[ \frac{1}{2}y^2 \right]^0_1 \\
&=& \frac{1}{2} + \frac{1}{2} = 1.0
\end{eqnarray}
$$

$(2)$ Xの周辺確率密度関数は？

周辺密度とは、同時確率において、対象の確率変数（ここでは$X$）以外の確率変数について積分した確率密度関数です。
$$
\begin{eqnarray}
f(x) &=& \int f(x, y) dy\\
&=& \int^1_0 (x+y) dy &=& \int^1_0 x dy + \int^1_0 y dy \\
&=& x \left[ y \right]^1_0 + \left[ \frac{1}{2}y^2 \right]^1_0 \\
&=& x+\frac{1}{2}
\end{eqnarray}
$$

$(3)$ $X$を与えた時の$Y$の条件付き確率密度関数は？

条件付き確率密度関数は同時確率密度関数と関係がありますので、その性質を使って導出します。

$$
\begin{eqnarray}
f(y | x) &=& \frac{f(x, y)}{f(x)} \\
&=& \frac{x+y}{\displaystyle x+\frac{1}{2}}
\end{eqnarray}
$$

問$2.2$

幾何分布に従う確率変数$X$について。幾何分布は以下で定義されている分布です。なお、$X$は非負の整数です。

$$
X \sim p(X=x) = \theta (1-\theta )^x
$$

確率母関数$G(s)$を求めよ

確率母関数$G(s)$は$E[s^X]$なので、定義通りに計算します。

$$
\begin{eqnarray}
G(s) &=& E[s^X] = \sum_x s^x p(x) \\
&=& \sum_x s^x \theta (1-\theta)^x \\
&=& \theta \sum_x (sq)^x
\end{eqnarray}
$$

ここで、$q = 1-\theta$としました。$x$は非負の整数なので、総和の範囲は$0$から$\infty$までです。
上記最後の式から、確率母関数$G(s)$は等比数列の無限和となります。公式を覚えていればそのまま当てはめることで以下のように解を導出できます。

$$
\begin{align}
G(s) &= \theta \sum_x (sq)^x \\
&= \theta \times \frac{1}{1-sq} \\
&= \frac{\theta}{1-s(1-\theta)}
\end{align}
$$

なお、等比数列の収束条件として、$\displaystyle s \lt \frac{1}{1-\theta}$との制限があります。

ここで、等比数列の無限和は以下のようにして導出ができます。

導出する等比数列を$\displaystyle \sum^{\infty}_{x=0} r^x$とします。この数列について、$0$から$(n-1)$までの和を$S_n$とします。

$$
\sum^{n-1}_{x=0} r^x = 1 + r + r^2 + \cdots + r^{n-1} = S_n
$$

続いて、$S_n$に$r$を掛けると次のようになります。

$$
rS_n = r + r^2 + \cdots + r^{n-1} + r^n = S_n
$$

これらの差から$S_n$を以下のように表現することができます。

$$
\begin{eqnarray}
S_n – rS_n = 1 – r^n \\
S_n = \frac{1-r^n}{1-r}
\end{eqnarray}
$$

これから、$n$を無限大まで増やして行った時に、$r \lt 1$と仮定すると以下の式になることがわかります。

$$
\lim_{n \rightarrow \infty} S_n = \frac{1}{1-r}
$$

なお、無限和が収束するためには、$r \lt 1$であることが条件となります。

確率母関数を利用して期待値と分散を求めよ

確率母関数の$1$階微分と$2$階微分から期待値と分散は以下のように導出できます。なお、$G^{(k)}(s)$とは、確率母関数$G(s)$の$k$階微分を表します。

$$
\begin{eqnarray}
G^{(1)}(s) &=& E[Xs^{X-1}] \\
G^{(1)}(1) &=& E[X] \\
G^{(2)}(s) &=& E[X(X-1) s^{X-2}] \\
G^{(2)}(1) &=& E[X(X-1)]
\end{eqnarray}
$$

ここから、

$$
\begin{eqnarray}
E[X] &=& G^{(1)}(1) \\
V[X] &=& E[X^2] – (E[X])^2 \\
&=& G^{(2)}(1) + G^{(1)}(1) – (G^{(1)}(1))^2
\end{eqnarray}
$$

なので、素直に当てはめることで導出できます。

$$
\begin{eqnarray}
G^{(1)}(s) &=& \frac{\theta (1 – \theta)}{( 1-s(1-\theta) )^2} \\
G^{(1)}(1) &=& \frac{\theta (1 – \theta)}{( 1-(1-\theta) )^2} \\
&=& \frac{1-\theta}{\theta} \\
G^{(2)}(s) &=& \frac{2 \theta (1-\theta)^2}{(1 – s(1-\theta))^3} \\
G^{(2)}(1) &=& \frac{2 \theta (1-\theta)^2}{(1 – (1-\theta))^3} \\
&=& \frac{2(1-\theta)^2}{\theta^2}
\end{eqnarray}
$$

これから、

$$
\begin{eqnarray}
E[X] &=& G^{(1)}(1) = \frac{1-\theta}{\theta}\\
V[X] &=& G^{(2)}(1) + G^{(1)}(1) – (G^{(1)}(1))^2 \\
&=& \frac{2(1-\theta)^2}{\theta^2} + \frac{1-\theta}{\theta} – \left( \frac{1-\theta}{\theta} \right) ^2 \\
&=& \frac{1-\theta}{\theta^2}
\end{eqnarray}
$$

なお、$s$の範囲には条件がありましたが、$\theta$は確率で$0$~$1$の範囲の整数なので、$s$の条件を満たします。

問$2.3$

指数分布に従う確率変数$X$について。指数分布は以下で定義されている分布です。なお、$X$は$0$以上の実数です($X \ge 0$)。

$$
p(X) = \lambda \exp\{ -\lambda x \}
$$

なお、$\lambda$は指数分布のパラメータで、正の実数です。

$X$のモーメント母関数$m(\theta)$を求めよ

モーメント母関数$m(\theta)$は以下の式で定義されています。

$$
m(\theta) = E[\exp\{\theta x\}] = \int^{\infty}_{-\infty} \exp\{\theta x\} p(x) dx
$$

この式に従って導出するだけです。$X$は指数分布に従う正の連続変数なので、積分範囲は$0$から$\infty$です。

$$
\begin{eqnarray}
m(\theta) &=& E[e^{\theta X}] = \int^{\infty}_{0} e^{\theta x} \lambda e^{-\lambda x } dx \\
&=& \lambda \int^{\infty}_{0} e^{(\theta -\lambda)x} dx \\
&=& \lambda \left[ \frac{1}{\theta – \lambda} e^{(\theta – \lambda)x} \right]^{\infty}_0 \\
&=& \lambda \left[ 0 – \frac{1}{\theta – \lambda} e^{0} \right] \\
&=& \frac{\lambda}{ \lambda – \theta}
\end{eqnarray}
$$

ここで、$3$行目から$4$行目の導出は$0 \lt \theta \lt \lambda$との条件を仮定しています。モーメント母関数の定義域に関しては「モーメント母関数とマクローリン展開」で詳しく取り扱いましたのでこちらも合わせてご確認ください。

モーメント母関数$m(\theta)$を利用して期待値と分散を求めよ

期待値と分散を導出するために、$1$階微分と$2$階微分を考えます。ここで、モーメント母関数の$k$階微分に対して、$\theta=0$とすると、$k$次のモーメント$E[X^k]$になることを利用します。（$m^{(k)}(0) = E[X^k]$）

$$
\begin{eqnarray}
E[X] &=& m^{(1)}(0) \\
V[X] &=& m^{(2)}(0) – (m^{(1)}(0))^2
\end{eqnarray}
$$

指数分布のモーメント母関数の$\theta$についての$1$階微分と$2$階微分は以下のようになります。

$$
\begin{eqnarray}
m^{(1)}(\theta) &=& \frac{\lambda}{(\lambda – \theta)^2} \\
m^{(2)}(\theta) &=& \frac{2\lambda}{(\lambda – \theta)^3}
\end{eqnarray}
$$

以上から、期待値と分散は次のように導出されます。

$$
\begin{eqnarray}
E[X] &=& m^{(1)}(0) = \frac{1}{\lambda} \\
V[X] &=& m^{(2)}(0) – (m^{(1)}(0))^2 = \frac{2 \lambda}{\lambda^3} – \frac{1}{\lambda^2} \\
&=& \frac{1}{\lambda^2}
\end{eqnarray}
$$

参考文献

・松原ら, 統計学入門, $1991$, 東京大学出版会

統計学入門 (基礎統計学Ⅰ)