python实现假设检验-t检验

一. 什么是t检验

设总体$X\sim N(\mu ,{\delta }^2)$，其中$\mu ,{\delta }^2$未知，统计量$t=\frac{\overline{X}-\mu }{S/\sqrt{n}}$服从标准正太分布，可以利用该统计量对总体均值$\mu$进行假设检验，这种检验方法称为$t$检验法。

适用条件：总体服从正态分布且方差未知；样本容量较小，一般小于30。

二. 常见t检验的实现

1. 单样本t检验

原假设$H_0$	备则假设$H_1$	拒绝域
$\mu \le {\mu }_0$	$\mu >{\mu }_0$	$t>t_{\alpha }(n-1)$
$\mu \ge {\mu }_0$	$\mu <{\mu }_0$	$t<-t_{\alpha }(n-1)$
$\mu ={\mu }_0$	$\mu \ne {\mu }_0$	$|t|>t_{\alpha /2}(n-1)$

某元件的寿命$X$服从参数为$N(\mu ,{\delta }^2)$的正太分布，$\mu ,\delta$未知，现抽取10个元件测得元件的寿命的平均值为240小时，标准差为50，问是否有理由认为元件总体的平均寿命大于225小时？

因为总体的方差$\delta$未知,且样本量小于10，因此需要使用$t$检验，设置信水平$\alpha =0.05$
$$H_0:u\ge 225H_1:u<225$$

代码实现：

from scipy.stats import t
import math

if __name__ == '__main__':
    n = 10
    sample_mean = 240
    sample_std = 50
    alpha = 0.05

    t_statistic = (sample_mean - 225) / (sample_std / math.sqrt(n))
    t_left = t(n - 1).ppf(alpha)
    print("t_left:", round(t_left, 3), "t_statistic:", round(t_statistic, 3))

    pval = t(n - 1).cdf(t_statistic)
    if t_statistic < t_left:
        print("reject null hypothesis, pval is", round(pval, 3))
    else:
        print("not reject null hypothesis, pval is", round(pval, 3))

运行结果：

t_left: -1.833 t_statistic: 0.949
not reject null hypothesis, pval is 0.816

2. 配对t检验

配对t检验一般用于比较同一研究对象处理前与处理后的效果比较，又称为重复测量设计下均值差异的检验。比如医学上药物效果的检验。

原假设$H_0$	备则假设$H_1$	拒绝域
${\mu }_D\le 0$	${\mu }_D>0$	$t>t_{\alpha }(n-1)$
${\mu }_D\ge 0$	${\mu }_D<0$	$t<-t_{\alpha }(n-1)$
${\mu }_D=0$	${\mu }_D\ne 0$	$|t|>t_{\alpha /2}(n-1)$

某研究者为了研究新型降压药的效果，现收集了8名高血压患者服药前和服药后的血压数据。

高血压患者	1	2	3	4	5	6	7	8
服药前	137	135	132	151	132	150	142	145
服药后	124	129	135	144	123	142	142	134
差值	13	6	-3	7	9	8	0	11

现在想知道该药物对高血压患者是否有明显效果？

要想验证药物是否有效，只需要验证服药前与服药后差异值的均值与0的关系，如果差值小于等于0则说明药物对高血压没有效果，反之则证明有效。设置信水平$\alpha =0.05$
$${\mu }_D:u\le 0H_1:{\mu }_D>0$$

代码实现：

from scipy.stats import t
import math
import numpy as np

if __name__ == '__main__':

    diff = [13, 6, -3, 7, 9, 8, 0, 11]
    diff_mean = np.mean(diff)
    diff_std = np.std(diff, ddof=1)
    diff_length = len(diff)
    alpha = 0.05

    t_statistic = diff_mean / (diff_std / math.sqrt(diff_length))
    t_right = t(diff_length - 1).ppf(1 - alpha)
    print("t_right:", round(t_right, 3))
    print("t_statistic:", round(t_statistic, 3))

    pval = t(diff_length - 1).sf(t_statistic)
    if t_statistic > t_right:
        print("reject null hypothesis, pval is", round(pval, 3))
    else:
        print("not reject null hypothesis, pval is", round(pval, 3))

运行结果：

t_right: 1.895
t_statistic: 3.341
reject null hypothesis, pval is 0.006

拒绝原假设，因此说明该药物对降压有明显效果。

3. 两样本t检验

设$X_1,X_2...X_{n1}$是来自正太总体$N({\mu }_1,{\delta }_1^2)$的样本,$Y_1,Y_2...Y_{n2}$是来自正太总体$N({\mu }_2,{\delta }_2^2)$的样本，且两样本独立，${\mu }_1,{\mu }_2,{\delta }_1^2,{\delta }_2^2$均未知。统计量$\frac{(\overline{X}-\overline{Y})-\delta }{S_w\sqrt{\frac{1}{n_1}+\frac{1}{n_2}}}$服从$t(n_1+n_2-2)$的分布。

其中$\overline{X},\overline{Y}$分别为两样本的均值，$n_1,n_2$为两样本的容量，$\delta ={\mu }_1-{\mu }_2$, $S_w^2=\frac{(n_1-1)S_1^2+(n_2-1)S_2^2}{n_1+n_2-2}$, $S_w=\sqrt{S_w^2}$

两样本t检验主要用来对两个总体的均值之间的关系做假设检验。

原假设$H_0$	备则假设$H_1$	拒绝域
${\mu }_1-{\mu }_2\le \delta$	${\mu }_1-{\mu }_2>\delta$	$t>t_{\alpha }(n_1+n_2-2)$
${\mu }_1-{\mu }_2\ge \delta$	${\mu }_1-{\mu }_2<\delta$	$t<t_{\alpha }(n_1+n_2-2)$
${\mu }_1-{\mu }_2=\delta$	${\mu }_1-{\mu }_2\ne \delta$	$|t|>t_{\alpha /2}(n_1+n_2-2)$

一家香水的制造商共有两条生产线，现在对两条生产线进行抽样，A生产线抽取18个样本，样本均值为80，样本标准差为5，B生产线抽取20个样本，样本均值为76，样本标准差为4。问这两条生产线所生产的香水的平均分量是否相同？

在$\alpha =0.05$的水平下做出如下假设：
$$H_0:{\mu }_A-{\mu }_B=0H_1:{\mu }_A-{\mu }_B\ne 0$$

代码实现：

from scipy.stats import t
import math

if __name__ == '__main__':

    n1 = 18
    n2 = 20
    sample_mean1 = 80
    sample_mean2 = 76
    sample_std1 = 5
    sample_std2 = 4
    alpha = 0.05

    sw = math.sqrt(((n1 - 1) * sample_std1**2 + (n2 - 1) * sample_std2**2) / (n1 + n2 - 2))
    t_statistic = (sample_mean1 - sample_mean2) / (sw * math.sqrt(1 / n1 + 1/ n2))
    t_left = t(n1 + n2 - 2).ppf(alpha/2)
    t_right = t(n1 + n2 - 2).ppf(1 - alpha/2)
    print("t_left:", round(t_left, 3), "t_right:", round(t_right, 3))
    print("t_statistic:", round(t_statistic, 3))

    if t_statistic >= t_right:
        pval = t(n1 + n2 -2).sf(t_statistic) * 2
    else:
        pval = 2 * t(n1 + n2 - 2)

    if t_statistic < t_left or t_statistic > t_right:
        print("reject null hypothesis, pval is", round(pval, 3))
    else:
        print("not reject null hypothesis, pval is", round(pval, 3))

运行结果：

t_left: -2.028 t_right: 2.028
t_statistic: 2.736
reject null hypothesis, pval is 0.01