のまとめやコードを実行してみた際のメモです．2群のデータに対するt検定
2つの変数の間で，平均値に差があるかどうかを判断してみます．以下の例では，薬を飲む前後の体温変化に関するデータを用いて，同じ対象を，異なった条件で2回測定して，その違いをみるという場合に対応のあるt検定を使っています．
読み込むデータ（.csvデータ）の中身は以下のようなものです．
person,medicine,body_temperature
A,before,36.2
B,before,36.2
C,before,35.3
D,before,36.1
E,before,36.1
A,after,36.8
B,after,36.1
C,after,36.8
D,after,37.1
E,after,36.9
薬を飲む前後における体温の差を計算します．まず，薬を飲む前後の標本平均を求めます．
この結果より，p値が0.05を下回っているので，帰無仮説は棄却され，薬を飲む前後で体温は有意に異なると主張できます．
ここまでは，対応のあるt検定でしたが，引き続いて対応のないt検定を行なってみます．対応のないt検定では，平均値の差に注目します．上記の対応のあるt検定では，データの差をとってから1郡のt検定を行なっていました．
平均値の差に基づいてt値を計算する場合は，t値の計算式が若干複雑になります．1軍のt検定におけるt値の計算式は以下のようになります．
読み込むデータ（.csvデータ）の中身は以下のようなものです．
person,medicine,body_temperature
A,before,36.2
B,before,36.2
C,before,35.3
D,before,36.1
E,before,36.1
A,after,36.8
B,after,36.1
C,after,36.8
D,after,37.1
E,after,36.9
このデータを対象として，対応のあるt検定を，以下の要領で行います．
帰無仮説：薬を飲む前と後で体温は変わらない．
対立仮説：薬を飲む前と後の体温が異なっている．
有意水準：5%（p値が0.05を下回れば，帰無仮説は棄却される）
t検定の実装は以下のようになります．
まずは，ライブラリの読み込みです．
>>> import numpy as np
>>> import pandas as pd
>>> import scipy as sp
>>> from scipy import stats
>>>
>>> from matplotlib import pyplot as plt
>>> import seaborn as sns
>>> sns.set()
データをデータフレーム形式で読み込みます．なお，サンプルサイズは10になります．
>>> paired_test_data = pd.read_csv("3-9-1-paired-t-test.csv")
>>> print(paired_test_data)
person medicine body_temperature
0 A before 36.2
1 B before 36.2
2 C before 35.3
3 D before 36.1
4 E before 36.1
5 A after 36.8
6 B after 36.1
7 C after 36.8
8 D after 37.1
9 E after 36.9
>>> before = paired_test_data.query('medicine == "before"')["body_temperature"]
>>> after = paired_test_data.query('medicine == "after"')["body_temperature"]
この時，シリーズ型のままだと計算がしにくいので，薬を飲む前後で抽出した後に，アレイ型に変換します．
>>> before = np.array(before)
>>> after = np.array(after)
差分を計算します．
>>> diff = after - before
>>> diff
array([ 0.6, -0.1, 1.5, 1. , 0.8])
この差分の平均値が0と異なるかどうかを，1軍のt検定で調べてみます．stats.ttest_rel関数を使うと検定が簡単です．
>>> stats.ttest_1samp(diff, 0)
Ttest_1sampResult(statistic=2.901693483620596, pvalue=0.044043109730074276)
>>> stats.ttest_rel(after, before)
Ttest_relResult(statistic=2.901693483620596, pvalue=0.044043109730074276)
ここまでは，対応のあるt検定でしたが，引き続いて対応のないt検定を行なってみます．対応のないt検定では，平均値の差に注目します．上記の対応のあるt検定では，データの差をとってから1郡のt検定を行なっていました．
平均値の差に基づいてt値を計算する場合は，t値の計算式が若干複雑になります．1軍のt検定におけるt値の計算式は以下のようになります．
$t値 = \cfrac{標本平均 - 比較対象値}{標準偏差 \div \sqrt{サンプルサイズ}} = \cfrac{標本平均 - 比較対象値}{標準誤差}$
変数$x$と$y$の平均値の差を検定するとします．対応の無いt検定のt値は以下のように計算されます．
$t = \cfrac{{\hat{\mu_{x}}} - {\hat{\mu_{y}}}}{\sqrt{{\hat{\sigma_{x}^{2}}}} / m + {\hat{\sigma_{y}^{2}}}}$
ここに，$\hat{\mu_{x}}$：$x$の標本平均，$\hat{\mu_{y}$：$y$の標本平均，$m$：$x$のサンプルサイズ，$n$：$y$のサンプルサイズ，${\hat{\sigma_{x}^{2}}}$：$x$の不偏分散，{\hat{\sigma_{y}^{2}}}：$y$の不偏分散です．
t値の計算は以下のようになります．
まず，平均値を求めます．
>>> mean_bef = sp.mean(before)
>>> mean_aft = sp.mean(after)
ここで，Scipyのmean関数を使っていますが，Scipyのmean関数はScipy2.0.0で廃止されるらしいので，Numpyのmean関数でも構いません．
分散を求めます．
>>> sigma_bef = sp.var(before, ddof = 1)
>>> sigma_aft = sp.var(after, ddof = 1)
ここで，Scipyのvar関数を使っていますが，Scipyのvar関数はいずれなくなるらしいので，Numpyのmean関数でも構いません．
サンプルサイズを求めます．
>>> m = len(before)
>>> n = len(after)
t値を計算します．
>>> t_value = (mean_aft - mean_bef) / \
... sp.sqrt((sigma_bef / m + sigma_aft / n))
>>> t_value
3.1557282344421034
ここで，Scipyのsqrt関数を使っていますが，Scipyのsqrt関数はScipy2.0.0で廃止されるらしいので，Numpyのsqrt関数でも構いません．
対応の無いt検定は，stats.ttest_ind関数を使えば計算できます．
>>> stats.ttest_ind(after, before, equal_var = False)
Ttest_indResult(statistic=3.1557282344421034, pvalue=0.013484775682079892)
p値が0.05を下回ったので，有意差があると判断できる結果となりました．ただし，先に行った対応のあるt検定で得られたp値（0.04...）とは異なっています．このことから，同じデータに対して同じ目的の検定を行なっても，検定手法が変わるとp値も変わるということが確認できます．
0 件のコメント :
コメントを投稿