Python3ではじめるシステムトレード：非復元単純無作為抽出と有限修正

大きさ$N$の有限母集団から大きさ$n$の標本を非復元単純無作為抽出で得た場合の標本平均の分散は母分散に有限修正を施すことで得られます。この有限修正項を導きたいと思ったのですが、期待値の意味がよくわかっていないので、導出に苦労しました。有限母集団の非復元抽出の際には、期待値は確率を用いて算出するよりも、順列と組み合わせ、場合の数を駆使して計算した方がいいかもしれません。基本的なところを調べメモしましたので、残しておきたいと思います。間違い等がありましたら、ご指摘いただけると幸いです。
コードはGitthubに掲載しました。
https://github.com/innovation1005/qiita_innovation1005

これは「データの分析方法を学ぶ会」の勉強会で悩んだ問題の１つです。日曜の朝9:30と火曜の午後7:00から勉強会を開いています。
https://study-data-analysis.connpass.com/
参考文献は統計学実践ワークブックです。今回の投稿は21章に関するものです。

非復元抽出 (Without Replacement)

非復元抽出の際の有限修正項を求めます。その際に期待値を用いますが、期待値は

$E(X) = \sum_{x} x \cdot P(X = x)$
$E(X) = \int x f(x)dx$

と統計学では表現します。しかし、非復元抽出のような際には、母集団の大きさは無限ではありません。そこで非復元抽出の場合、特定の確率変数に対する期待値を計算する際、各抽出が他の抽出に影響を与えることになります。したがって、確率が変化します。

そのために具体的な問題に応じた期待値の計算が可能です。その点に注意をしながら説明をしていきます。

非復元抽出の特徴:

• 一度選ばれた要素は再度選ばれません。
• 抽出が進むにつれて、選択される確率が変動します。(しかし、各々の出現確率は等確率です。)
• 母集団が小さい場合、非復元抽出の影響が顕著になります。

有限修正

有限修正は、母集団が有限で非復元抽出が行われる場合に、統計的推測の精度を高めるために使用される補正です。

• 標本が母集団に対して相対的に大きい場合、有限修正を使用することで、推測のバイアスを減らし精度を向上させます。
• 標本が母集団のごく一部である場合（例：標本が母集団の5%未満）、有限修正の影響はほとんどありません。

非復元抽出と有限修正は、実際のデータ収集や分析の際に考慮すべき重要な要因であり、母集団の大きさと標本の大きさがどのように相互作用するかを理解するのに役立ちます。特に、母集団が限られている場合や非復元抽出が行われる場合には、これらの概念が統計的推測の精度に重要な影響を及ぼす可能性があります。

非復元単純無作為抽出の標本平均と不偏推定量

大きさNの母集団から非復元単純無作為抽出されたn個のデータ点の標本平均が母平均の不偏推定量であることを示します。

例として、N=3,n=2を考えます。
$X_1,X_2,X_3$の結合分布は2個からなる組み合わせの数、$i_1,i_2$に対して

$$P_x(X_1=x_{i_1},X_2=x_{i_2})=\frac{1}{3!/(2!(3-2)!)}=\frac{1}{3} --------(1)$$

となります。

また、標本平均は$\bar{X}=\frac{1}{2}(x_{i_1}+x_{i_2})$です。
母集団を${1,2,3}$とすると可能な標本の組み合わせは

i1	i2
1	2
1	3
2	3

となります。
それぞれが起きる確率は$P_x(X_1=x_{i_1},X_2=x_{i_2})=\frac{1}{N!/(n!(N-n)!}=1/3$なので

$$E_x(\bar{X})=\frac{1}{3}\left[\frac{1+2}{2}+\frac{1+3}{2}+\frac{2+3}{2}\right]=\frac{12}{3\times 2}=2 --------(2)$$

これは母平均$\frac{1+2+3}{3}=2$に等しくなります。

一般化すると

$$E_x(\bar{X})=\frac{1}{\binom{N}{n}}\sum_{k=1}^{\binom{N}{n}} \frac{1}{n}(x_{i_1}+\cdots+x_{i_n}) -----------(3)$$
母集団から非復元でn個のデータ点を抽出する場合、特定のデータ点
$X_i$が選ばれる組み合わせの数は、残りの$n−1$個を残りの$N−1$個から選ぶ組み合わせの数である$\binom{N−1}{n-1}$と等しくなります。
$$=\frac{n!(N-n)!}{N!n}\frac{(N-1)!}{(N-n)!(n-1)!}\sum_{i=1}^NX_i=\frac{1}{N}\sum_{i=1}^NX_i=\mu ----(4)$$

非復元抽出における平均や分散の計算では、一般に順序は重要ではないため、組み合わせを使用することになります。しかし、特定の文脈で順序が意味を持ちます。3人の候補者から3人のメダリストを選ぶ場合、選ぶ順序が結果に影響するため、順列を用います。10人のプレイヤーから5人のチームメンバーを選ぶ場合、メンバー間の順序は通常は関係ないため、組み合わせを用います。

標本平均の標本分散と非復元単純無作為抽出

大きさNの母集団から非復元単純無作為抽出されたn個のデータ点の標本平均の分散を求めます。

この母分散は

$$\sigma^2=\frac{1}{N}\sum_{i=1}^N(X_i-\mu)^2 ----------(5)$$

で与えられます。

標本平均と総平均の差の二乗は

$$(\bar{X}-\mu)^2=\frac{1}{n^2}[(X_1-\mu)+\cdots+(X_n-\mu)]^2$$
$$=\frac{1}{n^2}\left[\sum_{i=1}^n(X_i-\mu)^2+\sum_{i=1}^n\sum_{j=1,i \neq j}^n(X_i-\mu)(X_j-\mu)\right] -----(6) $$

で与えられます。これは

$$=\frac{1}{n^2}\left[\sum_{i=1}^n(X_i-\mu)^2+\sum_{i}^n\sum_{j,i \neq j}^nX_iX_j-\mu\sum_{i}^n\sum_{j,i \neq j}^nX_i-\mu\sum_{i}^n\sum_{j,i \neq j}^nX_j+\sum_{i=1}^n\sum_{i \neq j}\mu^2\right]$$
$$=\frac{1}{n^2}\left[\sum_{i=1}^n(X_i-\mu)^2+\sum_{i}^n\sum_{j,i \neq j}^nX_iX_j-\mu\sum_{j,i \neq j}^nn\mu-\mu\sum_{j,i \neq j}^nn\mu+\sum_{i=1}^n\sum_{i \neq j}\mu^2\right]$$
$$=\frac{1}{n^2}\left[\sum_{i=1}^n(X_i-\mu)^2+\sum_{i}^n\sum_{j,i \neq j}^nX_iX_j-\sum_{j,i \neq j}^nn\mu^2\right]-----(7)$$

となります。つぎに標本平均と総平均の差の平方の期待値をとります。この際に要素の出現確率は等確率を仮定し、確率変数の独立であることは仮定していません。また、各試行の独立は仮定していません。ただし、母集団からデータ点は独立して選ばれます。期待値はすべての場合の数を取ることで得ています。

$$E(\bar{X}-\mu)^2=\frac{1}{n^2}E[(X_1-\mu)+\cdots+(X_n-\mu)]^2$$
$$=\frac{1}{n^2}\left[\sum_{i=1}^nE(X_i-\mu)^2\right]+\frac{1}{n^2}\left[\sum_{i=1}^n \sum_{i \neq j}^nE(X_i-\mu)(X_j-\mu)\right]$$
$$=\frac{1}{n}\sigma^2+\frac{n(n-1)}{n^2}\left[ E(X_i-\mu)(X_j-\mu)\right] ------------(8)$$

$E\sum_{i=1}^n(X_i-\mu)^2$は大きさ$N$の母集団から大きさ$n$の標本を無作為に取った時の期待値ですから、大きさ$N$の母集団について$(X_i-\mu)^2$を計算し、その場合の数で割る必要があります。したがって、

$E\sum_{i=1}^n(X_i-\mu)^2=\sum_{i=1}^nE(X_i-\mu)^2= \sum_{i=1}^nN^{-1}\sum_{j=1}^N(X_j-\mu)^2=n\sigma^2 --(9)$

これはそれぞれの $X_i$ は母集団から独立に選ばれるため、それぞれが同じ期待値を持つからです。

$E(X_i-\mu)(X_j-\mu)$

$$=\frac{1}{N(N-1)}\left[\sum_{i}^N\sum_{j,i \neq j}^NX_iX_j-\mu\sum_{i}^N\sum_{j,i \neq j}^NX_i-\mu\sum_{i}^N\sum_{j,i \neq j}^NX_j+\sum_{i=1}^N\sum_{i \neq j}^N\mu^2\right]$$
$$=E(X_iX_j)-\mu^2 -----(10)$$

この中で、$X_iX_j$の$X_j$を$X_j=N\mu-\sum_{k,k\neq i}X_k$に置き換えて展開すると、

$$\sum_{i=1}^N\sum_{j=1,i \neq j}^NX_iX_j=\sum_{i=1}^N\sum_{j=1,i \neq j}^NX_i(N\mu-\sum_{k,k\neq j}^NX_k)$$
$$=\sum_{j=1,i \neq j}^N\sum_{i=1}^NN\mu X_i-\sum_{i=1}^N\sum_{j=1,i \neq j}^N \sum_{k,k\neq j}^NX_iX_k$$
$$=\sum_{j=1,i \neq j}^N (N\mu)^2-\sum_{i=1}^N\sum_{j=1,i \neq j}^NX_i^2-\sum_{i=1}^N\sum_{j=1,i \neq j}^N \sum_{k,k\neq j,i}^NX_iX_k$$
$$=\sum_{j=1,i \neq j}^N (N\mu)^2-\sum_{i=1}^N\sum_{j=1,i \neq j}^NX_i^2-\sum_{i=1}^N\sum_{j=1,i \neq j}^N \sum_{k,k\neq j,i}^NX_iX_j$$

$$\sum_{i=1}^N\sum_{j=1,i \neq j}^NX_iX_j=\left(1+\sum_{k,k\neq j,i}^N\right)^{-1}\left[\sum_{j=1,i \neq j}^N (N\mu)^2-\sum_{i=1}^N\sum_{j=1,i \neq j}^NX_i^2\right]$$
$$=(N-1)^{-1}\left[\sum_{j=1,i \neq j}^N (N\mu)^2- \sum_{j=1,i \neq j}^N N(\sigma^2+\mu^2) \right]$$
$$=\left[(N\mu)^2- (N\sigma^2+N\mu^2) \right] -----(11)$$

よって、よって、非復元抽出の標本平均の分散は
$$E(X_iX_j)=\frac{1}{N(N-1)}\left((N\mu)^2- N\sigma^2-N\mu^2 \right)$$

$$E(X_iX_j)-\mu^2=\frac{1}{N(N-1)}\left[N^2\mu^2- N\sigma^2-N\mu^2-N(N-1)\mu^2\right]$$

$$E(\bar{X}-\mu)^2=\frac{\sigma^2}{n}-\frac{(n-1)}{(N-1)}\frac{\sigma^2}{n}=\frac{N-1-(n-1)}{(N-1)}\frac{\sigma^2}{n}=\frac{(N-n)}{(N-1)}\frac{\sigma^2}{n} --(12)$$
となります。有限修正が現れました。

例として、N=3,n=2を考えます。母集団を${1,2,3}$とすると

母分散は

$$\sigma^2=\frac{1}{3}[(1-2)^2+(2-2)^2+(3-2)^2]=2/3------(13)$$
となります。

表1

k	i1	i2	平均	平均偏差平方	偏差積	XiXj	muXi	muXj	mumu
1	1	2	1.5	(1.5-2)^2=0.25	(1-2)(2-2)=0	1x2	-2x2	-1x2	2x2
2	1	3	2.0	(2-2)^2=0	(1-2)(3-2)=-1	1x3	-2x3	-1x2	2x2
3	2	1	1.5	(1.5-2)^2=0.25	(2-2)(1-2)=0	2x1	-2x1	-2x2	2x2
4	2	3	2.5	(2.5-2)^2=0.25	(2-2)(3-2)=0	2x3	-2x3	-2x2	2x2
5	3	1	2.0	(2-2)^2=0	(3-2)(1-2)=-1	3x1	-2x1	-3x2	2x2
6	3	2	2.5	(2.5-2)^2=0.25	(3-2)(2-2)=0	3x2	-2x2	-3x2	2x2
						22	-24	-24	24

となります。この際の標本分散の期待値は
$$E(\bar{X}-\mu)^2=\frac{1}{3(3-1)}(0.25+0+0.25+0.25+0+0.25)=1/6 -----(14)$$

$$E(\bar{X}-\mu)^2=\frac{1}{n}\sigma^2+\frac{n(n-1)}{n^2}E\left[(X_i-\mu)(X_j-\mu)\right]=\frac{1}{2}\frac{2}{3}-\frac{2}{4}\frac{2}{6}=\frac{1}{6}-----(15)$$

表２

偏差積	XiXj	3muXi	XiXi	XiXk
(1-2)(2-2)=0	1x2	6x1	-1x1	1x3
(1-2)(3-2)=-1	1x3	6x1	-1x1	1x2
(2-2)(1-2)=0	2x1	6x2	-2x2	2x3
(2-2)(3-2)=0	2x3	6x2	-2x2	2x1
(3-2)(1-2)=-1	3x1	6x3	-3x3	3x2
(3-2)(2-2)=0	3x2	6x3	-3x3	3x1
	22	72	-28	-22

$$E\left[(X_i-\mu)(X_j-\mu)\right]$$
$$=E(X_iX_j)-3\mu E(X_i)- 3\mu E(X_j)+\mu\mu=22/6-24/6-24/6+24/6=-2/6-(16)$$

Pythonを用いた理解

さらにPythonを用いて非復元抽出と有限修正の問題を理解していきましょう。下記のプログラムではN=10,n=5としています。組み合わせを使うので、簡単に計算限界に達してしまいます。

このコードの主な目的は、非復元抽出の下で標本平均の分散を理解することです。特に、有限修正を用いると、この分散の計算がどう変化するかを示すことが重要です。この分析は、小さな母集団からの標本抽出において特に重要です。

以下は各部分の説明です。

母集団の設定:

大きさNの母集団をランダムに生成します。この母集団は0から1の間の連続一様分布から生成される値から成り立っています。

抽出の設定:

非復元抽出でn個のサンプルを取るすべての組み合わせを生成します。

母平均と母分散の計算:

母集団全体の平均（母平均）と分散（母分散）を計算します。

標本平均の分散の計算（有限修正なし）:

各標本の平均を計算し、その平均値の分散を求めます。これにより、標本平均の分散が計算されます。

標本平均の分散の計算（有限修正あり）:

標本平均の分散を有限修正を使用して計算します。有限修正は、標本が母集団全体の一部であるという事実を反映するために使用されます。

結果の出力と可視化:

上記の計算結果を出力し、標本平均のヒストグラムを表示します。

import numpy as np
from itertools import combinations
import matplotlib.pyplot as plt

def withoutreplacement(population,n=2):
    # 母集団の設定

    print('母集団の大きさ',N,'母集団',population)
    # 抽出の設定
    sample_indices = list(combinations(range(N), n))#非復元抽出のすべての組み合わせを生成
    #print('標本の大きさ',n,'標本',sample_indices)

    # 母平均と母分散の計算
    mu = np.mean(population)
    sigma_sq = np.var(population)

    # 標本平均の分散の計算
    sample_means = np.array([np.mean(population[np.array(indices)])\
         for indices in sample_indices])
    sample_var = np.var(sample_means)
    # 標本平均の分散の計算（有限修正あり）
    corrected_sample_var = (1/n)*sigma_sq*((N-n)/(N-1))

    print(f"Population mean: {mu}")
    print(f"Population variance: {sigma_sq}")
    print(f"Variance of sample means: {sample_var}")
    print(f"Variance of sample means \
    (with finite population correction): {corrected_sample_var}")
    plt.hist(sample_means)
    plt.show()
N=3
p=np.array([1,2,3])
withoutreplacement(population=p,n=2)

N=10
p = np.random.rand(N)  # 大きさNの母集団を生成します。
withoutreplacement(population=p,n=5)

# 標本が母集団に対して相対的に大きい場合、有限修正を使用することで、
# 推測のバイアスを減らし精度を向上させます。
withoutreplacement(population=p,n=9)

# 標本が母集団に対して相対的に大きい場合、有限修正を使用することで、
# 推測のバイアスを減らし精度を向上させます。
withoutreplacement(population=p,n=1)