静かなる名辞

pythonとプログラミングのこと

2019/03/22:TechAcademyがteratailの質問・回答を盗用していた件
2019/03/26:TechAcademy盗用事件 公式発表と深まる疑念



【python】statsmodelsでt検定する方法

はじめに

 この前はscipyでやる方法をまとめたわけですが、

www.haya-programming.com


 片側のオプションほしいなーと思ったのでstatsmodelsに浮気することにしました。

使い方の概要

 対応のないt検定はこれです。

statsmodels.stats.weightstats.ttest_ind — statsmodels 0.9.0 documentation

 引数は以下のようなものです。

statsmodels.stats.weightstats.ttest_ind(x1, x2,
alternative='two-sided', usevar='pooled', weights=(None, None), value=0)

 x1とx2が検定対象です。オプションで重要そうなのはalternativeとusevarです。

  • alternative

 "two-sided"と"larger"と"smaller"が渡せます。後ろの2つは片側のオプションです。

 こうしちゃうと、どっちがどっちに対して大きい/小さいのかわからん、とか考えないのかしら(x1 larger/smaller than x2(x1はx2より大きい/小さい)です。)。

  • usevar

 "pooled"と"unequal"が渡せます。スチューデントのt検定とウェルチのt検定に対応しているはずです。

 難しくはないけど、いちいちstrで重要なオプションを渡さないといけないのは少し面倒ですかね。間違った値にするとエラー出してくれるみたいなので、そういう意味では親切ですが。

実際にやってみる

 完全に前回と同じノリでやります。

 分布をとりあえず作る。

>>> from scipy import stats
>>> from statsmodels.stats.weightstats import ttest_ind
>>> d1 = stats.norm(loc=0, scale=1)
>>> d2 = stats.norm(loc=1, scale=1)

 やる。

>>> a = d1.rvs(3)
>>> b = d2.rvs(3)
>>> a
array([ 0.90831032, -0.88621515,  0.09060862])
>>> b
array([1.87384828, 2.92258811, 0.8539749 ])
>>> ttest_ind(a, b)
(-2.333627132285731, 0.07993392828229898, 4.0)

 t統計量、p値、自由度が返るみたいですね。

 ここで仮に有意水準0.05とすると、もう少しでいけそうだったけど切れなかったという残念な例です(そうなるまで何回か回しました)。そこで、片側検定にしてみます。

>>> ttest_ind(a, b, alternative="smaller")
(-2.333627132285731, 0.03996696414114949, 4.0)

 教科書通り半分のp値になって、めでたく「有意差」が出ました。教科書には「こういうことはやるな」と書いてあると思います。気をつけましょう。

 なお、片側の方向を間違えて指定した場合は、

>>> ttest_ind(a, b, alternative="larger")
(-2.333627132285731, 0.9600330358588506, 4.0)

 だいたい大きいp値になるので気づきます。でも、元のp値が0.42とかだったりすると案外気づかないかもしれない(どっちにしろ有意差ではないし、実害ないかもですが・・・)。

まとめ

 こっちの方が高機能だし、これでいいんじゃ? という感じもします。でもscipy入れててもstatsmodels入れてない人は多いと思うので、微妙っちゃ微妙ですね。