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

はじめに
記事リンク
4.3 データのもつ情報量
- 4.3.0 お気持ち
- 4.3.1 フィッシャー情報量
今日はここまで

はじめに

前回、承認欲求の権化と化した僕が承認欲求砲を大量に放ったお陰で、前回記事と前々回記事を多くの人に読んでもらうことができた。ありがとうございました。とってもうれしかったです。
読んだらちゃんと評価してくれよな。無反応は経済を停滞させ、クリエイターは減少の一途を辿り、そして人類は緩やかに死滅する。僕はまだ死滅したくないので、記事が面白いと思ったらツイッターでもっと拡散してほしい。記事がクソだと思ったらクソだクソだと叫びながら拡散してほしい。よろしゅうお願いします。

調子に乗ったのでクソどうでもいい話をするが、はてなブログには連続投稿日数なるものがあるようで、奴はここ3日間連続で記事を投稿した（めっちゃえらい）僕に無言の圧力を仕掛けてくる。
圧力と言っても、別に通知でうるさく喚いてきたりするわけではなく、記事を投稿した直後に「継続期間 3日」などと要らない情報を表示してくるのである。
ちょうどこんな感じだ。

f:id:kriver2135:20180501203154p:plain — Fig.1 余計なお世話

僕は昔から言われたことに唯々諾々と従うのが大好きなミスター唯々諾々だったので、こんなことをされたら継続期間を伸ばす以外一切のことが考えられなくなってしまう。なので今日も唯々諾々と記事を投稿していこうと思う。みんなも僕と一緒に唯々諾っていこう。

記事リンク

第一回
kriver-1.hatenablog.com

前回
kriver-1.hatenablog.com

次回
kriver-1.hatenablog.com

4.3 データのもつ情報量

4.3.0 お気持ち

前回までで、我々はデータの山 $\{(x_1, y_1), \cdots, (x_N, y_N)\}$ に対してモデルを仮定し、最尤推定量を計算するという一連の流れができるようになった。当初から予告していたとおり、今後は推定量の良さをいろいろと説明して、そのあとで最尤推定量がそれらをわりと満たすスーパーすごい推定量であることを示していこうと思う。

4.3.1項では、データにフィッシャー情報量という量を定義することについて話す。たぶん何の役に立つんだかいまいちわからないとは思うが、これは伏線を張っている段階だ。そしてネタバレになってしまうが、その伏線は4.3.2項で早々に回収される。すなわち、フィッシャー情報量の逆数が、クラメール・ラオの下限と呼ばれる不偏推定量の分散の下限に等しいという圧倒的事実が公開される。あまりに華麗な伏線回収に、読者諸氏は驚きに口を開け小便を漏らすことになるだろう。いまのうちにトイレに行っておくべきだ。

データに情報量がある、という話を、直感的な視点から説明しておこう。
我々は $x$ が与えられたときの $y$ の確率分布を知りたいのだった。例えばいま、 $x=1$ から $x=2$ までにおけるデータ点が既にたくさんあるとしよう。新しいデータ点はどこに欲しいだろうか。横軸が $x$ 、縦軸が $y$ であるようなグラフを考えて、 $1 \le x \le 2$ の範囲にたくさん点が打ってある状況を考えてほしい。足りないのは $x < 1$ と $2 < x$ の部分だ。それも $x=1$ や $x=2$ から十分離れている部分が圧倒的に足りていない。
この足りない部分についてのデータは、もう腐るほど余っている $1 \le x \le 2$ の部分に比べて多くのことを語っていると考えられる。つまりその部分のデータは情報量が多いとも言える。これを数式にしたのがフィッシャー情報量だ。
詳しい内容については次項で見ていくが、直感的にはそういう感じだと思ってほしい。

4.3.1 フィッシャー情報量

話を簡単にするために、この項ではとびきり単純なモデルを扱う。「原点を通る単回帰モデル」だ。
単回帰モデルというのは、説明変数 $x$ も目的変数 $y$ も1次元量なモデルのことだ。そして原点を通るというのは、 $x$ が0のとき $y$ の推定値も0になるということだ。要するに、次のようなモデルを考えるということになる。

$\displaystyle y_i = a x_i + \varepsilon_i$

考える必要があるのか怪しくなるレベルの単純さだ。ここで、 $\varepsilon_i$ は実データ $y_i$ と回帰値 $ax_i$ の差を表していて、この部分は独立で分散が一定な正規分布 $\mathcal{N}(0, \sigma^2)$ に従っていると仮定する。ここは最尤法のページであって線形モデルのページではないので詳しくは述べないが、（正規）線形モデルを仮定するときにはいっしょに誤差が独立一定な正規分布に従っていることも仮定しなければならないのである。
ここでは、簡単のため $\sigma$ についてももう値が分かっているものとしよう。知りたいのは $a$ だけだ。

これを線形モデルに対する最適化手法（最小二乗法）で解いたっていいのだが、せっかくなので読者が知ったばかりの最尤法でこの $a$ を最適化しよう。線形モデルは一般線形モデルの一種なので、最尤法は線形モデルにも使うことができる。ただし結果が同じになるとは限らない。最尤法はより一般性の高い手法なので、その分線形モデルに対する優位性は最小二乗法よりも少なくなっている（後述するように、最尤法はデータ数が少ないときは良さが保証されない）。常に上位互換というわけではないのである。

まずは対数尤度 $\log L(\theta; Y)$ を計算しよう。
尤度 $L$ は、そのパラメータで実際にサイコロを振ってみた時に、データ点がそっくりそのまま出てくる確率のことであった。いま、 $\varepsilon_i = y_i - ax_i$ が正規分布 $\mathcal{N}(0, \sigma^2)$ に従っていることを仮定しているのだから、それを使って確率を計算してあげればよい。すなわち、
$\displaystyle \begin{eqnarray} L(\theta; Y) &=& \prod_i^N (\text{正規分布から} y_i - ax_i \text{が飛び出してくる確率}) \\ &=& \prod_i^{N} Pr\left(\varepsilon_i = y_i - ax_i | \varepsilon_i \sim \mathcal{N}(0, \sigma^2) \right) \\ &=& \prod_i^N (2\pi \sigma^2)^{-\frac{N}{2}}\exp\left\{ - \frac{\left((y_i - ax_i) - 0\right)^2} {2\sigma^2} \right\} \\ \log L(\theta; Y) &=& \sum_i^N \left\{ -\frac{N}{2} \log (2\pi\sigma^2) \right\} -\sum_i^N \frac{(y_i - ax_i)^2}{2\sigma^2} \\ &=& -\frac{N^2}{2}\log(2\pi\sigma^2) - \frac{1}{2\sigma^2} \left\lVert \boldsymbol{y} - a\boldsymbol{x}\right\rVert^2 \end{eqnarray}$

のようになる。
一応言っておくと、 $\prod$ は $\sum$ の積バージョンで、たくさん掛けるときに便利な記号だ。大体 $\log$ に消される運命にある。
途中で「確率」を「確率密度関数」に変えてしまったり、二乗の和をベクトルでかっこよく表してみたりしているが、やっていることは単純だ。正規分布のキモい関数を頑張って思い出して、そこに値を代入しただけである。

ではこれを $\theta$ 、すなわち $a$ で微分することで、この尤度を最大化する $a$ を求めよう。
$\displaystyle \begin{eqnarray} \frac{\partial}{\partial a} \log L(a; Y) &=& \frac{\partial}{\partial a} \left\{ - \frac{1}{2\sigma^2} \left\lVert\boldsymbol{y} - a\boldsymbol{x}\right\rVert^2 \right\} \\ &=& - \frac{1}{2\sigma^2} \frac{\partial}{\partial a} \left\lVert\boldsymbol{y} - a\boldsymbol{x}\right\rVert^2 \\ &=& - \frac{1}{2\sigma^2} \left\{ -2 \boldsymbol{x}^T (\boldsymbol{y} - a\boldsymbol{x}) \right\} \\ &=& - \frac{2}{2\sigma^2} \left\{ a \boldsymbol{x}^T\boldsymbol{x} - \boldsymbol{x}^T\boldsymbol{y} \right\} &=& 0 \\ &\Leftrightarrow& a \boldsymbol{x}^T\boldsymbol{x} &=& \boldsymbol{x}^T\boldsymbol{y} \\ &\Leftrightarrow& a &=& \frac{\boldsymbol{x}^T\boldsymbol{y}}{\boldsymbol{x}^T\boldsymbol{x}} \end{eqnarray}$

とまあ、ざっとこんな感じである。ベクトルのスカラー微分に慣れていない人は読み飛ばしてもらっても構わない。なんやかんやあって、結局尤度を最大化する $a$ は $a=\hat{a} = \frac{\boldsymbol{x}^T\boldsymbol{y}}{\boldsymbol{x}^T\boldsymbol{x}}$ と求まった。なお、これは最小二乗法の結果と一致する。

さて、今回注目したいのは実はこの結果そのものではない。この $a$ が最適な $\hat{a}$ から少し外れたとき、どのくらい尤度が悪化するか知りたいのである。

もし仮に、 $a$ が $\hat{a}$ から少し外れただけで尤度が急速に悪化するというのなら、それは $\hat{a}$ がとても良い推定量になっているということを表している。逆に、 $a$ が $\hat{a}$ からいくら外れても尤度が全然変わらないとしたら、それは $\hat{a}$ が推定量として大した働きをしていないことになる。
この雰囲気がわかるだろうか。

わからなくてもいい。前述の通り、これは次項への伏線なのである。

さて、ある関数がある値の周りでどの程度急峻か知りたければ、2階微分してその値を見ればよい。1階微分は傾きを与えるが、2階微分は傾きの傾きを与えるからである。
実際に2階微分を計算してみよう。今回は対数尤度 $L(\theta; Y)$ を $\hat{\theta}$ の周りで2階微分してみる。

$\displaystyle \begin{eqnarray} \frac{\partial}{\partial a} \log L(a; Y) &=& - \frac{1}{\sigma^2} \left\{ a \boldsymbol{x}^T\boldsymbol{x} - \boldsymbol{x}^T\boldsymbol{y} \right\} \\ \frac{\partial^2}{\partial a^2} \log L(a; Y) &=& - \frac{1}{\sigma^2} \frac{\partial}{\partial a} \left\{ a \boldsymbol{x}^T\boldsymbol{x} - \boldsymbol{x}^T\boldsymbol{y} \right\} \\ &=& -\frac{ \boldsymbol{x}^T\boldsymbol{x} }{\sigma^2} \\ &=& -\frac{ \left\lVert\boldsymbol{x}\right\rVert^2 }{\sigma^2} \end{eqnarray}$

1階微分を既に計算してあったので、だいぶ楽をすることができた。 $\theta=\hat{\theta}$ は尤度の極大値になっているはずなので、2階微分は必ず負になっている（スカラーの場合）。急峻さを表す値になってほしいので、大きいほど急峻、という気持ちを表すため、-1倍してみよう。

$\displaystyle -\frac{\partial^2}{\partial a^2} \log L(a; Y) = \frac{ \left\lVert\boldsymbol{x}\right\rVert^2 }{\sigma^2}$

この値が、 $\hat{x}$ の良さを、ひいてはデータの良さを表している、というのである。まだ根拠について話していないので、ふーんそういう考え方もあるのね、くらいの気持ちで聞いてほしい。だが、後々これが有用になってくるので、これに名前を付けておく。フィッシャー情報量だ。

より一般的には、フィッシャー情報量 $I(\theta, \boldsymbol{x})$ は次のように定義される。説明データ $\boldsymbol{x}$ と推定量 $\theta$ に対して、
$\displaystyle I(\theta_0, \boldsymbol{x}) = -\mathbb{E}_{\boldsymbol{y} \sim f(\cdot; \theta_0, \boldsymbol{x})} \left[ \left. \frac{\partial^2}{\partial \theta^2} \log f(\boldsymbol{y}; \theta, \boldsymbol{x}) \right|_{\theta=\theta_0} \right]$
となる。

いきなりジャンプしすぎただろうか？ゆっくり見比べてみてほしい。
まず、尤度 $L(\theta; Y)$ の代わりに、 $\boldsymbol{y}$ の確率密度関数 $f(\boldsymbol{y}; \theta, \boldsymbol{x})$ を用いている。偶然表記が似通っているだけで、別にフィッシャー情報量は最尤法とは無関係だからだ。特に難しいことはなくて、 $L(\theta; Y) = f(\boldsymbol{y}; \theta, \boldsymbol{x})$ である。どちらも、 $\theta$ と $\boldsymbol{x}$ が与えられたときの $\boldsymbol{y}$ の出現確率、という点では変わらない。これまでは面倒だったので $\boldsymbol{x}$ を書いてこなかったが、せっかくなのでこれも明示している。
次に、 $\theta_0$ と $\theta$ が出てきているが、これは関数を微分したものがまた関数になる、ということを明示するためである。任意の $\theta$ に対して一つ値を返す $f(\boldsymbol{y}; \theta, \boldsymbol{x})$ という関数を、対数を取って $\theta$ で2階微分した結果（これもまた $\theta$ の関数である）に、 $\theta=\theta_0$ を代入したものが $\left. \frac{\partial^2}{\partial \theta^2} \log f(\boldsymbol{y}; \theta, \boldsymbol{x}) \right|_{\theta=\theta_0}$ である。
最後に、 $\boldsymbol{y} \sim f(\cdot; \theta_0, \boldsymbol{x})$ の下で期待値を取っているが、これは単に式の結果から $\boldsymbol{y}$ を消滅させたいからである。モデルと $\theta$ が決まれば $\boldsymbol{y}$ の確率分布は決まってしまうのだから、もはや変数として残す必要もないので、期待値を取ることで結果から消しておきたいのだ。

$\theta$ がベクトルなら、この式は
$\displaystyle I(\theta_0, \boldsymbol{x}) = -\mathbb{E}_{\boldsymbol{y} \sim f(\cdot; \theta_0, \boldsymbol{x})} \left[ \left. \frac{\partial^2}{\partial \theta_i \theta_j} \log f(\boldsymbol{y}; \theta, \boldsymbol{x}) \right|_{\theta=\theta_0} \right]_{ij}$
となる。 $I(\theta)$ の各 $i,j$ 要素が上式で表されるような行列になっているということだ。この場合、この行列をフィッシャー情報行列と呼ぶ。

フィッシャー情報量の意味が本格的に分かってくるのは、次項でクラメール・ラオの下限を知った後になる。今は簡単な単回帰モデルに対してフィッシャー情報量を実際に見ていくことでお茶を濁しておこう。

「原点を通る単回帰モデル」においては、フィッシャー情報量は $I(\theta, \boldsymbol{x}) = \frac{ \left\lVert\boldsymbol{x}\right\rVert^2 }{\sigma^2}$ のようになるのだった。このモデルにおいて、フィッシャー情報量がデータ点の持つ情報量を表していることを直感的に理解したい。
まず、モデルは原点を通るので、データ点として $x=0$ の点は無意味である。 $x=0$ に対応する点は $y=0$ に決まっているからだ。
逆に、データ点がどれだけ原点から離れていたとしても問題はない。今回のモデルでは単回帰モデル、つまりすべての点をいい感じに通るような巨大な直線を引きたいのであって、むしろデータ点は遠ければ遠いほど直線が引きやすくなる。データ点が原点付近に固まっていた場合、原点付近にしか綺麗な線は引けないのに、全くデータのない遥か遠くまで直線を引き伸ばさなければならなくなって、モデルの信頼性は落ちてしまう。
また、データの誤差の分散 $\sigma^2$ は、小さければ小さいほど嬉しい。これが大きいということは、データ自体に誤差が入りまくっているということで、データの情報量が小さくなるのはあたりまえである。

これらを踏まえて、フィッシャー情報量 $I(\theta, \boldsymbol{x}) = \frac{ \left\lVert\boldsymbol{x}\right\rVert^2 }{\sigma^2}$ をもう一度見てほしい。 $x$ が原点から離れるほど、データの誤差の分散 $\sigma^2$ が小さくなるほど、情報量が大きくなっていることがわかるだろうか。
これが「データの情報量」の直感的な理解である。

今日はここまで

ちょっと短いかもしれないが、今回はここで終わろうと思う。
というのも、次項で紹介するクラメール・ラオの下限はこの章の目玉商品の一つであり、そのぶん分量も大変なことになると想像できるからだ。ここで迂闊に書き続けてしまうと、この記事だけめちゃくちゃ長くなってしまうことが想定される。
フィッシャー情報量についてなんだかモヤモヤしているかもしれないが、それは次回の記事でスッキリ晴れ晴れすることになると思う。次回の僕にご期待下さい。

記事が短くて暇になってしまったという人は、ぜひ下の星ボタンを押すことで時間を潰してほしい。ログインや会員登録処理があって意外と暇つぶしになるはずだ。興が乗らなければ右下の青い鳥でも構わない。どうぞよろしくお願いします。

前回
kriver-1.hatenablog.com

次回
kriver-1.hatenablog.com