numpyを使って統計量「t検定」を計算する方法

Question

Python で作成したモデルに関する統計を生成したいと考えています。それに対して t 検定を生成したいのですが、numpy/scipy を使用してこれを行う簡単な方法があるかどうか疑問に思っていました。何か良い説明はありますか？

たとえば、次のような 3 つの関連データセットがあります。

[55.0, 55.0, 47.0, 47.0, 55.0, 55.0, 55.0, 63.0]

ここで、それらに対してスチューデントの t 検定を実行したいと思います。

Answer 1

のA scipy.stats のパッケージいくつかのttest_...機能があります。ここでのから例を参照してください。 >：

>>> print 't-statistic = %6.3f pvalue = %6.4f' %  stats.ttest_1samp(x, m)
t-statistic =  0.391 pvalue = 0.6955

Answer 2

scipyを使用したvanの答えは正確であり、 scipy.stats.ttest_* 機能がとても便利です。

しかし、私は純粋な解決策を探してこのページに来ました しこり, 、見出しに記載されているように、scipy 依存を避けるためです。この目的のために、ここで挙げた例を指摘しておきます。 https://docs.scipy.org/doc/numpy/reference/generated/numpy.random.standard_t.html

主な問題は、numpy には累積分布関数がないことです。したがって、私の結論は、実際には scipy を使用する必要があるということです。とにかく、numpy のみを使用することは可能です。

元の質問から、データセットを比較し、有意な偏差があるかどうかを t 検定で判断したいと推測します。さらに、サンプルはペアになっているのでしょうか?（見る https://en.wikipedia.org/wiki/Student%27s_t-test#Unpaired_and_paired_two-sample_t-tests ）その場合、次のようにt-とp値を計算できます。

import numpy as np
sample1 = np.array([55.0, 55.0, 47.0, 47.0, 55.0, 55.0, 55.0, 63.0])
sample2 = np.array([54.0, 56.0, 48.0, 46.0, 56.0, 56.0, 55.0, 62.0])
# paired sample -> the difference has mean 0
difference = sample1 - sample2
# the t-value is easily computed with numpy
t = (np.mean(difference))/(difference.std(ddof=1)/np.sqrt(len(difference)))
# unfortunately, numpy does not have a build in CDF
# here is a ridiculous work-around integrating by sampling
s = np.random.standard_t(len(difference), size=100000)
p = np.sum(s<t) / float(len(s))
# using a two-sided test
print("There is a {} % probability that the paired samples stem from distributions with the same means.".format(2 * min(p, 1 - p) * 100))

これで印刷されます There is a 73.028 % probability that the paired samples stem from distributions with the same means. これは正常な信頼区間 (たとえば 5%) をはるかに上回っているため、具体的なケースについては何も結論付ける必要はありません。

Answer 3

あなたはt値を取得したら、確率として、それをどのように解釈するかを不思議に思うかもしれません - 私がしました。ここで私はそれを助けるために書いた関数である。

これは、私が http://www.vassarstats.net/rsig.htmlと http://en.wikipedia.org/wiki/Student%27s_t_distribution に。

# Given (possibly random) variables, X and Y, and a correlation direction,
# returns:
#  (r, p),
# where r is the Pearson correlation coefficient, and p is the probability
# of getting the observed values if there is actually no correlation in the given
# direction.
#
# direction:
#  if positive, p is the probability of getting the observed result when there is no
#     positive correlation in the normally distributed full populations sampled by X
#     and Y
#  if negative, p is the probability of getting the observed result, when there is no
#     negative correlation
#  if 0, p is the probability of getting your result, if your hypothesis is true that
#    there is no correlation in either direction
def probabilityOfResult(X, Y, direction=0):
    x = len(X)
    if x != len(Y):
        raise ValueError("variables not same len: " + str(x) + ", and " + \
                         str(len(Y)))
    if x < 6:
        raise ValueError("must have at least 6 samples, but have " + str(x))
    (corr, prb_2_tail) = stats.pearsonr(X, Y)

    if not direction:
        return (corr, prb_2_tail)

    prb_1_tail = prb_2_tail / 2
    if corr * direction > 0:
        return (corr, prb_1_tail)

    return (corr, 1 - prb_1_tail)