検定を行う際によく出てくるのが「対応がある」「対応がない」という用語です。
自分も初めは何のことかわからなかったので、まず簡単に説明します。
まず「対応がある」というのは、同じ集団に違った処置を施し、その結果を検定する場合などに使われます。
例えば、ある新薬の効果を確かめる場合、同じ患者の集団に新薬Aと新薬Bを投与した結果に違いを検定する場合です。
逆に「対応がない」場合は、まったく別の集団の特性を検定する場合。
例えばA校の数学のテスト成績とB校のテスト成績の統計量を検定する場合などですね。
というわけで、さっそく対応のない２集団の統計量を検定してみましょう。
対応のない2集団の統計量検定
例によって問題は演習で身につける統計学入門の123ページの例題を使っていきましょう。
２種類のチーズAとBについて各種メーカーの製品をそれぞれ36個と34個集めその食塩濃度を測った時標本軽金と不偏標本分散はA：2.04と0.126。Bで1.75と0.113。
AとBの食塩濃度の平均に有意差はあるか有意水準5%で検定。
まず仮説を以下のように立てます。
- 帰無仮説 チーズAの平均食塩濃度 ＝ チーズBの平均食塩濃度
- 対立仮説 チーズAの平均食塩濃度 ≠ チーズBの平均食塩濃度
さっそく計算していきましょう。
データは、書籍のデータリンクからダウンロードしたものを読み込みます。
# 必要なライブラリーのインポート
import pandas as pd
import numpy as np
# エクセルファイルのインポート（書籍のデータリンクからダウンロードしたデータ）
c_salt = pd.read_excel('Ex6-1_z_test.xlsx', header = 0)
c_salt.head() # Aの最初の5行を表示
実行すると、以下のようなDataFrameができます。
次に読み込んだデータに不備がないかを確認していきます。
まずデータの数を確認してみます。
len(c_salt) # len()関数でデータの長さを確認する。
36 # データが36行あることがわかる。
データ数はAとBで36個あるようです。
データは36あるようですが、データがしっかりそろっているのか分かりませんので、欠損値を確認します。
c_salt.isna().any() # 欠損値があるかを確認する。
No. False # 欠損値無し
A False # 欠損値無し
B True # 欠損値あり
dtype: bool
結果を見ると、どうやらBに欠損値があるようですので、欠損値の数を確認してみましょう。
c_salt.isna().sum() # 欠損値の個数を確認する。
No. 0
A 0
B 2 # 欠損値2個
dtype: int64
Bに欠損値が2個あることがわかります。
この欠損値がある行をそのまま削除することもできますが、値を0で埋めておくことにします。
c_salt=c_salt.fillna(0) # 欠損値を0で埋める。
c_salt.isna().sum() # 再び欠損値の個数を確認する
No. 0
A 0
B 0
# 欠損値無し
ここから分析しやすいように、DataFrame の整形を行います。
参考サイト：
DataFrameの整形
ここから分析しやすいようにデータを整形していきます。
まず列Aの値を抜き出します。
c_saltA = c_salt.loc[:, "A"] # Aの列を抜き出す
c_saltA.head() # Aの最初の5行を表示
0 2.7
1 2.0
2 1.8
3 1.5
4 2.0
Name: A, dtype: float64
Bも同じように抜き出します。
c_saltB = c_salt.loc[:, "B"] # Bの列を抜き出す
c_saltB.head() # Bの最初の5行を表示
0 2.4
1 1.3
2 2.2
3 1.7
4 1.9
Name: B, dtype: float64
参考サイト：
ここで実際の検定を行うわけですが、検定方法にもいくつかの種類があるわけですね。
とりあえず以下の三つ。
- Ｚ検定（母分散既知）
- t検定（母分散未）
- ウェルチの検定（母分散に違いがある）
カイ二乗検定は質的データ（カテゴリカルデータ）に対しての検定なので、別の記事で解説。
今回のケースの場合、不偏標本分散は分かっていますが母分散は分からないので、Ｚ検定以外の検定を行います。
t検定かウェルチ検定かを決める場合、母分散に違いがあるかを確認する必要がありますで、母分散に差があるのかどうか確認するために、F検定を行います。
F検定
サンプルＡとＢの不偏標本分散は分かっているので、それをもとにＦ検定を行います。
F値を計算するモジュールはないのでプログラムを書くわけですが、ちょっとまだ書けないので以下の記事を参考にしました。
この場合の帰無仮説は「AとBの分散に差はない」であり、対立仮説は「AとBの分散には差がある」になります。
#有意水準を0.05、帰無仮説を「２群間の分散に差がないこと」とする。 A_var = 0.126 B_var = 0.113 A_df = len(c_saltA_n)-1 # A の自由度 B_df = len(c_saltB_n)-1 # B の自由度 #Fを求める場合は、数値が大きい方を分子にする。 bigger_var = np.max([A_var, B_var]) smaller_var = np.min([A_var, B_var]) f = bigger_var/smaller_var one_sided_pval1 = stats.f.cdf(f, A_df, B_df) one_sided_pval2 = stats.f.sf(f, A_df, B_df) two_sided_pval = min(one_sided_pval1, one_sided_pval2)*2 print('F: ', round(f, 3)) print('p-value: ', round(two_sided_pval, 3))
F: 1.115
p-value: 0.749
p>0.05なので、帰無仮説の「AとBの分散に差はない」が採用されます。
分散に差がないことがわかったので、さっそくt検定を行ってみましょう。
ここはt検定のモジュールがあるので、シンプルにそれを使って計算します。
statsmodelsを使ったstudentのt検定
分散に差がないことが分かったので、モジュールstatsmodelsでstudent のt検定を行います。
alternative で両側「two-sided」片側「smaller またはlarger 」を指定できます。
usevar では「pooled」でsutudentのt検定（分散に差がない）「unequal」でウェルチのt検定（分散に差がある）を指定できる。
# 必要なモジュールのインポート
from statsmodels.stats.weightstats import ttest_ind
# 対応無し両側t検定を行う場合
ttest_ind(c_saltA, c_saltB, alternative='two-sided', usevar='pooled')
(3.700279473352106, 0.0004252714975510123, 70.0)
検定の結果は t値が3.700、p値が 0.000425になります。
最後の70は自由度ということなのですが、なぜ70なのかよくわかりません。
この場合、p<0.05なので「AとBの塩分濃度に差がないとは言えない」という帰無仮説が採用されます。
参考サイト：
Pythonのt検定ではstatsmodelsを使う以外に、scipyのstatsを使っても求めることができます。
実際にやってみます。
import scipy.stats as stats
# 対応無しstudentのt検定を行う場合
stats.ttest_ind(c_saltTA_n, c_saltTB_n)
Ttest_indResult(statistic=3.700279473352104, pvalue=0.00042527149755101505)
ということで同じ結果を得ることができました。
statsmodelsを使ったウェルチの検定
次に2群間で分散に差がある場合の検定「ウェルチの検定」をやってみます。
例によって問題は演習で身につける統計学入門の130p、の例題12から。
ある男子高校で3年A組とB組の生徒からそれぞれ10人と12人を任意に選び、慎重を測定。A組の平均身長はB組の平均身長よりも低いか有意水準5% で検定。
まずはデータをダウンロード。
hight = pd.read_excel('Ex6-5 t test.xlsx', header = 3)
hight.head()
データの確認を行います。
len(hight['A']) # コラムAのデータの長さを確認。
12
len(hight['B']) # コラムBのデータの長さを確認。
12
hight.isna().sum() # 欠損値があるか確認。
No. 0
A 2
B 0
dtype: int64
# コラムAに欠損値が2個ある。
hight_comp=hight.fillna(hight.mean()) # 欠損値をコラムAの平均で埋める。
hight_comp.isna().sum()
No. 0
A 0
B 0
dtype: int64
次にDataFrameから、列ごとのデータを抜き出します。
まずはコラムAのデータ。
# コラムAのデータを抜き出す。
hight_compA = hight_comp.loc[:, "A"]
hight_compA.head()
0 167.0
1 171.0
2 160.0
3 172.0
4 165.0
Name: A, dtype: float64
次にコラムBのデータを抜き出します。
# コラムBのデータを抜き出す。
hight_compB = hight_comp.loc[:, "B"]
hight_compB.head()
0 177
1 175
2 165
3 190
4 171
Name: B, dtype: int64
F検定による分散の差の確認
では、F値を確認して2群間の分散に差があるかどうかを見てみましょう。
#有意水準を0.05、帰無仮説を「2群間の分散に差がないこと」とする。
A_var = np.var(hight_compA, ddof=1) # Aの分散を代入
B_var = np.var(hight_compB, ddof=1) # Bの分散を代入
A_df = len(hight_compA)-1 # A の自由度
B_df = len(hight_compB)-1 # B の自由度
#Fを求める場合は、数値が大きい方を分子にする。
bigger_var = np.max([A_var, B_var])
smaller_var = np.min([A_var, B_var])
f = bigger_var/smaller_var
one_sided_pval1 = stats.f.cdf(f, A_df, B_df)
one_sided_pval2 = stats.f.sf(f, A_df, B_df)
two_sided_pval = min(one_sided_pval1, one_sided_pval2)*2
print('F: ', round(f, 3))
print('p-value: ', round(two_sided_pval, 3))
F: 5.704
p-value: 0.008
結果、F値5.704、p値 0.008。
つまりp<0.05なので、帰無仮説は棄却され「2群間の分散に差がないとは言えない」が採用されます。
次にウェルチのt検定で、平均身長に差があるかを見てみましょう。
ウェルチの片側検定
大小関係を検定するので、片側検定を行うことにします。
参考サイト：
この場合使用するモジュールは、片側検定対応のstatsmodelsになり、片側右辺を使って検定を行います。
つまり「alternative=’smaller’」さらに、分散に差があるので「usevar=’unequal’」になります。
# 必要なモジュールのインポート
from statsmodels.stats.weightstats import ttest_ind
ttest_ind(hight_compA, hight_compB, alternative='smaller', usevar='unequal')
(-1.0042201353885785, 0.16573595640046906, 14.74210201305821)
結果はt値が-1.004、p値が0.165になります。
つまりp値>0.05なので、帰無仮説「A組の平均身長はB組の平均身長よりも低い」が棄却され、
「A組の平均身長はB組の平均身長よりも低いとは言えない」が採用されます。
コメント