ブログ

当まとめでは統計検定$2$級の公式テキストの副教材に用いることができるように、統計学入門に関して取り扱います。当記事では「統計検定$2$級対応 統計学基礎」の$4.4.2$節「対応ある$2$標本の場合」の内容を元に対応のある$2$標本に関する検定について確認を行います。
統計検定$2$級のテキストとの対応がわかりやすいように、目次を「統計検定$2$級対応 統計学基礎」と対応させました。学びやすさの観点からあえて目次を対応させましたが、当まとめは「統計の森」オリジナルのコンテンツであり、統計検定の公式とは一切関係ないことにご注意ください。

・統計検定$2$級対応・統計学入門まとめ
https://www.hello-statisticians.com/stat_basic

改訂版　日本統計学会公式認定　統計検定２級対応　統計学基礎

2,420円(11/04 19:31時点)

Amazon

「対応ある$2$標本の検定」の概要

概要

上記はダイエット処方の効果を確かめるにあたって、処方前の体重$d_1$、処方後の体重$d_2$と体重の変化$\Delta$を表にしたものです。以下、実際の効果を確かめるにあたって、帰無仮説$H_0: \, \Delta=0$と対立仮説$H_1: \, \Delta<0$に関して仮説検定を行います。

このときに$d_1$と$d_2$に相関があり「$4.4.1$ 母平均の差の検定」と同様な手法を用いることはできないので、当記事で以下取り扱う「対応ある$2$標本の取り扱い」の流れを抑えておく必要があります。

必要な数学

$\sqrt{x}$や$x^2$の取り扱いなど、基本的な計算を抑えておけば十分です。

対応ある$2$標本の検定

検定の流れ

前節の$d_1, d_2$の行の値は用いず、$\Delta$の行のみを用います。観測値を$\delta_1, \cdots , \delta_{10}$とおき、下記のように$\bar{\delta}, \hat{\sigma}_{\delta}^{2}$を定義します。
$$
\large
\begin{align}
\bar{\delta} &= \frac{1}{10} \sum_{i=1}^{10} \delta_{i} \\
\hat{\sigma}_{\delta}^{2} &= \frac{1}{9} \sum_{i=1}^{10} (\delta_{i}-\bar{\delta})^2
\end{align}
$$

このとき、下記のように検定統計量の実現値$t$を計算します。
$$
\large
\begin{align}
t &= \frac{\bar{\delta}}{\hat{\sigma}_{\delta}/\sqrt{n}} \\
&= \frac{\sqrt{n} \bar{\delta}}{\hat{\sigma}_{\delta}}
\end{align}
$$

上記に対応する統計量が自由度$10-1=9$の$t$分布$t(9)$に従うことを元に仮説検定を行うことができます。

具体的な計算例

対応ある$2$標本の検定法を用いる

前項の式に基づいて下記を実行することで検定を行うことができます。

import numpy as np
from scipy import stats

delta = np.array([0.1, -4.2, -6.4, -3.3, -3.0, -4.2, -2.8, 0.9, 3.5, -4.8])
sigma2 = np.sum((delta-np.mean(delta))**2)/(delta.shape[0]-1.)

t = np.sqrt(delta.shape[0]) * np.mean(delta) / np.sqrt(sigma2)

if stats.t.cdf(t,10-1)<0.05:
    print("t: {:.3f}, reject H_0.".format(t))
else:
    print("t: {:.3f}, accept H_0.".format(t))

・実行結果

t: -2.543, reject H_0.

通常の$2$標本の検定を用いる

「$4.4.1$ 母平均の差の検定」の内容に基づいて下記を実行することで検定を行うことができます。

import numpy as np
from scipy import stats

d = np.array([[65.0, 79.3, 55.9, 73.2, 58.0, 68.5, 68.1, 69.9, 71.7, 58.3], [65.1, 75.1, 49.5, 69.9, 55.0, 64.3, 65.3, 70.8, 75.2, 53.5]])
mean_d = np.mean(d,axis=1)
mean_d_mat = np.repeat(mean_d,10).reshape([2,10])

sum_sigma2 = np.sum((d-mean_d_mat)**2,axis=1)
sigma2 = np.sum(sum_sigma2)/(d.shape[1]+d.shape[1]-2.)

t = (mean_d[1]-mean_d[0]) / (np.sqrt(1./d.shape[1] + 1./d.shape[1]) * np.sqrt(sigma2))

if stats.t.cdf(t,10+10-2)<0.05:
    print("t: {:.3f}, reject H_0.".format(t))
else:
    print("t: {:.3f}, accept H_0.".format(t))

・実行結果

t: -0.652, accept H_0.

計算結果

「対応ある$2$標本」の検定法を用いた場合結果が有意であり、「通常の$2$標本の検定」を用いるとき有意でない結果が得られました。この解釈に関しては次項で取り扱います。

対応ある2標本の取り扱いの必要性の解釈

標本平均の差$\overline{X}-\overline{Y}$の分散$V[\overline{X}-\overline{Y}]$は$\overline{X}$と$\overline{Y}$に相関がある場合、下記のように表されます。
$$
\large
\begin{align}
V[\overline{X}-\overline{Y}] = V[\overline{X}] + V[\overline{Y}] – 2 \mathrm{Cov}(\overline{X},\overline{Y})
\end{align}
$$

「$4.4.1$ 母平均の差の検定」では$\mathrm{Cov}(\overline{X},\overline{Y})=0$を前提としているので、$V[\overline{X}-\overline{Y}] = V[\overline{X}] + V[\overline{Y}]$のように表されたことに注意が必要です。

このように「対応ある$2$標本」の場合は、$2$つのグループ間に正の相関があると$\overline{X}-\overline{Y}$の分散$V[\overline{X}-\overline{Y}]$が小さくなることに注意が必要です。冒頭の例に対し、$d_1$と$d_2$の散布図を描くと下記が得られます。

import numpy as np
import matplotlib.pyplot as plt

d = np.array([[65.0, 79.3, 55.9, 73.2, 58.0, 68.5, 68.1, 69.9, 71.7, 58.3], [65.1, 75.1, 49.5, 69.9, 55.0, 64.3, 65.3, 70.8, 75.2, 53.5]])

plt.scatter(d[0], d[1])
plt.show()

・実行結果

上記のように$2$つの標本には相関があることが確認できます。相関係数は下記のように計算できます。

import numpy as np

d = np.array([[65.0, 79.3, 55.9, 73.2, 58.0, 68.5, 68.1, 69.9, 71.7, 58.3], [65.1, 75.1, 49.5, 69.9, 55.0, 64.3, 65.3, 70.8, 75.2, 53.5]])
mean_d = np.mean(d,axis=1)
mean_d_mat = np.repeat(mean_d,10).reshape([2,10])

cov = np.mean((d[0]-np.mean(d[0]))*(d[1]-np.mean(d[0])))
sigma2 = np.mean((d-mean_d_mat)**2,axis=1)
r = cov / np.sqrt(sigma2[0]*sigma2[1]) 

print("r: {:.3f}".format(r))

・実行結果

r: 0.950

ここでの相関係数の値が$0.950$が得られるなど$2$標本間には強相関があるので、前項のような結果が得られたと解釈しておくと良いです。

発展事項①

「通常の$2$標本の検定を用いる」の計算で用いる分散から$2 \mathrm{Cov}(\overline{X},\overline{Y})$分を引いた際の検定は下記のように計算できます。

import numpy as np
from scipy import stats

d = np.array([[65.0, 79.3, 55.9, 73.2, 58.0, 68.5, 68.1, 69.9, 71.7, 58.3], [65.1, 75.1, 49.5, 69.9, 55.0, 64.3, 65.3, 70.8, 75.2, 53.5]])
mean_d = np.mean(d,axis=1)
mean_d_mat = np.repeat(mean_d,10).reshape([2,10])

sum_sigma2 = np.sum((d-mean_d_mat)**2,axis=1)
sum_cov = np.sum((d[0]-np.mean(d[0]))*(d[1]-np.mean(d[0])))
sigma2 = (np.sum(sum_sigma2)-2*sum_cov)/(d.shape[1]+d.shape[1]-2.)

t = (mean_d[1]-mean_d[0]) / (np.sqrt(1./d.shape[1] + 1./d.shape[1]) * np.sqrt(sigma2))

if stats.t.cdf(t,10+10-2)<0.05:
    print("t: {:.3f}, reject H_0.".format(t))
else:
    print("t: {:.3f}, accept H_0.".format(t))

・実行結果

t: -2.543, reject H_0.

この結果は「対応ある$2$標本の検定法を用いる」の結果と一致します。

発展事項②

$$
\large
\begin{align}
V[\overline{X}-\overline{Y}] = V[\overline{X}] + V[\overline{Y}] – 2 \mathrm{Cov}(\overline{X},\overline{Y})
\end{align}
$$

上記の式の導出は下記などで取り扱いました。

当まとめでは統計検定$2$級の公式テキストの副教材に用いることができるように、統計学入門に関して取り扱います。当記事では「統計検定$2$級対応 統計学基礎」の$4.4.1$節「母平均の差の検定」の内容を元に母分散既知・未知の場合の母平均の差の検定の方法について確認を行います。
統計検定$2$級のテキストとの対応がわかりやすいように、目次を「統計検定$2$級対応 統計学基礎」と対応させました。学びやすさの観点からあえて目次を対応させましたが、当まとめは「統計の森」オリジナルのコンテンツであり、統計検定の公式とは一切関係ないことにご注意ください。

・統計検定$2$級対応・統計学入門まとめ
https://www.hello-statisticians.com/stat_basic

改訂版　日本統計学会公式認定　統計検定２級対応　統計学基礎

2,420円(11/04 19:31時点)

Amazon

「母平均の差の検定」の概要

概要

上記のように$L$、$H$、$C$の$3$つのグループに関して観測値が得られたとき、$2$グループを選んでそれぞれのグループの母平均に差があるかを調べる場合を考えます。

このときの一連の手順が「母平均の差の検定」であり、次節で「母分散既知の場合」、「母分散未知かつ等しい場合」、「母分散未知かつ等しくない場合」の$3$つに場合分けし、それぞれについて検定の手順を確認します。

必要な数学

$\sqrt{x}$や$x^2$の取り扱いなど、基本的な計算を抑えておけば十分です。

母平均の検定

母分散既知の場合

$L$からの標本の実現値$x_1, \cdots , x_m$が$\mathcal{N}(\mu_1,\sigma^2)$に基づいて生成され、$H$からの標本の実現値$y_1, \cdots , y_n$が$\mathcal{N}(\mu_2,\sigma^2)$に基づいて生成されると仮定します。

このとき、標本平均$\overline{X}, \overline{Y}$を下記のように定めます。
$$
\large
\begin{align}
\overline{X} &= \frac{1}{m} \sum_{i=1}^{m} X_i \\
\overline{Y} &= \frac{1}{n} \sum_{j=1}^{n} Y_j \\
\overline{X} & \sim \mathcal{N} \left( \mu_1,\frac{\sigma^2}{m} \right), \quad \overline{Y} \sim \mathcal{N} \left( \mu_2,\frac{\sigma^2}{n} \right)
\end{align}
$$

ここで標本平均の差の$\overline{X}-\overline{Y}$を考えるとき、正規分布のモーメント母関数などを用いることで下記を示すことができます。
$$
\large
\begin{align}
\overline{X}-\overline{Y} \sim \mathcal{N} \left( \mu_1-\mu_2, \frac{\sigma^2}{m}+\frac{\sigma^2}{n} \right)
\end{align}
$$

上記より、$\displaystyle \overline{X}-\overline{Y} \sim \mathcal{N} \left( \mu_1-\mu_2, \left(\frac{1}{m}+\frac{1}{n}\right)\sigma^2 \right)$が成立します。よって、検定統計量$Z$を下記のように定義することができます。
$$
\large
\begin{align}
Z &= \frac{(\overline{X}-\overline{Y})-(\mu_1-\mu_2)}{\sqrt{\frac{1}{m}+\frac{1}{n}}\sigma} \\
Z & \sim \mathcal{N}(0,1)
\end{align}
$$

上記に対し$\mu_1-\mu_2$に対して帰無仮説$H_0: \, \mu_1-\mu_2=0$などを仮定し、$z_{\alpha=0.025}=1.96$などを用いることで両側$5$％検定を行うことができます。

発展事項

$\displaystyle \overline{X}-\overline{Y} \sim \mathcal{N} \left( \mu_1-\mu_2, \left(\frac{1}{m}+\frac{1}{n}\right)\sigma^2 \right)$が成立することは正規分布のモーメント母関数を用いることで示すことができます。詳しくは下記で取り扱いました。

上記の導出にはモーメント母関数の理解が必須なので、統計検定$2$級範囲では下記のように$E[\overline{X}-\overline{Y}], V[\overline{X}-\overline{Y}]$の計算で抑えておくでも十分です。
$$
\large
\begin{align}
E[\overline{X}-\overline{Y}] &= E[\overline{X}] – E[\overline{Y}] \\
V[\overline{X}-\overline{Y}] &= V[\overline{X}] + V[-\overline{Y}] \\
&= V[\overline{X}] + (-1)^2V[\overline{Y}] = V[\overline{X}] + V[\overline{Y}]
\end{align}
$$

母分散未知かつ等しい場合

母分散未知かつ等しくない場合

当まとめでは統計検定$2$級の公式テキストの副教材に用いることができるように、統計学入門に関して取り扱います。当記事では「統計検定$2$級対応 統計学基礎」の$1.2.2$節「幹葉図」の内容を元にヒストグラムと同様の目的で用いられる幹葉図に関して取り扱いました。
統計検定$2$級のテキストとの対応がわかりやすいように、目次を「統計検定$2$級対応 統計学基礎」と対応させました。学びやすさの観点からあえて目次を対応させましたが、当まとめは「統計の森」オリジナルのコンテンツであり、統計検定の公式とは一切関係ないことにご注意ください。

・統計検定$2$級対応・統計学入門まとめ
https://www.hello-statisticians.com/stat_basic

改訂版　日本統計学会公式認定　統計検定２級対応　統計学基礎

2,420円(11/04 19:31時点)

Amazon

幹葉図の概要

概要

幹葉図(stem and leaf plots)はヒストグラムと同様に、度数分布の可視化にあたってよく用いられます。幹葉図は下記のように表されます。

stem: 4, leaves: 5599
stem: 5, leaves: 22599
stem: 6, leaves: 0244445557777
stem: 7, leaves: 022224445555577999
stem: 8, leaves: 01359

上記は「人文・社会科学の統計学の演習$2.1$」の結果を用いました。図を元に$6$〜$7$の間のサンプルが多いことが確認できます。このように幹葉図はヒストグラムと同様な用途で用いることができます。

必要な数学

幹葉図の作成にあたって必要な数学知識は特にありませんが、プログラムを用いて作成する場合は割り算の商と余りの演算を重点的に抑えておくと良いと思います。

幹葉図

幹葉図のプログラム

幹葉図は下記のようなプログラムを実行することで作成できます。

import numpy as np

x = np.array([5.6, 8.8, 7.5, 6.2, 9.2, 7.2, 7.1, 7.7, 8.2, 7.5, 6.8, 8.5, 6.8, 8.8, 7.6])

stem = np.unique(x//1)
leaf = x%1

leaves = dict()

for i in range(stem.shape[0]):
    leaves[int(stem[i])] = list()

for i in range(x.shape[0]):
    leaves[int(x[i]//1)].append(int(np.floor(((x[i]+0.01)%1)*10)))

for i in range(stem.shape[0]):
    leaves[int(stem[i])] = np.sort(leaves[int(stem[i])])
    l = ""
    for j in range(len(leaves[int(stem[i])])):
        l += str(leaves[int(stem[i])][j])
    print("stem: {}, leaves: {}".format(int(stem[i]),l))

・実行結果

stem: 5, leaves: 6
stem: 6, leaves: 288
stem: 7, leaves: 125567
stem: 8, leaves: 2588
stem: 9, leaves: 2

9.2が9.1999999などのように表されたので、$14$行目でx[i]+0.01を計算しました。厳密な取り扱いではないので、プログラム自体は修正を行う可能性があります。

当まとめでは統計検定$2$級の公式テキストの副教材に用いることができるように、統計学入門に関して取り扱います。当記事では「統計検定$2$級対応 統計学基礎」の$1.2.2$節「ローレンツ曲線」の内容を元に経済などの分野でよく用いられるローレンツ曲線とジニ係数の概要に関して取り扱いました。
統計検定$2$級のテキストとの対応がわかりやすいように、目次を「統計検定$2$級対応 統計学基礎」と対応させました。学びやすさの観点からあえて目次を対応させましたが、当まとめは「統計の森」オリジナルのコンテンツであり、統計検定の公式とは一切関係ないことにご注意ください。

・統計検定$2$級対応・統計学入門まとめ
https://www.hello-statisticians.com/stat_basic

改訂版　日本統計学会公式認定　統計検定２級対応　統計学基礎

2,420円(11/04 19:31時点)

Amazon

ローレンツ曲線・ジニ係数の概要

概要

ローレンツ曲線(Lorenz curve)はM.O.ローレンツが考案した「所得などの量の集中度や格差を表すグラフ」です。ローレンツ曲線は下記のようなグラフで表されます。

上図は日本、アメリカ、スウェーデン、中国、ドイツの$5$か国の五分位所得割合に基づくローレンツ曲線を表したものです。ローレンツ曲線は上にいくほど「平等」で下にいくほど「不平等」であることを表します。図ではスウェーデン、ドイツ、日本、アメリカ、中国の順で平等であることが読み取れます。

ジニ係数(Gini’s coefficient)はグラフの$(0,0)$と$(5,100)$を通る直線とそれぞれのローレンツ曲線の間の面積の$2$倍で定義される指標であり、ジニ係数が大きいほど「不平等」であると考えられます。なお、グラフは$2011$〜$2014$年を元に作成されているので、その点には注意が必要です。

必要な数学

ジニ係数は「完全平等曲線」と「ローレンツ曲線」の間の面積の$2$倍で定義されるので、積分の定義や台形公式に基づく数値積分の基本的な計算法に関しては抑えておくと良いです。

ローレンツ曲線とジニ係数

ローレンツ曲線

ローレンツ曲線の描画にあたっては、まず所得の順に家計を並べ、世帯数を元に五等分します。たとえば$1$億世帯ある場合は所得の順に$1$億世帯をソートし、$2,000$万世帯ずつで五等分を行います。このときにそれぞれの$2,000$万世帯ずつが全体の所得に占める割合を計算します。計算結果は度数分布表と同様の要領で、表で表すことが可能です。前節で確認したローレンツ曲線に対応する表を作成すると、下記のように表せます。

上記の表を確認するにあたっては、どの行も「左から右にいくにつれて数字が必ず大きくなる」という点は抑えておく必要があります。このことは元々所得の順にソートを行なった上でグループ分けし、割合を計算したことを鑑みれば必然です。

次に「確率の累積」を左から順に計算します。左から「第$1$階級の数字」、「第$1$階級の数字＋第$2$階級の数字」〜「第$1$階級の数字＋第$2$階級の数字＋$\cdots$＋第$5$階級の数字」の数字を並べます。五カ国の例に基づいて計算すると下記のような結果が得られます。

上記を元にローレンツ曲線の描画を行うことができます。累積を計算するので一番右の列は$100$に近い値になりますが、有効数字や四捨五入を行う都合上、必ずしも$100$に一致しないことにも注意が必要です。一連の流れは下記を実行することでグラフを描くことができます。

import numpy as np
import matplotlib.pyplot as plt

x = np.array([[5.4, 10.7, 16.3, 24.1, 43.5], [5.1, 10.3, 15.4, 22.7, 46.4], [8.7, 14.3, 17.8, 23.0, 36.2], [5.2, 9.8, 14.9, 22.3, 47.9], [8.4, 13.1, 17.2, 22.7, 38.6]])

x_cum = np.zeros([x.shape[0],x.shape[0]+1])
x_cum[:,1:] = np.cumsum(x,axis=1)

nation_label = ["JPN", "USA", "SWE", "CHN", "DEU"]
for i in range(x_cum.shape[0]):
    plt.plot(np.arange(6),x_cum[i,:],label=nation_label[i])

plt.legend(loc="upper left")
plt.show()

・実行結果

ジニ係数

ジニ係数は全世帯が平等である場合の完全平等線(complete equality line)とローレンツ曲線の間の面積の$2$倍に対応します。下記を実行することで前項の例に完全平等線を追加することができます。

import numpy as np
import matplotlib.pyplot as plt

x = np.array([[5.4, 10.7, 16.3, 24.1, 43.5], [5.1, 10.3, 15.4, 22.7, 46.4], [8.7, 14.3, 17.8, 23.0, 36.2], [5.2, 9.8, 14.9, 22.3, 47.9], [8.4, 13.1, 17.2, 22.7, 38.6]])

x_cum = np.zeros([x.shape[0],x.shape[0]+1])
x_cum[:,1:] = np.cumsum(x,axis=1)

nation_label = ["JPN", "USA", "SWE", "CHN", "DEU"]
for i in range(x_cum.shape[0]):
    plt.plot(np.arange(6)/5.,x_cum[i,:]/100.,label=nation_label[i])

plt.plot(np.array([0.,1.]), np.array([0., 1.]), color="black",label="equality line")

plt.legend(loc="upper left")
plt.show()

・実行結果

上図のようなローレンツ曲線を元に「不平等」かどうかを論じることができますが、前項の表では第$1$階級では相対確率がアメリカ＜中国である一方で、それ以外ではアメリカ＞中国が成立するなど、取り扱いが難しい場合があります。

このような場合に面積に基づいて定量化を行うジニ係数のような指標があると用いやすいです。各国のジニ係数は台形公式に基づいて下記のような処理を行うことで計算することができます。

import numpy as np

x = np.array([[5.4, 10.7, 16.3, 24.1, 43.5], [5.1, 10.3, 15.4, 22.7, 46.4], [8.7, 14.3, 17.8, 23.0, 36.2], [5.2, 9.8, 14.9, 22.3, 47.9], [8.4, 13.1, 17.2, 22.7, 38.6]])

x_cum = np.zeros([x.shape[0],x.shape[0]+1])
x_cum[:,1:] = np.cumsum(x,axis=1)/100.

s = np.repeat(1.*1./2.,5)
for i in range(x.shape[1]):
    s -= 0.2*(x_cum[:,i]+x_cum[:,i+1])/2.

Gini = 2*s

print(Gini)

・実行結果

[ 0.3584  0.381   0.2548  0.3906  0.28  ]

上記より、ジニ係数を大きい順に並べると、「中国」、「アメリカ」、「日本」、「ドイツ」、「スウェーデン」であることが確認できます。