第２章 確率と統計的推論

内容

人間社会の活動の結果として生じ複雑に変動する事象は、運動方程式を計算して出てくる答えのような精密な予測をすることはできません。しかし統計的推論という強力な技術の蓄積のおかげで事象の背後にある母集団の特徴を理解する道が開かれています。
　この章では後続の章で必要となる確率・統計の知識を確認していいきます。特に確率の道具立てのところは、あれってどうなってたかな？と振り返ることがよくあるので補足資料もあわせ情報の整理を厚くしました。

1. 母集団と確率変数

本稿が想定する「母集団（Population）」は、物理や工学の世界と同じく調査や実験を通じて値を得ることができることを前提としています。例えば「日本の全ての小学６年生の身長からなる集合」は調査によってデータを収集できるので母集団の例になります。さらに、集団の要素の数が相当に大きくなおかつ観測や調査方法の制約のため要素がとりうる値には散らばりがある（ランダム性がある）、という類の対象を想定しています。

一般に母集団の要素の数が膨大になれば調査や実験の対象を一定数以下に限定せざるを得ません。この限定の仕方に人為的、作為的な偏りが生じないように抽出された部分集団のことを無作為抽出標本（以下では簡単に「標本」）と呼びます。標本のとりうる値全体が「標本空間」です。そして標本を得るための調査や実験のことを本稿では「試行」と呼ぶことにします。

母集団からｎ個の要素を観測する試行があったとしましょう。その要素が示す量が実数だとすれば標本空間は$\mathbb{R}^n$です。試行対象の要素のとりうる値はひとつの分布にしたがうという前提をおき、これを「母集団分布」といいます。「ひとつの分布」とわざわざいう理由は、観測対象の要素が二つの変量をもっている場合、二変量を規定する分布がひとつ存在している、という考え方をするからです。

母集団の特徴を要約した定量的数値が「母数（Modulus）」です。あるいは一般に母集団の特徴を決める量として「母集団のパラメータ」という言葉もよく使われます。母集団 $\{ x_1,..,x_N \}$ の平均や分散は代表的な母数です。式で表現すると、母平均 $\mu$ は $\frac{1}{N} \sum _ {i=1} ^ N x _ i$、母分散 $\sigma ^2$ は $\frac{1}{N} \sum _ {i=1} ^ N (x _ i - \mu) ^ 2$ です。母集団の重複しない要素の数が M 個であり、要素 $x _ i$ の出現頻度が $p _ i$ であったとするとき、一般に $\sum _ {i=1} ^ M (x _ i - \mu) ^ k \cdot p _ i$ を $\mu$ の周りの k 次のモーメント（Moment）と呼称します。母集団の要素が連続的に分布している場合であればモーメントは次のように積分として定義します。$\int _{-\infty }^{\infty }(x-\mu )^{k}p(x)\,dx$ 。母平均はゼロの周りの１次のモーメントであり、母分散は算術平均の周りの２次のモーメントと言い換えることができます。モーメントは統計的推論で使われる最も基本的な母数です。

一方、標本に対しても同様に特徴を要約した定量的数値を考えることができ、これを「統計量（Statistic）」といいます。例えばモーメントを標本に対して同様に定義した統計量を標本モーメントといいます。

1.1. 確率変数

標本の無作為性は標本、さらには標本を集めた集合である標本空間 $\Omega$ の部分集合の出現確率によって特徴づけられます。出現確率が定義された標本空間の部分集合のことを「事象（Event）」と呼称します。標本は根源事象と呼ばれることもあります。「確率変数（Random variable）」とは根源事象（標本）に対してその出現確率を対応づけた関数のことです。統計量を母集団分布に従う確率変数として扱うことには確率論の様々な知識を活用できるメリットがあります。

標本は独立同分布の確率変数の列 $X=(X_1,X_2,..,X_n )$ であるとします。列の各要素を標本変数と呼びます。また、１回の試行により実際に得られた値 $( x_1,x_2,..,x_n )$ を標本とは区別する意味でXの「実現値」とか「観測データ」と呼ぶことがあります。標本変数の算術平均（以後標本平均あるいは単に平均と呼ぶ）$\quad \overline{X}=(X_1+X_2+...+X_n)/n \quad$は統計量であり、$\quad \overline{x}=(x_1+x_2+..+x_n)/n \quad$は$\overline{X}$の実現値です。同様に標本変数から前述したモーメントを構成することができます。標本分散は標本平均の周りの２次のモーメントということができます。

一般に母数や母集団の特徴を決めるパラメータを直接計測することは困難です。なので、正しいと思われる「統計モデル」をもとに統計量からそれらの値を推定する試みが仮説検定の大きな目的であると言えます。統計モデルという用語は時系列分析や機械学習への応用では「予測モデル」と言われることもあるようです。

（少し長い注記）
ここでは標本という言葉から始めて次に確率変数を持ち出しましたが、現代的な確率の定義では、出現確率Pを定義できる事象の集まり $\mathscr{F}$ と根源事象の集まり $\Omega$ を考え三点セット {$\Omega,\mathscr{F},P$} を確率空間とよぶところからスタートします。そしてその空間の上で確率変数の概念を定義する流れになっています。ではあるのですが確率変数という言葉から浮かぶ常識的なイメージ、つまり「観測するたびに不規則に値が変動する量」という理解でも本稿ではあまり問題なさそうなので数学的に精密な説明は省きます。しっかりした説明を知りたい場合は確率論に関する教科書・文献は数多くあるので参照ください。参考書をひとつだけあげておきます（参考１）。確率の意味からはじまってファイナンス数理まで適度なページ数に収まっていて読み物としてもてうれしい構成の本です。
もう少しだけ追記します。上述した出現確率Pを定義できる事象の集まり$\mathscr{F}$はσ-加法族、あるいはσ-集合体と呼ばれ、以下のように確率がもつべき自然な性質を数学的に定式化したものです。$\mathscr{F}$の要素である事象は確率を測ることができる、という意味合いで「可測集合」であるといいます。

$$\begin{equation} \begin{aligned} \text{1.} & \;\emptyset \in \mathscr{F} \\ \text{2.} & \; A \in \mathscr{F} \Rightarrow A ^c \in \mathscr{F} \\ \text{3.} & \; \{ A _ n\} _ {n=1} ^ \infty \subset \mathscr{F} \Rightarrow \bigcup _ {n=1} ^{\infty} A _ n \in \mathscr{F} \end{aligned} \end{equation}$$

なかなか抽象的な概念なのですが、この定義から出発することで特殊な形の事象に対しても確率を積分計算により求めることできる大きなメリットがでてきます。また後の章で説明しますマルチンゲールという重要な概念の説明でも言及することになります。

一変量の確率変数

確率変数Xが離散的な場合（例えばサイコロの目の数）、値xをとる確率を関数として定めることで確率変数の性質は決まります。この関数を$f_X$と書き「確率関数（あるいは確率質量関数,probability mass function,pmf）」と呼びます。Xが従う「確率分布」は確率関数により以下のように表現できます。

$$\begin{equation} P(a \le X \le b) = \sum_{a \le x \le b}f_X(x) \end{equation}$$

なお「Xが確率分布Pに従う」ことを慣習的に以下のように表記します。

$$\begin{equation} X \sim P \end{equation}$$

$\; P(X \le b)\;$、すなわち X が b 以下となる確率を（累積）確率分布関数（cumulative distribution function , cdf）と呼ぶことが多いので少し紛らわしい感じがしますね。
累積分布関数の逆関数、言い換えれば確率 pを与えたときに, $\; P(X \le b)\;=p$ となる確率変数の値 b を返す関数のことをパーセントポイント関数（the Percent-Point Function,ppf）と呼びます。この関数は分位関数（quantile function）あるいは逆累積分布関数（inverse cumulative distribution function,icdf)とも呼ばれ統計分析でよく使われる概念です。分位数については「補足説明 分位数」にも少しまとめておきました。

一方、連続的な確率変数X（例えば身長や体重）に対しては、一点の出現確率はゼロになってしまうので代わりに区間の確率を積分形式で定めることで確率変数の性質は決まります。この被積分関数を$f_X$と書き「確率密度関数（ probability density function , pdf）」と呼びます。確率密度関数を使い、X が従う「確率分布」は以下のように表現できます。

$$\begin{equation} P(a \le X \le b) = \int_a^bf_X(\omega)d\omega \end{equation}$$

多変量の確率変数

サイコロとコインを同時に複数回投げる試行のように、標本の各要素を種類の異なる確率変数からなるベクトルとみなしたい場合があります。その場合標本の要素$\mathbb{X_t}= ( X_{1t},..,X_{nt} )$の同時確率（複数の変数の値を同時に決める意味で同時確率という）を与える確率質量関数$f_X(x_1,..,x_n)$によって性質は決まります。以後標本要素を区別する添え字は誤解の恐れが少ない場合は省略します。$\mathbb{X}$が従う確率分布は、 $A = A_1 \times..\times A_n:=\lbrace (a_1,..,a_n) \mid a_1 \in A_1,..,a_N \in A_n \rbrace$ を事象とすると、同時確率質量関数 $f_X(x_1,..,x_n)$ によって以下のように表されます。

$$\begin{equation} P(\mathbb{X} \in A) = \sum_{x_1 \in A_1}...\sum_{x_n \in A_n}f_X(x_1,..,x_n) \end{equation}$$

なお以後 $\lbrace \mathbb{X} \in A \rbrace$ を $\lbrace X_{1} \in A_1,..,X_{n} \in A_n \rbrace$ のように表すこともあります。

一つの変数 $X_i$ のみの確率分布をみるために他の変数について総和することを周辺化と呼び、その結果の分布を $X_i$ の周辺分布と呼びます。

$$\begin{equation} P(X_1 \in A_1) = \sum_{x_1 \in A_1}\sum_{x_2 \in \mathbb{R}}..\sum_{x_n \in \mathbb{R}}f_X(x_1,..,x_n) \end{equation}$$

工場ラインを流れる部品のサイズと重量を測定する試行のように標本要素が連続値の確率変数となる場合は、標本要素 $\mathbb{X_t}= ( X_1,...,X_n )$ の同時確率密度関数 $f_X(x_1,..,x_n)$ によって性質は決まります。$\mathbb{X_t}$ が従う確率分布は、Aを事象とすると以下のように表されます。

$$\begin{equation} P(\mathbb{X} \in A) = \int_{x_1 \in A_1}..\int_{x_n \in A_n}f_X(x_1,..,x_n)dx_1..dx_n \end{equation}$$

離散変数と同様周辺化、周辺分布を定義します。

$$\begin{equation} P(X_1 \in A_1) = \int_{x_1 \in A_1}\iint_{\mathbb{R} \times .. \times \mathbb{R}}f_X(x_1,..,x_n)dx_1..dx_n \end{equation}$$

条件付確率

確率変数 X 、確率変数 Y とするとき、$Y \in B$ のもとで $X \in A$ となる条件付確率は （同時確率 ÷ 周辺確率）として定義されます。

$$\begin{equation} P(X \in A| Y \in B) \;:= \; \frac{P(X \in A , Y \in B)}{P(Y \in B)} \end{equation}$$

これに対応する条件付確率質量関数（連続値のときは確率密度関数）に関する定義式を以下のように表記します。

$$\begin{equation} f_{X|Y}(x|y) \;:= \; \frac{f_{X,Y}(x,y)}{f_Y(y)} \end{equation}$$

また、$f$ の関数形がパラメータ $\theta$ によって特徴づけられていることを表すために $f(x;\theta)$のように表記することがあります。右辺の分子を同時確率質量関数（連続値のときは同時確率密度関数）、分母を周辺確率質量関数（連続値のときは周辺確率密度関数）と呼びます。

確率変数に関する基本的な定義

確率変数に対しては頻出する重要な量がいくつかあるので以下に定義を書いておきます。
（注）$\;:= \;$ は定義を意味する等号ですが、以降誤解の恐れがない場合には普通の等号を使う場合もあります。

【一変量確率変数】

確率変数Xを平均 $\mu$ の確率分布 P に従うとします。X の期待値 $E(X)$、分散 $V(X)$ 、事象 $\{Y=y\}$ が観測された条件下での確率変数 X の条件付き期待値 $E(X|Y=y)$ は、それぞれ以下のように定義されます。

Xが離散確率変数の場合：

$$\begin{equation} \begin{aligned} E(X) & := \sum_{x\in R:P(x)>0}P(x) \\ V(X) & := E((X-\mu)^2) := \sum_{x\in R:P(x)>0}(x-\mu)^2P(x) \\ E(X|Y=y)\; & := \;\sum_a a P(X=a \mid Y=y) \end{aligned} \end{equation}$$

Xが連続確率変数の場合：

$$\begin{equation} \begin{aligned} E(X) & := \int_{-\infty}^{+\infty}xf_X(x)dx \\ V(X) & := E((X-\mu)^2) := \int_{-\infty}^{+\infty}(x-\mu)^2f_{X}(x)dx \\ E(X|Y=y)\; & := \int_{-\infty}^{+\infty}xf_{X|Y}(x|y)dx \end{aligned} \end{equation}$$

【多変量確率変数】
変量数が2以上の場合の定義も同様なので二変量の場合を記述します。
確率変数のペア $\mathbb{X} = (X,Y)^T$ の各要素は平均 $(\mu_X,\mu_Y)^T$ の同時確率分布 P に従うとします。 $\mathbb{X}$ の期待値ベクトル $E(\mathbb{X})$、分散ベクトル $V(\mathbb{X})$、共分散行列 $Cov(\mathbb{X})$、共分散 $Cov(X,Y)$、相関係数 $\rho_{XY}$ は以下のように定義されます。

$$\begin{equation} \begin{aligned} E(\mathbb{X}) &:= (E(X),E(Y))^T \\ V(\mathbb{X}) &:= (V(X),V(Y))^T \\ Cov(\mathbb{X}) & := E[(\mathbb{X} - E[\mathbb{X}])(\mathbb{X} - E[\mathbb{X}])^T] \\ & := \begin{pmatrix} V(X) & Cov(X,Y) \\ Cov(Y,X) & V(Y) \end{pmatrix} \end{aligned} \end{equation}$$

共分散 $Cov(X,Y)$ は離散確率変数の場合、

$$\begin{equation} Cov(X,Y) := \sum_{x\in R:P(x)>0}\sum_{y\in R:P(y)>0}(x-\mu_X)(y-\mu_Y)P(x,y) \end{equation}$$

X,Y が連続確率変数の場合の共分散は、

$$\begin{equation} Cov(X,Y) := \iint_{-\infty}^{+\infty}(x-\mu_X)(y-\mu_Y)f_{XY}(x,y)dxdy \end{equation}$$

相関係数は、

$$\begin{equation} \rho_{XY} := \frac{Cov(X,Y)}{\sqrt{V(X)}\sqrt{V(Y)}} \end{equation}$$

役に立つ知識

確率変数に関して知っておくと後々確率に関わる計算や定理をフォローするときに役立つ基本知識を解説無しにあげておきます。詳しくは確率の教科書を参照してください。

• 互いに独立の確率変数 X、Y に対して

$$\begin{equation} \begin{aligned} E(XY) = E(X)E(Y) \\ V(X+Y) = V(X) + V(Y) \end{aligned} \end{equation}$$

• 上記加法性をつかうと、互いに独立でかつ同分布の確率変数列 $X_1,X_2,..$ の各変数の平均が $\mu$ 分散が $\sigma$ であるとき、算術平均 $\bar{X} = \frac{1}{n}(X_1+...+X_n)$ に対して

$$\begin{equation} \begin{aligned} E(\bar{X}) = \sigma \\ V(\bar{X}) = \frac{\sigma^2}{n} \end{aligned} \end{equation}$$

• 確率変数 X に対して

$$\begin{equation} \begin{aligned} V(X) = E(X^2) - E(X)^2 \end{aligned} \end{equation}$$

• 確率変数列$X_n,Y_n$に対して

$$\begin{equation} \begin{aligned} X _ n \rightarrow _ {P} X,Y _ n \rightarrow _ {P} Y \text{ならば} \quad X _ n + Y _ n \rightarrow _ {P} X+Y \\ X _ n \rightarrow _ {P} X,Y _ n \rightarrow _ {P} Y \text{ならば} \quad X _ n Y _ n \rightarrow _ {P} XY \\ X _ n \rightarrow _ {P} X,Y _ n \rightarrow _ {P} Y,Y \neq 0 \text{ならば} \quad X _ n / Y _ n \rightarrow _ {P} X / Y \end{aligned} \end{equation}$$

• （連続写像定理）確率変数列 $X_n$ 、連続関数 $g()$ に対して、

$$\begin{equation} \begin{aligned} X_n \rightarrow_{P} X \text{ならば} \quad g(X_n) \rightarrow_{P} g(X) \end{aligned} \end{equation}$$

• 条件付き期待値に関して

$$\begin{equation} \begin{aligned} E(Y) \;= \; \sum_{a\in A}E(\;Y\;|\;X = a \; )\;P(X = a ) \end{aligned} \end{equation}$$

• 事象の集合計算 確率変数 $X$ が値 $x$ を取る事象を $A$ 、$X$ の関数を $t(X)$ とし、$t(X)$ が値 $t(x)$ を取る事象を $B$ とする。このとき、$P(A \cap B) = P(A)$ である。 [証明] $X$ が値 $x$ を取る集合 $A$ と $t(X)$ が値 $t(x)$ を取る集合 $B$ について考える。

$$\begin{equation} \begin{aligned} A &= \{ \omega : X(\omega) = x \} \\ B &= \{ \omega : t(X(\omega)) = t(x) \} \end{aligned} \end{equation}$$

$B$ は関数 $t(X)$ によって定義されており、$X(\omega) = x$ であれば自動的に $t(X(\omega)) = t(x)$ になるので、$\omega \in A$ ならば $\omega \in B$。すなわち $A \subseteq B$ 。したがって、$A \cap B = A$ である。

補足説明 確率・統計のなかで二項分布、正規分布など基本的かつ重要な確率分布の事例を掲載していますので適宜参照ください。

1.2. 大数の法則と中心極限定理

標本 $X = (X_1,X_2,..,X_n )$ の各要素は独立同分布の確率変数でありその平均 $\overline{X}=(X_1+...+X_n)/n$ もまた確率変数です。標本平均 $\overline{X}$ の期待値 $E(\overline{X})$ を $\mu$ 、標本分散の期待値 $V(\overline{X})$ を $\sigma^2$ としましょう。 「大数の弱法則」（以後単に「大数の法則」とよぶ）は標本サイズ n が大きくなるにつれて標本平均 $\overline{X}$ の値と期待値 $\mu$ の差が極微になる確率が 1 に近づくという法則です。式で表すと、

$$\begin{equation} \begin{aligned} \lim_{n \to \infty} P(∣\overline{X}-\mu∣ < \epsilon) = 1 \quad \epsilon \ \text{は正の任意の数} \end{aligned} \end{equation}$$

あるいは用語的な表現では、$\overline{X}$は$\mu$へ確率収束する、と言い以下のように表します。

$$\begin{equation} \begin{aligned} \overline{X} \quad \rightarrow_P \quad \mu \end{aligned} \end{equation}$$

$\overline{X}$ のように正規分布に確率収束する分布に従う確率変数は「漸近正規性」を有すると表現されます。

$\overline{X}$ がしたがう分布が母集団の性質にかかわらず正規分布 $N(\mu,\sigma^2/n)$ へ近づいていくという非常に重要な性質があります。これは数学的に証明することができ古くから「中心極限定理」として知られています。標準化されたスケールで表すと

$$\begin{equation} \begin{aligned} Z = \frac{\overline{X}-\mu}{\sigma / \sqrt{n}}\quad \text{とおいたとき} \end{aligned} \end{equation}$$

$$\begin{equation} \begin{aligned} \lim_{n \to \infty} P(Z < z) = \int_{-\infty}^{z}\frac{1}{\sqrt{2\pi}}e^{-\frac{x^2}{2}} dx \end{aligned} \end{equation}$$

あるいは用語的な表現では、$Z$はN(0,1)へ分布収束する、と言い以下のように表します。

$$\begin{equation} \begin{aligned} Z \quad \rightarrow_d \quad N(0,1) \end{aligned} \end{equation}$$

分布収束は法則収束と呼ばれることもあります。元の集団が正規分布に従っていてもいなくても、標本サイズ（標本の要素数のこと） n が大きくなるにつれ標本平均の分布が正規分布に近づいていていくという事実はきわめて有用です。

大数の法則と中心極限定理について 補足説明 確率・統計 に補足していますので参照ください。

2. 統計的推論の基本知識

2.1. 推定の分類学

「推定（estimation）」とは、統計量をもとに、正しいと考えられる統計モデルに基づいて母集団のパラメータ値を推測する試みです。この推測のために標本において定義した統計量を「推定量（estimator）」とよびます。

「何を」推定するのか、そしてそれを「どのように」推定するのか、という観点で全然網羅的ではありませんが推定の技法の例を以下の表に整理してみました。

何を１	どのように１	どのように２	モデル例
点推定	パラメトリック	モーメント法
点推定	パラメトリック	最小二乗法	線形回帰モデル、ロジスティック回帰モデル、多項式回帰モデル、一般化線形モデル
点推定	パラメトリック	最尤法	線形回帰モデル、一般化線形モデル、一般化線形混合モデル
点推定	パラメトリック	ベイズ推定	階層ベイズモデル
点推定	パラメトリック	ロバスト推定	線形回帰モデル
点推定	パラメトリック	一般化モーメント法	資産価格モデル、動学的パネルデータモデル
点推定	ノンパラメトリック	カーネル密度推定	ガウスカーネル、エパネチニコフカーネル
区間推定	信頼区間		点推定の各種モデル
区間推定	ベイズ信頼区間		ベイズ推定のモデル

• 推定の対象を母数等、母集団のパラメータとするのが点推定、推定の確からしさを数値の幅をもって表現するのが区間推定です。この大別に基づいて以降の節は構成されます。

• パラメトリック推定は、母集団の分布関数を具体的に仮定し、関数形を特定するためのパラメータを推定する手法です。分布関数は正規分布、二項分布、ポアソン分布等、用途に応じて実に様々な関数形が考案されています。一方、ノンパラメトリック分布は特定の関数形を用いず観測データから近似的な分布を構成し、観測数を増やすことで近似の精度を高めていく手法です。これらの手法の特徴と差異について後節で説明します。

• ロバスト推定
  OLS などの従来の方法に比べて、観測結果に現れる「外れ値」外れ値の影響が少なく、極端な観測データが存在する場合でもより正確な推定値を提供できる手法です。本書では特に言及しません。

2.2. 一致性、不偏性、有効性

一般に母集団全体を計測できない限り母数は不明なので「推定」という作業が必要となるわけです。そうすると推定量を計算して得られる推定値がもっともらしいものであることを定量的に測り推定作業の良し悪しを判断する尺度が必要となります。それが表題に現れる三つのキーワードです。

標本サイズが増大するにつれて推定量が真のパラメータに収束することを「統計的推定の一致性」といい、そのような推定量を「一致推定量（Consistent estimator）」と呼称します。前述のように母集団と標本の性質から大数の法則が成立するので標本平均 $\overline{X}$ は母平均の一致推定量です。一致性が無い、とはすなわち真の値に収束しない推定です。そのような量が役に立つとは思えません。したがって推定に用いる統計量を採用する前提として一致性は必須の要件といえるでしょう。

標本平均はその期待値が母平均と一致しています。このように母集団の真のパラメータ値と標本の推定量の期待値が一致する場合、これを「不偏推定量（Unbiased estimator）」と呼びます。また、母集団の真のパラメータ値と推定量の差異を「バイアス（bias）」と呼びます。

一方、標本分散 $V( X ) = \frac{1}{n}\sum_i^n(X _ i -\overline{X})^2$ は一致推定量ですが、その期待値を計算してみると $\frac{(n-1)}{n}\sigma^2$ となります。これは $n \rightarrow \infty$ で母分散に収束するものの、不偏推定量にはなってはいません（証明はこちら）。このように、不偏推定量にはなっていないが標本サイズが大きくなるにつれ推定量の期待値が母集団の真のパラメータに収束する場合その推定量は「漸近不偏性（asymptotic unbiasedness）」を有するといいます。

まとめると、標本サイズを n、z を母集団のパラメータ、$\hat{Z}$をその推定量としたとき以下の定義が重要です。

• 一致推定量：

推定量 $\hat{\theta}$ が真のパラメータ $\theta_0$ に確率収束すること

$$\hat{\theta} \; \underset{p}{\operatorname{\rightarrow}} \theta_0 \;\;(n \rightarrow \infty )$$

• 不偏推定量：

推定量 $Z$ の標本分布による期待値が真のパラメータ $z$ に一致すること

$$E(\hat{Z}) = z$$

• 漸近不偏性：

推定量 $Z$ の標本分布による期待値が真のパラメータ $z$ に確率収束すること

$$E(\hat{Z}) \; \underset{p}{\operatorname{\rightarrow}} \;z$$

不偏推定量は分析のために非常に有用な性質ですが、不偏推定量Aに対して期待値ゼロの統計量Yを足したA＋Bも不偏推定量になります。つまり、ひとつの母数に対してユニークに不偏推定量が決まるわけではない点注意が必要です。ある母数に対して複数の不偏推定量があったときに、選択の指針になるのが「有効性（efficiency）」という概念です。有効性の尺度としては分散の大小を見ると良いといわれています。なおefficiencyは「効率性」と訳されることもあります。

不偏推定量 X の有効性が高い := 分散 V(X) が小さい

上述した標本平均は母平均に対する最小の不偏推定量になっていることを証明できます。このような不偏推定量を「一様最小分散不偏推定量（Uniformly Minimum Variance Unbiased Estimator,UMVUE）」、あるいは簡単に「有効推定量」と呼びます。有効推定量が存在する場合の下限値も知られています。

• クラメール－ラオ（Cramér-Rao）の不等式
  母集団のパラメータ$\theta$に対する不偏推定量を$\hat{\theta}$とするとき、フィッシャー情報量I(Θ)との間に以下の関係が成り立ちます。

$$\begin{equation} V(\hat{\theta}) \quad \ge \quad \frac{1}{I(\theta)} \end{equation}$$

ここでは下限が存在するという事実だけを記載しておきます。フィッシャー情報量など詳細は統計学の教科書を参照してください。

漸近不偏性を有する推定量の分散がクラメール－ラオの下限値に収束するとき「漸近有効性（Asymptotic efficiency）」を有するといいます。

• 漸近有効性の定義：

  $$\begin{equation} n \rightarrow \infty \;\text{then}\; V(\hat{\theta}) \quad \rightarrow \quad \frac{1}{I(\theta)} \end{equation}$$

漸近不偏性があれば大きな標本サイズでの推定の信頼性が高まると言えます。漸近有効性があれば大きな標本サイズでより小さな分散を得ることができます（つまり、より精度の高い推定ができる）。推定量の漸近性はそのような有用性を主張しています。

2.3. 点推定

「点推定（Point Estimation）」は、母集団のパラメータに対して単一の値（真の値）を推定する方法です。一般的な点推定方法には、最尤推定や最小二乗推定などがあります。以下主要な推定手法と具体的な使用例を見ていきましょう。

モーメント法

モーメント法（Method of Moments,MM ）は母集団の統計モデルにおいて成立する k 次モーメント（k=1,2,..）の関係式は標本においても成立する（あるいはほぼ成立する）はず、とのもっともらしい考え方に基づいます。この関係式をモーメント条件式とよび、標本モーメントで置き換えたモーメント条件式を解いて得られるパラメータを母数の推定量とする方法です。

確率変数 $X$ の 統計モデルを使い、母集団における $k$ 次モーメント $E(X^k)`$ が満たす関係式を以下のように表したものを理論的なモーメント条件式と呼ぶことにします。

$$\begin{equation} E(X^k) - m_k(\theta) = 0 \end{equation}$$

ここで、$m_k(\theta)$ は統計モデルに基づいて表現したパラメータ $\theta$ の関数です。 モーメント法による推定は以下のステップで行います。

step1. 理論モーメントの導出:

統計モデルに基づいて、パラメータ $\theta$ に関するモーメント条件式を定義します。たとえば、1次モーメント（平均）および2次モーメント（分散）などです。

$$\begin{equation} \begin{aligned} & E[X」- m_1(\theta) = 0\\ & E[X^2] - m_2(\theta) = 0 \end{aligned} \end{equation}$$

step2. 標本モーメントの計算:

観測データ $\{x_1, x_2, \ldots, x_n\}$ から標本モーメントを計算します。たとえば、標本平均および標本分散です。

$$\begin{equation} \begin{aligned} &\hat{s}_1 := \frac{1}{n} \sum_{i=1}^n x_i \\ &\hat{s}_2 := \frac{1}{n} \sum_{i=1}^n x_i^2 \end{aligned} \end{equation}$$

step3. モーメント条件式の設定:

モーメント条件式を標本モーメントにより置き換えた方程式の解がパラメータの推定値 $\hat{\theta}$ になると考えます。

$$\begin{equation} \begin{aligned} &\hat{s}_1 - m_1(\hat{\theta}) = 0\\ &\hat{s}_2 - m_2(\hat{\theta}) = 0 \end{aligned} \end{equation}$$

一般の次数では以下のように表現できます。

$$\begin{equation} \frac{1}{n} \sum_{i=1}^n x_i^k - m_k(\hat{\theta}) = 0 \end{equation}$$

step4. パラメータの推定:

モーメント条件式を解くことによってパラメータ $\theta$ を推定します。一般的に、方程式の数とパラメータの数が一致するようにモーメント条件式を設定します。

理論モーメントと標本モーメントを一致させるだけなので、モーメント法の実装は比較的簡単です。また、多くの異なる分布やモデルに適用可能です。

一方、モーメント法による推定量は、標本サイズに依存するため、後述する最尤法（MLE）と比較して効率性が劣ることがあります。特に、標本サイズが小さい場合バイアスが大きくなることは否めません。さらに、そもそもモーメント条件式の解が存在しない、または一意でない場合がありますので適用する標本サイズに対する配慮が必要です。

モーメント法による推定の例

例として母集団分布がパラメータ $a,b$ を持つ一様分布を考えてみましょう。この分布の定義式と平均 $\mu$ および分散 $\sigma^2$ は以下です:

$$\begin{equation} \begin{aligned} f(x) = \begin{cases} \frac{1}{b-a} & \text{if } a \leq x \leq b \\ 0 & \text{otherwise} \end{cases} \end{aligned} \end{equation}$$

$$\begin{equation} \mu = \frac{a+b}{2}, \quad \sigma^2 = \frac{(b-a)^2}{12} \end{equation}$$

標本モーメントは算術平均によりあらわされるので、標本を $\{ X _ 1,..,X _ n\}$、 平均、分散の推定量をぞれぞれ $\hat{s}_1,\hat{s}_2$ とすると、

$$\begin{equation} \begin{aligned} \hat{s}_1 &:= \bar{X} = \frac{1}{n}\sum _ {i=1} ^ n X _ i \\ \hat{s}_2 &:= \frac{1}{n}\sum _ {i=1} ^ n (X _ i - \hat{\mu})^2 \end{aligned} \end{equation}$$

標本モーメントに対してモーメント条件式が成立すると考えるので $a,b$ の推定値、$\hat{a},\hat{b}$ は以下の方程式を満たす解です。

$$\begin{equation} \begin{aligned} \frac{1}{n}\sum _ {i=1} ^ n X _ i - \frac{\hat{a}+\hat{b}}{2} = 0 \\ \frac{1}{n}\sum _ {i=1} ^ n (X _ i - \bar{X})^2 - \frac{(\hat{b}-\hat{a})^2}{12} = 0 \end{aligned} \end{equation}$$

簡単な式変形と$a \le \mu \le b$ であることから、以下の推定値を得ます。

$$\begin{equation} \begin{aligned} \hat{a} &= \bar{X} - \sqrt{ \frac{3}{n}\sum _ {i=1} ^ n ( X _ i - \bar{X} ^ 2)} \\ \hat{b} &= \bar{X} + \sqrt{ \frac{3}{n}\sum _ {i=1} ^ n ( X _ i - \bar{X} ^ 2)} \end{aligned} \end{equation}$$

モーメント法の考え方はシンプルで強力ですが、金融・経済分野で使われる統計モデルでは、モーメント条件式を構成したとき、パラメータ数よりもモーメント条件式の数が多くなることがしばしばあります。これは過剰適合（Over indentification）と呼ばれ、ここであげたモーメント法のようにパラメータ数とモーメント条件式の数が等しいことを前提とする手法ではうまく対処できません。

そこでハンセン（1982年）は以下のようにモーメント条件式の考え方を大きく見直すことで、非常に多くの応用に適用できる手法を提案しました。これは一般化モーメント法（Generalized Method of Moments,GMM）とよばれています。 本書ではこの手法は取り扱わないため詳しくは計量経済学の教科書などを参照ください。

一般化モーメント法

モーメント条件式を統計量$\ g _ k (X,\theta) := X^k - m_k$の期待値の形に書き換えてみましょう。

$$\begin{equation} E[g _ k (X,\theta)] = 0 \;\; where \; k=1,2 \end{equation}$$

古典的なモーメント法では各モーメント次数に対してこの方程式の右辺がゼロになることを要請するところが過剰適合の原因の根源でした。そこでこの対策としてモーメント条件式を以下のように緩和する見直しを行います。

モーメント条件式をつくる統計量 $g _ k ( X,\Theta )$ の数が K 個あり、$\Theta = \{\theta_1,..,\theta_L \}$ を推定の対象となる L 個のパラメータ、ただし$L \lt K$ であるとします。

$$\begin{equation} \begin{aligned} g _k ( \Theta ) & := E [ g _ k ( X,\Theta )] \;\; where \; k=1,2,..,K \end{aligned} \end{equation}$$

${g_k}$ に関する重み付きのコスト関数 $g(\Theta)$ を定義します。

$$\begin{equation} \begin{aligned} g( \Theta ) & := \sum _ {k=1} ^ K g _k( \Theta ) ^ T w _ k(\Theta) g _k( \Theta ) \end{aligned} \end{equation}$$

$w_k(\Theta)$ は重み係数として新に導入した正値の量です。これも新たなパラメータとして最適な値を推定する対象となります。

このもとで、$g(\Theta)$ を最小化する $\hat{\Theta}$ を推定量とします。

$$\begin{equation} \begin{aligned} \hat{ \Theta} & := \underset{\Theta}{\operatorname{argmin}} \; g( \Theta) \end{aligned} \end{equation}$$

前述の一様分布に対するモーメント法の例は、重み係数が 1 で、各モーメント条件式がゼロになること、 $g_1 = 0,g_2 = 0$　とみなすことができます。

最小二乗法推定

最小二乗法推定（Ordinary Least Squares Estimation, OLSEあるいは単にLSE）は、観測データである実現値とモデルの予測値との誤差の二乗和を最小化するパラメータを母集団のパラメータの推定値とする手法です。幅広い予測モデルへ適用できるポピュラーな手法です。以下に線形回帰モデルを例とし、モデル式の切片と回帰係数を母集団のパラメータと考え最小二乗法により推定するプログラムを statsmodels のサンプルコードを紹介します。

# 
# stat_estim_ols1.py
# テスト生成データを使い線形回帰モデルに対して最小二乗法推定を行う
import numpy as np
import matplotlib.pyplot as plt
import statsmodels.api as sm
import japanize_matplotlib

# サンプルデータ生成
np.random.seed(9876789)
X0 = np.random.normal(0, 0.9, size=100)  # 正規分布に従った数値を観測データとして100個生成
y = 2 * X0 + 5 + np.random.randn(100)  # 真の直線が傾き2、切片5になるよう係数を調整

# 線形回帰モデルクラス（OLS:Ordinary Least Squares）による最小二乗法推定
X = sm.add_constant(X0)  # 説明変数として定数項を第0列に挿入
res = sm.OLS(y, X).fit() # 予測を実行
print(res.summary()) # パラメータ推定結果の表示

# 予測直線のパラメータを取得
intercept, slope = res.params

# プロット用のデータ
x_values = np.linspace(-3, 3, 100)  # プロット用のX値（-2から2まで）
y_predicted = intercept + slope * x_values  # 予測直線の式
Y = 2 * np.linspace(-3, 3, 100) + 5  # 真の直線の式

# プロット
plt.scatter(X[:, 1], y, label='観測データ')  # サンプルデータ点の散布図
plt.plot(x_values, Y, color='blue', label='真の直線')  # 真の直線
plt.plot(x_values, y_predicted, color='red', label='予測直線')  # 予測直線
plt.xlabel('X')
plt.ylabel('y')
plt.legend()
plt.title('線形回帰')
plt.show()

これを実行したときの処理結果レポート

レポートの各項目の意味については 補足説明 ソフトウェアを準備する に説明をつけているので参照ください。

回帰モデルにおける最小二乗法推定が非常によい性質を持っていることが知られています。これについては第３章で補足します。

最尤推定法

母集団分布の形 $f(x;θ)$ はわかっているがパラメータ $\theta$ は不明であるようなケースにおいて、標本の統計量から $\theta$ を推定したいと考えたとします。$f$ が正規分布ならば $\theta$ として平均 $\mu$ や分散 $\sigma^2$ を想像してください。「最尤推定法（Maximum Likelihood Estimation, MLE）」では、試行により標本 X の観測データ $x=(x_1,...,x_n)$が得られたとしたとき、x を得る確率 $L(\theta;x)$ を定義しこれがが最大となる $\theta$ を母数の推定値と考えます。これを「最尤推定値」（＝もっともらしい推定値）と呼称します。式で表すと

$$\begin{equation} \theta_{MLE}=\underset{\theta}{\operatorname{argmax}} \; L(\theta;x) \end{equation}$$

ここで $\underset{\theta}{\operatorname{\argmax}}$: A は、Aを$\theta$の関数としたとき値が最大になる$\theta$のこと

母集団のパラメータが一つではなく複数の場合は多変量関数$L(\theta_1,...,\theta_m)$に対して同様に考えます。点推定は計算が比較的容易であり、1つの値を得ることができるため直感的です。ただし、推定の不確実性についての情報を提供しないので、不確実性の度合いを知りたい場合は後述の区間推定を使うことになります。

最尤推定法の基本的なステップを示します。

step1. 尤度関数の設定：

まず、母集団のパラメータ $\theta$ を変数とみなして、確率変数 X の観測データ x が得られる確率を $\theta$ の関数 $L(x;\theta)$ として定義します。これを「尤度関数（Likelihood Function）」と呼びます。この関数は、与えられたパラメータに対して、観測データが得られる確率を与えるものであるため、母集団の確率密度関数（pdf）、あるいは離散型確率変数の時は確率質量関数、（pmf）の積として構成することができます。なお記述を短くするために以下の議論は連続型確率変数を前提とします。観測データを$x=\{x_1,...,x_n\}$、確率密度関数を $f_X(x_i;\theta)$ とすると

$$\begin{equation} L(\theta;x) = f _ X (x_1;θ)f _ X (x_2;θ)...f _ X (x_n;θ) = \prod_{i=1}^{n}f _ X (x_i;θ) \end{equation}$$

例えば $X$ が正規分布に従う場合には

$$\begin{equation} L(\mu,\sigma^2;x)= \frac{1}{\sqrt{2\pi\sigma^2}}\exp(-\frac{\sum_{i=1}^n(x_i-\mu)^2}{2\sigma^2}) \end{equation}$$

step2. 対数尤度関数への変換：

尤度関数を対数変換して計算を簡単にし、数値的な見通しをよくすることができます。このステップにより、母集団分布関数の積が対数和に変換されます。

$$\begin{equation} \log L(\mu,\sigma^2;x) = -\frac{n}{2}\log( 2\pi ) -\frac{n}{2}\log( \sigma^2 ) -\frac{1}{2\sigma^2}\sum_{i=1}^n(x_i-\mu)^2 \end{equation}$$

step3. argmaxの計算：

対数尤度関数を最大化するパラメータ $\theta$ の値を見つける。このためには微分や勾配降下法などの最適化アルゴリズムを活用します。最尤推定では、尤度関数の対数尤度をパラメータ（上記例は $\mu$ と $σ^2$ ）について偏微分し、それをゼロと等置した方程式を解くことになります。
実際に尤度関数の $\mu$ に関する偏微分方程式を解いてみると $\mu$ の最尤推定値 $\hat{\mu}$ は観測データの平均に一致することがわかります。

$$\begin{equation} \frac{ \partial }{ \partial \mu} \log ( L ( \mu,\sigma^2;x) ) = 0 \end{equation}$$

$\mu$ に関する定数項は消えるので

$$\begin{equation} -\frac{1}{\sigma^2}\sum_{i=1}^n(x_i-\mu) = 0 \end{equation}$$

これを変形することで最尤推定値 $\hat{\mu}$ が求まりました。

$$\begin{equation} \hat{\mu} = \frac{1}{n}\sum_{i=1}^n(x_i) \end{equation}$$

同じように $\sigma^2$ に関して計算すると $\sigma^2$ の最尤推定値 $\hat{\sigma^2}$ は観測データの分散に一致することがわかります。

$$\begin{equation} \frac{\partial }{\partial \sigma^2}\log (L(\mu,\sigma^2;x)) = 0 \end{equation}$$

$\log(p)$ の微分が1/pであることを思い出すと

$$\begin{equation} -\frac{n}{2}\frac{1}{\sigma^2}+\frac{1}{(\sigma^2)^2}\sum_{i=1}^n(x_i-\mu) = 0 \end{equation}$$

これを変形することで最尤推定値 $\hat{\sigma^2}$ が求まりました。

$$\begin{equation} \hat{\sigma^2} = \frac{1}{n}\sum_{i=1}^n(x_i-\mu)^2 \end{equation}$$

step4. 推定されたパラメータの利用：

最尤推定法によって得られたパラメータの値を使用して、統計モデルの予測、推論、モデル適合度の評価などを行います。

科学技術計算用ライブラリ scipy を使って実際に最尤推定を行ってみましょう。まずscipyをインストールします。

pip install scipy

以下は scipy の minimize 関数を活用した最尤推定のサンプルコードです。

# stat_estim_mle1.py
# 正規分布の平均を最尤推定する
#
import numpy as np
from scipy.stats import norm
from scipy.optimize import minimize

# 標本の生成
loc = 5 # 分布の平均値
scale = 2 # 分布の標準偏差
np.random.seed(0)
data = np.random.normal(loc=loc, scale=scale, size=100)

# 尤度関数から対数尤度関数へ変換
def negative_log_likelihood(params, data):
    mu = params[0]  # 推定したい平均
    # norm.pdf 確率密度関数
    # loc 平均値 scale 標準偏差
    likelihood = norm.pdf(data, loc=mu,scale=1)
    log_likelihood = np.log(likelihood)
    return -np.sum(log_likelihood)

# 最適化を用いて最尤推定を行う
initial_guess = [0.0]  # 初期値
result = minimize(negative_log_likelihood, initial_guess, args=(data,))
mu_mle = result.x[0]  # 最尤推定された平均

print(f"ランダム生成した標本の平均値: {loc} 最尤推定された平均: {mu_mle:.2f}")

最尤推定が示す重要な性質

最尤推定は、先にあげた推定における重要な三つの条件（一致性、不偏性、有効性）を含め、以下のように優れた特徴を示すことが知られています。

一致性

最尤推定量 $\hat{\theta}$ は Wald(1949) によって正則性条件としてまとめられたいくつかの条件のもとで一致性を示すことが知られています。すなわち $\hat{\theta}$ は真のパラメータ $\theta_0$ に確率収束します。

$$\hat{\theta} \; \underset{p}{\operatorname{\rightarrow}} \theta_0 \;\;(n \rightarrow \infty )$$

一致性のための正則性条件

• パラメータ空間がコンパクトであり、最尤推定量は有限となること。
• 尤度関数が識別可能である、言い換えれば異なるパラメータが異なる確率分布を生成すること。
• 標本対数尤度関数 $\frac{1}{n} \log L_n(\theta)$ が、期待対数尤度関数 $E[\log L(\theta)]$ に一様に収束すること。
• 真のパラメータで連続的に最大化されること。$\lim_{\theta \to \theta_0} \frac{1}{n} \log L_n(\theta) = \frac{1}{n} \log L_n$
• 尤度関数 $L_n(\theta)$ が、パラメータ $\theta$ の連続関数であること

漸近正規性

漸近正規性（Asymptotic Normality）とは、以下に記す正則性条件のもとで、標本サイズが大きくなるにつれて、最尤推定量の分布が正規分布に近づく性質を指します。

$$\begin{equation} \sqrt{n} (\hat{\theta} - \theta_0) \xrightarrow{d} N(0, I^{-1}(\theta_0)) \end{equation}$$

漸近正規性のための正則性条件 - 対数尤度関数が一階および二階微分可能であり、真のパラメータ周りで適切な近似が可能であること。
- 標本サイズが増えると、中心極限定理によりスコア関数（対数尤度関数の一階微分）が平均ゼロの正規分布に収束すること。
- フィッシャー情報量行列が正定値であること。これは情報量行列が逆行列を持つことを意味するとともに、最尤推定量が、真のパラメータ $0_0$ の周りで平均 0、分散 $𝐼(\theta_0)^{-1}$ の正規分布へ漸近するための条件となります。 - 標本が独立同分布（i.i.d.）であること。

この正則性条件は、真のパラメータの周りで正規分布に従うための条件であり、尤度関数の滑らかさやスコア関数の正則性が必要となります。

この事実は、最尤推定量が大きな標本サイズにおいて、平均 $\theta$、分散 $I^{-1}(\theta) / n$ の正規分布に従うことを意味します。漸近正規性は、標本サイズが十分に大きい場合に適用されるため、推定値の信頼区間や仮説検定に利用されます。

漸近不偏性

漸近不偏性は、十分大きな標本サイズで推定量が不偏性を示すとみなせることを主張します。つまり、

$$\lim_{n \to \infty} E[\hat{\theta}] = \theta_0$$

一般に最尤推定量は、不偏推定量にはなりませんが、正則性条件の下で漸近不偏性を満足することを証明することができます。

この事実は、前出の漸近正規性からの帰結です。

$$\sqrt{n} (\hat{\theta} - \theta_0) \xrightarrow{d} N\left(0, I(\theta_0)^{-1}\right)$$

従って、

$$\hat{\theta} \simeq \theta_0 + \frac{1}{\sqrt{n}} Z$$

ここで、 $Z \sim N(0, I(\theta_0)^{-1})$ です。

期待値をとれば、

$$E[\hat{\theta}] \simeq E\left[\theta_0 + \frac{1}{\sqrt{n}} Z\right] = \theta_0 + \frac{1}{\sqrt{n}} E[Z]$$

$n \rightarrow \infty$ で右辺の２項目はゼロへ収束するので、

$$E[\hat{\theta}] \rightarrow \theta_0$$

漸近有効性

Cramér-Raoの不等式は、不偏推定量の分散には下限値が存在することを示すものです。

$$\begin{equation} V(\hat{\theta}) \quad \ge \quad \frac{1}{I(\theta_0)} \end{equation}$$

先に示した漸近不偏性は、最尤推定量の分散が下限値へ収束する、漸近有効性を有することを主張しています。

$$\begin{equation} n \rightarrow \infty \;\text{then}\; V(\hat{\theta}) \quad \rightarrow \quad \frac{1}{I(\theta_0)} \end{equation}$$

十分統計量

まず最尤推定に関わる重要な統計量である「十分統計量」、そして十分統計量に関する重要な定理を紹介しておきます。

【十分統計量の定義】

標本 X の統計量 $T=t(X)$ に関して、$T$ の実現値 $t(x)$ のもとで $X$ が実現値 $x$ を得る条件付確率 $P _ \theta (X=x \mid t(X)=t(x))$ がパラメータ $\theta$ に依存しない時、そしてその時に限り、$T$ はパラメータ $\theta$ の十分統計量（Sufficient statistic）であるという。あるいは、$T$ が $\theta$ のために十分であるという。

Tが十分統計量であることを式で表現すると、

$$\begin{equation} \begin{aligned} & \Longleftrightarrow P_ \theta (X=x\mid t(X)=t(x))=P(X=x\mid t(X)=t(x)) \\ & \Longleftrightarrow f_{X|T}(x|t(x);\theta) = f_{X|T}(x|t(x)) \end{aligned} \end{equation}$$

十分統計量は観測データから母数 $\theta$ を推定するために必要な情報を保持しており、他の情報を追加で利用する必要がない統計量であることを示すものです。

例として母集団分布がベルヌーイ分布の場合の十分統計量を考えてみましょう。 独立したベルヌーイ試行を n 回繰り返す試行の標本を $X=(X_1, X_2, \ldots, X_n)$ とします。 ベルヌーイ試行の確率質量関数は次の通りです。

$$\begin{equation} f_X(x_i; p) = p^{x_i} (1-p)^{1-x_i} \end{equation}$$

標本 $X$ に対して統計量 $t(X) = \sum_{i=1}^n X_i$ を導入します。標本Xが観測データxを示す確率はベルヌーイ分布の結合分布として表現でき、その結果 $X$ と $t(X)$ の同時確率の形に書き直すことができます。これはパラメータ $p$ の関数と考えれば尤度関数 $L(p;X)$ と同じです。

$$\begin{equation} \begin{aligned} P _ p (X=x) &= \prod_{i=1}^n p^{x_i} (1-p)^{1-x_i} \\ &= p^{\sum_{i=1}^n x_i} (1-p)^{n - \sum_{i=1}^n x_i} \\ &= p^{t(x)} (1-p)^{n - t(x)} \\ &= p^{t(x)} (1-p)^{n - t(x)} \mathbf{1}(t(X)=t(x)) \\ &= P _ p (X=x,t(X)=t(x)) \end{aligned} \end{equation}$$

上記で $\mathbf{1}(t(X)=t(x))$ は $t(X) = t(x)$ のとき１, それ以外のときは ０ の値をとるインディケーター関数です。このように結合分布が標本Xの個々観測データ $x_i$ に依存せず統計量 $t(X)$ の関数として表現できているところがポイントです。

条件付分布の定義から、

$$\begin{equation} \begin{aligned} P _ p (X = x | t(X) = t(x)) &= \frac{P _ p (X = x, t(X) = t(x))}{P _ p (t(X) = t(x))} \\ &= \frac{P _ p (X = x )}{P _ p (t(X) = t(x))} \end{aligned} \end{equation}$$

分母の $P _ p (t(X) = t(x) )$ は、ベルヌーイ分布の性質により次のように求められます。

$$\begin{equation} P _ p (t(X) = t(x)) = \binom{n}{t(x)} p^{t(x)} (1-p)^{n-t(x)} \end{equation}$$

これを用いると、条件付き分布は次のように書けます。

$$\begin{equation} \begin{aligned} P_p(X = x | t(X) = t(x)) &= \frac{p^{t(x)} (1-p)^{n - t(x)} \mathbf{1}(t(X)=t(x)) }{\binom{n}{t(x)} p^{t(x)} (1-p)^{n-t(x)}} \\ &= \frac{\mathbf{1}(t(X)=t(x))}{\binom{n}{t(x)}} \end{aligned} \end{equation}$$

これは、パラメータ $p$ に依存していないことがわかります。従って、 $t(X) = \sum_{i=1}^n X_i$ は $p$ に対して十分統計量です。

次の因子分解定理は、十分統計量を見つけるための方法を与えてくれます。

Fisher-Neymanの因子分解定理

X の確率密度関数（離散的な場合には確率質量関数）を $f_X(x ;\theta)$ とすると、ある関数 g と h が次の関係にある場合、かつその場合に限り、$T$ は $\theta$ に対する十分統計量である。

$$\begin{equation} f_X(x_1,..,x_n;\theta) = g(t(x); \theta) h(x) \end{equation}$$

証明は補足説明　確率・統計を参照ください。

パラメータ $\theta$ のための十分統計量はいくつも存在します。ではどのような基準で十分統計量を決めればよいでしょうか？　ひとつの答えは次の最小十分統計量というアイデアです。

最小十分統計量

統計量 T を パラメータ $\theta$ のための十分統計量とする．T が最小十分統計量であるとは，パラメータ $\theta$ のための任意の十分統計量 U に対し，１対１写像 h が存在して $T = h(U)$　のように変換できることをいう。

このような都合のよい T が存在するのであれば明らかに　いろいろな十分統計量を探さずとも最初から T を使うことが効率的です。

【最小十分統計量を見つける方法】

$f_X (x;θ)$ を標本 X が従う確率密度関数（離散的な場合には確率質量関数）とします。統計量 $t(X)$ は以下の性質をもつとき最小十分統計量となります。

任意の二つの観測データ x と y について，$T (x) = T (y)$ （標本空間の T による分割において x と y が同値クラス）のときかつそのときに限り比率 $f_X(x;θ) / f _X (y;θ)$ が $\theta$ の関数として一定となる。すなわち、

$$\begin{equation} \begin{aligned} & T (x) = T (y) \\ \Longleftrightarrow & f_X(x;θ) / f _X (y;θ) \;\; \text{がθに依存しない} \\ \Longleftrightarrow & L(θ;x) / L (θ;y) \;\; \text{がθに依存しない} \end{aligned} \end{equation}$$

証明は専門の確率論の教科書を参照ください。この事実を用いて最小十分統計量を見つけることができます。

尤度比を使用して最小十分統計量を見つける手順：

1. 観測データ $x$ に基づいて尤度関数 $L(\theta; x)$ を構築します。
2. 2つの異なる値 $x_1$ と $x_2$ に対して、尤度比を計算します：

$$\begin{equation} \frac{L(\theta; x_1)}{L(\theta; x_2)} \end{equation}$$

1. 尤度比が $\theta$ の関数として一定となる場合を特定します。これは、パラメータ $\theta$ に関する同じ情報を提供する観測データ $x$ の同値クラスを見つけることを意味します。
2. この同値クラスを定義する統計量が最小十分統計量です。

例：母集団分布が正規分布の場合

仮に $X_1, X_2, \ldots, X_n$ が平均 $\mu$ と分散 $\sigma^2$ の正規分布 $N(\mu, \sigma^2)$ からの独立同分布（i.i.d.）のサンプルであるとします。尤度関数は次のようになります：

$$\begin{equation} L(\mu, \sigma^2; x) = \left(\frac{1}{\sqrt{2 \pi \sigma^2}}\right)^n \exp\left(-\frac{1}{2\sigma^2} \sum_{i=1}^n (X_i - \mu)^2\right) \end{equation}$$

2つのサンプル $x_1$ と $x_2$ に対して：

$$\begin{equation} \begin{aligned} \frac{L(\mu, \sigma^2; x_1)}{L(\mu, \sigma^2; x_2)} &= \exp\left(-\frac{1}{2\sigma^2} \left[\sum_{i=1}^n (x_{1i} - \mu)^2 - \sum_{i=1}^n (x_{2i} - \mu)^2\right]\right) \\ &= \exp\left(-\frac{1}{2\sigma^2} \left[ \sum_{i=1}^n x_{1i}^2 - \sum_{i=1}^n x_{2i}^2 + 2 \mu( \sum_{i=1}^n x_{1i} - \sum_{i=1}^n x_{2i} )\right]\right) \end{aligned} \end{equation}$$

尤度比がパラメータ $\mu$ の関数とみたとき一定となるのは次の場合に限られます：

$$\begin{equation} \sum_{i=1}^n x_{1i} = \sum_{i=1}^n x_{2i} \end{equation}$$

従って、$\mu$ に対する最小十分統計量は１次モーメント $\sum_{i=1}^n X_{i}$ です。

また、尤度比がパラメータ $\sigma^2$ の関数とみたとき一定となるのは次の場合に限られます：

$$\begin{equation} \sum_{i=1}^n x_{1i} = \sum_{i=1}^n x_{2i} \;\;\text{and}\;\; \sum_{i=1}^n x_{1i}^2 = \sum_{i=1}^n x_{2i}^2 \end{equation}$$

従って、$\sigma^2$ に対する最小十分統計量は １次モーメントと２次モーメントのペア $(\sum_{i=1}^n X_{i},\sum_{i=1}^n X_{i}^2)$ です。

十分統計量 $t(X)$ が与えられた場合、その条件付確率分布 $P(X|T)$ はパラメータ $\theta$ に関する情報をすべて含むため、$t(X)$ 自体から尤度関数をつくることでを用いてパラメータ $\theta$ を推定できます。これは $t(X)$ に関する条件付確率分布の特性を利用して推定する方法です。

$$\begin{equation} \theta_{MLE}= \underset{\theta}{\operatorname{argmax}} \; L_T(\theta;x) = \underset{\theta}{\operatorname{argmax}} \; h(x)\,g(t(X),\theta ) \end{equation}$$

最小十分統計量を使い最尤推定する

母集団分布が前出のベルヌーイ分布の場合を例とします。 ベルヌーイ分布の尤度関数は次の通りです：

$$\begin{equation} L(\theta; X) = \prod_{i=1}^n \theta^{X_i} (1 - \theta)^{1 - X_i} \end{equation}$$

因子分解定理により、尤度関数は次の形に因子分解できます：

$$\begin{equation} L(\theta; X) = \theta^{\sum_{i=1}^n X_i} (1 - \theta)^{n - \sum_{i=1}^n X_i} \end{equation}$$

尤度比を作ります。

$$\begin{equation} \begin{aligned} \frac{L(\mu, \sigma^2; x_1)}{L(\mu, \sigma^2; x_2)} &= \frac{\theta^{\sum_{i=1}^n x_{1i}} (1 - \theta)^{n - \sum_{i=1}^n x_{1i}}}{\theta^{\sum_{i=1}^n x_{2i}} (1 - \theta)^{n - \sum_{i=1}^n x_{2i}}} \\ &= (\frac{\theta}{1-\theta})^{\sum_{i=1}^n x_{1i}-\sum_{i=1}^n x_{2i}} \end{aligned} \end{equation}$$

尤度比がパラメータ $\theta$ の関数とみたとき一定となるのは次の場合に限られます：

$$\begin{equation} \sum_{i=1}^n x_{1i} = \sum_{i=1}^n x_{2i} \end{equation}$$

従って、$\theta$ に対する最小十分統計量は１次モーメント $\sum_{i=1}^n X_{i}$ です。

最小十分統計量 $t(X) = \sum_{i=1}^n X_i$ を用いると、尤度関数は次のように簡略化できます：

$$\begin{equation} L(\theta; t(X)) = \theta^{t(X)} (1 - \theta)^{n - t(X)} \end{equation}$$

尤度関数を最大化するために、対数尤度関数を考えます：

$$\begin{equation} \log L(\theta; t(X)) = t(X) \log \theta + (n - t(X)) \log (1 - \theta) \end{equation}$$

対数尤度関数の $\theta$ に関する一階導関数を求めて 0 と置き、解を求めます：

$$\begin{equation} \frac{d}{d\theta} \log L(\theta; t(X)) = \frac{t(X)}{\theta} - \frac{n - t(X)}{1 - \theta} = 0 \end{equation}$$

これを解くと、最尤推定量 $\hat{\theta}$ は次のようになります：

$$\begin{equation} \hat{\theta} = \frac{t(X)}{n} = \frac{\sum_{i=1}^n X_i}{n} \end{equation}$$

小標本サイズでの留意点

最尤推定法（MLE）が大きな標本サイズにおいて「漸近的に良い性質を持つ」ことがわかりました。これは逆に言えば、標本サイズが小さいケースに対して最尤推定を適用する際には、バイアスや分散の影響等を考慮しておく必要があることを意味しています。

これらの影響に対しては、パラメータに制約を課すことで過度な推定バイアスを防ぐ「正則化」、事前情報を用いてパラメータの不確実性を軽減する「ベイズ推定法」、小標本での分散や信頼区間評価のための「ブートストラップ法」、等の対策を活用することが考えられます。

ベイズ推定法

最尤推定法をはじめこれまでに述べた点推定法は「頻度主義」という前提にたっています。頻度主義の立場では母集団の特徴を定めるパラメータはある定められた値をもち、観測データは確率的に変動すると考えます。推定値は観測データの統計量から導出される（最小二乗法推定）、あるいは観測データがあらわれる確率を最大化するパラメータとして算出されます（最尤推定法）。
これに対してベイズ推定法（Bayesian Estimation）は逆に母集団のパラメータは確率的にしか値を知ることができない不可知性（Uncertainty）、不確実性をもったものと考えます。そして今パラメータに対してもっている不確実性をパラメータの「事前分布（Prior distribution）」として表現し、これを新しい観測データによってアップデートすることでパラメータの不確実性を減じていくというアプローチをとります。このアップデートされた分布をパラメータの「事後分布（Posteriori distribution）」と呼びます。ポイントは事前分布をどのように定めるのかとアップデートの方法です。なお分布と名のつく用語がいくつも出てきましたが、事前分布とは母集団の分布関数とは異なる概念でありパラメータ自体の分布を決定する確率分布であることを再度強調しておきます。

事前分布の決定手法:

事前分布を適切に決定することは、ベイズ推定の結果に大きな影響を与えるため、慎重に選択する必要があります。以下に事前分布を決定する一般的な方法を示します

1. 事前知識の利用: 事前分布はパラメータについての既存の知識や信念に基づいて設定されることがあります。ドメインの専門知識や関連する研究から得られた情報を利用して、事前分布を決定します。
   注）この決め方が、確率＝頻度で慣らされた頭にはとっつきにくい感じがしてします。

2. 非情報事前分布（Non-Informative Prior）:特定の情報を反映しない分布を利用します。例えば、一様分布はわかりやすい例。

3. 経験的ベイズ法（Empirical Bayes）: 過去のデータからパラメータの事前分布を決める方法です。対象となる統計モデルに対する尤度関数を選択する、そして最尤推定法を使用して観測データから事前分布のパラメータを推定します。事前分布は通常は特定の確率分布（例: ベータ分布、指数分布）を仮定し、その分布のパラメータを観測データから推定します。

4. 感度分析: ベイズ推定では事前分布に対する感度分析を行うことが一般的です。事前分布のパラメータを変化させて結果の感受性を評価し、選択した事前分布が結果にどれだけ影響を与えるかを調査します。

事前分布のアップデート:

ベイズの更新式は、ベイズの定理に基づいて尤度関数と事前分布から事後分布を計算する式であり、次のように表されます：

事後分布  =  尤度関数  × 事前分布  / 正規化定数 

これを式で書くと以下のように表現できます:

$$\begin{equation} P(θ∣x)= \frac{P(x∣θ)⋅P(θ)}{P(x)} \end{equation}$$

​ここで、各項の意味は以下の通り：
$P(θ∣x)$: 事後分布。パラメータθ の条件付き確率分布で、データx が与えられた条件下での分布。
$P(x∣θ)$: 尤度関数（Likelihood　function）。パラメータ θ のもとでデータ x が得られる条件付確率を表します。これは実際の観測データに対するモデルの適合度を示します。
$P(θ)$: 事前分布。パラメータθ についての不確実性を表す確率分布です。事前分布はデータを考慮する前の分布です。
$P(x)$: 正規化定数（Normalization Constant）、事後分布を確率密度として正確に表現するための正規化定数です。通常は尤度関数と事前分布の積を全パラメータ空間にわたって積分することによって求められる（周辺化：marginalize）ことから「周辺尤度(marginal likelihood)」とも呼ばれています。

MAP推定法により事後分布の統計量をえる

MAP（Maximum A Posteriori）推定法は、計算量の問題で事後分布自体を求めることが難しい場合、かわりに事後分布内で最も確率密度が高い（最も尤もらしい）パラメータの値を求める手法です。これは事後分布の「ピーク」に対応します。MAP 推定値は事後分布から得られる一つの統計値であり、パラメータの不確実性を考慮する点推定の手法といえます。

MAP 推定値を以下の式で表します：

$$\begin{equation} \theta_{MAP}=\underset{\theta}{\operatorname{argmax}} \; P(\theta∣x) \end{equation}$$

ここで、各記号の意味は次の通り：
$\theta_{MAP}$: MAP推定されたパラメータの値
$\theta$: 推定の対象となるパラメータ
x: 観測データ
$P(θ∣x)$: パラメータ$\theta$の事後分布

右辺をベイズの更新式で置き換えたとき正規化定数の部分は$\theta$に依存していないので事後分布のargmaxは尤度関数×事前分布のargmaxに等しい。

$$\begin{align} \underset{\theta}{\operatorname{argmax}} \; P(θ∣x) &= \underset{\theta}{\operatorname{argmax}} \; \frac{P(x∣θ)⋅P(θ)}{P(x)} \\ &= \underset{\theta}{\operatorname{argmax}} \; P(x∣θ)⋅P(θ) \end{align}$$

例えば事前分布が$\theta$によらないものを考えれた場合には（例えば一様分布）、これは最尤推定そのものです。この式の形から一般的には尤度関数に対する事前分布をうまく決めることがポイントになることがわかります。

うまい決め方のひとつとして「共役事前分布（Conjugate　Prior）」をあげておきます。 これは尤度関数と事後分布が同じ関数の形になるように事前分布を決める方法です。

例えば母集団の観測データ $X = \{x_1, x_2, \ldots, x_n\}$ の分布関数がパラメータ $\theta$を平均にもつ正規分布$N(\theta,\sigma^2)$ であり、$\theta$ の MAP 推定値を計算したいとします。このとき $\theta$ の事前分布として平均 $\mu _ 0$ と分散 $\sigma_0^2$ の正規分布を設定してみましょう。すると、ベイズ更新式から

$$\begin{equation} P \left( \theta | X \right) \propto \exp \left( -\frac{1}{2 \sigma ^ 2} \sum _ {i=1} ^ n (x _ i - \theta) ^ 2 \right) \times \exp\left( -\frac{1}{2\sigma _ 0 ^ 2} (\theta - \mu _ 0) ^ 2 \right) \end{equation}$$

両辺の対数をとることで、

$$\begin{equation} \log P \left( \theta | X \right) = - \frac {1} { 2 \sigma ^ 2 } \sum _ {i=1} ^ n ( x _ i - \theta ) ^ 2 - \frac {1}{ 2 \sigma _ 0 ^ 2 } ( \theta - \mu _ 0 ) ^ 2 + C \end{equation}$$

ここで $C$ は $\theta$ に依存しない定数です。 $\log P(\theta | X)$ を最大化する$\theta$はこれに $-1$ を乗じたものの最小値をもとめることと同じなので以下のように書き換えます。

$$\begin{equation} -\log P(\theta | X) = \frac{1}{2\sigma^2} \sum_{i=1}^n (x_i - \theta)^2 + \frac{1}{2\sigma_0^2} (\theta - \mu_0)^2 + C' \end{equation}$$

ここで $C'$ は $\theta$ に依存しない定数です。 右辺を $\theta$ で微分した式がゼロとなる値を求めます。

$$\begin{equation} \frac{d}{d\theta} \left( -\log p(\theta | X) \right) = \frac{1}{\sigma^2} \sum_{i=1}^n (x_i - \theta) + \frac{1}{\sigma_0^2} (\theta - \mu_0) \end{equation}$$

式を変形していけば以下のようになります。

$$\begin{equation} \frac{1}{\sigma^2} \sum_{i=1}^n (x_i - \theta) + \frac{1}{\sigma_0^2} (\theta - \mu_0) = 0 \end{equation}$$

$$\begin{equation} \frac{1}{\sigma^2} \sum_{i=1}^n x_i - \frac{n}{\sigma^2} \theta + \frac{1}{\sigma_0^2} \theta - \frac{\mu_0}{\sigma_0^2} = 0 \end{equation}$$

$$\begin{equation} \left( \frac{n}{\sigma^2} + \frac{1}{\sigma_0^2} \right) \theta = \frac{1}{\sigma^2} \sum_{i=1}^n x_i + \frac{\mu_0}{\sigma_0^2} \end{equation} 従って、 \begin{equation} \theta _ {MAP}= \frac{\frac{1}{\sigma^2} \sum_{i=1}^n x_i + \frac{\mu_0}{\sigma_0^2}}{\frac{n}{\sigma^2} + \frac{1}{\sigma_0^2}} \end{equation}$$

正規分布のように解析的な計算により MAP 推定値を得ることができる母数の分布としては他に二項分布（事前分布はベータ分布）、ポアソン分布（事前分布はガンマ分布）等の指数型分布族に属する分布が知られています。

このほかにマルコフ連鎖モンテカルロ法（MCMC）と呼ばれる機械学習や予測で使われる非常に興味深い技法がありますがここでは踏み込みません。ここでは医療統計の記載が充実した書籍を参考としてあげておきます（参考２）。

2.4. 区間推定

不確実性の認識と信頼度：
母集団のパラメータは確定した真の値をもつという前提は点推定と共通しつつ、区間推定（Interval Estimation）では、我々が行う統計的推定行為には不確かさがあるという前提にたちます。その不確かさを測る指標として「信頼度（信頼水準）」を事前に設定します。一般的には95％や99％などが使われるようです。推定対象パラメータを$Param$とし、信頼度を95%に設定したとしましょう。このとき$\hspace{0.3em} |Param|<v \hspace{0.3em}$ で決まる幅$2v$の区間に95%の確率で真のパラメータを発見できる、と判定した場合この区間のことを「信頼区間」と呼びます。信頼区間を決定することを一般に区間推定といいます。

信頼区間の計算： 以下に主要な区間推定の方法をあげておきます。

Zスコアによる推定：
Z スコアは、正規分布に従う観測データ $X$ が母集団の平均からどれだけ標準偏差の数値分離れているかを理解しやすいよう変換したものです。

$$\begin{equation} Z = \frac{X-\mu}{\sigma} \end{equation}$$

ここで X は観測データ、μ は母集団の平均、σ は母集団の標準偏差です。
標本平均 $\overline{X}$ もまた正規分布に従い、そのZスコアは以下の式となります。

$$\begin{equation} Z = \frac{\overline{X}-\mu}{\sigma/\sqrt{n}} \end{equation}$$

この式の分母に現れる量を標準誤差（Standard Error, SE）とよびます。

$$\begin{equation} SE = \frac{\sigma}{\sqrt{n}} \end{equation}$$

観測データ $X$ が正規分布に従わない場合であっても中心極限定理が示すように標本平均の Z スコアは n が十分大きければ近似的に正規分布に従います。そこで信頼度を α 、これを実現する区間幅を $z_{conf}$ とすると、信頼区間の定義から

$$\begin{equation} P(-z_{conf} \ \le\ \frac{\overline{X}-\mu}{SE} \ \le\ +z_{conf}) \ =\ \alpha \end{equation}$$

Pの括弧内の論理式は以下と同値

$$\begin{equation} P(\overline{X}-z_{conf}\cdot SE \ \le\ \mu \ \le\ \overline{X}+z_{conf}\cdot SE) \ =\ \alpha \end{equation}$$

つまり信頼区間の推定値は標本平均からずれの下限、上限として以下のように求められます。

下限 = $\overline{X}-z_{conf}\cdot SE$
上限 = $\overline{X}+z_{conf}\cdot SE$

具体的な数値は公開されているZ分布表（例：産総研のZ分布表）を利用することで容易に信頼区間の推定値を得ることができます。

上図で左側の赤丸はZスコアの小数点1位までの数字が「1.9」であることを示します。上方の赤丸は小数点２位の数字が「6」であることを示します。1.9 の横軸と6の縦軸の交点の数字「0.024998」が信頼度を示します。すなわちZ<=1.96となる確率が約 0.25% であることを示しています。信頼度の設定を正あるいは負どちらかに限定することを片側推定といい、正負両方設定した場合を両側推定と呼ぶ習慣があるので信頼度5%での両側推定の信頼区間は |Z|<=1.96 ということになります。

言葉の使い方
- 信頼度95%の信頼区間は |Z| <= 1.96 である。 - |Z| <= 2 であるとき 「 Z は ２ シグマ以下である」という。信頼度 95.45% の信頼区間を示す。

Tスコアによる推定：

Z スコアによる推定は母分散が既知であることを前提としています。またnが小さくなるほど推定の誤差は増えていきます。Tスコアによる推定はこのようなZスコアが適用できないケースにも対応できるため実務ではよく使われています。

T スコアの背後にある考え方は、母集団の平均に関する統計的な推定を行う際に、標本サイズの大きさが小さい場合に生じる不確かさ（標本誤差）を補正することです。 不明な母分散 $\sigma^2$ の代わりとなる量が必要となります。答えは「標本不偏分散」 $s^2$ と呼ばれる確率変数を定義し置き換えることです。標本不偏分散の平方根を「標本標準偏差」と呼ばれることがあります。少し紛らわしい言葉の使い方のように感ずるのですが。

$$\begin{equation} s := \sqrt{\frac{1}{n-1}\sum_i^n(x_i-\overline{X})^2} \end{equation}$$

日次の株価 $P_t$ からなる時系列から算出した標本標準偏差は株価分析の世界ではボリンジャーバンドとしてよく知られています。また、株価収益率 $R_{t+1} := (P_{t+1} - P_t)/P_t$ に対する標本標準偏差はヒストリカルボラタリティと呼ばれ、株価の値動きの大きさを示す指標としてよく使われています。

$$\begin{equation} H_v := \sqrt{\frac{1}{n-1}\sum_{i=t-n}^t(R_i-\overline{R})^2}*\sqrt{market\_days}*100 \end{equation}$$

上記で、market_days は市場の年間営業日数を表し、例えば250日等の数字を指定します。

標準誤差の式が $\sigma/\sqrt{n}$ であったのでこれを置き換えるのであれば式 $\sum_i^n(x_i-\overline{X})^2/n$（標本分散と呼ばれる確率変数）を使ったほうが自然のようにも思えますが、これはうまくいきません。実際にこの確率変数の期待値を計算してみると得られるのは $(n-1)\sigma^2/n$ という母数 $\sigma^2$ よりも小さな数字であり不偏性が損なわれているからです。そこで式 (97) のように n を n-1 で置き換えた式を使えば、$E(s^2)$ は $\sigma^2$ に一致する不偏推定量となっていることが少し計算するとわかります。

標本標準偏差（$s/\sqrt{n}$）により Z スコア式を修正した結果が以下のＴスコアです。

$$\begin{equation} T = \frac{\overline{X}-\mu}{s/\sqrt{n}} \end{equation}$$

ここで$\overline{X}$は標本の平均、μ は母集団の平均、s は標本標準偏差（＝標本不偏分散の平方根）、n は標本サイズ。

Tスコアが従うT分布の信頼区間の推定値は標本平均からずれの下限、上限として以下のように求められます。信頼度をα、これを実現する区間幅を$t_{conf}$として、

下限 = $\overline{X}-t_{conf}\cdot s/\sqrt{n}$
上限 = $\overline{X}+t_{conf}\cdot s/\sqrt{n}$

上図は比較のため T 分布の自由度１、２のグラフと標準正規分布のグラフを excel を使って書いてみました。T 分布は T.DIST 関数、正規分布は NORM.DIST 関数で計算したものを散布図で容易にグラフ化にできます。自由度は標本サイズ（＝標本あたりの観測データ数）から１を引いた数として定義されるため、自由度１のグラフは標本サイズ = 2 のときに適用されます。
信頼度 αを決める（例えば95%）と、T 分布表（参考：産総研の T 分布表）から信頼区間を確認することができます。T 分布の密度関数の式形は少し複雑なので割愛します。考案者（ゴセット）がどのような過程で T 分布を発見したのか筆者は調べていませんが本稿の目的にはT分布の公式を知らなくてもよさそうなのでこれ以上立ち入らないことにします。

参考 - Z分布表　産総研有機基準物研究グループ - T分布表　産総研有機基準物研究グループ

3. 統計的に有意であるとは

統計の世界ではある命題の真偽を判定するために計画された実験や調査があったとして、その命題を肯定する結果が偶然に出現したとは考えにくい場合、命題が真であると考え「統計的に有意な」結果が得られたと言明します。ここに出てきた「偶然に出現したとは考えにくい」を判定する基準を数値化し「有意水準（significance_level）」あるいは「危険率」と呼称します。これを指すのにしばしば記号αを用いられます。有意水準は実験や調査の内容に応じて決められるべき数値ですが以下のようなとらえ方がよく行われているようです。
- 5% 以下の出現率 → 小さい - 1% 以下の出現率 → 非常に小さい - 0.1 %以下の出現率 → 極端に小さい

「 100% 正しい」、も「 100% 間違っている」も無い世界で真偽を定義するアイデアとして、便宜的な感じがしつつも便利には違いないと思います。

3.1 帰無仮説

目的とする母集団の特徴を知るということを統計の手法では真偽を判定できる命題の形で表現します。仮説した統計モデルの正しさを検証する等、真偽を証明したい命題を対立仮設と呼びます。これに対して、その命題が偽であると仮定することを「帰無仮説（Null Hypothesis）」を設定するといいます。そして試行から得られたデータをもとに帰無仮説を肯定する結果が出現する確率を計算します。それがあらかじめ定めた有意水準（例えば 5% ）以下であれば、帰無仮説を棄却（否定）する,つまり帰無仮説の反対を主張するもともとの命題（「対立仮説」と呼ばれる）が真であると判定します。

少し言い換えると帰無仮説は対象とする命題が偶然に成立することはさほど珍しいことではないと主張します。しかしその前提で計算してみると極めて稀にしか成立しない事象であることが判明したとします。するともっともらしい結論とは、帰無仮説は誤りであり対象とする命題が成り立つことは偶然ではなく必然だろうとする判断です。

まわりくどい感じもするのですが、「仮説」が正しくないことを証明するためには反例がひとつでも見つかればOKなので、背理法に従ったこの考え方は合理的であるといえます。ただし帰無仮説が棄却されなかった場合何がわかったのかというと、「対立仮説が正しいとはいえない」であって「対立仮説は誤っている」ではない、ということに注意が必要です。

このようにして得られた試行結果をもとに帰無仮説が真であることを証明する作業を「検定(Test)」と呼称します．具体的な手順として帰無仮説を証明するための定量的指標をまず定めます．これを「検定統計量（test statistic）」とよびます．帰無仮説が成立すると仮定したうえで標本に対する検定統計量を算出し、その出現確率が有意水準を下回るのであれば、帰無仮説を棄却する。そうでなければ対立仮説を受容します。

3.2. 超能力は存在するか？

参考３の「推計学のすすめ」という本に記載されていた題材がおもしろっかたのでこれを補強して帰無仮説を説明してみましょう。

唐突だけれどサイコロ投げという行為を通じて「超能力は存在する」という命題をたてその真偽を判定したいと私は考えました。これを証明するための実験としてサイコロを1000回振って念じて偶数を出現させる試行を行おうと思います。

母集団は無限に続くサイコロ投げ、標本は一定回数投げて出た目の集合とする。帰無仮説（H0）は「超能力は存在しない」、対立仮説（H1）は「超能力は存在する」。有意水準を 5% と設定します。 検定統計量は、念じて偶数を出現させた回数（例えば、偶数の出現回数）を用います。帰無仮説の下では超能力が存在しないと仮定されているため、偶数が出現する確率は通常のサイコロの確率（1/2）と同じであると仮定します。

検定のステップ:
1. サイコロを1000回振る試行を実行し、偶数が出現する回数を検定統計量と決め記録する。 2. 試行で得られた検定統計量の値を与える確率を算出し有意水準と比較する。 3. 計算した確率が有意水準以下であれば帰無仮説を棄却し「超能力は存在する」と判定する。そうでない場合ば「超能力が存在するとはいえない」と判定する

試行の結果：
- 偶数の出現数 → 542 - 奇数の出現数 → 458

帰無仮説の判定：
均等に目が出るサイコロを振ったときに指定回数k以上の偶数が出る確率は二項分布を使って求めることができます。

$$\begin{equation} P(X \ge k) = 1- \sum_{i=0}^{k-1} \left( \begin{array}{c} n \\ i \end{array} \right) \cdot p^i \cdot (1-p)^{n-i} \end{equation}$$

$P(X \ge k)$ は、k 回以上成功（偶数が出る回数）する確率。
n は１回の試行に含まれる標本数（この試行では 1000 ）。
k は成功回数（ = 542 ）。
p は1回の試行で成功する確率（偶数が出る確率）で、1/2 。
$\left(\begin{array}{c}n \\ k \\ \end{array} \right)$はn回中k回成功する組み合わせの数 $_nC_k= \frac{n!}{k!(n-k)!}$ であり二項係数と呼びます。

ExcelのBINOM.DIST関数あるいはscipyのbinom.cdf関数など、二項分布の確率分布関数をサポートするツールで計算すると約0.358%を得ます。これは有意水準5%をかなり下回るので偶数の目が出た回数542回は偶然に出現するとは考えられず、結論として帰無仮説は棄却され、「超能力は存在する」と推定できます。 有意水準との比較に使われた$P(X \ge k)$をp値と呼びます。

蛇足ですが有意水準に対応する偶数出現数を計算すると526という数字が出ます。この程度のばらつきは何回か試行するとなんとなく容易に出現しそうな気もするのですが実際にはありえないほどレアということです。

3.3. 手法の一般化

手法の一般化: Z検定

上述した Excel の BINOM.DIST 関数が二項分布をどういうアルゴリズムで出力しているのか筆者は知りませんが、一般に検定統計量が現実的な時間では計算できない場合もあるでしょう。そのような場合であっても中心極限定理が成り立つ条件さえ満たせていれば検定統計量が従う分布として正規分布を使うことができます。二項分布の場合では、分布の期待値 E 、分散 V はそれぞれ np、np(1-p) であること、さらに二項分布は標本サイズ n が大きい場合期待値 np、分散 np(1-p) の正規分布の確率密度関数へ近似することが古くからド・モアブル=ラプラスの極限定理として知られています。

$$\begin{equation} X \quad \sim \quad \frac{1}{\sqrt{2\pi np(1-p)}}\exp(-\frac{(k-np)^2}{2np(1-p)}) \end{equation}$$

したがって二項分布に従う統計量をXとしたとき標本サイズが大きければ

$$\begin{equation} Z \quad = \quad \frac{X-np}{\sqrt{np(1-p)}} \end{equation}$$

は近似的に標準正規分布に従います。この式と標準正規分布表をつかって前述のP(X>=k)を算出してみると0.40%というけっこう近い値が出てきます。離散分布の二項分布を連続分布の正規分布で近似する際に現る差異を補正するために定数0.5を上記Zの定義式分子に加えることを連続補正といいます。連続補正によりZの定義式分子をX+0.5-npで置き換えると近似はもっと良くなり約0.362%という値がでてきます。

平均と分散があらかじめわかっている正規分布に従うZスコアを統計量として使用する検定手法はZ検定（Z-test）と呼ばれます．この事例では母集団（＝サイコロを振って出た目の列）の平均と分散が既知であったためうまい具合にZ検定を適用することができました． $\overline{X}$ を標本平均、$\mu$ を母集団の平均、n を標本サイズ、$\sigma$ を母集団の標準偏差（分散の平方根のこと）としたときのZスコアを再掲します

$$\begin{equation} Z \quad = \quad \frac{\overline{X}-\mu}{\sigma / \sqrt{n}} \end{equation}$$

Z 検定の適用シーン
Z 検定は、特に母集団の平均や比率に関する仮説を検定するために広く使用される手法です。以下に、Z 検定がよく使われる適用シーンをいくつかあげておきます。

• 標本平均の母集団平均との比較（１標本の検定）:
  単一の標本から得られた平均値 $\mu$ が、ある既知の母集団平均 $\mu_0$ と統計的に異なるかどうかを判定するために Z 検定を利用できます。例えば、ある製品の平均寿命が既知の平均寿命と等しいかどうかを検証する場合です。$\mu - \mu_0=0$ を帰無仮説とし、これが棄却される結果を得たならば標本の平均と母集団の平均には有意差があると判定できます。

• 母集団の平均の比較（２標本の検定）:
  Z 検定は、2 つの独立した標本それぞれの平均値 $\mu_0,\mu_1$ が統計的に異なるかどうかを判定するためにも利用できます。例えば、新しい薬の治療効果を評価するために、治療群と対照群の平均体温を比較する場合などです。$\mu_0-\mu_1=0$ を帰無仮説とし、これが棄却される結果を得たならば二つの母集団の平均には有意差があると判定できます。

• 母集団の比率の比較:
  カテゴリカルデータの比率を比較する場合、Z検定は非常に有用です。例えば、2つの異なる広告キャンペーンのクリック率を比較して、どちらが優れているかを判断する場合など。

• 標本サイズが大きいケースへの適用:
  母集団が正規分布にしたがっていないとしても標本サイズが十分に大きい場合、中心極限定理に基づいてZ検定を適用することができます。十分に大きいの目安はn=30とされているようです。

手法の一般化: T 検定

残念ながら現実の応用場面では母集団の平均、分散が既知であることはまれです．そのため母集団の分布は正規分布であるが平均、分散は不明であるという前提での検定手法が一般によく用いられる「T検定」です．T 検定では前出の T スコアを検定統計量として使用する手法です．

$$\begin{equation} T \quad = \quad \frac{\overline{X}-\mu}{S / \sqrt{n}} \end{equation}$$

ここで$S^2$は標本不偏分散、$S / \sqrt{n}$ は一般に標本標準誤差と呼ばれます。

T検定の適用シーンはZ検定と同様です。

１標本検定の場合は母集団の平均$\mu_0$は既知として、$\mu - \mu_0=0$を帰無仮説として設定します。これが棄却された場合標本の平均と母集団の平均には有意差があると判定できます。

２標本検定の場合は、ふたつの母集団の平均の差を検定します。母集団の平均の値自体は問わず、分散は不明という前提となります。Ｚ検定のときと同じく標本平均の差$\mu_0 - \mu_1=0$を帰無仮説とし、これが棄却される結果を得たならば母集団の平均には有意差があると判定できます。

3.4. 棄却域とP値:

棄却域（Rejection Region）:
棄却域とは、帰無仮説が成立すると仮定した場合に、検定統計量の出現確率が有意水準と同じか下回る領域を意味します（下図の黄色でハッチした部分）。棄却域に入らない領域を信頼区間と呼ぶ。棄却域と信頼区間の境界の検定統計量の値を「臨界値（critical value）」とよびます。

P 値（P-value）:
P 値は、帰無仮説が成立すると仮定した場合、標本から計算した検定統計量の出現確率以下になる領域に検定統計量が入る確率を示すものです。従ってP値が有意水準以下となる場合には帰無仮説は棄却されます。一方P値が有意水準以上の場合、帰無仮説は棄却されません。
棄却域として分布の右側あるいは左側の一方だけを使う検定を片側検定といい、両方使う場合を両側検定といいます。両側検定では有意水準を 5% とした場合右側、左側それぞれ 2.5% の有意水準として計算します。

4. 推定の性能を測る指標

仮説検定では命題の真偽が便宜的であるだけに、帰無仮説 $H_0$ を判定した結果偽であるとし棄却した場合であっても実際のところは真であったという場合が当然起こりえます。この判定と真実の関係は以下の混同行列と呼ばれる四象限の表により表すことができます。この行列表現は機械学習におけるクラス分類モデル（分類アルゴリズム）の予測結果と実際のクラス（正解ラベル）との一致や不一致を表現する際にも頻出する重要な概念です。
検定の適用シーンによってはよく使われる用語（陽性、陰性）があるので表に付記しました。感染症の陽性を検出する機器による検査を想定し「($H_0$)検査機器は無効であり、正確な結果を提供しない」を帰無仮説とし、「($H_1$)検査機器は有効である」を対立仮説とする検定を想定すると理解しやすいと思います。

表中4つの主要な要素が含まれます：

帰無仮説を棄却する：
- 真陽性（True Positive, TP）：対立仮説が正しいと検定し、真実も対立仮説が正しかった。
- 偽陽性（False Positive, FP）：対立仮説が正しいと検定したが、真実は対立仮説が誤っていた。

帰無仮説を棄却しない：
- 真陰性（True Negative, TN）：帰無仮説が正しいと検定し、真実も帰無仮説が正しかった。
- 偽陰性（False Negative, FN）：帰無仮説が正しいと検定したが、真実は帰無仮説が誤っていた。

試行により得られた観測データの構成比率は推定作業の性能（成績）を評価する指標となります。以下は主要な指標です。

正解率（Accuracy）：真の陽性、真の陰性が正しく判定された観測データの割合

$$\begin{equation} Accuracy = \frac{TP + TN} {TP + TN + FP + FN} \end{equation}$$

適合率（精度）（Precision）：陽性と判定された観測データのうち、真は陽性である割合

$$\begin{equation} Precision = \frac{TP} {TP + FP} \end{equation}$$

再現率（感度）（Recall）：真は陽性の観測データのなかで、正しく陽性と判定された割合

$$\begin{equation} Recall = \frac{TP} {TP + FN} \end{equation}$$

F1スコア（F1 Score）：適合率と再現率の調和平均で、対立仮説のバランスを示す指標。

$$\begin{equation} F1-score = 2 * (\text{適合率} * \text{再現率}) / (\text{適合率} + \text{再現率}) \end{equation}$$

第一種過誤と第二種過誤は表中トレードオフの位置になっているため同時に減らすことはできません。統計の標準的な考え方では検定の誤りを減らすためにはまず適合率（精度）を向上させます。第一種の過誤となる確率 $\alpha$ （有意水準）を一定以下にすることを優先したうえで、再現率（感度）を向上させることで第二種の過誤となる確率βを減らす検討をすることが一般的のようです。$1-\beta$ は検出力とも呼ばれます。ただし重要な指標は検定の目的によって決めるべきものであるため必ずしもこの指針が絶対というわけではありません。実際、新型インフルエンザ流行が激しかった頃の報道を思い出すと真の陽性を見逃さないために検査機器の感度の高さ（第二種過誤を減らすこと）を求める声が大きかったように思います。

次に株価予測モデルの有効性を検定する作業を考えてみましょう。有効な株価モデルとは「上昇すると予測した株価が実際に上昇する」、上昇的中率が高いモデルであると考えることにします。

「($H_0$)予測モデルは無効であり、正確な結果を提供しない」を帰無仮説とし、「($H_1$)予測モデルは有効である」を対立仮説とする検定を想定します。

上昇予測数に対する上昇的中数の割合の分布を測定し検証した結果帰無仮説が棄却されたなら予測モデルは有効だ！と判定する、という流れです。下表でいえば観測データのなかで $(TP+FP)$ に占めるTPの割合が非常に高い場合が予測モデルが有効である考える根拠になると言えます。筆者はロスカット発生を恐れる性分なのでこの予測モデルで重視したいのは、実際には値下がりする株を値上がりすると予測する誤りが少ないことです。そうすると第一種の過誤を減らす指標として偽陽性率 (False Positive Rate, FPR) $FP/(FP+TN)$ の低減を指標として予測モデルを改良していくことが目的に合っていると結論できます。

補足説明

補足説明 確率・統計
補足説明 ソフトウェアを準備する