統計 10 · 点推定 — 良い推定量とは何か

Chapter 10

点推定 — 良い推定量とは何か

この章がなぜ必要なのか——「もっともらしい値」を原理から決める.

母平均を推定するのに標本平均を使う。これは自然に思える。 だがなぜ標本平均なのか。中央値ではだめなのか。上位 10% を捨てた平均ではだめなのか。

実はどれも「推定量」として使える。ならばどれが良いのか、良いとは何かを決める必要がある。 この章は、その基準(不偏性・一致性・有効性)を定め、 そして推定量を機械的に作り出す方法(最尤法)を与える。

最尤法は 18 章の回帰、21 章のロジスティック回帰、 さらには現代の機械学習(分類モデルの学習で最小化する「クロスエントロピー」は、負の対数尤度と同じものである)まで、すべての基礎になっている。

この章で使う既出の用語(定義は各リンク先). 母集団・母数・標本(01 章 2 節)、標本平均 \(\bar x\)・不偏分散 \(s^2\)(02 章)、外れ値(03 章 1 節)、 オッズ \(a/b\)(03 章 8 節。起きる確率を \(p\) とすれば \(p/(1-p)\) と同じもの)、コーシー・シュワルツの不等式(03 章 3 節)、 確率変数・確率質量関数・確率密度関数(06 章 1〜2 節)、期待値・分散・共分散(06 章 3〜5 節)、モーメント \(E[X^k]\)(06 章 7 節)、 二項・ポアソン分布(07 章)、正規・指数・一様・ガンマ分布(08 章)、 標本分布と標準誤差(09 章 1 節, 5 節)、大数の法則(09 章 2 節)、中心極限定理と分布収束 \(\xrightarrow{d}\)(09 章 3 節)。 母平均 \(\mu\)・母分散 \(\sigma^2\)・成功確率 \(p\)・発生率 \(\lambda\) はいずれも母集団の未知の値(母数、パラメータ)であり、本章ではこれらを総称して \(\theta\) と書く。

1. 推定量と推定値

用語意味例
推定量 (estimator)データから推定値を計算する規則・関数。確率変数\(\bar{X}=\frac{1}{n}\sum X_i\)
推定値 (estimate)実際のデータを入れて出た数値170.2

推定したい母数を \(\theta\)、その推定量を \(\hat\theta\) と書く(ハットは「推定」の印。\(\theta\) が \(\mu\) なら \(\hat\mu\)、\(p\) なら \(\hat p\))。推定量は確率変数なので分布を持つ—— これが 09 章の標本分布であり、良し悪しを議論できる理由である。

2. 良い推定量の 3 条件

(a) 不偏性 (unbiasedness)

\[ \boxed{E[\hat\theta] = \theta} \]

平均的に当たること。系統的なズレ(バイアス)がない。

\[ \text{Bias}(\hat\theta) = E[\hat\theta]-\theta \]

例 1: \(E[\bar{X}]=\mu\) なので標本平均は不偏(06 章)。

例 2: 不偏分散 \(s^2 = \frac{1}{n-1}\sum(X_i-\bar{X})^2\) は \(E[s^2]=\sigma^2\) で不偏(02 章で導出)。 一方 \(\frac{1}{n}\sum(X_i-\bar X)^2\) は \(\frac{n-1}{n}\sigma^2\) となり、過小評価する。

不偏性の落とし穴 (1): 関数を通すと壊れる. \(s^2\) は \(\sigma^2\) の不偏推定量だが、\(s\) は \(\sigma\) の不偏推定量ではない。 イェンセンの不等式(上に凸な関数 \(g\)——グラフ上の 2 点を結ぶ弦がグラフの下に来る関数——について \(E[g(Y)] \le g(E[Y])\)。ばらつきがあれば等号は成り立たない)より、 \(\sqrt{\ }\) は上に凸なので \(E[\sqrt{s^2}] < \sqrt{E[s^2]} = \sigma\)。 標本標準偏差は平均的に真値をわずかに下回る。

一般に \(E[g(\hat\theta)] \neq g(E[\hat\theta])\) である。不偏性は変換で保存されない。

不偏性の落とし穴 (2): 不偏なら良いとは限らない. 極端な例として「最初の 1 個 \(X_1\) だけを使う」推定量も \(E[X_1]=\mu\) で不偏である。 しかし分散は \(\sigma^2\) で、\(\bar X\) の \(\sigma^2/n\) よりずっと大きい。不偏だが役に立たない。

不偏性は「当たりやすさ」の一側面にすぎず、ばらつきも見なければならない。

(b) 一致性 (consistency)

\[ \hat\theta_n \xrightarrow{P} \theta \qquad (n\to\infty) \]

データを増やせば真値に収束すること。\(\xrightarrow{P}\) は確率収束の記号で、「どんな \(\varepsilon>0\) についても \(P(\lvert\hat\theta_n-\theta\rvert>\varepsilon)\to 0\)」という意味。 これは 09 章 2 節 の大数の法則(弱法則)の主張そのものなので、\(\bar X\) は \(\mu\) の一致推定量である。

不偏性より一致性の方が大事だと考える立場が強い. 不偏でなくても、\(n\) を増やせば正しい値に行くなら実用上は問題ない。 逆に、不偏でも収束しないなら使えない。 実際、最尤推定量は一般に不偏ではないが、一致性は持つ。

(c) 有効性 (efficiency)

不偏推定量の中で分散が最小のものが良い。

バイアスのある推定量も含めて「外れ具合」を 1 つの数で測るには、真値からのずれの二乗の期待値——平均二乗誤差 (MSE)——を使う。

\[ \text{MSE}(\hat\theta) = E[(\hat\theta-\theta)^2] = V[\hat\theta] + \left(\text{Bias}(\hat\theta)\right)^2 \]

導出: \(m = E[\hat\theta]\) と置く。\(m\) も \(\theta\) も定数(確率変数ではない)なので、期待値の外に出せる。

\[ E[(\hat\theta-\theta)^2] = E[(\hat\theta - m + m - \theta)^2] = E[(\hat\theta-m)^2] + 2(m-\theta)\underbrace{E[\hat\theta-m]}_{=0} + (m-\theta)^2 \]
\[ = V[\hat\theta] + \text{Bias}^2 \]

∎

\[ \boxed{\text{平均二乗誤差} = \text{分散} + \text{バイアス}^2} \]

不偏推定量どうしならバイアスは 0 なので MSE = 分散であり、「分散が最小」と「MSE が最小」は同じことである。有効性はその意味での最良を指す。

これが「バイアス-バリアンス分解」であり、22 章のモデル選択の中心概念である. 誤差には 2 種類ある——的の中心からずれている(バイアス)か、 ばらばらに散らばっている(バリアンス)か。 そして重要なことに、バイアスを少し許す代わりに分散を大きく減らせば、総合誤差は小さくなる。 20 章の Ridge 回帰、22 章の正則化はすべてこの取引を利用している。

クラメール・ラオの下界

不偏推定量の分散には、原理的な下限がある。母集団の確率質量関数または密度関数を、母数 \(\theta\) を含むことを明示して \(f(x;\theta)\) と書く。 \(X_1,\dots,X_n\) を \(f(x;\theta)\) からの i.i.d. 標本とすると、任意の不偏推定量 \(\hat\theta\) について

\[ V[\hat\theta] \geq \frac{1}{nI(\theta)}, \qquad I(\theta)=E\left[\left(\frac{\partial}{\partial\theta}\log f(X;\theta)\right)^2\right] \]

\(I(\theta)\) をフィッシャー情報量と呼ぶ。この下限を達成する推定量を有効推定量という。 (\(\partial/\partial\theta\) は「複数の変数のうち \(\theta\) だけで微分する」偏微分の記号。ここでは \(f\) が \(x\) と \(\theta\) の 2 変数関数なのでこう書く。\(\theta\) が 1 個の母数なら \(d/d\theta\) と同じ計算である。)

フィッシャー情報量の意味. \(\log f(X;\theta)\) を \(\theta\) の関数と見たもの(3 節 で対数尤度と呼ぶ)の傾きの、二乗の期待値。 傾きが急なら「\(\theta\) が少し変わるだけで \(f\) が大きく変わる」—— つまりデータが \(\theta\) について多くを語っている。だから情報量と呼ぶ。 情報が多いほど分散の下限が小さくなる、という自然な関係になっている。

導出: \(L(\theta) = \prod_i f(X_i;\theta)\)(標本全体の同時密度)、\(S = \frac{\partial}{\partial\theta}\log L(\theta) = \sum_i \frac{\partial}{\partial\theta}\log f(X_i;\theta)\) と置く(スコアと呼ぶ)。 以下、微分と積分(離散なら和)の順序が交換できることを仮定する(これが後述の「正則条件」の中身である)。

  1. \(E[S]=0\): \(\int L\,dx = 1\) を \(\theta\) で微分すると \(0 = \int \frac{\partial L}{\partial\theta}dx = \int \frac{\partial\log L}{\partial\theta}L\,dx = E[S]\)(\(\partial L/\partial\theta = L\cdot\partial\log L/\partial\theta\) を使った)。
  2. \(V[S] = nI(\theta)\): \(S\) は独立な \(n\) 項の和で、各項は平均 0、分散 \(I(\theta)\)(定義そのもの)だから。
  3. \(\hat\theta\) が不偏なら \(\theta = E[\hat\theta] = \int\hat\theta L\,dx\)。両辺を \(\theta\) で微分すると \(1 = \int\hat\theta\frac{\partial L}{\partial\theta}dx = \int \hat\theta S L\,dx = E[\hat\theta S] = \text{Cov}[\hat\theta, S]\)(最後は \(E[S]=0\) より)。
  4. 共分散は \(\text{Cov}^2 \le V[\hat\theta]\,V[S]\) を満たす(相関係数が \(\pm1\) を超えないこと、すなわちコーシー・シュワルツの不等式)。よって \(1 \le V[\hat\theta]\cdot nI(\theta)\)。∎

例: 正規分布 \(N(\mu,\sigma^2)\) の \(\mu\) について、\(\frac{\partial}{\partial\mu}\log f = (x-\mu)/\sigma^2\) なので \(I(\mu) = E[(X-\mu)^2]/\sigma^4 = 1/\sigma^2\)。下限は \(\sigma^2/n\) で、これは \(V[\bar X]\) に等しい。つまり標本平均は有効推定量である。

3. 最尤推定法 (Maximum Likelihood Estimation, MLE)

考え方

「観測されたデータが最も起こりやすくなるようなパラメータを選ぶ」

サイコロを 10 回振って 1 が 8 回出た。このサイコロの「1 が出る確率 \(p\)」はいくつだと思うか。 \(p=1/6\) なら 8 回も出るのは極めて不自然だ。\(p=0.8\) なら自然である。 観測されたことを最もうまく説明する \(p\) を選ぶ——これが最尤法である。

尤度関数

観測データ \(x_1,\dots,x_n\) を固定して、\(\theta\) の関数と見る。

\[ L(\theta) = \prod_{i=1}^n f(x_i;\theta) \]

確率と尤度の違い(重要). 同じ式 \(f(x;\theta)\) を見ているが、

だから「尤度が 0.3」と言っても「確率 30%」ではない。尤度の絶対値に意味はなく、比だけに意味がある。

対数尤度

積は扱いにくい(微分すると積の微分法で \(n\) 個の項が並び、また小さな確率を何百個も掛けると値が極端に小さくなる)ので対数を取り、積を和に変える。\(\log\) は単調増加なので最大化点は変わらない。

\[ \ell(\theta) = \log L(\theta) = \sum_{i=1}^n \log f(x_i;\theta) \]
\[ \hat\theta_{\text{MLE}} = \arg\max_\theta \ell(\theta) \]

通常は \(\dfrac{\partial\ell}{\partial\theta}=0\)(尤度方程式)を解く。

4. 最尤推定の実例

(a) ベルヌーイ(コインの表確率)

\(n\) 回中 \(k\) 回成功。各回の結果 \(x_i\in\{0,1\}\) の確率は \(p^{x_i}(1-p)^{1-x_i}\)(07 章 2 節)なので、観測された並びの確率は

\[ L(p) = \prod_i p^{x_i}(1-p)^{1-x_i} = p^k(1-p)^{n-k}, \qquad \ell(p) = k\log p + (n-k)\log(1-p) \]

(「\(k\) 回成功」という事象の確率なら二項係数 \(\binom{n}{k}\) が掛かるが、それは \(p\) によらない定数なので、掛けても最大化点は変わらない。)

\[ \frac{d\ell}{dp} = \frac{k}{p}-\frac{n-k}{1-p} = 0 \]
\[ k(1-p) = (n-k)p \quad\Longrightarrow\quad k = np \quad\Longrightarrow\quad \boxed{\hat p = \frac{k}{n}} \]

標本比率が出てきた。直感どおりである。

(b) 正規分布(\(\mu\) と \(\sigma^2\))

密度は \(f(x;\mu,\sigma^2) = \dfrac{1}{\sqrt{2\pi\sigma^2}}\exp\left(-\dfrac{(x-\mu)^2}{2\sigma^2}\right)\)(08 章 4 節)なので \(\log f = -\frac12\log(2\pi) - \frac12\log\sigma^2 - \frac{(x-\mu)^2}{2\sigma^2}\)。\(n\) 個足して

\[ \ell(\mu,\sigma^2) = -\frac{n}{2}\log(2\pi) - \frac{n}{2}\log\sigma^2 - \frac{1}{2\sigma^2}\sum(x_i-\mu)^2 \]

\(\mu\) について:

\[ \frac{\partial\ell}{\partial\mu} = \frac{1}{\sigma^2}\sum(x_i-\mu) = 0 \quad\Longrightarrow\quad \boxed{\hat\mu = \bar{x}} \]

\(\sigma^2\) について(\(\sigma\) ではなく \(v=\sigma^2\) を 1 つの変数と見て微分する。\(\log v\) の微分は \(1/v\)、\(1/(2v)\) の微分は \(-1/(2v^2)\)):

\[ \frac{\partial\ell}{\partial\sigma^2} = -\frac{n}{2\sigma^2}+\frac{1}{2\sigma^4}\sum(x_i-\mu)^2 = 0 \]

両辺に \(2\sigma^4\) を掛けて \(\sigma^2 = \frac1n\sum(x_i-\mu)^2\)、\(\mu\) に \(\hat\mu=\bar x\) を入れて

\[ \boxed{\hat\sigma^2 = \frac{1}{n}\sum(x_i-\bar{x})^2} \]

\(\sigma^2\) の最尤推定量は \(n\) で割る形になり、不偏ではない! \(E[\hat\sigma^2] = \frac{n-1}{n}\sigma^2 < \sigma^2\) で、わずかに過小評価する。

つまり最尤法は不偏性を保証しない。それでも最尤法が標準的なのは、 \(n\to\infty\) でバイアスが消える(一致性)ことと、次に述べる漸近的な良い性質のためである。

実務では不偏分散(\(n-1\))を使うのが慣例だが、 「\(n\) で割るのが最尤、\(n-1\) で割るのが不偏」と整理しておけば混乱しない。

(c) ポアソン分布

\(P(X=x) = \lambda^x e^{-\lambda}/x!\)(07 章 4 節)の対数は \(x\log\lambda - \lambda - \log x!\) なので

\[ \ell(\lambda) = \sum\left(x_i\log\lambda - \lambda - \log x_i!\right) = \left(\sum x_i\right)\log\lambda - n\lambda + \text{const} \]
\[ \frac{d\ell}{d\lambda} = \frac{\sum x_i}{\lambda}-n = 0 \quad\Longrightarrow\quad \boxed{\hat\lambda = \bar{x}} \]

(d) 指数分布

密度 \(f(x;\lambda) = \lambda e^{-\lambda x}\)(08 章 5 節)の対数は \(\log\lambda - \lambda x\) なので

\[ \ell(\lambda) = n\log\lambda - \lambda\sum x_i, \qquad \frac{d\ell}{d\lambda}=\frac{n}{\lambda}-\sum x_i = 0 \]
\[ \boxed{\hat\lambda = \frac{1}{\bar x}} \]

08 章 5 節 で示した \(E[X]=1/\lambda\) の関係をそのまま反転した形になっている。

(e) 一様分布 \(\text{Unif}(0,\theta)\) — 微分では解けない例

密度は \(f(x;\theta) = 1/\theta\)(\(0\le x\le\theta\))、それ以外で 0(08 章 3 節)。すべての \(x_i\) が \([0,\theta]\) に入っていないと積が 0 になるので

\[ L(\theta) = \begin{cases}\dfrac{1}{\theta^n} & \theta \geq \max_i x_i \\ 0 & \text{それ以外}\end{cases} \]

\(1/\theta^n\) は \(\theta\) の減少関数なので、制約 \(\theta\ge\max x_i\) を満たす最小の \(\theta\) が最尤。

\[ \boxed{\hat\theta = \max_i x_i} \]

微分して 0 と置く方法が通用しない. 尤度が滑らかでなく、最大が境界にあるからだ。 必ず尤度関数の形を確認することという教訓が得られる。

なお \(M = \max_i X_i\) は不偏でない。\(M\le m\) となるのは全部の \(X_i\) が \(m\) 以下のときなので \(P(M\le m) = (m/\theta)^n\)、微分して密度は \(nm^{n-1}/\theta^n\)。よって \(E[M] = \int_0^\theta m\cdot\frac{nm^{n-1}}{\theta^n}dm = \frac{n}{\theta^n}\cdot\frac{\theta^{n+1}}{n+1} = \frac{n}{n+1}\theta < \theta\)。 不偏にするには \(\frac{n+1}{n}\max x_i\) とすればよい。

5. 最尤推定量の性質

正則条件——母数の取りうる範囲が開区間で、\(f(x;\theta)>0\) となる \(x\) の範囲が \(\theta\) によらず、\(\log f\) が \(\theta\) について 2 回微分でき、微分と積分の順序が交換できる——のもとで、\(n\to\infty\) のとき次が成り立つ。 ((e) の一様分布は「範囲が \(\theta\) に依存する」ので正則条件を破る例である。)

性質内容
一致性\(\hat\theta_{\text{MLE}} \xrightarrow{P} \theta\)
漸近正規性\(\sqrt{n}(\hat\theta-\theta) \xrightarrow{d} N\left(0, \dfrac{1}{I(\theta)}\right)\)(\(\xrightarrow{d}\) は分布収束。09 章 3 節)
漸近有効性分散が \(n\to\infty\) でクラメール・ラオの下界 \(1/(nI(\theta))\) に一致する(最良)
不変性\(g(\theta)\) の最尤推定量(\(\widehat{g(\theta)}\) と書く)は \(g(\hat\theta)\) に等しい

漸近正規性の導出のあらすじ: スコア \(\ell'(\theta) = \sum_i \frac{\partial}{\partial\theta}\log f(X_i;\theta)\) は \(\hat\theta\) で 0 になる(尤度方程式)。\(\theta\) のまわりでテイラー展開して 1 次で止めると

\[ 0 = \ell'(\hat\theta) \approx \ell'(\theta) + \ell''(\theta)(\hat\theta-\theta) \quad\Longrightarrow\quad \sqrt{n}(\hat\theta-\theta) \approx \frac{\ell'(\theta)/\sqrt{n}}{-\ell''(\theta)/n} \]

分子は「平均 0、分散 \(I(\theta)\) の i.i.d. 項の和を \(\sqrt n\) で割ったもの」なので中心極限定理で \(N(0, I(\theta))\) に近づく。 分母は i.i.d. 項の平均なので大数の法則で \(E[-\frac{\partial^2}{\partial\theta^2}\log f(X;\theta)]\) に収束し、これは \(I(\theta)\) に等しい (クラメール・ラオの導出の 1. の式 \(\int \frac{\partial\log f}{\partial\theta} f\,dx = 0\) をもう一度 \(\theta\) で微分すると \(\int\frac{\partial^2\log f}{\partial\theta^2}f\,dx + \int\left(\frac{\partial\log f}{\partial\theta}\right)^2 f\,dx = 0\) となるから)。 よって全体は \(N(0, I(\theta))\) を \(I(\theta)\) で割った \(N(0, 1/I(\theta))\) に近づく。一致性は、この式で \(\hat\theta-\theta\) が \(1/\sqrt n\) の速さで 0 に行くことから従う。

不変性の理由: \(g\) が 1 対 1 なら、\(\eta = g(\theta)\) と置き直した尤度 \(L(g^{-1}(\eta))\) を最大にする \(\eta\) は、\(L(\theta)\) を最大にする \(\theta\) を \(g\) で写したもの、すなわち \(g(\hat\theta)\) である(1 対 1 でない場合も、同じ \(\eta\) を与える \(\theta\) の中の最大値で定義すれば同じ結論になる)。

不変性は実務で非常に便利である. たとえば \(\hat p = 0.3\) なら、 オッズ \(p/(1-p)\) の最尤推定は \(\hat p/(1-\hat p) = 0.43\) とそのまま計算してよい。 不偏性はこの性質を持たない(前述)ので、この点でも最尤法は扱いやすい。

漸近正規性が意味すること. \(n\) が大きければ \(\hat\theta \approx N(\theta, 1/(nI(\theta)))\) となる。 つまり最尤推定量の標準誤差が計算でき、信頼区間が作れる。 統計ソフトが出す係数の標準誤差は、たいていこの理論から来ている。

6. モーメント法

もっと古典的で単純な方法。\(k\) 次の理論モーメント \(E[X^k]\)(母数の式。06 章 7 節)と、データから計算した標本モーメント \(\frac1n\sum x_i^k\) を等しいと置き、母数の個数だけ方程式を作って解く。

\[ E[X] = \bar{x}, \qquad E[X^2] = \frac{1}{n}\sum x_i^2, \qquad \dots \]

例(ガンマ分布 \(\text{Ga}(\alpha,\beta)\)): 08 章 6 節 より \(E[X]=\alpha/\beta\)、\(V[X]=\alpha/\beta^2\)。1 次モーメントと分散(2 次の中心モーメント)を標本の値 \(\bar x, s^2\) と等しいと置く。

\[ \frac{\alpha}{\beta}=\bar{x}, \qquad \frac{\alpha}{\beta^2}=s^2 \]

1 本目を 2 本目で割ると \(\beta = \bar x/s^2\)、これを 1 本目に戻して \(\alpha = \beta\bar x = \bar x^2/s^2\)。

\[ \hat\beta = \frac{\bar{x}}{s^2}, \qquad \hat\alpha = \frac{\bar{x}^2}{s^2} \]
モーメント法最尤法
計算簡単(連立方程式)数値最適化が必要なことも多い
効率劣ることが多い漸近的に最良
使いどころ最尤法の初期値を作る本命

7. その他の推定原理

方法考え方主な用途
最小二乗法残差(観測値とモデルの予測値の差)の平方和を最小化回帰(18 章)。誤差が正規なら最尤と一致
ベイズ推定事後分布(23 章)の平均や最頻値事前情報がある場合
M 推定一般の損失関数(推定値の外れ具合に付ける罰点。二乗なら最小二乗)を最小化ロバスト推定(外れ値の影響を受けにくい推定)
モーメント法モーメントを合わせる簡便・初期値

正規誤差なら最小二乗 = 最尤である. \(y_i = \mu + \varepsilon_i\)、\(\varepsilon\sim N(0,\sigma^2)\) の対数尤度は \(-\frac{1}{2\sigma^2}\sum(y_i-\mu)^2 + \text{const}\) なので、 尤度の最大化と残差平方和の最小化はまったく同じ計算になる。 18 章で最小二乗法を使うとき、裏では正規分布を仮定した最尤推定をしていることになる。

8. まとめ

基準意味
不偏性\(E[\hat\theta]=\theta\)。平均的に当たる
一致性\(n\to\infty\) で真値に収束
有効性不偏の中で分散最小
MSE\(V[\hat\theta] + \text{Bias}^2\)
分布最尤推定量
ベルヌーイ\(\hat p = k/n\)
正規\(\hat\mu=\bar x\), \(\hat\sigma^2=\frac{1}{n}\sum(x_i-\bar x)^2\)
ポアソン\(\hat\lambda=\bar x\)
指数\(\hat\lambda = 1/\bar x\)
一様 \((0,\theta)\)\(\hat\theta=\max x_i\)