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

はじめに

前回の記事は全然PVが伸びなかった。さすがに誰も気付いていないとは思うが、実は僕は第一回から第三回までの記事をそれぞれ21時、15時、18時に投稿している。ここまで僕の記事を読んできた読者諸氏であればすぐにピンときたであろうが、ここにはPV数の確率分布を推定するためのフィッシャー情報量を最大化しようという意図がある。その結果、記事は15時に投稿するのがよいということが示唆された。今日はその結果を活かして15時に投稿しようと思ったのだが、寝て起きたら15時を過ぎていたので残念ながら叶わなかった。実に残念だ。

4.3.2 クラメール・ラオの下限

さて、前回の記事で、我々はフィッシャー情報量というなんだかよくわからないがデータの情報としての質を表すらしい値を定義した。その定義式をここに書いておこう。というのも、この式が今回の記事でもう一度出てくるからだ。

 \displaystyle
 I(\theta', \boldsymbol{x}) = -\mathbb{E}_{\boldsymbol{y} \sim f(\cdot; \theta', \boldsymbol{x})} \left[ \left. \frac{\partial^2}{\partial \theta^2} \log f(\boldsymbol{y}; \theta, \boldsymbol{x}) \right|_{\theta=\theta'} \right]

 \theta_0がこの項で使う別の変数と名前被りを起こしていたので、 \theta'に変更させてもらった。ただの変数名なので、僕の気まぐれで勝手に変えてもなんの問題もないはずだ。

この項の目的は、不偏推定量の分散には下界があること、その下界の逆数がこのフィッシャー情報量に等しいこと、の2つを示すことにある。なので、まず不偏推定量というものが何者なのかについて説明する。その後に、その分散を実際に計算してみて(これは任意の不偏推定量について計算する)、下界があることを確かめる。最後に、その下界をわちゃわちゃ式変形して、最終的にフィッシャー情報量の逆数と等しくなることを示す。
これらの話は、最尤法や最尤推定量とは一切無関係であることに注意してほしい。これらの話は任意の推定量について成り立つものであって、この話と最尤法の実際の関係については次節で述べることになる。
この話を先にしておくのは、次節で最尤法の素晴らしさを語るための伏線なのである。

4.3.2.0 不偏推定量

前述した通り、この項で話す内容は最尤推定量だけでなく任意の推定量について成り立つ話である。
定量とは、その他の条件から確率変数を推定する関数(あるいは、その関数によって推定された予測値)のことであった。今回はデータの山からパラメータを当ててほしいので、これを t(\boldsymbol{x}, \boldsymbol{y})と表すことにしよう。 \boldsymbol{x} \boldsymbol{y}を受け取って、 \thetaの予測値を返す、任意の関数ということだ。

定量は、 \boldsymbol{x} \boldsymbol{y}を受け取って \thetaを予測してさえくれればなんだっていい。例えば、どんな \boldsymbol{x} \boldsymbol{y}が来ても無視して常に 0を返す関数も推定量だ。
もちろん、これが無意味な推定量なのはわかると思う。推定量の定義は満たしているだろうが、そんなやる気のない推定量に興味はない。我々が対象にしているのは、もっと真面目に推定するつもりのある推定量だけだ。

そこで、推定量として当然持っていてほしい性質を一つ紹介しよう。
その性質とは、「データを何度も取り直したときの予測値の期待値が、パラメータの真の値に一致する」というものである。

いま、パラメータ \thetaの真の値を \theta_0としよう。もう忘れかけているかもしれないが、パラメータは 1. データの山に対して 2. モデルを決めて 3. 最後に推定するものである。3. の工程ではモデルは無限に正しいものとしてパラメータを最適化することになる。モデルは無限に正しいのだから、当然パラメータをうまく決めればもとの確率分布を完全に再現できるはずだ。なので、その神のパラメータを \theta_0とするのである。我々の目的は \theta_0がいくらなのかを当てることだ。

真の値を \theta_0とおいたので、データ \boldsymbol{y}の確率分布は、モデル f(\boldsymbol{y};  \theta, \boldsymbol{x})と説明変数 \boldsymbol{x}を用いて、
 \displaystyle
 \boldsymbol{y} \sim f(\cdot; \theta_0, \boldsymbol{x})
と表されるはずだ。これは真の値 \theta_0の定義と言ってしまってもいい。

これを用いて、「データを何度も取り直したときの予測値の期待値が、パラメータの真の値に一致する」を数式で表そう。特に難しいことはなくて、言われたとおりに数式にしてみればよい。

 \displaystyle
 \mathbb{E}_{\boldsymbol{y} \sim f(\cdot; \theta_0, \boldsymbol{x})} \left[ t(\boldsymbol{x}, \boldsymbol{y}) \right] = \theta_0

あるいは、期待値を積分で表すのが好きな人は、

 \displaystyle
 \int f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) \; t(\boldsymbol{x}, \boldsymbol{y}) \; d\boldsymbol{y} = \theta_0

としてもよい。

この式は、まともな推定量なら当然満たしていてほしい式である。この式を満たすようなまともな推定量のことを、不偏推定量と言う。真の値から見て偏っていない推定量、というような意味だ。英語の unbiased estimator という言い方のほうが意味が伝わりやすいかもしれない。

4.3.2.1 分散の下界の存在

さて、役者は揃った。いよいよこの章の目玉の一つ、クラメール・ラオの下限の説明に入ろう。

前目において、我々は不偏推定量という推定量に関する強力な性質を知った。不偏推定量でないような推定量、例えば前述したような「どんな \boldsymbol{x} \boldsymbol{y}が来ても無視して常に 0を返す」推定量は、もし奇跡的に真の値が \theta_0=0であれば、分散0で真の値を当てたことになる。
しかしながら、これはもちろん不偏推定量の条件を満たしていない。不偏推定量という強力な制限を掛けてしまえば、もう分散0で真の値を当てるのは不可能だろう。このことをきちんと数式で示すのがこの目の目的だ。


まずは、推定量 t(\boldsymbol{x}, \boldsymbol{y})について我々が唯一仮定している情報、すなわち t(\boldsymbol{x}, \boldsymbol{y})が不偏推定量である、というところから出発しよう。

 \displaystyle
 \int f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) \; t(\boldsymbol{x}, \boldsymbol{y}) \; d\boldsymbol{y} = \theta_0

真の値をわざわざ \theta_0などと明示しているが、実際には真の値が \theta_1だろうが \theta_2だろうが上式は成り立つはずである。すなわち、これは両辺が \theta_0の関数であると見ても問題ないはずだ。

僕はこの部分で詰まってしまって2時間ほど消し飛ばし、ここがわからないというブログ記事まで書いた(下書き一覧に今も残っている)。分かる人には分かる書き方をすると、 \theta_0は一見すると「真の値」という文脈で束縛された束縛変数に思えるが、式中の f(\boldsymbol{y}; \theta_0, \boldsymbol{x})が「 \theta_0が真の値である」ということ、すなわち「 \boldsymbol{y} \theta_0をパラメータとするモデルに従って生起すること」を完全に表現しているので、 \theta_0はただの自由変数、つまり \theta_0として別の値を持ってきてもこの式はそのまま成り立つのである。

両辺が \theta_0の関数になっているのだから、両辺を \theta_0微分したっていいはずだ。やってみよう。

 \displaystyle
\begin{eqnarray}
 \frac{\partial}{\partial \theta_0} \int f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) \; t(\boldsymbol{x}, \boldsymbol{y}) \; d\boldsymbol{y} &\equiv& \frac{\partial}{\partial \theta_0} \theta_0 \\
 \Leftrightarrow \int  t(\boldsymbol{x}, \boldsymbol{y}) \; \frac{\partial}{\partial \theta_0} f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) \; d\boldsymbol{y} &\equiv& 1
\end{eqnarray}

微分した後の式はどんな \thetaを代入しても成り立つので、 =ではなく \equivで表してみた。真面目に書くと、

 \displaystyle
 \forall \theta', \int  t(\boldsymbol{x}, \boldsymbol{y}) \; \left. \frac{\partial}{\partial \theta} f(\boldsymbol{y}; \theta, \boldsymbol{x}) \right|_{\theta=\theta'} \; d\boldsymbol{y} = 1

ということだ。いちいちこんなことを書くのは面倒なので、以降は特に断らない限りこういう恒等式も普通の =で表すことにする。

また、さっきの式変形では何も言わずにこっそり積分微分の順序を入れ替えてしまったが、実はこういうことをしていいのはある正則条件を満たしているときだけだ。ここでは詳しくは触れないが、この正則条件が要求されるのはモデル fの話で、モデルは最適化を行う前に決まっているのだから、ここでは微積の入れ替え可能条件は満たしているものとしよう。
主に \thetaの定義域の端っこで \boldsymbol{y}についての積分が変なことになっていないことが要求されるのだが、そういう変なモデルを考える場合は分散の下限はもっと気持ち悪い値になってしまう(あるいは、存在しない?)という話だ。

さて、この両辺を微分するというワザはかなり使えそうだ。というのも、我々は f(\boldsymbol{y}; \theta, \boldsymbol{x})についての積分で表される恒等式をもう一つ知っているからである。
確率密度関数の総和は1」という式だ。

 \displaystyle
 \int f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) \; d\boldsymbol{y} = 1

これも両辺を \theta_0微分してみる。

 \displaystyle
 \int \frac{\partial}{\partial \theta_0} f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) \; d\boldsymbol{y} = 0

 f微分の総和は0ということがわかった。0は何倍しても0なのだから、これを有効活用しない手はない。とりあえず、不偏推定量 tからその平均 \theta_0を引いておこう。

 \displaystyle
\begin{eqnarray}
 \int  t(\boldsymbol{x}, \boldsymbol{y}) \; \frac{\partial}{\partial \theta_0} f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) \; d\boldsymbol{y} &=& 1 \\
 \Leftrightarrow \int  t(\boldsymbol{x}, \boldsymbol{y}) \; \frac{\partial}{\partial \theta_0} f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) \; d\boldsymbol{y} - 0 \theta_0 &=& 1 - 0 \theta_0 \\
 \Leftrightarrow \int  t(\boldsymbol{x}, \boldsymbol{y}) \; \frac{\partial}{\partial \theta_0} f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) \; d\boldsymbol{y} - \theta_0 \int \frac{\partial}{\partial \theta_0} f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) \; d\boldsymbol{y} &=& 1 - 0 \\
 \Leftrightarrow \int  \left\{ t(\boldsymbol{x}, \boldsymbol{y}) - \theta_0 \right\} \left\{ \frac{\partial}{\partial \theta_0} f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) \right\} d\boldsymbol{y} &=& 1
\end{eqnarray}

ここで、少しトリッキーなことをする。我々はこの積分に意味を見出したいので、なんとか期待値 \int g(\boldsymbol{y}) f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) d\boldsymbol{y}の形に持ち込みたいのだが、ここで対数についての微分の式
 \displaystyle
\begin{eqnarray}
 \left. \frac{d}{dx} \log g(x) \right|_{x=x_0} &=& \frac{1}{g(x_0)} \left. \frac{d}{dx} g(x) \right|_{x=x_0}  \\
 \Leftrightarrow \left. \frac{d}{dx} g(x) \right|_{x=x_0} &=& g(x_0)  \left. \frac{d}{dx} \log g(x) \right|_{x=x_0} 
\end{eqnarray}
を使えば、微分の中に \logが生えてくることを犠牲に、外側に fを引っ張り出してくることができるのだ。早速やってみよう。
この変形は確率の微分を含む積分が出てきたときの常套手段なので、ぜひ覚えて帰ってほしい。

 \displaystyle
\begin{eqnarray}
&&\int  \left\{ t(\boldsymbol{x}, \boldsymbol{y}) - \theta_0 \right\} \left\{ \frac{\partial}{\partial \theta_0} f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) \right\} d\boldsymbol{y} \\
&=& \int  \left\{ t(\boldsymbol{x}, \boldsymbol{y}) - \theta_0 \right\} \left\{ \frac{\partial}{\partial \theta_0} \log f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) \right\} f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) d\boldsymbol{y} \\
&=& \mathbb{E}_{\boldsymbol{y}} \left[ 
 \left\{ 
  t(\boldsymbol{x}, \boldsymbol{y}) - \theta_0
 \right\}
 \left\{
  \frac{\partial}{\partial \theta_0} \log f(\boldsymbol{y}; \theta_0, \boldsymbol{x})
 \right\}
 \right] 
 = 1
\end{eqnarray}

同様のことが、先ほどの \int \frac{\partial}{\partial \theta_0} f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) \; d\boldsymbol{y} = 0についても行える。

 \displaystyle
\begin{eqnarray}
&&\int \frac{\partial}{\partial \theta_0} f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) \; d\boldsymbol{y}\\
&=& \int  \frac{\partial}{\partial \theta_0} \log f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) d\boldsymbol{y} \\
&=&  \mathbb{E}_{\boldsymbol{y}} \left[ 
  \frac{\partial}{\partial \theta_0} \log f(\boldsymbol{y}; \theta_0, \boldsymbol{x})
 \right] 
 = 0
\end{eqnarray}

これらを組み合わせれば、以下のような式が導けるだろう。

 \displaystyle
\begin{eqnarray}
&& 
\mathbb{E}_{\boldsymbol{y}} \left[ 
 \left\{ 
  t(\boldsymbol{x}, \boldsymbol{y}) - \theta_0
 \right\}
 \left\{
  \frac{\partial}{\partial \theta_0} \log f(\boldsymbol{y}; \theta_0, \boldsymbol{x})
 \right\}
 \right]  \\
&=& 
\mathbb{E}_{\boldsymbol{y}} \left[ 
 \left\{ 
  t(\boldsymbol{x}, \boldsymbol{y}) - \theta_0
 \right\}
 \left\{
  \frac{\partial}{\partial \theta_0} \log f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) - 0
 \right\}
 \right]  \\
&=& \mathbb{E}_{\boldsymbol{y}} \left[ 
 \left\{ 
  t(\boldsymbol{x}, \boldsymbol{y}) 
  - \mathbb{E}_{\boldsymbol{y}} \left[
   t(\boldsymbol{x}, \boldsymbol{y})
  \right] 
 \right\}
 \left\{
  \frac{\partial}{\partial \theta_0} \log f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) 
  - \mathbb{E}_{\boldsymbol{y}} \left[
   \frac{\partial}{\partial \theta_0} \log f(\boldsymbol{y}; \theta_0, \boldsymbol{x})
  \right]
 \right\}
 \right] \\
&=& \text{Cov} \left( 
  t(\boldsymbol{x}, \boldsymbol{y}), 
  \frac{\partial}{\partial \theta_0} \log f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) 
 \right) = 1
\end{eqnarray}

期待値というか、これは共分散だ。いま、長いので t(\boldsymbol{x}, \boldsymbol{y})を単に t \frac{\partial}{\partial \theta_0} \log f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) を単に Fと書くことにしてみると、 t Fの共分散が1という式が出てきたことになる。知りたいのは tの分散 V(t)なのだが、共分散と分散の関係性といえば相関係数 r(X, Y) = \text{Cov}(X, Y) / \left(V(X) V(Y)\right)が有名なので、これを用いて式変形してみよう。

 \displaystyle
\begin{eqnarray}
 \text{Cov}(t, F) = r(t, F) V(t) V(F) &=& 1 \\
 \Leftrightarrow r(t, F) &=& \frac{1}{V(t)V(F)}
\end{eqnarray}

相関係数が、不偏推定量 tと対数尤度の微分 Fの分散の積の逆数になるという式が導かれた。
ところが、一般に相関係数 rは絶対値が必ず1以下になることが知られている。これは、内積の二乗は二乗の積より小さくなるというコーシー・シュワルツの不等式から導かれることなのだが、ここではその証明については触れず道具として用いることにする。すなわち、
 \displaystyle
\begin{eqnarray}
 \left| r(t, F) \right| = \frac{1}{V(t)V(F)} &\le& 1 \\
\Leftrightarrow V(t) &\ge& \frac{1}{V(F)}
\end{eqnarray}

となるのである。これは不偏推定量の分散に下界が存在して、その値が V(F)であると主張していることにほかならない。

以上より、任意の不偏推定量の分散には下界が存在し、 V(F)がその下界となることが示された。
めでたしめでたし。

4.3.2.2 分散の下界とフィッシャー情報量

いや、全然めでたくない。分散の下界はフィッシャー情報量の逆数になるという触れ込みだったのに、実際に計算してみたら対数尤度の微分の分散の逆数になってしまった。
これはこれで綺麗な結果ではあるが、求めていたものではない。式変形を続けてフィッシャー情報量を導き出そう。

まず、 Fの期待値 \mathbb{E}_{\boldsymbol{y}} \left[  \frac{\partial}{\partial \theta_0} \log f(\boldsymbol{y}; \theta_0, \boldsymbol{x}) \right] は、上の方で計算したとおり0になるのだった。分散と期待値の関係式
 \displaystyle
V(X) = \mathbb{E}(X^2) - \mathbb{E}(X)^2
を用いれば、
 \displaystyle
\begin{eqnarray}
V(F) &=& \mathbb{E}(F^2) - \mathbb{E}(F)^2 \\
&=& \mathbb{E}(F^2)
\end{eqnarray}
となることがただちにわかる。また、フィッシャー情報量 I(\theta)を変形すれば、
 \displaystyle
\begin{eqnarray}
I(\theta', \boldsymbol{x}) &=& -\mathbb{E}_{\boldsymbol{y} \sim f(\cdot; \theta', \boldsymbol{x})} \left[ \left. \frac{\partial^2}{\partial \theta^2} \log f(\boldsymbol{y}; \theta, \boldsymbol{x}) \right|_{\theta=\theta'} \right] \\
&=& -\mathbb{E}_{\boldsymbol{y} \sim f } \left[ 
  \frac{\partial^2}{\partial \theta^2} \log f 
\right] \\
&=& -\mathbb{E}_{\boldsymbol{y} \sim f } \left[ 
  \frac{\partial}{\partial \theta} \left\{ \frac{\partial}{\partial \theta}  \log f  \right\}
\right] \\
&=& -\mathbb{E}_{\boldsymbol{y} \sim f } \left[ 
  \frac{\partial}{\partial \theta} \left\{ \frac{1}{f} \frac{\partial f}{\partial \theta}  \right\}
\right] \\
&=& -\mathbb{E}_{\boldsymbol{y} \sim f } \left[ 
  \frac{-1}{f^2} \left( \frac{\partial f}{\partial \theta} \right)^2 + \frac{1}{f} \frac{\partial^2 f}{\partial \theta^2}  
\right] \\
&=& \mathbb{E}_{\boldsymbol{y} \sim f } \left[ 
  \left( \frac{1}{f} \frac{\partial f}{\partial \theta} \right)^2
\right] 
 - \mathbb{E}_{\boldsymbol{y} \sim f } \left[ 
  \frac{1}{f} \frac{\partial^2 f}{\partial \theta^2}
\right] \\
&=& \mathbb{E}_{\boldsymbol{y} \sim f } \left[ 
  \left( \frac{\partial }{\partial \theta} \log f \right)^2
\right] 
 - \int
  \frac{1}{f} \frac{\partial^2 f}{\partial \theta^2}
 \; f\; d\boldsymbol{y} \\
&=& \mathbb{E}_{\boldsymbol{y} \sim f } (
  F^2
)
 - \int
  \frac{\partial^2 f}{\partial \theta^2}
 \; d\boldsymbol{y} 
\end{eqnarray}

となるが、第一項 \mathbb{E}(F^2)は先ほどの計算から V(F)に等しく、第二項 -\int \frac{\partial^2 f}{\partial \theta^2} \; d\boldsymbol{y}  -\int \frac{\partial f}{\partial \theta} \; d\boldsymbol{y} = 0の両辺を \theta微分することで簡単に0だとわかる。したがって、
 \displaystyle
\begin{eqnarray}
I(\theta', \boldsymbol{x}) &=& \mathbb{E}_{\boldsymbol{y} \sim f } (
  F^2
)
 - \int
  \frac{\partial^2 f}{\partial \theta^2}
 \; d\boldsymbol{y} \\
&=& \mathbb{E}(F^2) - 0 \\
&=& V(F)
\end{eqnarray}

となって、見事フィッシャー情報量の逆数が不偏推定量の分散の下界になることがわかったのだった。

今日はここまで

途中重い式変形をはさみつつ、駆け足でここまで来てしまった。僕も疲れたし、読者諸氏もかなり体力を使っただろうから、今日はこのあたりで終わっておこうと思う。次項ではクラメール・ラオの下限の具体例を見たり、フィッシャー情報量の独立性を見たりする予定だ。あるいは、いきなり最尤法の良さの話に突入してしまうかもしれない。

今回の記事は(主にmathjaxのせいだが)13800文字もあるらしい。100文字打つのに1分としたって2時間以上掛かっている計算だ。実際はもっと掛かっている。反響がないとパッタリやめてしまうので、続きが読みたいと少しでも思った人はぜひ宣伝をしてほしい。ツイッターリツイートしたりするのが手軽でおすすめだ。下にある星や鳥を押すのだってかまわない。どうぞよろしくお願いします。



前回
kriver-1.hatenablog.com

次回
kriver-1.hatenablog.com