【数学溢れ話】【Token】確率密度空間と累積分布空間①記述統計との狭間

【Python演算処理】行列演算の基本④大源流における記述統計学との密接な関連性?
こちらの投稿の内容をまとめるうちに頭がこんがらがってきたので、別投稿にまとめる事にした次第です。大学時代に覚えた筈の領域の話なのに、いざ再学習に挑戦すると、いやはや…

統計科学事典，清水良一訳(https://www.ism.ac.jp/~ayaka/2017_gairon_1.pdf)

記述統計学とはデータのもっている主要な特性をより鮮明に表現するために，データを要約したり作表をしたりすること一般を指す。

実際の統計学は「とりあえず手元データ自体を母集合と考える」記述統計学(Descriptive Statistics)=古典統計学(Classical Statistics)と「手元データは母集合から抽出した部分標本に過ぎない」と考える推計統計学(Inferential Statistics)の端境期に「大数の法則」概念より発展した正規分布の概念が絶対視され、実質上記述統計の一部と考えられていた「パラメトリック・モデル絶対主義」の時代が挟まり、ベイズ統計の再評価を発端にMCMC技法や統計モデリングの時代に推移するので話がややこしくなってくる訳です。
記述統計学と推計統計学の違い
古典統計学・ベイズ統計・統計モデリングの関係について

#　確率密度空間(Probability Distribution Space)

Wikipediaでは以下の「バクテリア死亡率モデル」が紹介されています。
確率密度関数 - Wikipedia

寿命が4〜6時間程度のバクテリアがいると仮定する。

• 特定のバクテリアが丁度5時間で死亡する確率はどれ位だろうか？ 答えは0%である。およそ5時間で寿命を迎えるバクテリアはたくさん居るが、正確に5.0000000000…時間で死ぬ確率は限りなく小さい。
• 一方、5〜5.01時間で死亡する確率はどうだろうか？ 例えばこれを2%とし、その$\frac{1}{10}$の範囲の5〜5.001時間である確率は？ 答えはおよそ2%$\frac{1}{10}$=0.2%となる。
• さらにその$\frac{1}{10}$の範囲の5〜5.0001時間である確率は、およそ0.02%となる。

上記の3例において「特定の時間範囲内に死亡する確率」を「その範囲の長さ」で割った値に着目すると、1時間につき2に定まる。

• 例えば、5〜5.01時間の0.01時間の範囲でバクテリアが死亡する確率は0.02であり、$\frac{確率0.02}{0.01時間}=2時間^{−1}$である。この$2時間^{−1}$（毎時200%）という量を、5時間時点での確率密度と呼ぶ。
• 従って、「バクテリアの寿命が5時間である確率」を問われた時、真の答えは0%であるが、より実用的には、$2時間^{−1}dt$であると言える。これは、無限小の時間範囲dt内で、バクテリアが死亡する確率である。例えば、丁度5時間〜5時間+1ナノ秒($10^{−9}$)の寿命である確率は、$2時間^{−1}(6x10^{-4})×1ナノ秒≈6×10^{−13}$となる。

これを確率密度関数fを用いて$f（5時間）=2時間^{−1}$と表現するならfを任意の（微小に限らない）時間範囲で積分することで、当該時間範囲内でバクテリアの寿命が尽きる確率を求める事ができる。

とりあえずこれが正規分布に従う場合を考えてみましょう。

f(x)=(x \in \mathbb{R}| \frac{1}{\sqrt{2πσ^2}}exp(\frac{-(x-μ)^2}{２σ^2}))\\
μ=0,σ^2=1の時\\
f(x)=(x \in \mathbb{R}| \frac{1}{\sqrt{2π}}exp(\frac{-x^2}{２}))\\
累積分布関数\\
F(x)=\int_{-\infty}^{x}f(x)dx

この式をこのままSympyに入力すると思わぬ形に整理されてしまう上、定積分した場合の累積分布関数のグラフが出力出来ません。

import numpy as np
import sympy as sp
σ,μ,x = sp.symbols('σ,μ,x')
a=1/(sp.sqrt(2*sp.pi*σ**2))
b=sp.exp(-(x-μ)**2/(2*σ**2))
c=a*b
d=c.subs([(μ,0),(σ,1)])
e=sp.integrate(c,x)
f=sp.integrate(c,(x,sp.oo,x))
e2=e.subs([(μ,0),(σ,1)])
f2=f.subs([(μ,0),(σ,1)])
sp.init_printing()
display(a) 
print(sp.latex(a))
display(b) 
print(sp.latex(b))
display(c) 
print(sp.latex(c))
display(d) 
print(sp.latex(d))
display(e) 
print(sp.latex(e))
display(f) 
print(sp.latex(f))
display(e2) 
print(sp.latex(e2))
display(f2) 
print(sp.latex(f2))
sp.plot(d, (x, -4, 4))
sp.plot(e2, (x, -4, 4))

確率密度関数\\
f(x)=\frac{\sqrt{2}}{2 \sqrt{\pi} \sqrt{σ^{2}}}e^{- \frac{\left(x - μ\right)^{2}}{2 σ^{2}}}=\frac{\sqrt{2} e^{- \frac{\left(x - μ\right)^{2}}{2 σ^{2}}}}{2 \sqrt{\pi} \sqrt{σ^{2}}}\\
μ(平均)=0,σ^2=1の時\\
\frac{\sqrt{2} e^{- \frac{x^{2}}{2}}}{2 \sqrt{\pi}}\\
累積分布関数\\
不定積分F(x)=\int f(x)dx=\frac{σ \operatorname{erf}{\left(\frac{\sqrt{2} \left(x - μ\right)}{2 σ} \right)}}{2 \sqrt{σ^{2}}}\\
μ(平均)=0,σ^2=1の時\\
F(x)=\frac{\operatorname{erf}{\left(\frac{\sqrt{2} x}{2} \right)}}{2}\\
定積分F(x)=\int_{-\infty}^{x}f(x)dx=\frac{\sqrt{2} e^{- \frac{μ^{2}}{2 σ^{2}}} \int\limits_{\infty}^{x} e^{- \frac{x^{2}}{2 σ^{2}}} e^{\frac{x μ}{σ^{2}}}\, dx}{2 \sqrt{\pi} \sqrt{σ^{2}}}\\
μ(平均)=0,σ^2=1の時\\
F(x)=\frac{\sqrt{2} \int\limits_{\infty}^{x} e^{- \frac{x^{2}}{2}}\, dx}{2 \sqrt{\pi}}

関数化しても結果は変わりません。

import numpy as np
import sympy as sp
σ,μ,x = sp.symbols('σ,μ,x')

def a(x):
    return 1/(sp.sqrt(2*sp.pi*σ**2))

def b(x):
    return sp.exp(-(x-μ)**2/(2*σ**2))

def c(x):
    return a(x)*b(x)

def d(x):
    return c(x).subs([(μ,0),(σ,1)])

def e(x):
    return sp.integrate(c(x),x)

def f(x):
    return sp.integrate(c(x),(x,sp.oo,x))

def e2(x):
    return e(x).subs([(μ,0),(σ,1)])

def f2(x):
    return f(x).subs([(μ,0),(σ,1)])

sp.init_printing()
display(a(x)) 
print(sp.latex(a(x)))
display(b(x)) 
print(sp.latex(b(x)))
display(c(x)) 
print(sp.latex(c(x)))
display(d(x)) 
print(sp.latex(d(x)))
display(e(x)) 
print(sp.latex(e(x)))
display(f(x)) 
print(sp.latex(f(x)))
display(e2(x)) 
print(sp.latex(e2(x)))
display(f2(x)) 
print(sp.latex(f2(x)))
sp.plot(d(x), (x, -4, 4))
sp.plot(e2(x), (x, -4, 4))

確率密度関数\\
f(x)=\frac{\sqrt{2}}{2 \sqrt{\pi} \sqrt{σ^{2}}}e^{- \frac{\left(x - μ\right)^{2}}{2 σ^{2}}}=\frac{\sqrt{2} e^{- \frac{\left(x - μ\right)^{2}}{2 σ^{2}}}}{2 \sqrt{\pi} \sqrt{σ^{2}}}\\
μ(平均)=0,σ^2=1の時\\
\frac{\sqrt{2} e^{- \frac{x^{2}}{2}}}{2 \sqrt{\pi}}\\
累積分布関数\\
不定積分F(x)=\int f(x)dx=\frac{σ \operatorname{erf}{\left(\frac{\sqrt{2} \left(x - μ\right)}{2 σ} \right)}}{2 \sqrt{σ^{2}}}\\
μ(平均)=0,σ^2=1の時\\
F(x)=\frac{\operatorname{erf}{\left(\frac{\sqrt{2} x}{2} \right)}}{2}\\
定積分F(x)=\int_{-\infty}^{x}f(x)dx=\frac{\sqrt{2} e^{- \frac{μ^{2}}{2 σ^{2}}} \int\limits_{\infty}^{x} e^{- \frac{x^{2}}{2 σ^{2}}} e^{\frac{x μ}{σ^{2}}}\, dx}{2 \sqrt{\pi} \sqrt{σ^{2}}}\\
μ(平均)=0,σ^2=1の時\\
F(x)=\frac{\sqrt{2} \int\limits_{\infty}^{x} e^{- \frac{x^{2}}{2}}\, dx}{2 \sqrt{\pi}}

これでは使い物にならないので、scipy.stats&NumPy&matplotlibを使った方法に切り替えます。

import sympy as sp
import numpy as np
import pandas as pd
from scipy import stats
import matplotlib.pyplot as plt

# 等差数列を生成
X = np.arange(start=-4, stop=4, step=0.1)

# pdfで確率密度関数を生成
norm_pdf = stats.norm.pdf(x=X, loc=0, scale=1) # 期待値(平均)=0, 標準偏差(分散)=1 
norm_cdf = stats.norm.cdf(x=X, loc=0, scale=1) # 期待値(平均)=0, 標準偏差(分散)=1 

plt.style.use('default')
fig = plt.figure()

# 可視化
plt.plot(X, norm_pdf,color="blue", label="Probability Density Function(PDF)")
plt.plot(X, norm_cdf,color="red", label="Cumulative Distribution Function(CDF)")
plt.ylim([-0.1,1.1])
plt.xlim([-4,4])
plt.xlabel("Random Variable")
plt.ylabel("Probability Density")
plt.title("Standard Normal Distribution")
ax = fig.add_subplot(111)
ax.legend(loc='upper left')
plt.show()

import sympy as sp
import numpy as np
import pandas as pd
from scipy import stats
import matplotlib.pyplot as plt

# 等差数列を生成
X = np.arange(start=-6, stop=6, step=0.1)

# pdfで確率密度関数を生成
norm_pdf1 = stats.norm.pdf(x=X, loc=0, scale=0.4) # 期待値(平均)=0, 標準偏差(分散)=1 
norm_pdf2 = stats.norm.pdf(x=X, loc=0, scale=1) # 期待値(平均)=0, 標準偏差(分散)=1 
norm_pdf3 = stats.norm.pdf(x=X, loc=0, scale=4.0) # 期待値(平均)=0, 標準偏差(分散)=1 
norm_pdf4 = stats.norm.pdf(x=X, loc=-4, scale=0.7) # 期待値(平均)=0, 標準偏差(分散)=1 

plt.style.use('default')
fig = plt.figure()

# 可視化
plt.plot(X, norm_pdf1,color="blue", label="μ=0,σ^2=0.4")
plt.plot(X, norm_pdf2,color="red", label="μ=0,σ^2=1.0")
plt.plot(X, norm_pdf3,color="yellow", label="μ=0,σ^2=4.0")
plt.plot(X, norm_pdf4,color="green", label="μ=-4,σ^2=0.7")
plt.ylim([-0.1,1.1])
plt.xlim([-6,6])
plt.xlabel("Random Variable(x)")
plt.ylabel("Probability Density(y)")
plt.title("Probability Density Functions(PDF)")
ax = fig.add_subplot(111)
ax.legend(loc='upper left')
plt.show()

import sympy as sp
import numpy as np
import pandas as pd
from scipy import stats
import matplotlib.pyplot as plt

# 等差数列を生成
X = np.arange(start=-6, stop=6, step=0.1)

# pdfで確率密度関数を生成
norm_cdf1 = stats.norm.cdf(x=X, loc=0, scale=0.4) # 期待値(平均)=0, 標準偏差(分散)=1 
norm_cdf2 = stats.norm.cdf(x=X, loc=0, scale=1) # 期待値(平均)=0, 標準偏差(分散)=1 
norm_cdf3 = stats.norm.cdf(x=X, loc=0, scale=4.0) # 期待値(平均)=0, 標準偏差(分散)=1 
norm_cdf4 = stats.norm.cdf(x=X, loc=-4, scale=0.7) # 期待値(平均)=0, 標準偏差(分散)=1 

plt.style.use('default')
fig = plt.figure()

# 可視化
plt.plot(X, norm_cdf1,color="blue", label="μ=0,σ^2=0.4")
plt.plot(X, norm_cdf2,color="red", label="μ=0,σ^2=1.0")
plt.plot(X, norm_cdf3,color="yellow", label="μ=0,σ^2=4.0")
plt.plot(X, norm_cdf4,color="green", label="μ=-4,σ^2=0.7")
plt.ylim([-0.1,1.1])
plt.xlim([-6,6])
plt.xlabel("Random Variable(x)")
plt.ylabel("Probability Density(y)")
plt.title("Cumulative Distribution Function(CDF)")
ax = fig.add_subplot(111)
ax.legend(loc='lower right')
plt.show()

とりあえず上掲のバクテリア死亡率が正規分布に従うと考えてみましょう。まずはこの場合確率変数xが時間(t)に対応すると考え、平均μを5(時間)に設定します。

import sympy as sp
import numpy as np
import pandas as pd
from scipy import stats
import matplotlib.pyplot as plt

#タイトル
TTL="Probability of Bacterial survival"
#とりあえずの有意最低限
Min=0
#とりあえずの有意最低限
Max=10
# mean(平均)
MN = 5
# standard deviation(分散)
SCL = 1

# 等差数列を生成
X = np.arange(start=Min, stop=Max, step=0.1)
#下側(左側)確率変数
LFT=4 
#上側(右側)確率変数
RGT=6 
#下側(左側)確率密度
LFTq=stats.norm.cdf(LFT, loc=MN, scale=SCL) 
#上側(右側)確率密度
RGTq=stats.norm.sf(RGT, loc=MN, scale=SCL) 

# pdfで確率密度関数を生成
norm_pdf = stats.norm.pdf(x=X, loc=MN, scale=SCL) # 期待値(平均)=0, 標準偏差(分散)=0.5 
norm_cdf = stats.norm.cdf(x=X, loc=MN, scale=SCL) # 期待値(平均)=0, 標準偏差(分散)=0.5 

x0=np.arange(start=LFT,stop=RGT,step=0.1)
xp=np.concatenate((
    [LFT],
    x0,
    [RGT,LFT])) 
yp=np.concatenate((
    [0],
    stats.norm.pdf(x=x0, loc=MN, scale=SCL),
    [0,0])) 

xc=np.concatenate((
    [LFT],
    x0,
    [RGT,LFT])) 
yc=np.concatenate((
    [0],
    stats.norm.cdf(x=x0, loc=MN, scale=SCL),
    [0,0])) 

plt.style.use('default')
fig = plt.figure()

# 可視化
plt.plot(X, norm_pdf,color="blue", label="Probability Density Function(PDF)")
plt.plot(X, norm_cdf,color="red", label="Cumulative Distribution Function(CDF)")
plt.fill(xp,yp,color="skyblue",alpha=0.5)
plt.fill(xc,yc,color="orange",alpha=0.5)
plt.ylim([-0.1,1.1])
plt.xlim([Min,Max])
plt.xlabel("Random Variable(g)")
plt.ylabel("Probability Density")
plt.title(TTL)
ax = fig.add_subplot(111)
ax.legend(loc='upper left')
plt.show()

print("4時間以前に死亡している確率="+str(np.round(LFTq*100,decimals=2))+"%")
print("6時間以降に死亡する確率="+str(np.round(RGTq*100,decimals=2))+"%")
print("4時間から6時間にかけて死亡する確率="+str(np.round((1-LFTq-RGTq)*100,decimals=2))+"%")

4時間以前に死亡している確率=15.87%
6時間以降に死亡する確率=15.87%
4時間から6時間にかけて死亡する確率=68.27%

ううむ「平均寿命4時間~6時間」といったら、この区間内で95%くらいは死んで欲しいものです。そういう方向で分散$σ^2$の最適値を探すと、大体0.5程度に落ち着く様ですね。

import sympy as sp
import numpy as np
import pandas as pd
from scipy import stats
import matplotlib.pyplot as plt

#タイトル
TTL="Probability of Bacterial survival"
#とりあえずの有意最低限
Min=0
#とりあえずの有意最低限
Max=10
# mean(平均)
MN = 5
# standard deviation(分散)
SCL = 0.5

# 等差数列を生成
X = np.arange(start=Min, stop=Max, step=0.1)
#下側(左側)確率変数
LFT=4 
#上側(右側)確率変数
RGT=6 
#下側(左側)確率密度
LFTq=stats.norm.cdf(LFT, loc=MN, scale=SCL) 
#上側(右側)確率密度
RGTq=stats.norm.sf(RGT, loc=MN, scale=SCL) 

# pdfで確率密度関数を生成
norm_pdf = stats.norm.pdf(x=X, loc=MN, scale=SCL) # 期待値(平均)=0, 標準偏差(分散)=0.5 
norm_cdf = stats.norm.cdf(x=X, loc=MN, scale=SCL) # 期待値(平均)=0, 標準偏差(分散)=0.5 

x0=np.arange(start=LFT,stop=RGT,step=0.1)
xp=np.concatenate((
    [LFT],
    x0,
    [RGT,LFT])) 
yp=np.concatenate((
    [0],
    stats.norm.pdf(x=x0, loc=MN, scale=SCL),
    [0,0])) 

xc=np.concatenate((
    [LFT],
    x0,
    [RGT,LFT])) 
yc=np.concatenate((
    [0],
    stats.norm.cdf(x=x0, loc=MN, scale=SCL),
    [0,0])) 

plt.style.use('default')
fig = plt.figure()

# 可視化
plt.plot(X, norm_pdf,color="blue", label="Probability Density Function(PDF)")
plt.plot(X, norm_cdf,color="red", label="Cumulative Distribution Function(CDF)")
plt.fill(xp,yp,color="skyblue",alpha=0.5)
plt.fill(xc,yc,color="orange",alpha=0.5)
plt.ylim([-0.1,1.1])
plt.xlim([Min,Max])
plt.xlabel("Random Variable(g)")
plt.ylabel("Probability Density")
plt.title(TTL)
ax = fig.add_subplot(111)
ax.legend(loc='upper left')
plt.show()

print("4時間以前に死亡している確率="+str(np.round(LFTq*100,decimals=2))+"%")
print("6時間以降に死亡する確率="+str(np.round(RGTq*100,decimals=2))+"%")
print("4時間から6時間にかけて死亡する確率="+str(np.round((1-LFTq-RGTq)*100,decimals=2))+"%")

4時間以前に死亡している確率=2.28%
6時間以降に死亡する確率=2.28%
4時間から6時間にかけて死亡する確率=95.45%

条件を緩めて4時間から6時間の間に死亡する確率が90%程度で良いとするなら分散は0.6前後に落ち着きます。

import sympy as sp
import numpy as np
import pandas as pd
from scipy import stats
import matplotlib.pyplot as plt

#タイトル
TTL="Probability of Bacterial survival"
#とりあえずの有意最低限
Min=0
#とりあえずの有意最低限
Max=10
# mean(平均)
MN = 5
# standard deviation(分散)
SCL = 0.6

# 等差数列を生成
X = np.arange(start=Min, stop=Max, step=0.1)
#下側(左側)確率変数
LFT=4 
#上側(右側)確率変数
RGT=6 
#下側(左側)確率密度
LFTq=stats.norm.cdf(LFT, loc=MN, scale=SCL) 
#上側(右側)確率密度
RGTq=stats.norm.sf(RGT, loc=MN, scale=SCL) 

# pdfで確率密度関数を生成
norm_pdf = stats.norm.pdf(x=X, loc=MN, scale=SCL) # 期待値(平均)=0, 標準偏差(分散)=0.5 
norm_cdf = stats.norm.cdf(x=X, loc=MN, scale=SCL) # 期待値(平均)=0, 標準偏差(分散)=0.5 

x0=np.arange(start=LFT,stop=RGT,step=0.1)
xp=np.concatenate((
    [LFT],
    x0,
    [RGT,LFT])) 
yp=np.concatenate((
    [0],
    stats.norm.pdf(x=x0, loc=MN, scale=SCL),
    [0,0])) 

xc=np.concatenate((
    [LFT],
    x0,
    [RGT,LFT])) 
yc=np.concatenate((
    [0],
    stats.norm.cdf(x=x0, loc=MN, scale=SCL),
    [0,0])) 

plt.style.use('default')
fig = plt.figure()

# 可視化
plt.plot(X, norm_pdf,color="blue", label="Probability Density Function(PDF)")
plt.plot(X, norm_cdf,color="red", label="Cumulative Distribution Function(CDF)")
plt.fill(xp,yp,color="skyblue",alpha=0.5)
plt.fill(xc,yc,color="orange",alpha=0.5)
plt.ylim([-0.1,1.1])
plt.xlim([Min,Max])
plt.xlabel("Random Variable(g)")
plt.ylabel("Probability Density")
plt.title(TTL)
ax = fig.add_subplot(111)
ax.legend(loc='upper left')
plt.show()

print("4時間以前に死亡している確率="+str(np.round(LFTq*100,decimals=2))+"%")
print("6時間以降に死亡する確率="+str(np.round(RGTq*100,decimals=2))+"%")
print("4時間から6時間にかけて死亡する確率="+str(np.round((1-LFTq-RGTq)*100,decimals=2))+"%")

4時間以前に死亡している確率=4.78%
6時間以降に死亡する確率=4.78%
4時間から6時間にかけて死亡する確率=90.44%

確率変数Xと確率密度Yに平均μや分散$σ^2$が与える影響は大体こういうイメージ。サンプル数が十分大きい場合には、大数の法則や中心極限定理を根拠にこの様な「母集合が正規分布に従う」前提の推計も記述統計的に扱われるケースが多く、この数理を前提とするなら、こんな荒っぽい手口も使えるという話なんですね。かくして記述統計の概念は大幅な拡張を受ける展開を迎えた訳です。

補助的関数の追加

確率密度空間と累積分布空間の実際の計算では、さらに複数の補助的関数が使われます。

例題その1(パーセント点関数)

逆に分布中のパーセント点から確率密度を求める手口も存在します。
パーセント点の求め方

データ集合が特定の分布に従うとき、そのデータの範囲を例えば「パーセント点」という形で確率的に評価することが可能となります。

例題
「農園Aのみかんの重量が以下の正規分布に従うことがわかっているとします。

• 平均：100g
• 標準偏差：5g

この場合、農園Aの適当なみかんはおおよそ90%の確率で92g以上、108g以下に収まると推定されます。」

import sympy as sp
import numpy as np
import pandas as pd
from scipy import stats

#農園Aのみかんの重量

# mean(平均)
loc = 100
# standard deviation(分散)
scale = 5

#下側(左側)5%点
print("下側(左側)5%点") 
print(stats.norm.ppf(q=0.05, loc=loc, scale=scale) )

#上側(右側)5%点
print("上側(右側)5%点") 
print(stats.norm.ppf(q=0.95, loc=loc, scale=scale) )

下側(左側)5%点
91.77573186524263
上側(右側)5%点
108.22426813475737

確率密度空間と確率分布空間に対応する全体像

import sympy as sp
import numpy as np
import pandas as pd
from scipy import stats
import matplotlib.pyplot as plt

#タイトル
TTL="Harvested Vegetable Weight"
#とりあえずの有意最低限
Min=80
#とりあえずの有意最低限
Max=120
# mean(平均)
MN = 100
# standard deviation(分散)
SCL = 5

# 等差数列を生成
X = np.arange(start=Min, stop=Max, step=0.1)
#下側(左側)確率密度
LFTq=0.05 
#上側(右側)確率密度
RGTq=0.95 
#下側(左側)確率変数
LFT=stats.norm.ppf(LFTq, loc=MN, scale=SCL) 
#上側(右側)確率変数
RGT=stats.norm.ppf(RGTq, loc=MN, scale=SCL) 

# pdfで確率密度関数を生成
norm_pdf = stats.norm.pdf(x=X, loc=MN, scale=SCL) # 期待値(平均)=0, 標準偏差(分散)=0.5 
norm_cdf = stats.norm.cdf(x=X, loc=MN, scale=SCL) # 期待値(平均)=0, 標準偏差(分散)=0.5 

x0=np.arange(start=LFT,stop=RGT,step=0.1)
xp=np.concatenate((
    [LFT],
    x0,
    [RGT,LFT])) 
yp=np.concatenate((
    [0],
    stats.norm.pdf(x=x0, loc=MN, scale=SCL),
    [0,0])) 

xc=np.concatenate((
    [LFT],
    x0,
    [RGT,LFT])) 
yc=np.concatenate((
    [0],
    stats.norm.cdf(x=x0, loc=MN, scale=SCL),
    [0,0])) 

plt.style.use('default')
fig = plt.figure()

# 可視化
plt.plot(X, norm_pdf,color="blue", label="Probability Density Function(PDF)")
plt.plot(X, norm_cdf,color="red", label="Cumulative Distribution Function(CDF)")
plt.fill(xp,yp,color="skyblue",alpha=0.5)
plt.fill(xc,yc,color="orange",alpha=0.5)
plt.ylim([-0.1,1.1])
plt.xlim([Min,Max])
plt.xlabel("Random Variable(g)")
plt.ylabel("Probability Density")
plt.title(TTL)
ax = fig.add_subplot(111)
ax.legend(loc='upper left')
plt.show()

print("5%以下="+str(np.round(LFT,decimals=2))+"g以下")
print("95%以上="+str(np.round(RGT,decimals=2))+"g以上")

5%以下=91.78g以下
95%以上=108.22g以上

ここでq値に対応する確率変数を導出するのに用いたパーセント点関数(PPF=Percent Point Function)は累積分布関数CDFの逆関数で、累積分布関数をq％と指定するとその値をとる変数を返します。
SciPy.stats

import numpy as np
from scipy import stats
import matplotlib.pyplot as plt

p = np.linspace(0,1,100)
y_ppf = stats.norm.ppf(p)  # 標準正規分布
plt.plot(p,y_ppf)
plt.xlabel('P')
plt.ylabel('Random Variable(x)')
plt.title('Percent Point Function')
plt.show()

###例題その2(CDF関数そのもの)
パーセント点の求め方

また逆に、正規分布を前提としてX軸のある点が何％点となるのか？を求める場合もあります。

例題
「ある工場では製品に品質管理用の製品評価指数を設けており、この指数は以下の正規分布に従うことがわかっているとします。

• 平均：100
• 標準偏差：3

この製品評価指数が90を下回る製品は不良品扱いとして出荷しないことにしました。不良品の発生は割合は1万個に4個程度と推定されます。」

import sympy as sp
import numpy as np
import pandas as pd
from scipy import stats

#工場の不良品発生率

# mean
loc = 100
# standard deviation
scale = 3

#不良品発生率
print("不良品発生率") 
print(stats.norm.cdf(x=90, loc=loc, scale=scale) )

不良品発生率
0.0004290603331968372

確率密度空間と確率分布空間に対応する全体像

import sympy as sp
import numpy as np
import pandas as pd
from scipy import stats
import matplotlib.pyplot as plt

#タイトル
TTL="Defective Product Rate"
#とりあえずの有意最低限
Min=88
#とりあえずの有意最低限
Max=90.5
# mean(平均)
MN = 100
# standard deviation(分散)
SCL = 3

# 等差数列を生成
X = np.arange(start=Min, stop=Max, step=0.001)
#下側(左側)確率変数(限りなくゼロに近い)
LFT=0.0000001 
#上側(右側)確率変数
RGT=90
#下側(左側)確率密度(限りなくゼロに近い)
LFTq=LFT 
#上側(右側)確率密度
RGTq=stats.norm.cdf(RGT, loc=MN, scale=SCL) 

# pdfで確率密度関数を生成
norm_pdf = stats.norm.pdf(x=X, loc=MN, scale=SCL) # 期待値(平均)=0, 標準偏差(分散)=0.5 
norm_cdf = stats.norm.cdf(x=X, loc=MN, scale=SCL) # 期待値(平均)=0, 標準偏差(分散)=0.5 

x0=np.arange(start=LFT,stop=RGT,step=0.001)
xp=np.concatenate((
    [LFT],
    x0,
    [RGT,LFT])) 
yp=np.concatenate((
    [0],
    stats.norm.pdf(x=x0, loc=MN, scale=SCL),
    [0,0])) 

xc=np.concatenate((
    [LFT],
    x0,
    [RGT,LFT])) 
yc=np.concatenate((
    [0],
    stats.norm.cdf(x=x0, loc=MN, scale=SCL),
    [0,0])) 

plt.style.use('default')
fig = plt.figure()

# 可視化
plt.plot(X, norm_pdf,color="blue", label="Probability Density Function(PDF)")
plt.plot(X, norm_cdf,color="red", label="Cumulative Distribution Function(CDF)")
plt.fill(xp,yp,color="skyblue",alpha=0.5)
plt.fill(xc,yc,color="orange",alpha=0.5)
plt.ylim([0.000,0.0005])
plt.xlim([Min,Max])
plt.xlabel("Random Variable(g)")
plt.ylabel("Probability Density")
plt.title(TTL)
ax = fig.add_subplot(111)
ax.legend(loc='upper left')
plt.show()

print("不良品発生率="+str(np.round(RGTq*10000,decimals=2))+"個/10000個")

不良品発生率=4.29個/10000個

想像していたより苦戦させられました。

• とりあえず「下側数値と上側数値の2点表示」を流用すべく下側数値を0,上側数値を算出値と置いたらStat関数が0%や確率密度0を処理出来ない事が判明(確かにそれは分布上無限遠点の彼方に存在する訳である)。とりあえず疑似的に0.0000001を置く。
• 全体像が俯瞰可能な表示スケールの模索。もはや人間の主観の領域？
• CDFから下ろす法線が斜めに傾く。x軸の分割数を増大させる事で緩和に向かうが、処理が加速度的に重くなっていく。さらに範囲関係で工夫を重ねないと実用にならない感じ。

そういえばこの領域ではPDFの値がCDFの累積結果を上回るのですね。実際に目にすると感慨ひとしおです。

このパーセント点を求める操作が統計学の入門部分の一つの肝で、ここができれば今後説明する区間推定や帰無仮説検定の理解がスムーズに進むかと思います。また、今回は正規分布を利用しましたが、正規分布以外でもパーセント点を求められる場合が多々ありますので様々に活用することができるようになるはずです。

例題その3(生存関数)

標準正規分布表の使い方1

例題

「あるクラスのテスト結果は平均72.8点、標準偏差15点の正規分布に従っています。この時、88点以上の人は15.6%と推定されます。」

import sympy as sp
import numpy as np
import pandas as pd
from scipy import stats

#あるクラスのテスト結果は平均72.8点、標準偏差15点
#の正規分布に従っています。この時、88点以上の人は
#何%いるでしょうか。(答え15.6%)

z =  (88- 72.8) /15
print("標準化データZ")
print(z)

#この場合、求めるのは上側確率
print("上側確率") 
print(stats.norm.sf(z))

#ちなみに下側確率
print("下側確率") 
print(stats.norm.cdf(z))

#上側確率+下側確率=1
print("上側確率+下側確率") 
print(stats.norm.sf(z)+stats.norm.cdf(z))

#どちらか覚えれば十分ですが、ただし
#1からsfを引くよりは直接cdfを使った方が
#精度は良くなる可能性があります。

標準化データZ
1.0133333333333334
上側確率
0.15545048547780277
下側確率
0.8445495145221973
上側確率+下側確率
1.0

確率密度空間と確率分布空間に対応する全体像

import sympy as sp
import numpy as np
import pandas as pd
from scipy import stats
import matplotlib.pyplot as plt

#タイトル
TTL="Test Results"
#とりあえずの有意最低限
Min=0
#とりあえずの有意最低限
Max=100
# mean(平均)
MN = 72.8
# standard deviation(分散)
SCL = 15

#データ標準化
z =  (88- MN) /SCL

# 等差数列を生成
X = np.arange(start=Min, stop=Max, step=0.1)
#下側(左側)確率密度
LFTq=stats.norm.sf(88,loc=MN, scale=SCL) 
#上側(右側)確率密度(100点満点)
RGTq=1.0
#下側(左側)確率変数
LFT=88
#上側(右側)確率変数(100点満点)
RGT=100

# pdfで確率密度関数を生成
norm_pdf = stats.norm.pdf(x=X, loc=MN, scale=SCL) # 期待値(平均)=0, 標準偏差(分散)=0.5 
norm_cdf = stats.norm.cdf(x=X, loc=MN, scale=SCL) # 期待値(平均)=0, 標準偏差(分散)=0.5 

x0=np.arange(start=LFT,stop=RGT,step=0.1)
xp=np.concatenate((
    [LFT],
    x0,
    [RGT,LFT])) 
yp=np.concatenate((
    [0],
    stats.norm.pdf(x=x0, loc=MN, scale=SCL),
    [0,0])) 

xc=np.concatenate((
    [LFT],
    x0,
    [RGT,LFT])) 
yc=np.concatenate((
    [0],
    stats.norm.cdf(x=x0, loc=MN, scale=SCL),
    [0,0])) 

plt.style.use('default')
fig = plt.figure()

# 可視化
plt.plot(X, norm_pdf,color="blue", label="Probability Density Function(PDF)")
plt.plot(X, norm_cdf,color="red", label="Cumulative Distribution Function(CDF)")
plt.fill(xp,yp,color="skyblue",alpha=0.5)
plt.fill(xc,yc,color="orange",alpha=0.5)
plt.ylim([-0.1,1.0])
plt.xlim([Min,Max])
plt.xlabel("Random Variable(g)")
plt.ylabel("Probability Density")
plt.title(TTL)
ax = fig.add_subplot(111)
ax.legend(loc='upper left')
plt.show()

print("標準化データ"+str(np.round(z,decimals=2)))
print("88点以上の得点者率="+str(np.round(LFTq*100,decimals=2))+"%")

標準化データ1.01
88点以上の得点者率=15.55%

この作図過程ではこんな事を学びました。

• 結局、標準化データは計算で使わなかった。おそらくこの様な計算が必要だったのは紙の標準正規分布表を検索する都合上であって、コンピューター時代には平均μや標準偏差$σ^2$の値ごとライブラリに丸投げしてしまう方が手取り早いのである。この辺りの都合、コンピューターの普及によって常用対数表の検索技量が不要となっていった過程に対応しそうである。
  【数理考古学】常用対数表(Table of Common Logarithms)を使った計算
• 上掲の不良品率の計算に際して「下側確率が無限遠点の彼方の0に設定出来ない」問題との格闘があった様に、元来上側確率の計算に際しても、それが「無限遠点の彼方の0に設定出来ない」問題が存在する訳だが、テストの成績評価に関しては0点以下も100点以上も存在しないので、この問題を意識する必要がない。ただしこの事はおそらく「テストの成績によって評価可能な事」の限界とも密接に関係してくる。
  【無限遠点を巡る数理】無限遠点(Infinity)としての正規分布と分散概念の歴史

またここではさりげなく累積分布関数CDFとは逆に、確率変数Xがある値x以上となる確率を計算する生存関数(Survival Function)が導入されています。まぁ実は上掲のバクテリア生存率計算でも既にこっそり使っていた訳ですが。

1. Pythonで学ぶ統計学　2. 確率分布[scipy.stats徹底理解]

import numpy as np
from scipy import stats
import matplotlib.pyplot as plt

X = np.linspace(start=1, stop=7, num=100)
y_ppf = stats.norm.sf(x=X, loc=4, scale=0.8)  # 標準正規分布
plt.plot(x,y_ppf)
plt.xlabel('Random Variable(x)')
plt.ylabel('p')
plt.title('Survival Function')
plt.show()

｢1-cdf｣と定義されることもありますが、sfの方がより正確な場合もあるとの指摘も見られます。

例題その4(逆生存関数)

逆に累積確率から横軸の値(確率変数x)の値を取り出したい時には逆生存関数(ISF=Inverse Survival Function)を使います。
1. Pythonで学ぶ統計学　2. 確率分布[scipy.stats徹底理解]

逆生存関数は、生存関数sfの逆関数ですが、やることはパーセント点関数ppfと逆のことです。生存関数の残余をq％と指定するとその値をとる変数を返します。

import numpy as np
from scipy import stats
import matplotlib.pyplot as plt

p = np.arange(start=0.001, stop=1, step=0.001)
pdf0=[]
for num in range(len(p)):
    pdf0.append(stats.norm.isf(q=p[num], loc=4, scale=0.8))
y_ppf=np.array(pdf0)
plt.plot(p,y_ppf)
plt.xlabel('p')
plt.ylabel('Random Variable(x)')
plt.title('Inverse survival function')
plt.show()

Python​/標準正規分布

import sympy as sp
import numpy as np
import pandas as pd
from scipy import stats

#累積確率から横軸の値を取り出す。
print("累積確率から横軸の値を取り出す") 
print(stats.norm.ppf(q=0.9))

#上側確率から横軸の値を取り出す
print("上側確率から横軸の値を取り出す") 
print(stats.norm.isf(q=0.9))

#標準正規分布の性質から、累積確率の結果とは符号が逆転します。

#上側確率から横軸の値を求める
print("95%信頼区間") 
print(stats.norm.isf(q=0.025))

#90%信頼区間の係数
print("90%信頼区間") 
print(stats.norm.isf(q=0.05))

#一般に100(1-a)%の信頼区間の係数を求める
#ならnorm.isf(q=a/2)と書けます。

累積確率から横軸の値を取り出す
1.2815515655446004
上側確率から横軸の値を取り出す
-1.2815515655446004
95%信頼区間
1.9599639845400545
90%信頼区間
1.6448536269514729

符号反転するし使い道が良くわかりません。とりあえずメモがてら…

信頼区間推定(Confidence Interval Estimation)

繰り返しになりますが、この様に推定統計の結果をあたかも記述統計の結果の様に断定的に扱えるのは母集団の平均(母平均)と分散(母分散)が明らかになっている場合、あるいはサンプル数が十分大きくて大数の法則や中心極限定理による「母集合の正規分布に従う」条件が使える場合だけに限ります。この辺りが明確でない場合にまず使われるのが信頼区間推定(Confidence Interval Estimation)となります。
19-1. 区間推定とは

母平均の区間推定では、母分散が分かっている場合と分からない場合とで、その算出方法が異なります。

• 母分散が分かっている場合（母分散既知）…母分散の値を使い、標準正規分布を用いて信頼区間を算出する。
  19-2. 母平均の信頼区間の求め方（母分散既知）

• 母分散が分からない場合（母分散未知）…次の式から算出される不偏分散$S^2$の値を使い、t分布を用いて信頼区間を算出する。
  20-1. 標本とt分布

s^2=\frac{1}{n-1}\sum_{i=1}^{n}(x_i-\bar{x})^2

ただし、母平均が分からないのに母分散だけは分かっているという状況は現実にはほとんどないので、通常母平均の区間推定を行う場合にはt分布を用いた方法が使われます。

t分布は自由度によって分布の形が変わり(正規分布に近づく)、サンプルサイズが30以上であれば、ほぼ正規分布と同一視可能となります。要するにこれまで述べてきた「サンプル数が十分大きくて大数の法則や中心極限定理による母集合の正規分布近似が使える場合」とはこのケースを指していたとも言えそうなのです。

t分布の確率密度関数

import numpy as np
from scipy import stats
import matplotlib.pyplot as plt

x = np.linspace(-4,4,100)
y1 = stats.t.pdf(x, df=1)
y2 = stats.t.pdf(x, df=2)
y3 = stats.t.pdf(x, df=5)
y4 = stats.t.pdf(x, df=30)

plt.style.use('default')
fig = plt.figure()

plt.plot(X, y1,color="blue", label="df=1")
plt.plot(X, y2,color="red", label="df=2")
plt.plot(X, y3,color="yellow", label="df=5")
plt.plot(X, y4,color="green", label="df=30")

plt.xlabel('Random Variable(x)')
plt.ylabel('p')
plt.title("Student's t-distribution")

ax = fig.add_subplot(111)
ax.legend(loc='upper left')
plt.show()

t分布の累積分布関数

import numpy as np
from scipy import stats
import matplotlib.pyplot as plt

x = np.linspace(-4,4,100)
y1 = stats.t.cdf(x, df=1)
y2 = stats.t.cdf(x, df=2)
y3 = stats.t.cdf(x, df=5)
y4 = stats.t.cdf(x, df=30)

plt.style.use('default')
fig = plt.figure()

plt.plot(X, y1,color="blue", label="df=1")
plt.plot(X, y2,color="red", label="df=2")
plt.plot(X, y3,color="yellow", label="df=5")
plt.plot(X, y4,color="green", label="df=30")

plt.xlabel('Random Variable(x)')
plt.ylabel('p')
plt.title("Student's t-distribution")

ax = fig.add_subplot(111)
ax.legend(loc='upper left')
plt.show()

Scipy.statには信頼区間推定(Confidence Interval Estimation)を直接行うinterval関数が搭載されています。

1. Pythonで学ぶ統計学　2. 確率分布[scipy.stats徹底理解]

import sympy as sp
import numpy as np
import pandas as pd
from scipy import stats
import matplotlib.pyplot as plt

#タイトル
TTL="Confidence Interval Estimation"
#とりあえずの有意最低限
Min=1
#とりあえずの有意最低限
Max=7
# mean(平均)
MN = 4
# standard deviation(分散)
SCL = 0.8

# intervalで信頼区間95％に当たる変数を取得
# 信頼係数=0.95, 期待値=4, 標準偏差=0.8
lower, upper = stats.norm.interval(alpha=0.95, loc=MN, scale=SCL) 

# 等差数列を生成
X = np.arange(start=Min, stop=Max, step=0.1)
#下側(左側)確率変数
LFT=lower 
#上側(右側)確率変数
RGT=upper
#下側(左側)確率密度
LFTq=stats.norm.cdf(LFT, loc=MN, scale=SCL) 
#上側(右側)確率密度
RGTq=stats.norm.sf(RGT, loc=MN, scale=SCL) 


# pdfで確率密度関数を生成
norm_pdf = stats.norm.pdf(x=X, loc=MN, scale=SCL) # 期待値(平均)=0, 標準偏差(分散)=0.5 
norm_cdf = stats.norm.cdf(x=X, loc=MN, scale=SCL) # 期待値(平均)=0, 標準偏差(分散)=0.5 

x0=np.arange(start=LFT,stop=RGT,step=0.1)
xp=np.concatenate((
    [LFT],
    x0,
    [RGT,LFT])) 
yp=np.concatenate((
    [0],
    stats.norm.pdf(x=x0, loc=MN, scale=SCL),
    [0,0])) 

xc=np.concatenate((
    [LFT],
    x0,
    [RGT,LFT])) 
yc=np.concatenate((
    [0],
    stats.norm.cdf(x=x0, loc=MN, scale=SCL),
    [0,0])) 

plt.style.use('default')
fig = plt.figure()

# 可視化
plt.plot(X, norm_pdf,color="blue", label="Probability Density Function(PDF)")
plt.plot(X, norm_cdf,color="red", label="Cumulative Distribution Function(CDF)")
plt.fill(xp,yp,color="skyblue",alpha=0.5)
plt.fill(xc,yc,color="orange",alpha=0.5)
plt.ylim([-0.1,1.1])
plt.xlim([Min,Max])
plt.xlabel("Random Variable(g)")
plt.ylabel("Probability Density")
plt.title(TTL)
ax = fig.add_subplot(111)
ax.legend(loc='upper left')
plt.show()

print('信頼区間95％の下限：', lower)
print('信頼区間95％の上限：', upper)

信頼区間95％の下限： 2.4320288123679568
信頼区間95％の上限： 5.567971187632043

• 区間推定では母数θを含む区間の確率を指定します。この例ではalpha=0.95がそれで、これを信頼係数と呼びます。多くの場合、99％、95％、90％が使われますが、これは任意です。
• その指定された確率のもとで母数θθを含む区間の下限と上限が求められます。この例では95％の確率で、母数θが確率変数2.43～5.57の間に含まれていることを示しています。

19-1. 区間推定とは

母平均の区間推定では「95%信頼区間（95％CI）」を求めることが多いですが、これは「母集団から標本を取ってきて、その標本平均から母平均の95%信頼区間を求める、という作業を100回やったときに、95回はその区間の中に母平均が含まれる」ということを意味します。

また、95％信頼区間以外に「99%信頼区間」や「90%信頼区間」といった区間を求めることもあります。この95％や99%、90%のような、ある区間に母数が含まれる確率のことを「信頼係数」あるいは「信頼度」といいます。

この辺りまで来ると次第に「記述統計の一部」とはいえなくなってくる様です。

統計的仮説検定(Testing of Statistical Hypothes)

実は以下の理解でつまづいたのが今回の投稿の発端だったのです。

とある記述(大学時代に学習した内容との再邂逅に過ぎないにも関わらず、最初読んだ時は上手く頭に入ってこなかった)
Pythonで統計学を学ぶ(4)

以下の手順に従う。

①母集団に関する帰無仮説と対立仮説を設定する。

「帰無仮説」とは提案する手法が従来の手法と「差がない」、すなわち提案する手法は「効果がない」という本来主張したいこととは逆の仮説である。 この仮説が棄却されることを目標として仮説検定を行う。具体的には、母平均μ=0(母平均は0である),母相関係数ρ=0(相関がない),母平均の差μ1−μ2=0(差がない)というような仮説を指す。

一方「対立仮説」とは帰無仮説が棄却されたときに採択される逆の仮説であり、実験などで示したい・主張したいことを表したもの。具体的には、母平均μ≠0(母平均は0でない),母相関係数ρ≠0(相関がある),母平均の差μ1−μ2≠0(差がある)というような仮説を指す。

対立仮説の設定により、検定は次のどちらかで行う(両側検定の方がより厳しい条件であり、普通は両側検定で行う)。

• 片側検定：対立仮説が、母平均μ＞0(もしくはμ＜0) 、母相関係数ρ＞0(もしくはρ＜0)、母平均の差μ1＞μ2(もしくは μ1＜μ2)の場合
• 両側検定：対立仮説が、母平均μ≠0、母相関係数ρ≠0、母平均の差μ1-μ2≠0の場合。要するに、両側検定では例えば母平均μ≠0を調べるには、母平均μ>0とμ<0の両方を調べなければならない。

帰無仮説が正しいものとして分析を行う。 実際に得られたデータから計算された検定統計量の値によって採択を判断する。帰無仮説が正しいとしたとき、検定統計量が、ほぼ起こり得ない値(それほど極端な値)であれば、帰無仮説を棄却する(つまり、本来の主張を表す対立仮説が採択される)。そうでなければ(確率的に十分起こりうるような値であれば、帰無仮説を採択する(この場合は、本来主張したかった対立仮説が棄却されてしまう)。

②検定統計量を選択する。

「検定統計量」とは統計的仮説検定のために用いられる標本統計量のこと。代表的な検定統計量の例:t,χ2、F

③有意水準αの値を決定する。

「有意水準」とは対立仮説を採択するか決定するときに基準である(αで表されます)。

• 有意水準は5％または1％(α=0.05またはα=0.01)に設定することが多い(つまり、標本が100回に5回(5%の場合)以下にしか現れないデータであった---こんなことは偶然では起こりえない---、だから帰無仮説が成り立たないと考えて良いのではないか、という判断基準)
• 帰無仮説が正しいものとして考えた時の標本分布を「帰無分布」という---帰無分布に基づいて確率が計算される。
• 棄却域以外の部分すなわち「帰無仮説が採択される領域」を「採択域」という。
• 棄却域と採択域の境目の値を「臨界値」という。
• 帰無仮説のもとで、非常に生じにくい検定統計量の値の範囲を「棄却域」という---帰無仮説が棄却される領域(だから、この範囲に入るのが『望ましい』)。

④データを収集した後、データから検定統計量の実現値を求める。

「検定統計量の実現値」とは実際のデータ（手に入った標本）を基に計算してえられる具体的な値のことで、対立仮説に合うほど 0から離れた値を示す。

⑤検定統計量の実現値が棄却域に入れば帰無仮説を棄却し、対立仮説を採択する。そうでなければ、帰無仮説を採択する。

帰無仮説が正しいという仮定のもとで、標本から計算した検定統計量の実現値以上の値が得られる確率を「p値」という。p値が有意水準より小さい時に帰無仮説を棄却する。

• p値の大きさこそが対立仮説を採択する(帰無仮説を棄却する)決め手となる。p値が小さいということは 『帰無仮説が正しいとすると』確率的にほとんど起こりえないことが起きた(有意水準が5%なら100回中5回以下、1%なら100回中1回以下)ということを意味する。逆にp値が大きいということは、確率的にはよくあることが起きた(だから、この結果では差があるとはいえない)ということになる。

検定統計量の実現値が棄却域に入った場合、「差がない」という帰無仮説を棄却し、 「差がある」という対立仮説を採択する。この時は以下の様に記述する。

• [検定結果は5%(または1%)水準で有意である」
• 「p<.05(またはp<.01)で有意差が見られた」 。

帰無仮説が棄却できない場合は、「検定の結果、差が有意でなかった」または「有意差が認められなかった」と書く。

##また別の記述(直接関係ない内容にも関わらず、再理解の突破口となる)
[６]　母集団と標本

「推測統計」では母集団からとられた比較的小さな標本を手がかりとして母集団のあるべき姿を推測する。標本に基づいて計算された量（統計量）から、母集団の数理的特性（母数）を推定する。

①標本抽出の方法は主に以下。

• ランダムサンプリング…母集団の中の要素を、すべての要素が等しい確率で、無作為に抽出する。
• 層化無作為抽出法…母集団を予めいくつかの層に分け、その層の中で無作為抽出する。

②この結果「標本抽出分布=ある母集団から、特定サイズの標本の抽出を繰り返したときに、それらの標本から得られる統計量（平均など）の分布」が得られる。その特徴は以下。

• 理論的標本抽出分布の平均は、母平均と同じ…標本抽出を多数繰り返すと、標本抽出分布の平均値や標準偏差は理論値に近づく。
• 理論的標本抽出分布の標準偏差（標準誤差）の算出方法は以下で、その値は標本サイズｎが大きくなるほど、小さくなっていく。

理論的標準誤差＝\frac{母標準偏差σ}{\sqrt{標本サイズ}}

• 母標準偏差σが不明の時は、標本標準偏差を代用して、標準誤差を推定する。

推定標準誤差＝\frac{標本標準偏差S}{\sqrt{標本サイズ}}

標本抽出分布と正規性（無作為抽出の場合

• 母集団分布が正規(母集団の分布が正規分布)であれば、平均値の標本抽出分布も正規分布となる。
• 母集団が正規分布でなくても、標本サイズが大きくなると、平均値の標本抽出分布は正規分布に近づく。（中心極限定理）

標本誤差自体については、以下でも触れてますね。
【Python演算処理】行列演算の基本④大源流における記述統計学との密接な関連性?

両記述の合算

まずどちらも「正規分布絶対主義」時代独特の思考様式に基づく記述である事から出発しないといけません。機械学習理論が深層学習概念に到達するまでしばらく停滞期を迎えた様に、統計学の世界もまた「単一ガウス・モデル」から主成分分析（PCA=Principal Component Analysis）やEMアルゴリズムが前提とする「混合ガウス・モデル」に移行するまでに相応の停滞期を迎えていたのです。
ディープラーニング - Wikipedia
主成分分析 - Wikipedia
EMアルゴリズム - Wikipedia

【初心者向け】正規分布(Normal Distribution)とは何か？
混合ガウス分布（GMM=Gaussian Mixture Model)

• とりあえず前者記述には「主語=操作対象」に関する記述が大幅に欠落しており、それが後者記述によって補われた感じとなります。こうして全体像を俯瞰する立場から、前者記述を改めて読み直すと、その発想の本質自体は十分「混合ガウス・モデルの時代」にも通用する事が明らかとなるのです。
• 意外だったのが、多変量解析では活躍する「Z値計算」の出番が結局全然なかった事。
  【Python演算処理】行列演算の基本④大源流における記述統計学との密接な関連性?
  逆を言えばここまでは「(各評価軸同士の対応を擦り合わせる必要のない)単変量解析の世界」だったという結論になりそうですね。
  多変量解析の活用

そんな感じで以下続報…