$ \newcommand\diff[2]{\frac{d{#1}}{d{#2}}} $
 

順序統計量が従う分布

順序統計量が従う分布.

$[0, 1]$ の一様分布から10個サンプルをとってくる試行を, 100セット繰り返す. 各セットで3番目に大きい値のヒストグラムを描くと, どんな形になるだろうか? そのような問題について, 考えてみる.


非専門家が書いています. 十分に批判的に読んで頂くようお願いいたします. 間違い・疑問点などあれば, ぜひコンタクトフォームへ連絡いただけると幸いです.

順序統計量

確率変数 $X_1, X_2, \ldots, X_n \overset{iid}{\sim} F$ とする. つまり, 確率変数 $X_1, X_2, \ldots, X_n$ は独立に同一の分布に従うとする. これを, 小さい順から並べ替えたものを $X_{(1)}, X_{(2)}, \ldots, X_{(n)}$ としたとき, これらを 順序統計量という. 換言すると, $i$ 番目の順序統計量 $X_{(i)}$ は, 小さい方から $i$ 番目の値である.

たとえば $F$ を $[0, 1]$ 上の一様分布 $\mathrm{Uni}[0, 1]$ として, $X_1, X_2, \ldots, X_5 \overset{iid}{\sim} \mathrm{Uni}[0, 1]$ とする. 実現値を得るために, MATLABで乱数を発生させてみる.

MATLAB
>> d = rand(1, 5)

d =

    0.0975    0.2785    0.5469    0.9575    0.9649

これが, $X_1, X_2, \ldots, X_5$ の実現値である. ソートすると,

MATLAB
>> sort(d)

ans =

    0.0975    0.2785    0.5469    0.9575    0.9649

これが, $X_{(1)}, X_{(2)}, \ldots, X_{(5)}$ の実現値である.

順序統計量が従う分布に関する一般理論

方針

では, 順序統計量 $X_{(i)}$ が従う密度函数を求めよう. いきなり密度函数を求めるのは難しいので, まず累積分布函数を求めよう. つまり,

\begin{align} \large F_{(i)}(x) \equiv \Pr{(X_{(i)} \leq x)} \end{align}

を求める.

累積分布函数を求める

$X_{(i)} \lt x$ となるのは,

  • $n$ 個の実現値うち $i$ 個が $x$ 以下,
  • $n$ 個の実現値うち $i+1$ 個が $x$ 以下,
  • $n$ 個の実現値うち $i+2$ 個が $x$ 以下,
    $\qquad \qquad \vdots$
  • $n$ 個の実現値すべてが $x$ 以下,
という事象の和である ( わかりにくければ, $i=3$ などと置いてみればよい ) .

すなわち,

$B_k = \{ X_1, X_2, \ldots, X_n$ のうち $k$ 個が $x$ 以下 $\}$
と定義すると,

\begin{align} \large F_{(i)}(x) = \sum_{k=i}^{n} \Pr{(B_k)} \end{align}

となる.

では, $\Pr{(B_k)}$ をもとめてみよう. これは, 高校レベルの問題である. $X_1, X_2, \ldots, X_n \overset{iid}{\sim} F$ であることを踏まえると, 各 $X_i$ が $x$ 以下である確率は $F(x)$, $x$ より大きい確率は $1 - F(x)$ である. よって, $\Pr{(B_k)}$ は,

\begin{align}\large{ \Pr{(B_k)} = \binom{n}{k} p^k (1 - p)^{n-k}, \quad p = F(x) }\end{align}

である. よって,

\begin{align}\large{ F_{(i)}(x) = \sum_{k=i}^{n} \binom{n}{k} p^k (1 - p)^{n-k} }\end{align}

と, 累積分布が求まった.

密度函数を求める

$X_{(i)}$ の密度函数 $f_{(i)}(x)$ は, 累積分布函数 $F_{(i)}(x)$ を $x$ で微分すればよい. しかし, 直接微分するのは難しいので, 合成函数の微分を用いて,

\begin{align} \large f_{(i)}(x) &\large= \diff{F_{(i)}}{x} \label{aaa} \\[8pt] &\large= \diff{F_{(i)}}{p} \diff{p}{x} \notag \\[8pt] &\large= \color{#cf201f}{\diff{F_{(i)}}{p}} f(x) \notag \\[8pt] \end{align}

である. ここで, $F$ の密度函数を $f$ とおいた.

では, 赤字部分を求めよう.

\begin{align} \large \color{#cf201f}{\diff{F_{(i)}}{p}} &\large= \diff{}{p} \sum_{k=i}^{n} \binom{n}{k} p^k (1 - p)^{n-k} \\[8pt] &\large= \sum_{k=i}^{n} \binom{n}{k} \diff{}{p} p^k (1 - p)^{n-k} \notag \\[8pt] &\large= \sum_{k=i}^{n} \frac{n!}{(n-k)!k!} \{ k p^{k-1} (1-p)^{n-k} - p^k (n-k) (1-p)^{n-k-1} \} \notag \\[8pt] &\large= \sum_{k=i}^{n} \left[ \frac{n!}{(n-k)!(k-1)!} p^{k-1} (1-p)^{n-k} - \frac{n!}{(n-k-1)!k!} p^k (1-p)^{n-k-1} \right] \notag \\[8pt] &\large= n \sum_{k=i}^{n} \left[ \frac{(n-1)!}{(n-k)!(k-1)!} p^{k-1} (1-p)^{n-k} - \frac{(n-1)!}{(n-k-1)!k!} p^k (1-p)^{n-k-1} \right] \notag \\[8pt]. \end{align}

ここで, 二項分布に従う確率変数が値 $k$ を取る確率を $\mathrm{Bin}(k; n, p)$ と表すと,

\begin{align}\large{ \mathrm{Bin}(k; n, p) = \binom{n}{k}p^k (1-p)^{n-k} }\end{align}

であるので, この記号を用いると,

\begin{align} \large \color{#cf201f}{\diff{F_{(i)}}{p}} &\large= n \sum_{k=i}^{n} \left[ \mathrm{Bin}(k-1; n-1, p) - \mathrm{Bin}(k; n-1, p) \right] \\[8pt] &\large= n \left[ \mathrm{Bin}(i-1; n-1, p) - \mathrm{Bin}(n; n-1, p) \right] \notag \\[8pt] &\large= n \mathrm{Bin}(i-1; n-1, p) \notag \\[8pt] &\large= n \binom{n-1}{i-1}p^{i-1} (1-p)^{n-i} \notag \\[8pt] &\large= n \frac{(n-1)!}{(n-i)!(i-1)!} p^{i-1} (1-p)^{n-i} \notag \\[8pt] &\large= \frac{n!}{(n-i)!(i-1)!} p^{i-1} (1-p)^{n-i} \notag \\[8pt] \end{align}

これを式\eqref{aaa}に $p=F(x)$ とともにもどして,

\begin{align}\large{ f_{(i)}(x) = f(x) \frac{n!}{(n-i)!(i-1)!} F(x)^{i-1} (1-F(x))^{n-i}. \label{eqmain} }\end{align}

これで, $i$ 番目の順序統計量 $X_{(i)}$ の確率密度函数が求まった.

例と数値計算

$F$ が一様分布の場合

$F$ が $[0, 1]$ 上の一様分布 $\mathrm{Uni}[0, 1](x)$ の場合を考えてみよう. この場合,

\begin{align}\large{ F(x) = x ,\quad f(x) = 1 }\end{align}

なので, これを式\eqref{eqmain}に代入すると,

\begin{align}\large{ f_{(i)}(x) = \frac{n!}{(n-i)!(i-1)!} x^{i-1} (1-x)^{n-i} }\end{align}

となる. ここで, ガンマ函数の性質

\begin{align}\large{ \Gamma(n+1) = n! }\end{align}

を思い出すと, 上の式は,

\begin{align}\large{ f_{(i)}(x) = \frac{\Gamma(n+1)}{\Gamma(n-i+1)\Gamma(i)} x^{i-1} (1-x)^{n-i} }\end{align}

さらに, ベータ函数

\begin{align}\large{ \mathrm{B}(\alpha, \beta) = \frac{\Gamma(\alpha) \Gamma(\beta)}{\Gamma(\alpha + \beta)} }\end{align}

を用いると,

\begin{align}\large{ f_{(i)}(x) = \frac{1}{\mathrm{B}(\color{#cf201f}{i}, \color{#cf201f}{n-i+1})} x^{\color{#cf201f}{i}-1} (1-x)^{(\color{#cf201f}{n-i+1})-1} }\end{align}

となり, これはベータ分布 $\mathrm{Beta}(\alpha, \beta)$, $\alpha = i$, $\beta = n-i+1$ にほかならない.

これを数値的に検証してみよう.

$[0, 1]$ 上の一様分布に従う乱数を10個サンプリングするという試行を10000回繰り返す. そして, 各試行の中で下から3番目の大きさのデータ ( 10000個 ) を抽出し, ヒストグラムを描く.

MATLAB
rng(1234)
N    = 10;
ITER = 10000;
I    = 3;

alpha = I;
beta  = N - I + 1;

d = rand(N, ITER);
d = sort(d);
ithData = d(I, :);
histogram(ithData)

結果.

name

さて, これはベータ分布に従うだろうか? $\alpha = 3$, $\beta = 10 - 3 + 1$ のパラメータを持つベータ分布をプロットしてみよう.

MATLAB
X = 0:0.01:1;
y = betapdf(X, alpha, beta);
plot(X, y)

結果. したがってそうである.

name

$F$ が正規分布の場合

正規分布の場合はどうだろうか. 正規分布の確率密度函数を $\phi(x)$, 累積分布函数を $\Phi(x)$ として, これらを式\eqref{eqmain}に代入すると,

\begin{align} \large f_{(i)}(x) &\large= \phi(x) \frac{n!}{(n-i)!(i-1)!} \Phi(x)^{i-1} (1-\Phi(x))^{n-i} \\[8pt] &\large= n \phi(x) \binom{n-1}{i-1} \Phi(x)^{i-1} (1-\Phi(x))^{n-i} \notag \\[8pt] \end{align}

これがどんな形状か皆目検討もつかないが, 数値実験はできる.

標準正規分布に従う乱数を10個サンプリングするという試行を10000回繰り返す. そして, 各試行の中で下から2番目の大きさのデータ ( 10000個 ) を抽出し, ヒストグラムを描く.

MATLAB
t = tiledlayout(1, 2); t.Padding = 'compact'; t.TileSpacing = 'compact';
nexttile
rng(1234)
N    = 10;
ITER = 10000;
I    = 2;

% サンプルする
d = randn(N, ITER);
d = sort(d);
ithData = d(I, :);
histogram(ithData, 'EdgeColor', 'none')
xlim([-5, 5])

% 順序統計量の確率密度分布のプロット
nexttile
X = -5:0.01:5;
y = N * normpdf(X) .* nchoosek(N-1, I-1) .* normcdf(X).^(I-1) .* (1 - normcdf(X)).^(N-I); 
plot(X, y)

結果. したがってそうである.

name

感想・参考文献

感想

いまひとつ, 順序統計量には興味が持てなかったので, いっそのこと記事にすることで勉強してしまえという思いで書いた. 数値実験してみると面白い.

参考文献

竹村 新装改訂版 現代数理統計学

この記事のTOP    BACK    TOP