統計学 学習ノート

はじめに

「最尤推定量を求めよって、どこまで書けば満点なの？」

準1級では、対数尤度を微分して0とおき、推定値を計算しました。ところが、論述の答案で「微分して0とおくと $\hat{\theta} = \dots$」とだけ書くと、それでは足りないことがあります。その点が本当に最大なのかを、確かめていないからです。

霧の中の山登りを想像してください。足元が平らな場所を見つけても、そこが山頂とは限りません。尾根の途中かもしれないし、谷底かもしれませんよね。尤度の山でも、同じことが起きます。

実は、最尤推定の理論は、対数尤度を微分した1つの関数に集約されます。この関数の性質を押さえれば、推定量の計算も、次の記事以降の理論も、すべて見通しよく書けるようになります。

そこで登場するのが、今回の主役 有効スコア関数 です。有効スコア関数とは、対数尤度関数を母数で微分した関数で、単にスコア関数とも呼びます。

前回の記事では、情報を失わない要約である十分統計量を学びました。今回の記事では、尤度関数の書き下ろし方、有効スコア関数の性質の証明、そして最尤推定を答案として書き切る方法について解説します。

尤度関数：モデルから式を書き下ろす

はじめに、尤度関数の定義と書き方について解説します。

準1級の復習です。データ $\boldsymbol{x}$ を固定して、同時確率（密度）を母数 $\theta$ の関数と見たものを 尤度関数 と呼び、その対数を 対数尤度関数 と呼びます。

$$L(\theta) = f(\boldsymbol{x};\theta) = \prod_{i=1}^{n} f(x_i;\theta), \qquad \ell(\theta) = \log L(\theta) = \sum_{i=1}^{n} \log f(x_i;\theta)$$

尤度は、オーディションの点数のようなものです。審査員はデータで、候補者は $\theta$ の値です。「その値が真なら、このデータはどのくらい起きやすいか」で、候補者に点数がつきます。

ここで注意が2つあります。

1つ目は、尤度は $\theta$ の確率ではないということです。候補者全員の点数を足しても、1になる保証はありません。例えば、12回中3回当たる二項分布の尤度を $p$ で積分すると、$\int_0^1 \binom{12}{3} p^3 (1-p)^9\,dp = 1/13$ です。

2つ目は、尤度の定数倍には意味がないということです。点数は比べるためのものなので、全員に同じ数を掛けても順位は変わりません。例えば、次の2つの実験を比べてみましょう。

【同じ尤度をもつ2つの実験】

• 12回引くと決めて、3回当たった： $L(p) = \binom{12}{3} p^3 (1-p)^9$（二項分布）
• 3回当たるまで引いて、12回かかった： $L(p) = \binom{11}{2} p^3 (1-p)^9$（負の二項分布）

実験の設計は違いますが、尤度はともに $p^3(1-p)^9$ に比例します。したがって、最尤推定値はどちらも $\hat{p} = 1/4$ です。

1級の腕の見せどころは、複雑な観測の尤度を自分で書くことです。原則は1つだけで、「観測された事象の確率を、そのまま掛ける」です。

例えば、寿命が率 $\lambda$ の指数分布に従う部品5個を試験し、2, 3, 5 の時点で3個が故障して、残り2個は6の時点で動いたまま試験を打ち切ったとします。故障した部品の情報は密度 $\lambda e^{-\lambda t}$ で、打ち切った部品の情報は「6より長くもった」という確率 $e^{-6\lambda}$ です。

$$L(\lambda) = \prod_{\text{故障}} \lambda e^{-\lambda t_i} \prod_{\text{打ち切り}} e^{-\lambda c_j} = \lambda^{d}\, e^{-\lambda \tau}$$

$d = 3$ は故障数で、$\tau = 2 + 3 + 5 + 6 + 6 = 22$ は総試験時間です。微分すると $\hat{\lambda} = d/\tau = 3/22 \approx 0.136$ で、平均寿命の推定値は約7.3です。打ち切りの2個を無視して故障の3個だけで平均すると3.3になり、長生きした部品を捨てたぶん、大きく過小評価してしまいます。

以上が、尤度関数についての解説です。

有効スコア関数：2つの等式を証明する

次に、最尤法の理論の中心である有効スコア関数について説明します。

有効スコア関数（efficient score function）は、対数尤度関数の傾きです。

$$U(\theta) = \frac{\partial \ell(\theta)}{\partial\theta} = \sum_{i=1}^{n} u(X_i;\theta), \qquad u(x;\theta) = \frac{\partial}{\partial\theta}\log f(x;\theta)$$

$U(\theta) = 0$ が 尤度方程式 です。$\boldsymbol{X}$ が確率変数なので、$U(\theta)$ も確率変数です。ここからが1級の話です。$U$ は次の2つの等式を満たします。

$$E_\theta[U(\theta)] = 0, \qquad V_\theta[U(\theta)] = -E_\theta\left[\frac{\partial U(\theta)}{\partial\theta}\right]$$

イメージはシーソーです。真の $\theta$ の位置で対数尤度の傾きを測ると、データしだいで右上がりにも右下がりにもなります。でも、平均するとちょうど釣り合って、傾きは0です。

等式には、次の 正則条件 が必要です。

【正則条件】

1. $f(x;\theta) > 0$ となる $x$ の範囲（台）が $\theta$ によらない
2. $f$ は $\theta$ について2回微分できる
3. $\int f\,dx$ を $\theta$ で2回微分するとき、微分と積分の順序を交換できる

【証明】

1. $\int f(x;\theta)\,dx = 1$ を $\theta$ で微分する。$\partial f/\partial\theta = u f$ なので、$\int u f\,dx = 0$、すなわち $E_\theta[u] = 0$ である
2. $\int u f\,dx = 0$ をもう一度 $\theta$ で微分すると、$\int (\partial u/\partial\theta) f\,dx + \int u\,(\partial f/\partial\theta)\,dx = 0$ である
3. 第2項は $\int u^2 f\,dx$ なので、$E_\theta[u^2] = -E_\theta[\partial u/\partial\theta]$ を得る
4. $U$ は独立な $u(X_i;\theta)$ の和なので、$E_\theta[U] = 0$ で、$V_\theta[U] = n E_\theta[u^2] = -E_\theta[\partial U/\partial\theta]$ である

この $E_\theta[u^2]$ が、準1級で学んだフィッシャー情報量 $I(\theta)$ です。スコアの揺れの大きさが、データのもつ情報の量になっています。

具体的に見てみましょう。ポアソン分布では $u = x/\lambda - 1$ です。$E[u] = \lambda/\lambda - 1 = 0$ で、$V[u] = V[X]/\lambda^2 = 1/\lambda$ です。もう一方も $-E[\partial u/\partial\lambda] = E[X]/\lambda^2 = 1/\lambda$ で、確かに一致します。

別の例も見てみましょう。正則条件が崩れる $U(0, \theta)$ では、台の中で $u = -1/\theta$ です。つまり $E[u] = -1/\theta \ne 0$ で、等式が成り立ちません。台が $\theta$ とともに動くので、手順1の微分で積分の端の項が余分に出るからです。

例題：スコア関数の性質を確かめる

密度 $f(x;\theta) = (x/\theta)\,e^{-x^2/(2\theta)}$（$x > 0$、$\theta > 0$）からの大きさ $n$ の標本について、有効スコア関数 $U(\theta)$ を求め、$E_\theta[U] = 0$ を確かめよ。さらに $V_\theta[U]$ を求めよ。

【答案】

1. $\log f = \log x - \log\theta - x^2/(2\theta)$ より、$u = -\dfrac{1}{\theta} + \dfrac{x^2}{2\theta^2}$ で、$U(\theta) = -\dfrac{n}{\theta} + \dfrac{\sum_i X_i^2}{2\theta^2}$ である
2. $Y = X^2$ と変換すると、$Y$ の密度は $\frac{1}{2\theta}e^{-y/(2\theta)}$ で、平均 $2\theta$ の指数分布である
3. $E[X^2] = 2\theta$ より、$E[u] = -1/\theta + 2\theta/(2\theta^2) = 0$ である。よって $E_\theta[U] = 0$ である
4. $V[X^2] = 4\theta^2$ より、$V[u] = 4\theta^2/(4\theta^4) = 1/\theta^2$ で、$V_\theta[U] = n/\theta^2$ である

答えは $V_\theta[U] = n/\theta^2$ です。$-E[\partial u/\partial\theta] = -1/\theta^2 + 2\theta/\theta^3 = 1/\theta^2$ でも、同じ値になります。

最後に、母数が複数ある場合です。関心のある母数 $\theta$ のほかに、じゃまな母数 $\eta$（局外母数）があるとします。このとき、$\theta$ のスコアには $\eta$ の不確かさの影響が混ざるので、その分を取り除いた次の量を考えます。

$$\tilde{U}_\theta = U_\theta - I_{\theta\eta}\, I_{\eta\eta}^{-1}\, U_\eta$$

$I_{\theta\eta}$ などは、スコアの共分散行列の成分です。これは、$U_\theta$ を $U_\eta$ に回帰したときの残差にあたります。

例えば、$y_i = \alpha + \beta x_i + \varepsilon_i$（$\sigma^2$ 既知）で、$\beta$ に関心があるとしましょう。$U_\alpha = \sum_i \varepsilon_i/\sigma^2$、$U_\beta = \sum_i x_i\varepsilon_i/\sigma^2$ で、$I_{\beta\alpha}/I_{\alpha\alpha} = \bar{x}$ です。代入すると、$\tilde{U}_\beta = \sum_i (x_i - \bar{x})\varepsilon_i/\sigma^2$ となります。その分散は $S_{xx}/\sigma^2$ で、逆数は最小二乗推定量の分散 $\sigma^2/S_{xx}$ に一致します。切片が分からないぶん、$\beta$ の情報が $\sum_i x_i^2/\sigma^2$ から $S_{xx}/\sigma^2$ に減ったのです。

補足（ここは読み飛ばしてOK）： 「有効スコア」は efficient score の訳で、局外母数がないときは対数尤度の微分そのものです。局外母数があるときは、上の残差を指す本もあります。試験では問題文の定義を確かめましょう。

以上が、有効スコア関数についての説明です。

最尤推定：答案を書き切る

続いて、最尤推定 を答案として書き切る手順を解説します。

【最尤推定の答案の手順】

1. 尤度関数と対数尤度関数を、母数の範囲とともに書く
2. 台が母数で決まるなら微分を使わず、尤度の形を直接見る
3. 尤度方程式 $U(\theta) = 0$ を解く
4. 求めた点が最大であることを示す

手順4の方法は、主に3つあります。全域で $\ell''(\theta) < 0$ を示す（凹関数なら、尤度方程式の解は最大点です）。$U(\theta)$ の符号の変化を調べる。解が複数あるなら、尤度の値を比べる。

さきほどの例題の分布で書いてみましょう。$U(\theta) = (n/\theta^2)(\hat{\theta} - \theta)$（$\hat{\theta} = \sum_i x_i^2/(2n)$）と変形すると、$\theta < \hat{\theta}$ で $U > 0$、$\theta > \hat{\theta}$ で $U < 0$ です。ここで $\ell$ は増加から減少に変わるので、$\hat{\theta}$ は最大点です。$\ell''(\hat{\theta}) < 0$ を1点で確かめただけでは、局所的な最大しか言えない点に注意しましょう。

では、確認をしないと何が起きるのでしょうか？ 位置母数 $\theta$ のコーシー分布 $f = 1/[\pi\{1 + (x-\theta)^2\}]$ で、データが $-3, 3$ の2個だとします。尤度方程式を解くと、解は $\theta = 0$ と $\theta = \pm\sqrt{8} \approx \pm 2.83$ の3つです。ところが、$\ell(0) \approx -4.61$ に対して、$\ell(\pm\sqrt{8}) \approx -3.58$（定数を除く）です。$\theta = 0$ は山頂ではなく谷底でした。冒頭の霧の山が、本当に現れたわけです。

例題：遺伝子型の頻度を推定する

3つの型の確率が $\theta^2$、$2\theta(1-\theta)$、$(1-\theta)^2$（$0 < \theta < 1$）で、度数が $n_1, n_2, n_3$（$n = n_1 + n_2 + n_3$）である。$\theta$ の最尤推定量を導け。$(n_1, n_2, n_3) = (30, 50, 20)$ のときの値も求めよ。

【答案】

1. 定数を除いて、$\ell(\theta) = a\log\theta + b\log(1-\theta)$ である。ここで $a = 2n_1 + n_2$、$b = n_2 + 2n_3$ である
2. $U(\theta) = a/\theta - b/(1-\theta) = 0$ を解くと、$\hat{\theta} = a/(a + b) = (2n_1 + n_2)/(2n)$ である
3. $\ell''(\theta) = -a/\theta^2 - b/(1-\theta)^2 < 0$ なので $\ell$ は凹関数で、$\hat{\theta}$ は最大点である（$a, b > 0$ のとき）
4. 値は $(60 + 50)/200 = 0.55$ である

答えは $\hat{\theta} = (2n_1 + n_2)/(2n)$ で、値は0.55です。数え方の意味は素直で、全部で $2n$ 個ある遺伝子のうちの片方の型の割合です。

次に、微分が効かない場合です。1つずつ見ていきましょう。

1. 尤度に折れ目がある場合

ラプラス分布 $f = \frac{1}{2}e^{-|x-\theta|}$ では、$\ell(\theta) = -\sum_i |x_i - \theta| + \text{定数}$ です。これはデータ点で折れる折れ線で、$\theta$ より大きいデータの個数が小さいデータの個数を上回る間は増加します。したがって、最尤推定量は標本中央値です。

データが 1, 3, 7 なら $\hat{\theta} = 3$ です。データが偶数個のときは、中央の2つの値の間で傾きが0になるので、その間のすべての値が最尤推定値です。最尤推定値は1つとは限りません。

2. 台が母数で決まる場合

密度 $f(x;\alpha, c) = \alpha c^\alpha x^{-\alpha-1}$（$x \ge c$）のパレート分布を考えます。$L = \alpha^n c^{n\alpha} (\prod_i x_i)^{-\alpha-1}\,\mathbf{1}\{c \le x_{(1)}\}$ は $c$ について増加なので、まず $\hat{c} = x_{(1)}$ です。これを代入して $\alpha$ を微分すると、$\hat{\alpha} = n/\sum_i \log(x_i/\hat{c})$ です。

データが 2, 3, 4, 8 なら、$\hat{c} = 2$ で、$\hat{\alpha} = 4/(\log 1.5 + \log 2 + \log 4) \approx 1.61$ です。片方を先に最大化してから、もう片方を解く方法を覚えておきましょう。

例題：尤度方程式を数値で解く

ゼロが観測されないポアソン分布 $P(X = x) = e^{-\lambda}\lambda^x/\{x!\,(1 - e^{-\lambda})\}$（$x = 1, 2, \dots$）がある。(1) $\bar{x}$ を使って尤度方程式を書け。(2) $\bar{x} = 2$ のとき、$\lambda_0 = 2$ から始めて、ニュートン法で $\hat{\lambda}$ を小数第2位まで求めよ。

【答案】

1. $\ell(\lambda)/n = -\lambda + \bar{x}\log\lambda - \log(1 - e^{-\lambda}) + \text{定数}$ を微分して整理すると、次の式を得る

$$\frac{U(\lambda)}{n} = \frac{\bar{x}}{\lambda} - \frac{1}{1 - e^{-\lambda}} = 0$$

2. $U'(\lambda)/n = -\bar{x}/\lambda^2 + e^{-\lambda}/(1 - e^{-\lambda})^2$ を使い、$\lambda_{\text{new}} = \lambda - U(\lambda)/U'(\lambda)$ を繰り返す
3. $\lambda$ の値は、$2 \to 1.509 \to 1.590 \to 1.594 \to 1.594$ と進む

答えは約1.59です。$\lambda/(1 - e^{-\lambda})$ は増加関数なので、$\bar{x} > 1$ なら解はただ1つです。$\bar{x} = \lambda/(1 - e^{-\lambda})$ は、標本平均と母平均をそろえた式でもあります。

試験のポイント：解が陽に書けない問題でも、「尤度方程式を書く」「解の一意性を示す」までは必ず答案にできます。

以上が、最尤推定の答案の書き方についての解説です。

不変性と十分統計量との関係

最後に、最尤推定量がもつ2つの構造を見ていきましょう。

1つ目は 不変性 です。$\hat{\theta}$ が最尤推定量なら、$g(\theta)$ の最尤推定量は $g(\hat{\theta})$ です。円で測って一番高い買い物は、ドルで測っても一番高いですよね。1対1の変換では山頂の位置が動かないので、これは当然です。

$g$ が1対1でないときは、$\tau = g(\theta)$ の尤度を、$L^*(\tau) = \sup_{\theta: g(\theta) = \tau} L(\theta)$ と定めます。

【不変性の証明】

1. $\hat{\tau} = g(\hat{\theta})$ とおく。$\hat{\theta}$ は $g(\theta) = \hat{\tau}$ の範囲に含まれるので、$L^*(\hat{\tau}) \ge L(\hat{\theta})$ である
2. 任意の $\tau$ について、$L^*(\tau) = \sup_{\theta: g(\theta) = \tau} L(\theta) \le \sup_\theta L(\theta) = L(\hat{\theta})$ である
3. したがって $L^*(\tau) \le L^*(\hat{\tau})$ で、$\hat{\tau}$ は $L^*$ の最大点である

具体的に見てみましょう。身長が正規分布に従い、$\hat{\mu} = 170$、$\hat{\sigma} = 5$ だったとします。180を超える確率は2つの母数の関数で、1対1ではありません。それでも、その最尤推定値は $1 - \Phi((180 - 170)/5) = 1 - \Phi(2) \approx 0.023$ です。

ただし、不変性は不偏性を保ちません。例えば、$S^2$ が不偏でも $S$ は不偏でないように、関数をとると期待値はずれます。「最尤推定量は不偏」とは言えないのです。

2つ目は、十分統計量との関係です。分解定理より $\ell(\theta) = \log g(T;\theta) + \log h(\boldsymbol{x})$ で、第2項は $\theta$ によりません。つまり、対数尤度の山の形は $T$ だけで決まります。

$$\hat{\theta} = \arg\max_\theta\, \log g(T;\theta)$$

したがって、最尤推定量がただ1つなら、それは十分統計量の関数です。先ほどの打ち切りの例では、$(d, \tau)$ が十分統計量で、$\hat{\lambda} = d/\tau$ はその関数になっています。

なお、最尤推定量が1つでないと、話は少し変わります。$U(\theta, \theta + 1)$ では、尤度は $x_{(n)} - 1 \le \theta \le x_{(1)}$ で一定の1です。この区間の中のどの点も最尤推定値で、十分統計量と無関係な選び方もできます。その場合は、十分統計量の関数である最尤推定量が「選べる」という言い方が正確です。

以上が、不変性と十分統計量との関係についての解説です。

細かい例を忘れてしまっても問題ありません。ここでは、「スコアの平均は0、分散は情報量」「最大の確認まで書いて答案は完成」という点だけ押さえておきましょう。

今回のまとめ

それでは、今回の内容を振り返ります。

今回は、尤度関数を書き下ろし、有効スコア関数の性質を証明して、最尤推定を答案として書き切る方法を紹介しました。尤度方程式の解は、霧の山の平らな場所にすぎません。山頂だと確かめるところまでが、最尤推定です。

今回のポイントを振り返りましょう。

振り返り

• 尤度は確率ではなく、定数倍に意味がない。複雑な観測でも、観測された事象の確率をそのまま掛ければ書ける。
• 有効スコア関数は対数尤度の傾き。正則条件のもとで、平均は0で、分散は2階微分の期待値の符号を変えたものに等しい。
• 最尤推定の答案は、尤度方程式を解いたあとに最大の確認まで書く。解が複数のこともある。
• 折れ目のある尤度や、台が母数で決まる尤度では、微分に頼らず形を直接見る。
• 不変性は成り立つが不偏性は保たれない。最尤推定量がただ1つなら十分統計量の関数である。

キーワードのおさらい

• 尤度関数: データを固定し、同時確率（密度）を母数の関数と見たもの
• 対数尤度関数: 尤度関数の対数。独立な標本では各観測の対数密度の和
• 有効スコア関数: 対数尤度関数を母数で微分した関数
• 尤度方程式: 有効スコア関数を0とおいた方程式
• 正則条件: 台が母数によらず、微分と積分を交換できるなどの条件
• 局外母数: 関心はないが、モデルに含まれていて未知の母数
• 最尤推定: 尤度関数を最大にする値を推定値とする方法

さて、最尤法は強力ですが、推定量の作り方はこれだけではありません。尤度が書けないとき、分布を仮定したくないときは、どうすればよいのでしょうか？

次回は、モーメント法、最小二乗法、線形推定（BLUE）など、最尤法以外の推定法について紹介します。

それでは、次の記事へ進んでいきましょう。