【統計】Pythonでt検定をやってみる

Python

前回はt検定について、中身と流れをまとめました。

https://www.learning-nao.com/?p=2589

今回は、そのt検定をPythonでやってみようと思います。

目次

使用するデータ

今回は定番のアイリスデータセットを使っていきます。

import pandas as pd
from sklearn.datasets import load_iris
iris=load_iris()
df = pd.DataFrame(iris.data, columns=iris.feature_names)

さらに種類を追加します。iris.targetには種類が格納されており、0(=setosa), 1(=versicolor), 2(=virginica)の3値が入っています。

df['species']=iris.target
df
df

今回は2標本のt検定を実施するので、上のdfから2つの標本を作成します。species=0の標本と、species=1の標本を作成し、「petal width（cm）」で検定をします。

この両標本は、「対応のないデータ」として扱います。

petal_len0 = df_iris[df_iris["species"]==0]["petal length (cm)"]
petal_len1 = df_iris[df_iris["species"]==1]["petal length (cm)"]

検定の前に

さて、早速t検定を…と行きたいところですが、その前にデータの特徴を確認しておきます。t検定には前提条件があるのでしたね。その前提とは

  • 母集団が正規分布に従う
  • 母集団の分散が等しい

です。この前提を満たし対応のない標本同士のt検定は「Studentのt検定」と呼ばれます。まずは、標本がこの前提を満たすかを確認します）。と言っても母分散はわからないので、標本の分散（標準偏差）と正規性あたりを確認します。

平均と標準偏差

各標本の平均と標準偏差を確認します。

#平均
print(petal_len0.mean())
print(petal_len1.mean())
>>1.4620000000000002
  4.26 

#標準偏差
print(petal_len0.std())
print(petal_len1.std())
>>0.1736639964801841
  0.46991097723995806

平均は3cm程度異なるようです。また、標準偏差も2.5倍くらい異なっていそうです。標本数は同じのため母分散も等しくはなさそうです。

ヒストグラム

データのばらつきを見るために、それぞれの標本でヒストグラムを描いてみます。

import matplotlib.pyplot as plt

plt.hist(petal_len0,bins=10)
plt.show()
plt.hist(petal_len1,bins=10)
plt.show()
plt.hist(petal_len0,bins=10)
plt.hist(petal_len0,bins=10)
plt.hist(petal_len1,bins=10)
plt.hist(petal_len1,bins=10)

ほんまに正規分布に従ってるんか？って感じです。

Q-Qプロット

グラフで正規性を確認する方法を紹介しました。としてQ-Qプロットも有効です。正規分布に近い場合、きれいな右肩上がりの直線になります。逆に、正規分布に従わない場合はその線が歪んで出てきます。

Q-Qプロットはscipyのstatsモジュールから描画可能です。

from scipy import stats

stats.probplot(petal_len0, plot=plt)
plt.show()

stats.probplot(petal_len1, plot=plt)
plt.show()
stats.probplot(petal_len0, plot=plt)
stats.probplot(petal_len0, plot=plt)
stats.probplot(petal_len1, plot=plt)

どちらも端の方でやや歪んでいるように見えます。きれいな正規分布とはいえなさそうな気がします。

シャピロ・ウィルク検定

シャピロ・ウィルク検定はデータの正規性を検定する手法の1つです。シャピロ・ウィルク検定の帰無仮説は「分布は正規分布に従う」です。今回は有意水準5%(p<=0.05)で検定を行ってみます。p値が0.05以下であれば分布は正規分布であるという帰無仮説が棄却されます。

この検定はscipyのstats.shapiro()で実行可能です。

from scipy import stats

#2つ目の戻り値がp値
p0 = stats.shapiro(petal_len0)[1]
p1 = stats.shapiro(petal_len1)[1]

print("sepal_len_0のp値：",round(p0,3))
print("sepal_len_1のp値：",round(p1,3))

>> sepal_len_0のp値： 0.055
   sepal_len_1のp値： 0.158

scipy.stats.shapiro()の引数には、検定したいデータを渡します。検定結果はp値にて判断しますが、そのp値は2つ目の戻り値にて返されます（なので最後に[1]がついています）。

結果を見ると、両分布とも帰無仮説が採択されました（sepal_len_0はギリですが）。なので、一応正規分布に従うとして検定を進めていきます。

t検定

ここまででデータの分散および正規性を確認しました。ここまでの前提をまとめると、

  • 対応のないデータ
  • 分散は等しくない
  • 正規分布に従う

です。この条件では「Welchのt検定」という手法を用います。Pythonでのt検定は、scipy.statsのttest()で行います。以下で実行してみましょう。

Welchのt検定

先述の通り、Welchのt検定は等分散ではないと仮定した対応の無いt検定です。等分散でない場合は引数equal_varをFalseに指定します。

# Welchのｔ検定
stats.ttest_ind(petal_len0,petal_len1, equal_var=False)

>> Ttest_indResult(statistic=-39.492719391538095, pvalue=9.934432957587695e-46)

p値はpvalue（2つ目の引数）です。9.93で帰無仮説を棄却できないのかと思いがちですが、注意してください。e-46は小数点が右に46個ずれるという意味です。実際はかなり小さい値で、5%水準を軽く下回ります。

よって、2つの標本平均に差がないとは言えないという結論に至りました。

その他t検定

今回の場合、Welchのt検定を用いるのが適切でしたが、他の手法のt検定をPythonではどのように実行するのかも紹介しておきます。

studentのt検定

この検定は対応なし、等分散の前提のもとで実施する検定です。

# studentのｔ検定
stats.ttest_ind(petal_len0,petal_len1)

>> Ttest_indResult(statistic=-39.492719391538095, pvalue=5.404910513441677e-62)

対応のあるt検定

標本に対応がある場合は、scipy.stats.ttest_rel()を用います。

# 対応のあるｔ検定
stats.ttest_rel(petal_len0,petal_len1)

>> Ttest_relResult(statistic=-37.24112591165691, pvalue=1.3246949994286501e-37)

まとめ

Pythonでt検定を行う方法を紹介しました。

検定自体は1行で簡単にできますが、適切な手法を選ぶためにも、その前にデータの特性を確認しておく必要はあると思います。

2標本の比較は実際の場面でもよくあるケースだと思うので、t検定を活用して根拠を持って判断ができるといいですね。

ではでは👋

Pythonの最新記事8件