はじめに
つい最近、統計を用いた分析の講義的なものを聴きました。 統計分析は目的に応じて手法を適宜使い分けなければなりません。
その講義では、分散の計算には標本分散ではなく不偏分散を使いましょうという具合の話がありました。きっとその分野では不偏分散を扱ったほうが良い分析ができるということなのでしょう。ところが、多くの人は単にそれを鵜呑みにしているだけだったりしないでしょうか。
不偏分散がなぜ良いのか、あるいは本当に良いのであろうか、ということは一考の余地があります。今回は推定量の良し悪しの基準の1つである不偏性について説明し、分散の推定を例にその不偏性をもたらす推定が必ずしも目的には合わないこともあるということを見ていきます。
不偏性
推定量
まず、仮定した確率分布 $p(X \mid \theta)$ の点推定量 $\hat \theta$ は $N$ 個の確率変数 $\mathbf X = \{X _ 1, \cdots, X _ N \}$ の関数です。すなわちデータを $N$ 個採取することによって $\bf X$ が実現値として値が定まり、その定まった値をゴニョゴニョ計算して推定量を得るというわけです。
もちろん計算の仕方なんで無数にあります。それが役に立つのかどうかはともかくとして、 $\hat \theta = \frac{1}{2} (X _ 1 + X _ N)$ と、最初と最後のデータだけを使って計算しても、それは推定量と呼ばれる資格があります。このように無数に考えうるからこそ、どんな計算の仕方が良いのかを評価する必要があり、その評価の仕方の1つとして "不偏性を有すると望ましい" としているのです。
不偏推定量
さて、$\bf X$ の関数である推定量 $\hat \theta (\mathbf X)$ が不偏性を有するというのは一体どういうことでしょうか。まず $\hat \theta(\mathbf X)$ が $\bf X$ の関数であるので、$\bf X$ の採取をやり直したら、同じ計算の仕方でも結果が少し変わってしまうことに着目します。
$\hat \theta (\mathbf X)$ 自体がバラツキを持っているわけです。もしも $\hat \theta$ がバラツキを持っているとしても、平均的には、推定したかった $\theta$ の値になっていてほしいというのは当然の願いです。すなわち
$$ \mathbb E _ \theta [\hat \theta(\mathbf X)] = \theta $$
となっていてほしいということになります。このような性質を持つ推定量 $\hat \theta (\mathbf X)$ のことを不偏推定量と呼びます。
不偏性を有さない推定
さて、不偏性を有さない推定というのはどういうものでしょうか。これは単純に
$$ \mathbf E _ \theta [\hat \theta(\mathbf X)] = \theta $$
を満たさない推定量 $\hat \theta$ のことを指します。従って、不偏推定量であれば
$$ 0 = \mathbf E _ \theta [\hat \theta(\mathbf X)] - \theta $$
となっているのですが、不偏性を有さない推定量はこの左辺が $0$ ではなく何らかの値が残ってずれてしまっているわけです。従って
$$ {\rm Bias}(\hat \theta) = \mathbf E _ \theta [\hat \theta(\mathbf X)] - \theta $$
と記述し直してバイアスと呼ぶことにします。またバイアスを使えば ${\rm Bias} (\cdot) = 0$ となるような推定量を不偏推定量と呼ぶことになります。このバイアスは小さければ小さいほど良いと考えられ、その中でもバイアスが $0$ のものを特別に不偏推定量呼んでいるわけですね。
正規分布の分散の推定
さて、正規分布の分散パラメータ $\sigma ^ 2$ の推定を実施してみましょう。ここでは既に天下り的に、不偏推定量として
$$ \hat \sigma ^ 2 _ U = \frac{1}{ N - 1} \sum _ {i = 1} ^ N (x _ i - \hat x ) ^ 2 $$
を用いてみます。ここで $\hat x =(1/N) \sum _ {i=1} ^ N x _ i $ という標本平均です。下記のシミュレーションを見てください。
まず真の分布の例として
$$ {\rm Normal}(X \mid \mu=1, \sigma ^ 2 = 4) $$
を利用します。TensorFlow Probabilityではスケールパラメータ(標準偏差)を与える仕様なので下記のようにします。
import tensorflow as tf import tensorflow_probability as tfp import matplotlib.pyplot as plt tfd = tfp.distributions var_true = 4.0 normal = tfd.Normal(loc=1.0, scale=tf.sqrt(var_true))
さて、このような分布から $N$ 個のデータを獲得し、不偏分散を計算するコードを書いてみます。
def calc_var(N): X = normal.sample(100) mean = tf.reduce_mean(X) var = 1/(N-1) * tf.reduce_sum((X - mean)**2) return var
これによって、$100$ 個のデータを採取した場合の不偏分散を見てみましょう。
N = 100 var1 = calc_var(N) var1.numpy() # 4.4412413
この計算は calc_var(N) を実行する度に結果が変わります(データを採取し直しているから当然)。今の計算は
$$ \hat \theta(\mathbf X) = \frac{1}{ N - 1} \sum _ {i = 1} ^ N (x _ i - \hat x ) ^ 2 $$
を1回計算してみたということになります。 さて、不偏分散を計算するという試みを $M$ 回繰り返してみます。その $M$ 回の結果を平均することで
$$ \mathbf E _ \theta [\hat \theta(\mathbf X)] \simeq \frac{1}{M} \sum _ {j = 1} ^ M \hat \theta _j (\mathbf X) $$
を評価することができるというわけです。上記の計算は
def reduce_mean_calc_var(M, N): var_list = [] for _ in range(M): var_list.append(calc_var(N)) return tf.reduce_mean(var_list)
によって実施できます。1回のデータの採取は $100$ 個で、繰り返しを10回行ってみましょう。
N = 100 M = 10 expected_var = reduce_mean_calc_var(M, N) expected_var.numpy() # 3.825538
さて、より正確に期待値を近似するべく繰り返しを増やしてみます。
N = 100 M = 50 expected_var = reduce_mean_calc_var(M, N) expected_var.numpy() # 4.085039
さて、真の値が $4$ なのでいい感じに近づいていますね。下記に$M = 2, ..., 300$ と増やしていたときの真の値との剥離、すなわちバイアスをプロットしてみます。
expected_var_list = [] num_M = range(2, 300) for m in num_M: expected_var_list.append(reduce_mean_calc_var(m, N) - var_true) plt.plot(num_M, expected_var_list) plt.grid()
結果として正規分布に対して、バイアスを持たない推定が不偏分散の計算でできていそうだということが確認できました。
不偏推定量は必ずしも良い推定量ではない
バリアンス
さて、実際、我々が100個のデータを集めて推定するという試みを300回も行うことなんてありえません。 実際には $N$ 個のデータを1回だけ採取して、推定を行わなければならず、不偏性がもたらす性質は、単にその"推定1回1回は期待値の意味でバイアスがありません"というだけなのです。
無論1回1回の推定は真の値を必ずしも捉えることはできないのです。そこで下記のような素朴な指標を考えたくなります。
$$ \hat \theta(\mathbf X) - \theta $$
これは単に、1回の推定量 $\hat \theta$ が $\theta$ とどれくらいずれているかを計算しているものです。こんなものが計算できるなら最初から困らないわけですから、$\theta$ を使わずに似たようなことをやりたくなります。そこで推定量の期待値を使ってやることにします。(先程見てきたとおり、不偏推定量であれば推定量の期待値は $\theta$ に一致するし、そうでなくとも今まで考えてきたバイアスというものでズレを評価できているので素性が良い)
すなわち
$$ \hat \theta(\mathbf X) - \mathbb E [ \hat \theta (\mathbf X) ] $$
という指標を考えてみたくなるのです。これは1回の推定 $\hat \theta$ が推定の期待値からどれくらいずれていると言えるかという指標になります。このズレの大きさに興味がある場合は2乗を取ってやり、
$$ (\hat \theta(\mathbf X) - \mathbb E [ \hat \theta (\mathbf X) ]) ^ 2 $$
を評価することになります。 そしてこの指標について更に期待値を取れば、推定量 $\hat \theta$ がどれくらいバラつくのかを見ることができます。この指標を下記のように書き
$$ {\rm Var} (\hat \theta) = \mathbb E [ (\hat \theta(\mathbf X) - \mathbb E [ \hat \theta (\mathbf X) ] ) ^ 2] $$
これをバリアンスと呼びます。このバリアンスという指標は、1回1回の推定が推定の期待値からどれくらいバラツキを持ちうるのかを評価しており、バイアスと並んで重要な指標になります。
バイアスとバリアンス
さて、バイアスは推定量の期待値が真の値からどれくらいずれているかという指標でした。実際には我々は推定量の期待値を計算することはなく、1回の推定を行うのみであることに注意が必要です。一方でバリアンスは1回1回の推定量が推定量の期待値からどれくらいバラついているかを見る指標でした。こちらはあくまで推定量の期待値からのバラツキなので真の値との比較はできませんが、1回1回の推定がやたらめったら違う値にならないことを望むのならば重要な指標です。
そして、バイアスとバリアンスの2つを合わさえれば、総合的な良さを測れそうです。仮にバイアスが小さくバリアンスも小さければ、推定の期待値が真の値に近く、1回1回の推定は推定の期待値に近いということになり、結果1回1回の推定が真の値に近いことを望めるという最高の性質を持ちます。
ところが、そんな都合の良い推定量は普通は作れません。バイアスとバリアンスは多くの場合トレードオフの関係にあります。
もしもバイアスは小さいがバリアンスが大きいということになれば、推定の期待値は真の値に近いかもしれないが、推定の期待値からのバラツキが大きく1回1回の推定は信用がならんということになりかねません。一方でバリアンスは小さいがバイアスが大きいという場合には、1回1回の推定はその推定の期待値の付近に固まっており推定のバラツキは小さいが、肝心な推定の期待値が真の値から遠いということになります。
無論、不偏推定量とはバイアスが $0$ と言っているだけの話です。
両者を考慮した平均二乗誤差
そこで、この両者を考慮して天下り的ですがパラメータ $\theta$ を推定量 $\hat \theta$ で推定する際の評価指標として ${\rm MSE} (\theta, \hat \theta)$ なる指標を下記で表します。
$$ {\rm MSE}(\theta, \hat \theta) = {\rm Var}(\hat \theta) + \{ {\rm Bias}(\hat \theta) \} ^ 2 $$
これが小さい推定量 $\hat \theta$ を選ぶとバリアンスとバイアスの両者をバランス良く小さくすることができます。ところで、こんな単純な足し合わせ方で良いのでしょうか。実はこれは1回の推定の真の値からのずれ
$$ \hat \theta(\mathbf X) - \theta $$
を考慮したものになっています(バリアンスを見ていく中で、こんなものを評価できるなら最初から苦労しない、と切り捨てたものである)。この値の二乗の期待値を取ることにしてしまえば
$$ \begin{align} \mathbb E _ \theta [ (\hat \theta(\mathbf X) - \theta) ^ 2 ] &= \mathbb E _ \theta [ \{ \hat \theta - ( \mathbb E _ \theta [\hat \theta] - \mathbb E _ \theta [\hat \theta] ) - \theta) \} ^ 2] \\ &= \mathbb E _ \theta [\{ (\hat \theta - \mathbb E _ \theta [\hat \theta] )+ (\mathbb E _ \theta [\hat \theta] - \theta) \} ^ 2]\\ &= \mathbb E _ \theta [(\hat \theta - \mathbb E _ \theta [\hat \theta] ) ^ 2] + \mathbb E _ \theta [(\mathbb E _ \theta [\hat \theta] - \theta) ^ 2] \\ & = {\rm Var}(\hat \theta) + \{ {\rm Bias}(\hat \theta) \} ^ 2 \end{align} $$
となります。
仮に不偏推定量を選んだ場合、$\rm MSE$ の第二項であるバイアスの二乗は $0$ となるため、$\rm MSE = Var$ が成り立ちます。実はパラメータ $\theta$ に対する不偏推定量の作り方は複数ありえます。なので、不偏推定量の中で更に優劣を付けたい場合は $\rm MSE$ の観点で言うと、$\rm Var$ が最も小さいものを選びたいということになります。
正規分布の不偏分散と標本分散の比較
正規分布に対する不偏分散の計算は、"不偏推定量の中では"最も分散が小さいので、バイアスの無い推定がしたい場合には確かに最適なのです。しかし、バイアスがあっても良いから、更にバリアンスの小さな推定量が作れないのか？ということを考えたときには、もっと選択肢があることになります。
今良く知られている不偏分散推定量
$$ \sigma ^ 2 _ U = \frac{1}{N - 1} \sum _ i (x _ i - \bar x) $$
と標本分散
$$ \sigma ^ 2 _ {ML} = \frac{1}{N } \sum _ i (x _ i - \bar x) $$
について検討します。補足として標本分散は正規分布の分散パラメータ $\sigma ^ 2$ の最尤推定量 (Maximum Likelihood estimator) となります。
バイアスとバリアンスの両方を考慮している $\rm MSE$ の観点で評価したときには、正規分布の分散パラメータ $\sigma ^2 $に対する不偏分散推定量 $\sigma ^ 2 _ U$ の $\rm MSE$ を計算をすると
$$ {\rm MSE} (\sigma ^ 2, \sigma ^ 2 _ U) = \frac{2\sigma ^ 4}{N - 1} $$
となります。一方で $\sigma ^ 2 _ {ML}$ の場合は
$$ {\rm MSE} (\sigma ^ 2, \sigma ^ 2 _ {ML}) = \frac{(2N - 1)\sigma ^ 4}{N ^ 2} $$
となります。これらの大小を比較すると、
$$ \begin{align} & {\rm MSE} (\sigma ^ 2, \sigma ^ 2 _ U) - {\rm MSE} (\sigma ^ 2, \sigma ^ 2 _ {ML}) \\ =& \frac{2\sigma ^ 4}{N - 1} - \frac{(2N - 1)\sigma ^ 4}{N ^ 2} \\ =& \frac{ 2n ^ 2 - (2n - 1)(n - 1) }{n ^ 2 (n-1) } \sigma ^ 4 \\ =& \frac{ 3n - 1 }{n ^ 2 (n-1) } \sigma ^ 4 > 0 \end{align} $$
となっており、実は常に $\rm MSE$ の観点では最尤推定量である標本分散の方が良いのです。
では不偏分散を使う意義はないのでしょうか。そんなことはありません。まず $\rm MSE$ という指標が最も良いかが分からないからです。この指標に基づけば上記の議論ができますが、果たして大本となった $(\mathbb E [\hat \theta] - \theta ) ^2$ という測り方が良いのかは議論の余地があります。
また、不偏分散を積極的に使うのは、バイアスを徹底的になくしておき、バリアンスが大きくても良いと思っているときになります。ところでバリアンスは推定量によって決まってくるものですが、そもそも推定がばらつくのはなぜでしょう。もしもサンプルを100個取る際に、毎回正確に同じようなサンプリングが実施できる整った環境である場合には、推定値のバラツキも小さくなることが想定でき、そもそもバリアンスをそんな気にする必要がない可能性もあります。
いずれにしても統計的推定は、分析の目的に応じて使い分けなければならないということになります。
