堅牢な統計

Gelmanらは著書『Bayesian Data Analysis』（2004年）の中で、サイモン・ニューカムが行った光速測定に関するデータセットを取り上げています。この書籍のデータセットは「Classic data sets」ページから入手でき、書籍のウェブサイトにはデータに関する詳細情報が掲載されています。

データの大部分はほぼ正規分布に従っているように見えるが、明らかに外れ値が2つ存在する。これらの外れ値は平均値に大きな影響を与え、平均値を外れ値に向かって、そしてデータの中心から遠ざける。したがって、平均値がデータの中心位置を示す指標として用いられる場合、外れ値が存在すると、ある意味で偏りが生じることになる。

また、平均値の分布は中心極限定理により漸近的に正規分布になることが知られています。しかし、外れ値が存在すると、かなり大きなデータセットであっても平均値の分布は非正規分布になる可能性があります。このような非正規性に加えて、平均値は外れ値が存在する場合には非効率的であり、より変動の少ない位置尺度が利用可能です。

下の図は、光速データの密度プロットとラグプロット（パネル（a））を示しています。また、通常のQ-Qプロット（パネル（b））も示しています。これらのプロットには外れ値が明確に示されています。

図のパネル（c）と（d）は、平均値（c）と10%トリミング平均値（d）のブートストラップ分布を示しています。トリミング平均値は、データの両端から一定割合（ここでは10%）の観測値を削除し、通常の方法で平均値を計算する、シンプルで頑健な位置推定値です。解析はRで行われ、生の平均値とトリミング平均値それぞれについて10,000個のブートストラップサンプルが使用されました。

平均値の分布は、10%トリム平均の分布よりも明らかに広くなっています（グラフは同じスケールです）。また、トリム平均の分布は正規分布に近いように見えますが、生の平均値の分布は左に大きく偏っています。したがって、この66個の観測値のサンプルでは、​​外れ値が2つあるため、中心極限定理は適用できません。

頑健な統計的手法（トリム平均はその簡単な例である）は、外れ値が存在する場合、あるいはより一般的には、基礎となるパラメトリックな仮定が必ずしも正しくない場合に、古典的な統計的手法よりも優れた性能を発揮することを目指している。

この例では、トリミング平均は平均値と比較して良好な結果を示していますが、より頑健な推定値も存在します。実際、平均値、中央値、トリミング平均はすべてM推定量の特殊なケースです。詳細は以下のセクションをご覧ください。

光速データにおける外れ値は、平均値に悪影響を与えるだけでなく、それ以上の影響を及ぼします。通常、尺度の推定値として用いられる標準偏差は、平均値からの偏差の二乗が計算に含まれるため、外れ値の影響をさらに大きく受けます。つまり、外れ値の影響が増幅されるのです。

以下の図は、標準偏差、中央絶対偏差(MAD)、および Rousseeuw–Croux (Qn)スケール推定量のブートストラップ分布を示しています。^{[ 4 ]}これらの図は、各推定量について 10,000 個のブートストラップサンプルに基づいており、再サンプリングされたデータにガウスノイズが追加されています (平滑化ブートストラップ)。パネル (a) は標準偏差の分布、(b) は MAD の分布、(c) は Qn の分布を示しています。

標準偏差の分布は外れ値の影響で不規則かつ広範囲に及んでいます。MADはMADよりも安定した挙動を示し、QnはMADよりもわずかに効率的です。この簡単な例から、外れ値が存在する場合、標準偏差を尺度の推定値として推奨できないことがわかります。

従来、統計学者は手作業でデータから外れ値を選別し、それらを除去し、通常はデータの出所を調べて外れ値が誤って記録されたかどうかを確認していました。実際、上記の光速の例では、さらなる分析に進む前に2つの外れ値を見つけて除去することは容易です。しかし、現代では、データセットは多数の実験単位で測定された多数の変数から構成されることが多く、そのため、外れ値を手作業で選別することは非現実的な場合が多いのです。

外れ値は互いに影響し合い、互いを覆い隠してしまうことがよくあります。簡単な例として、小さな単変量データセットに、中程度の外れ値と大きな外れ値がそれぞれ1つずつ含まれている場合を考えてみましょう。推定標準偏差は、大きな外れ値によって著しく膨らまされます。その結果、中程度の外れ値は比較的正常に見えます。大きな外れ値が除去されると、推定標準偏差は縮小し、今度は中程度の外れ値が異常に見えます。

データの複雑さが増すにつれて、このマスキングの問題は悪化します。例えば、回帰分析では、外れ値を特定するために診断プロットが使用されます。しかし、いくつかの外れ値を除去すると、他の外れ値が目立つようになることがよくあります。この問題は、高次元になるとさらに深刻になります。

堅牢な手法は、外れ値を自動的に検出、重み付けを下げる（または除去する）、およびフラグ付けする方法を提供し、手動によるスクリーニングの必要性を大幅に排除します。注意が必要です。南極上空にオゾンホールが最初に現れたことを示す初期データは、非人間によるスクリーニングによって外れ値として却下されました。^{[ 5 ]}

本稿では単変量統計手法の一般原則について述べるが、回帰問題、一般化線形モデル、および様々な分布のパラメータ推定についても、より堅牢な手法が存在する。

直感的に言えば、推定量の破綻点とは、推定量が誤った（例えば、任意に大きい）結果を出す前に処理できる誤った観測値（例えば、任意に大きい観測値）の割合のことです。通常、漸近的（無限標本）限界が破綻点として引用されますが、有限標本破綻点の方がより有用な場合もあります。^{[ 6 ]}例えば、${\displaystyle n}$独立な確率変数${\displaystyle (X_{1},\dots ,X_{n})}$そしてそれに対応する実現${\displaystyle x_{1},\dots ,x_{n}}$我々は、${\displaystyle {\overline {X_{n}}}:={\frac {X_{1}+\cdots +X_{n}}{n}}}$平均を推定する。このような推定量のブレークダウンポイントは 0 (または有限標本ブレークダウンポイントは${\displaystyle 1/n}$なぜなら私たちは${\displaystyle {\overline {x}}}$任意の値を大きくする${\displaystyle x_{1},\dots ,x_{n}}$。

推定量のブレークダウンポイントが高いほど、その推定量はより頑健になります。直感的に、観測値の半分以上が汚染されている場合、基礎となる分布と汚染された分布を区別することができないため、ブレークダウンポイントは 50% を超えることはできないと理解できます( Rousseeuw & Leroy (1987))。したがって、最大ブレークダウンポイントは 0.5 であり、そのようなブレークダウンポイントを達成する推定量があります。たとえば、中央値のブレークダウンポイントは 0.5 です。X% トリム平均のブレークダウンポイントは、選択された X のレベルに対して X% です。Huber (1981)およびMaronna et al. (2019)には、より詳細な情報が記載されています。検定のレベルと検出力のブレークダウンポイントについては、He, Simpson & Portnoy (1990)で調査されています。

ブレークダウンポイントが高い統計は、耐性統計と呼ばれることがある。^{[ 7 ]}

経験的影響関数は、推定量が標本中のいずれかの点の値にどの程度依存しているかを示す尺度です。これは、異なる標本を用いて推定量を再度計算するだけで済むという意味で、モデルに依存しない尺度です。右側に示されているのは、テューキーのバイウェイト関数です。これは、後述するように、「良い」（定義は後述します）経験的影響関数がどのようなものであるべきかを示す一例です。

数学的に言えば、影響関数は推定量の空間におけるベクトルとして定義され、推定量は母集団の部分集合である標本に対して定義される。

1. ${\displaystyle (\Omega ,{\mathcal {A}},P)}$確率空間である。
2. ${\displaystyle ({\mathcal {X}},\Sigma )}$は測定可能な空間（状態空間）であり、
3. ${\displaystyle \Theta }$次元のパラメータ空間${\displaystyle p\in \mathbb {N} ^{*}}$、
4. ${\displaystyle (\Gamma ,S)}$測定可能な空間であり、

例えば、

1. ${\displaystyle (\Omega ,{\mathcal {A}},P)}$は任意の確率空間であり、
2. ${\displaystyle ({\mathcal {X}},\Sigma )=(\mathbb {R} ,{\mathcal {B}})}$、
3. $\displaystyle \Theta =\mathbb {R} \times \mathbb {R} ^{+}}$
4. ${\displaystyle (\Gamma ,S)=(\mathbb {R} ,{\mathcal {B}})}$、

経験的影響関数は次のように定義される。

させて$\displaystyle n\in \mathbb {N} ^{*}}$そして${\displaystyle X_{1},\dots ,X_{n}:(\Omega ,{\mathcal {A}})\to ({\mathcal {X}},\Sigma )}$はiidであり、${\displaystyle (x_{1},\dots ,x_{n})}$これは、これらの変数からのサンプルです。${\displaystyle T_{n}:({\mathcal {X}}^{n},\Sigma ^{n})\to (\Gamma ,S)}$は推定量である。${\displaystyle i\in \{1,\dots ,n\}}$経験的影響関数${\displaystyle EIF_{i}}$観察時${\displaystyle i}$定義は以下のとおりです。

$${\displaystyle EIF_{i}:x\in {\mathcal {X}}\mapsto n\cdot (T_{n}(x_{1},\dots ,x_{i-1},x,x_{i+1},\dots ,x_{n})-T_{n}(x_{1},\dots ,x_{i-1},x_{i},x_{i+1},\dots ,x_{n}))}$$

これはつまり、サンプル中のi番目の値を任意の値に置き換えて、推定量の出力を調べているということです。あるいは、EIF は、n ではなく n+1 でスケーリングした、点を追加したときの推定量への影響として定義されます。${\displaystyle x}$サンプルへ。

データだけに頼るのではなく、確率変数の分布を利用することもできます。このアプローチは前の段落とはかなり異なります。ここで私たちがしようとしているのは、データの分布を少し変更したときに推定量がどうなるかを確認することです。推定量は分布を仮定し、この分布の変化に対する感度を測定します。対照的に、経験的影響はサンプルセットを仮定し、サンプルの変化に対する感度を測定します。^{[ 8 ]}

させて${\displaystyle A}$上のすべての有限符号付き測度の集合の凸部分集合である。${\displaystyle \Sigma }$パラメータを推定したい${\displaystyle \theta \in \Theta }$分布の${\displaystyle F}$で${\displaystyle A}$関数を${\displaystyle T:A\to \Gamma }$ある推定量列の漸近値とする${\displaystyle (T_{n})_{n\in \mathbb {N} }}$この関数はフィッシャー整合的であると仮定します。${\displaystyle \forall \theta \in \Theta ,T(F_{\theta })=\theta }$これは、モデルにおいて${\displaystyle F}$推定器のシーケンスは漸近的に正しい量を測定する。

させて${\displaystyle G}$何らかの分布がある${\displaystyle A}$データがモデルに従わない場合はどうなるのか？${\displaystyle F}$まさにその通りだが、少し異なる「向かっている」${\displaystyle G}$？

私たちは以下の点を検討しています。 $${\displaystyle dT_{G-F}(F)=\lim _{t\to 0^{+}}{\frac {T(tG+(1-t)F)-T(F)}{t}},}$$

これは、 ​${\displaystyle T}$で${\displaystyle F}$の方向に${\displaystyle G-F}$。

させて${\displaystyle x\in {\mathcal {X}}}$。${\displaystyle \Delta _{x}}$は、質量 1 を与える確率測度です。${\displaystyle \{x\}}$私たちは選びます${\displaystyle G=\Delta _{x}}$影響関数は次のように定義されます。

$${\displaystyle IF(x;T;F):=\lim _{t\to 0^{+}}{\frac {T(t\Delta _{x}+(1-t)F)-T(F)}{t}}.}$$

これは、その地点での微小汚染の影響を説明しています。${\displaystyle x}$我々が求めている推定値について、質量で標準化${\displaystyle t}$汚染（観測値の汚染によって生じる漸近的なバイアス）について。ロバストな推定量の場合、影響関数は有界である必要があり、つまり、x が任意に大きくなっても無限大にならない関数が必要です。

経験的影響関数は経験分布関数を使用する${\displaystyle {\hat {F}}}$分布関数の代わりに${\displaystyle F}$ドロップインの原則を利用して。

影響関数が望ましい性能を発揮するために必要な特性は以下のとおりです。

1. 有限の拒否点${\displaystyle \rho ^{*}}$、
2. 小さな大誤差に対する感度${\displaystyle \gamma ^{*}}$、
3. 局所的なシフトに対する感度が小さい${\displaystyle \lambda ^{*}}$。

M推定量は必ずしも確率密度関数と関連しているわけではありません。したがって、尤度理論から派生した既成の推論手法は、一般的には使用できません。

M推定量は漸近的に正規分布に従うことが示されており、標準誤差が計算できる限り、近似的な推論手法が利用可能である。

M推定量は漸近的にのみ正規分布に従うため、サンプルサイズが小さい場合は、ブートストラップなどの別の推論手法を用いるのが適切かもしれません。しかし、M推定量は必ずしも一意ではありません（つまり、方程式を満たす解が複数存在する可能性があります）。また、特定のブートストラップサンプルには、推定量のブレークダウンポイントよりも多くの外れ値が含まれる可能性もあります。したがって、ブートストラップ手法を設計する際には、注意が必要です。

もちろん、光速の例で見たように、平均値は漸近的に正規分布に従うだけであり、外れ値が存在する場合は、かなり大きなサンプルであっても近似精度が非常に悪くなる可能性があります。しかし、平均値に基づく検定を含む古典的な統計的検定は、通常、検定の公称サイズによって上限が定められています。M推定量の場合はそうではなく、第一種過誤率は公称レベルを大幅に上回る可能性があります。

これらの考察は、M推定をいかなる意味でも「無効にする」ものではありません。単に、他のあらゆる推定方法と同様に、M推定を用いる際には注意が必要であることを明確にしているにすぎません。

M推定量の影響関数は、${\displaystyle T}$に比例する${\displaystyle \psi }$[ ^{9 ]}これは、推定量の特性（棄却点、大誤差感度、局所シフト感度など）は、その推定量が分かっている場合に導出できることを意味します^。${\displaystyle \psi }$関数。

$${\displaystyle IF(x;T,F)=M^{-1}\psi (x,T(F))}$$

と共に${\displaystyle p\times p}$提供元:

$${\displaystyle M=-\int _{\mathcal {X}}\left({\frac {\partial \psi (x,\theta )}{\partial \theta }}\right)_{T(F)}\,dF(x).}$$

多くの実際的な状況では、${\displaystyle \psi }$関数は、良好で頑健な推定値を得る上で重要ではなく、多くの選択肢が同様の結果をもたらし、外れ値が存在する場合の従来の推定値よりも効率とバイアスの点で大幅な改善をもたらします。^{[ 10 ]}

理論的には、${\displaystyle \psi }$関数が推奨され、Tukeyのバイウェイト（バイスクエアとも呼ばれる）関数が一般的な選択肢です。Maronnaら^{[ 11 ]}は、正規分布での効率を85%に設定したバイウェイト関数を推奨しています。

光速データの場合、尖度パラメータを変化させ、尤度を最大化すると、次のようになります。

$${\displaystyle {\hat {\mu }}=27.40,\quad {\hat {\sigma }}=3.81,\quad {\hat {\nu }}=2.13.}$$

修正${\displaystyle \nu =4}$そして、可能性を最大化すると

$${\displaystyle {\hat {\mu }}=27.49,\quad {\hat {\sigma }}=4.51.}$$