母分散の区間推定（Python）
今回の設定
母集団のμとσ^2が未知のとき、標本を用いて母集団の分散を区間推定する。
ラインに流れる商品のばらつきを調べるときなどに使われる。
分散σ^2の正規分布に従う独立な確率変数（X1.X2...Xn）があるとする。
このとき、T＝（n-1）U^2/σ^2は自由度n-1のカイ二乗分布に従う
※U^2は不偏分散を示す
Pythonで計算するには、scipyのstats.chi2.intervalを用いる。
例題
内容量が100gとされているお菓子を10個選択して、重さを測定した。
お菓子の内容量の分布は正規分布に従うものとし、このデータから母分散の95%信頼区間を求める。
import scipy as sp
import numpy as np
from scipy import statsdata = [102, 99, 100, 100, 105, 107, 103, 100, 99, 100]
mu = sp.mean(data) #標本平均
n = len(data) #サンプルサイズ
df = n - 1 #自由度
alpha = 0.95 #信頼係数
sigma2 = np.var(data)#不偏分散bottom, up= stats.chi2.interval(alpha=0.95, df=df)#下側2.5%点と上側2.5%点
#下側2.5%点 <= 自由度＊不偏分散/σ^2 <= 上側2.5%点
#自由度＊不偏分散/上側2.5%点 <= σ^2 <= 自由度*不偏分散/下側2.5%点(f"{round(df*sigma2/up, 3)}<= σ^2 <={round(df*sigma2/bottom, 3)}")
>'3.146<= σ^2 <=22.163'
よって、母分散の95%信頼区間は3.146<= σ^2 <=22.163となる。
標本平均 = 101.5
参考
母比率の推定（Python）
母比率の推定とは？
例えば集団の中で、○✗クイズをしたときに、○と答えた人の確率を母比率pという。
一方、標本の○率を標本比率Rとすると、Rは母比率ｐに対応する。
よって、標本比率Rを用いてｐを推定することが可能になります。
推定には中心極限定理を利用するため、サンプルサイズは十分に取る必要がある。
Pythonで計算するには、scipyのstats.binom.intervalを用いて推定が可能になります。
例題
サイコロを400回投げたところ、6の目が80回出た。このサイコロで6の目が出る母比率の95％信頼区間を求めよ。
実際に計算してみます。
import scipy as sp
from scipy import stats
alpha = 0.95
n = 400
p = 80/400 #標本確率
sita, ue = stats.binom.interval(alpha=alpha, n=n, p=p, loc=0)
bottom = sita/400 #下側確率
up = ue/400 #上側確率print(f"{bottom} <= p <= {up}")
>0.1625 <= p <= 0.24 #答え
よって、 このサイコロで6の確率が出る信頼確率95%は、0.1625 <= p <= 0.24となる。
参考
母集団分布が未知の区間推定（Python）
母集団分布が未知の区間推定
これの続きです
母集団分布が未知の場合で、以下のパラメータがわかっているとき母平均の推定を行ってみる。
N = 100
標本平均¯X = 25000
不偏標準偏差U=5000
ここで中心極限定理より、
Nが十分に大きいとき、母集団の平均がμ、分散σ^2から取り出した標本平均¯Xは近似的に平均μ、分散σ^2/nの正規分布に従う。
つまり、nが十分に大きいとき¯Xは正規分布に従う。
stats.t.intervalを使って計算を行う。
n = 100 #N
mu = 25000 #標本平均
sigma = 5000 #不偏標準偏差
df = 99 #自由度（N-1）
alpha = 0.95 #信頼係数stats.t.interval(alpha = alpha, df = df, loc = mu, scale = sigma/sp.sqrt(n))
>(24007.891524245657, 25992.108475754343)
だいたい95%信頼区間は、24007.89<= μ <=25992.11であることがわかる。
中心極限定理を利用した場合の計算式は、
25000-1.96*5000/√100<=μ<=25000+-1.96*5000/√100
24020<=μ<=25980となる。
参考
Pythonで学ぶあたらしい統計学の教科書（馬場 真哉）｜翔泳社の本
区間推定（Python）
区間推定とは？
例えば母集団のμとσ^2がわかっていないとき（片方がわかっている場合もある）、標本から母幅をもたせて母集団を推定することである。
区間推定には、幅における信頼度として信頼区間という言葉が使われる。
以下、定理
平均μ、分散σ^2の正規分布に従う独立な確率変数X1、X2...Xnがあるとするとき、T=¯Xーμ/U/√nは自由度n-1のt分布に従う。
まずはPythonのライブラリを用いて簡易的に区間推定を行う。(有効桁数は適当です)
95%信頼区間として、μを求める。
df = [94,99,86,101] #N=4のサンプル
ライブラリは、stats.t.intervalを用いる。
引数(alpha = 信頼係数、df = 自由度、loc = 標本平均、scale = 標準誤差)
※標準誤差は、σ/√Nで計算できる。
stats.t.interval(alpha = 0.95, df = len(df)-1, loc = mu, scale = se)
>(84.365, 105.634)
この結果によると、95%信頼区間は(84.365 <= μ <=105.634)であることがわかる。
次にT分布を用いてμの区間推定を行う。
t値の計算式は、
T = ¯X(標本平均)-μ(母平均)/U(不偏標準偏差)/√n
※ここでの分母「U(不偏標準偏差)/√n」を標準誤差(standart error)という。
＞ T分布は P(-t0.975 <= T <= t0.975)=95%（T = ¯X-μ/U/√n）となる。
上記式を μについて解くと
＞¯X-t0.975*U/√n <= μ <= ¯X+t0.975*U/√n となる。
ここで、T分布の97.5%点を求める。
t_975 = stats.t.ppf(0.975, len(df)-1)
>3.182
¯X-t0.975*U/√n <= μ <= ¯X+t0.975*U/√nにそれぞれ代入。
＞mu - 3.182*se <= μ <= mu + 3.182*se
mu - 3.182*se
>84.367
mu + 3.182*se
>105.63
stats.t.intervalで計算した結果とほぼ一緒になることがわかる。
結論、Pythonライブラリ使っとけばいい。
参考
Pythonで学ぶあたらしい統計学の教科書（馬場 真哉）｜翔泳社の本
点推定（Python）
点推定とは？
母集団の平均μと、分散σ^2が未知のときに、標本の値から母集団をピンポイント推定することである。標本平均を推定量として使用する。
今回は適当な母集団から標本を抽出したことにする。(N=6)
※母集団は正規分布に従うとする。
データは、[9270, 33662, 18298, 36092, 17745, 7865]を標本抽出した。
まずは必要なライブラリを読み込む。
import pandas as pd
import numpy as np
import scipy as sp
from scipy import stats#グラフ作成
from matplotlib import pyplot as plt
import seaborn as sns
sns.set()%matplotlib inline
%precision 3
numpyを用いてリスト化。
df = np.array([9270, 33662, 18298, 36092, 17745, 7865])
実際に点推定を行う。母平均の推定値は標本平均で求まる。
mu = sp.mean(df)
>20488.66
母分散の点推定では不偏分散を用いる。（標本分散は点推定に適さない）
bunsan = sp.var(df, ddof=1)#不偏分散
>142908578.26
点推定の結果、μ=20488.66、σ^2=142908578.26と母分散を推定することができる。
↓以下は点推定に関係ありません。
続きを読む二項分布（Python）
二項分布のパラメーターは、
試行回数：N
一つの試行で一つのイベントが発生する確率：p
二項分布はscipy.statsに関数が用意されている。
確率質量変数：pmf（引数には成功確率・試行回数・成功確率を入れる（k, n, p））
■では例として、表が出る確率50%のコインを10回投げて、そのうちの3回が表である確率を計算する。
import matplotlib.pyplot as plt
%matplotlib inlineimport seaborn as sns
sns.set()n = 10
p = 0.5m = np.arange(0, int(n+1), 1)#0回からn+1回当たる確率に注目
pmf_binomial = sp.stats.binom.pmf(k = m, n = n, p = p)
plt.plot(m, pmf_binomial)
■次に12人でレースを行い、30レース行うと何回1位を取ることができるかを調べます。
試行回数はN=30、確率はp=1/12となる
import matplotlib.pyplot as plt
%matplotlib inlineimport seaborn as sns
sns.set()n = 30
p = 1/12m = np.arange(0, int(n+1), 1)#0回からn+1回当たる確率に注目
pmf_binomial = sp.stats.binom.pmf(k = m, n = n, p = p)
plt.plot(m, pmf_binomial)#グラフ出力
pmf_binomial.argmax() #30レースした結果、何回1位を取る確率が高いのか
>2
つまり、12人で30レース走ると2回一位を取る確率が高いらしい。
参考
