自然科学の統計学4章最尤法を読んだその6

はじめに
記事リンク
4.4 最尤推定量の最適性
4.4.1 一致性
- 4.4.2 漸近有効性
おわりに

はじめに

5月からちまちま書き続けているこのブログだが、完全にハマってしまってすっかり困っている。
ハマるという単語はかなり文脈依存の多義語なのだが、この場合は be addicted to の意味だ。毎日統計や数学の記事が書きたくて仕方がない。

僕は現状金を支払う代わりに良い環境で自由に時間を使える立場にいるので（ニートに近似できそうだ）、生活を律して研究を律しなければならない立ち位置にいる。なのだが、こういう時間を吸う楽しい事象があると簡単に堕落してしまう。ゲームとかと違って自分のためになっている錯覚に陥るのもよくないポイントだ。

早々に切り上げて主活動に戻らなければならないのだが、身体は数学を求める。難しい問題だ。

記事リンク

第一回
kriver-1.hatenablog.com

前回
kriver-1.hatenablog.com

4.4 最尤推定量の最適性

前節までで、我々は推定量に関するいくつかの重要な指標を手に入れた。すなわち、

不偏性。期待値が真のパラメータに等しいような推定量を不偏推定量という。
有効性。不偏推定量の中でも、分散がフィッシャー情報量の逆数、クラメール・ラオの下界に等しいような不偏推定量を有効推定量という。
最小分散性。分散があらゆる不偏推定量の中で最小になるような推定量を最小分散不偏推定量という。有効推定量は必ず最小分散不偏推定量である。

この節では、最尤推定量が上記の指標を漸近的にすべて満たすこと、すなわち、データ数 $N$ が増加するにつれて最尤推定量の期待値は真のパラメータにいくらでも近づき、しかもその分布の分散はクラメール・ラオの下界にどんどん近づいていくことを示す。
さらに、その漸近分布が正規分布になっていることも示す。

これまでの議論では、最尤推定量はあくまで尤度を最大化しただけのものであって、パラメータの真の値とどの程度近いかは一切わからない、ただの推定量でしかなかった。この節では最尤推定量の良さをきちんと示すことで、最尤推定量が（ $N$ 大で）ちゃんと良い推定量になっていることを理解していこうと思う。

4.4.1 一致性

漸近的に不偏推定量となる推定量のことを、一致推定量と言う。
もう少し真面目に書くと、推定量 $t_n(\theta)$ がパラメータ $\theta$ の一致推定量であるとは、 $t_n(\theta)$ が $\theta_0$ に確率収束することを言う。

確率収束や分布収束、概収束などの違いについては、別途記事を書いたのでそちらを参照してほしい。

kriver-1.hatenablog.com

分布収束は累積分布関数が各点で収束することで、確率収束はある値以外の出る確率がどんどん小さくなっていくこと、概収束はその値以外が最終的に一切出なくなることを表している。
今回は確率収束なので、真の値以外が出る可能性はゼロというわけではないが、いくらでも小さくできることを表している。

さて、では最尤推定量が一致推定量であることを簡単に示していこう。
一致推定量であることを示すには、最尤推定量 $\hat\theta$ が真の値 $\theta_0$ に確率収束することを示せばよい。そこで、まずは最尤推定量の満たす方程式から出発しよう。

$\displaystyle \begin{eqnarray} \left. \frac{\partial }{\partial \theta} \log f(y; \theta) \right|_{\theta=\hat\theta} = \left. \sum_{i=1}^N \frac{\partial }{\partial \theta} \log f(y_i; \theta) \right|_{\theta=\hat\theta} = 0 \end{eqnarray}$

説明変数 $x$ は、どうでもいいので省略させてもらった。また、各 $y_i$ が独立であることも用いている。

さて、この式は「独立な $N$ 個の確率変数 $\frac{\partial }{\partial \theta} \log f(y_i; \theta)$ の総和」についての式になっている。せっかくなので、両辺に $1/N$ を掛けて総和を平均に変換してみよう。

$\displaystyle \begin{eqnarray} \frac{1}{N} \left. \sum_{i=1}^N \frac{\partial }{\partial \theta} \log f(y_i; \theta) \right|_{\theta=\hat\theta} = 0 \end{eqnarray}$

$N$ が十分大きい世界では、左辺は（ $\theta=\hat\theta$ でなくとも）期待値に確率収束する。大数の弱法則だ。

$\displaystyle \begin{eqnarray} \frac{1}{N} \sum_{i=1}^N \frac{\partial }{\partial \theta} \log f(y_i; \theta) \xrightarrow{p} \mathbb{E}_{y \sim f(y; \theta_0)} \left[ \frac{\partial }{\partial \theta} \log f(y; \theta) \right] \end{eqnarray}$

期待値を取るときの確率分布が $\theta$ でも $\hat\theta$ でもなく $\theta_0$ になっている点についてきちんと理解しておいてほしい。いま我々はモデルが無限に正しいと仮定し、モデルに真の値 $\theta_0$ があると仮定している。データはこの真の値からサンプルされている（と仮定している）のだから、期待値の確率分布は $\theta_0$ から取るべきなのである。

いま、 $\text{（右辺）}=0$ を満たすような $\theta$ は1つしかない、という仮定をおこう。これは背理法の仮定ではなく、最後まで正しいと思って話を進める仮定だ。
もちろん実際には、こんな仮定を置かなくとも最尤推定量の一致性を示すことはできるのだが、自然科学の統計学に合わせて、ここではこの仮定を置くことで証明を簡単にすることにする。きちんとした証明が見たい人はこのpdfとかがわかりやすいんじゃないだろうか。イェンセンの不等式と大数の強法則を用いて、最尤推定量が漸近的に真の値に概収束することを示している。

左辺は右辺に確率収束するのだから、 $\text{（左辺）}=0$ の解も $\text{（右辺）}=0$ の解に確率収束するはずだ（この理屈もわりと怪しい気がするが、略証なので許してほしい）。
ここで、もし $\theta=\theta_0$ が $\text{（右辺）}=0$ を満たすとしたら、仮定より $\theta=\theta_0$ は $\text{（右辺）}=0$ の唯一の解であるから、 $\text{（左辺）}=0$ の解の一つである $\hat\theta$ が $\theta_0$ に確率収束することになって、証明が完了したことになる。

実際、 $\theta=\theta_0$ は $\text{（右辺）}=0$ を満たす。というのも、

$\displaystyle \begin{eqnarray} \left. \text{（右辺）} \right|_{\theta=\theta_0} &=& \left. \mathbb{E}_{y \sim f(y; \theta_0)} \left[ \frac{\partial }{\partial \theta} \log f(y; \theta) \right] \right|_{\theta=\theta_0} \\ &=& \left. \int f(y; \theta_0) \frac{\partial }{\partial \theta} \log f(y; \theta) dy \right|_{\theta=\theta_0} \\ &=& \left. \int f(y; \theta_0) \frac{1}{f(y; \theta)} \frac{\partial }{\partial \theta} f(y; \theta) dy \right|_{\theta=\theta_0} \\ &=& \int \left. \frac{\partial }{\partial \theta} f(y; \theta) \right|_{\theta=\theta_0} dy \\ &=& \left. \frac{\partial }{\partial \theta} \int f(y; \theta) dy \right|_{\theta=\theta_0} \\ &=& \left. \frac{\partial }{\partial \theta} 1 \right|_{\theta=\theta_0} = 0 \end{eqnarray}$

となるからである。対数の微分を商と微分に分解する式変形や、確率分布の総和が1であることを用いた式変形は、前回までにも頻出していたので、式変形を真面目に追いかけている読者にとってはいつもどおりの式変形だったと思う。

以上から、最尤推定量が一致推定量であること、すなわち、データ数 $N$ が十分大きければ最尤推定量は漸近的に不偏推定量になることが示された。

4.4.2 漸近有効性

最尤推定量が漸近的な不偏推定量であることがわかったところで、次は最尤推定量が漸近的な有効推定量であることを示そう。
ここでも、証明を簡単にするべくいくつかの雑な近似を行う。雑な近似が見たくない潔癖症の方は、先程と同様このpdfを読んでほしい。近似の有無でだいぶ遠回りをすることになることがわかるはずだ。

まずは先程と同様、最尤推定量の満たすべき方程式から出発する。

雑な近似を使っていいとのことなので、早速中辺を雑にテイラー展開していこう。 $\hat\theta$ が $\theta_0$ のまわりでどういう挙動をするかが知りたいので、 $\theta_0$ のまわりでテイラー展開する。 $N$ が大きければ $\hat\theta$ は $\theta_0$ にいくらでも近づくのだから、項は1次くらいまで展開すれば十分だろう。

$\displaystyle \begin{eqnarray} \left. \sum_{i=1}^N \frac{\partial }{\partial \theta} \log f(y_i; \theta) \right|_{\theta=\hat\theta} \simeq \left. \sum_{i=1}^N \frac{\partial }{\partial \theta} \log f(y_i; \theta) \right|_{\theta=\theta_0} + (\hat\theta - \theta_0) \left. \sum_{i=1}^N \frac{\partial^2 }{\partial \theta^2} \log f(y_i; \theta) \right|_{\theta=\theta_0} \end{eqnarray}$

$N$ が大きいときの式変形としては、大数の法則や中心極限定理が有名だ。なので、 $N$ 個の独立な確率変数を見つけ出してあげよう。
ここでは、以下の2つを新たに確率変数として抜き出してみる。

$\displaystyle \begin{eqnarray} z &=& \left. \frac{\partial }{\partial \theta} \log f(y; \theta) \right|_{\theta=\theta_0} \\ w &=& \left. \frac{\partial^2 }{\partial \theta^2} \log f(y; \theta) \right|_{\theta=\theta_0} \end{eqnarray}$

これを用いて、上式はこのように変形できる。

$\displaystyle \begin{eqnarray} \left. \sum_{i=1}^N \frac{\partial }{\partial \theta} \log f(y_i; \theta) \right|_{\theta=\hat\theta} \simeq \sum_{i=1}^N z_i \; + \: (\hat\theta - \theta_0) \sum_{i=1}^N w_i = 0 \end{eqnarray}$

さて、気になるのは $z$ と $w$ の期待値と分散だ。これがわかれば、大数の法則や中心極限定理などを適切に用いることで、肝心の $\hat\theta-\theta_0$ の挙動がわかるだろう。

まず、 $z$ の期待値だが、これは先程も言ったように0になる。一応式変形も載せておこう。

$\displaystyle \begin{eqnarray} \mathbb{E}_{y \sim f(y; \theta_0)} \left[ z \right] &=& \mathbb{E}_{y \sim f(y; \theta_0)} \left[ \left. \frac{\partial }{\partial \theta} \log f(y; \theta) \right|_{\theta=\theta_0} \right] \\ &=& \left. \int f(y; \theta_0) \frac{\partial }{\partial \theta} \log f(y; \theta) \right|_{\theta=\theta_0} dy \\ &=& \left. \int f(y; \theta_0) \frac{1}{f(y; \theta_0)} \frac{\partial }{\partial \theta} f(y; \theta) \right|_{\theta=\theta_0} dy \\ &=& \int \left. \frac{\partial }{\partial \theta} f(y; \theta) \right|_{\theta=\theta_0} dy \\ &=& \left. \frac{\partial }{\partial \theta} \int f(y; \theta) dy \right|_{\theta=\theta_0} \\ &=& \left. \frac{\partial }{\partial \theta} 1 \right|_{\theta=\theta_0} = 0 \end{eqnarray}$

次に $z$ の分散だが、これは第四回、第五回で説明したクラメール・ラオの下界を用いることできれいに表せる。
クラメール・ラオの下界 $I(\theta)$ の定義を再掲しよう。

$\displaystyle I_1(\theta) = -\mathbb{E} \left[ \frac{\partial^2}{\partial \theta^2} \log f(y; \theta) \right] = \text{Var} \left[ \frac{\partial}{\partial \theta} \log f(y; \theta) \right]$

第四回の計算で出てきた同値表現も書いておいた。
いま、 $z=\left. \frac{\partial }{\partial \theta} \log f(y; \theta) \right|_{\theta=\theta_0}$ なのだから、この分散 $\text{Var}(z)$ は明らかに $I_1(\theta_0)$ に他ならない。
したがって、 $z$ の期待値は $0$ 、分散は $I_1(\theta_0)$ となる。

では、今度は $w$ について考えよう。
$w$ は定義から

$\displaystyle \begin{eqnarray} w &=& \left. \frac{\partial^2 }{\partial \theta^2} \log f(y; \theta) \right|_{\theta=\theta_0} \end{eqnarray}$

である。
この期待値を考えたいのだが、ここで上にあるクラメール・ラオの下界の式を見てほしい。なんと、中辺に全く同じ表現がある！
そういうわけで、 $w$ の期待値は $-I_1(\theta_0)$ となる。

さて、 $w$ の分散がまだ求まっていないが、実はこれだけあればもう十分だ。少し前の式

$\displaystyle \begin{eqnarray} \sum_{i=1}^N z_i \; + \: (\hat\theta - \theta_0) \sum_{i=1}^N w_i = 0 \end{eqnarray}$

の両辺を $N$ で割って、平均を出す。

$\displaystyle \begin{eqnarray} \frac{1}{N} \sum_{i=1}^N z_i \; + \: \frac{1}{N} (\hat\theta - \theta_0) \sum_{i=1}^N w_i = 0 \end{eqnarray}$

$z$ は期待値も分散もわかっているので、中心極限定理を使って正規分布に分布収束させてやろう。
中心極限定理は「期待値 $\mu$ 、分散 $\sigma^2$ の独立同一分布に従う確率変数 $N$ 個の平均は、期待値 $\mu$ 、分散 $\frac{\sigma^2}{N}$ の正規分布に分布収束する」という定理だった。したがって、

$\displaystyle \begin{eqnarray} \frac{1}{N} \sum_{i=1}^N z_i &\xrightarrow{d}& \mathcal{N}\left(0, \frac{I_1(\theta_0)}{N}\right) \end{eqnarray}$

となる。

一方、 $w$ は期待値しか求めていないので、大数の弱法則で期待値に確率収束させてやろう。
大数の弱法則は「期待値 $\mu$ の独立同一分布に従う確率変数 $N$ 個の平均は、期待値 $\mu$ に確率収束する」という定理だった。これを用いて、

$\displaystyle \begin{eqnarray} \frac{1}{N} \sum_{i=1}^N w_i &\xrightarrow{p}& - I_1(\theta_0) \end{eqnarray}$

となる。これらを代入して、 $N$ が十分大きいときには

$\displaystyle \begin{eqnarray} && \frac{1}{N} \sum_{i=1}^N z_i \; + \: \frac{1}{N} (\hat\theta - \theta_0) \sum_{i=1}^N w_i \\ &\simeq& \mathcal{N}\left(0, \frac{I_1(\theta_0)}{N}\right) \; - \; \left(\hat\theta - \theta_0\right) I_1(\theta_0) = 0 \\ &\Leftrightarrow& \hat\theta - \theta_0 \sim \frac{1}{I_1\left(\theta_0\right)} \mathcal{N}\left(0, \frac{I_1(\theta_0)}{N}\right) \\ &\Leftrightarrow& \hat\theta - \theta_0 \sim \mathcal{N}\left(0, \frac{1}{N I_1\left(\theta_0\right)} \right) \\ &\Leftrightarrow& \hat\theta \sim \mathcal{N}\left(\theta_0, \frac{1}{I\left(\theta_0\right)} \right) \\ \end{eqnarray}$

となる。これはまさに、 $N$ 大のとき、最尤推定量が真の値を期待値とし、かつ分散最小を達成するような正規分布に従うことを表している。

以上から、冒頭で言っていた

この節では、最尤推定量が上記の指標を漸近的にすべて満たすこと、すなわち、データ数 $N$ が増加するにつれて最尤推定量の期待値は真のパラメータにいくらでも近づき、しかもその分布の分散はクラメール・ラオの下界にどんどん近づいていくことを示す。
さらに、その漸近分布が正規分布になっていることも示す。

が、すべて示された。

おわりに

全6回というとんでもない長さになってしまったが、無事最尤推定量が様々な意味で漸近的に「良い」推定量であることを示すことができた。
最尤推定を行う上でこれらの知識が必須というわけではないが、最尤推定の意義や標本数の重要性、推定量の分散下界の存在など、重要な要素がいくつも詰まった面白い章だったと思う。

実は、自然科学の統計学の4章にはまだ続きがあるのだが（4.5節、最尤推定を利用した検定とエフィシエントスコアについての話）、これを含めると「最尤推定量の良さを調べることが目標です」という道筋が見えにくくなると感じたので、今回は含めないことにした。
むしろこの話は6章、検定力の話のあとに持ってくるべきなように思える。この節についてもいずれ紹介しようと思うが、この連載記事とは別枠にした方が良さそうだ。

なにはともあれ、最尤法についての記事はこれで最終回となる。ここまで読んでくださった皆さんには本当に感謝したい。
もし意見や感想があったら積極的に受け入れていこうと思うので、今後もぜひコメントを付けていってほしい。
よろしくお願いします。

ありがとうございました。

前回
kriver-1.hatenablog.com

第一回
kriver-1.hatenablog.com