python 实现方差检验

前面介绍的$z$检验和$t$检验主要用于对总体的均值进行假设检验，下面主要介绍对总体方差进行假设检验的方法。

一. 单个总体方差检验

设总体$X\sim N(\mu ,{\delta }^2)$，$\mu ,{\delta }^2$均未知，$X_1,X_2...,X_n$是来自$X$的样本，则
$$\frac{(n-1)S^2}{{\delta }^2}\sim {\chi }^2(n-1)$$

其中S为样本标准差，$\delta$为总体标准差，n为样本容量

根据以上的统计量，可以对单个总体的方差$\delta$进行假设检验。

原假设$H_0$	备则假设$H_1$	拒绝域
${\delta }^2={\delta }_0^2$	${\delta }^2\ne {\delta }_0^2$	${\chi }^2>{\chi }_{\alpha /2}^2(n-1)$ 或 ${\chi }^2<{\chi }_{1-\alpha /2}^2(n-1)$
${\delta }^2\le {\delta }_0^2$	${\delta }^2>{\delta }_0^2$	${\chi }^2>{\chi }_{\alpha }^2(n-1)$
${\delta }^2\ge {\delta }_0^2$	${\delta }^2<{\delta }_0^2$	${\chi }^2<{\chi }_{1-\alpha }^2(n-1)$

某厂生产某种型号的电池，其寿命长期以来服从方差${\delta }^2=5000$的正太分布，现有一批这种电池，从它的生产情况看，寿命的波动性有所改变。现随机抽26只电池，测出其寿命的样本方差$S^2=9200$，请问根据这一数据能否推断这批电池寿命的波动性较以前有显著的改变？

1. 双侧检验

在$\alpha =0.05$的水平下做出如下假设：
$$H_0:{\delta }^2=5000H_1:{\delta }^2\ne 5000$$

从备则假设可知，该假设检验是一个双边检验。

代码实现：

from scipy.stats import chi2

if __name__ == '__main__':

    alpha = 0.05
    n = 26

    # 总体方差
    var = 5000
    # 样本方差
    sample_var = 10000
    # 卡方统计量
    chi2stats = (n - 1) * var / sample_var

    rv = chi2(df=n - 1)
    chi2left = rv.ppf(alpha / 2)
    chi2right = rv.ppf(1 - alpha / 2)

    print("chi2left: {}, chi2right: {}".format(round(chi2left, 3), round(chi2right, 3)))
    print("chi2stats: ", round(chi2stats, 3))

    # 计算p值
    if chi2stats > chi2right:
        pval = 2 * (1 - rv.cdf(chi2stats))
    else:
        pval = 2 * (rv.cdf(chi2stats))

    if chi2stats > chi2right or chi2stats < chi2left:
        print("reject null hypothesis, pval: ", round(pval, 3))
    else:
        print("not reject null hypothesis, pval: ", round(pval, 3))

运行结果：

chi2left: 13.12, chi2right: 40.646
chi2stats:  12.5
reject null hypothesis, pval:  0.036

2. 单侧检验

保持$\alpha =0.05$不变，修改原假设和备则假设$$H_0:{\delta }^2\ge 5000H_1:{\delta }^2<5000$$，此时就变成了一个左侧检验。

代码实现：

from scipy.stats import chi2

if __name__ == '__main__':

    alpha = 0.05
    n = 26

    # 总体方差
    var = 5000
    # 样本方差
    sample_var = 10000
    # 卡方统计量
    chi2stats = (n - 1) * var / sample_var

    rv = chi2(df=n - 1)
    chi2left = rv.ppf(alpha)

    print("chi2left: {}".format(round(chi2left, 3)))
    print("chi2stats: ", round(chi2stats, 3))

    # 计算p值
    pval = rv.cdf(chi2stats)

    if chi2stats < chi2left:
        print("reject null hypothesis, pval: ", round(pval, 3))
    else:
        print("not reject null hypothesis, pval: ", round(pval, 3))

运行结果：

chi2left: 14.611
chi2stats:  12.5
reject null hypothesis, pval:  0.018

二. 两个总体方差检验

设$X_1,X_2...,X_n$是来自总体$N({\mu }_1,{\delta }_1^2)$的样本，$Y_1,Y_2...,Y_n$是来自总体$N({\mu }_2,{\delta }_2^2)$的样本，且两样本独立，其样本方差分别为$S_1^2,S_2^2$，且${\mu }_1,{\mu }_2,{\delta }_1^2,{\delta }_2^2$均未知，则
$$\frac{S_1^2/S_2^2}{{\delta }_1^2/{\delta }_2^2}\sim F(n_1-1,n_2-1)$$

${\delta }_1,{\delta }_2,$为两个总体的标准差，$n_1,n_2$为两个样本的容量

根据以上的统计量，可以对两个总体的方差之间的大小关系${\delta }_1^2/{\delta }_2^2$进行假设检验。

原假设$H_0$	备则假设$H_1$	拒绝域
${\delta }_1^2={\delta }_2^2$	${\delta }_1^2\ne {\delta }_2^2$	$F>F_{\alpha /2}(n_1-1,n_2-1)$ 或 $F<F_{1-\alpha /2}(n_1-1,n_2-1)$
${\delta }_1^2\le {\delta }_2^2$	${\delta }_1^2>{\delta }_2^2$	$F>F_{\alpha }(n_1-1,n_2-1)$
${\delta }_1^2\ge {\delta }_2^2$	${\delta }_1^2<{\delta }_2^2$	$F<F_{1-\alpha }(n_1-1,n_2-1)$

两工厂生产同一份零件且零件尺寸服从正态分布，A生产的16个零件中，均值${\mu }_A$为3.1，样本方差${\delta }_A^2$为0.04；B生产的25个零件中，均值${\mu }_B$为3.08，样本方差${\delta }_B^2$为0.09。问能否在$\alpha =0.05$的显著性水平下认为两工厂生产零件的方差相同？

1. 双边检验

在$\alpha =0.05$的水平下做出如下假设：
$$H_0:{\delta }_A^2={\delta }_B^2{H}_1:{\delta }_A^2\ne {\delta }_B^2$$
从备则假设可知，该假设检验是一个双边检验。

代码实现：

from scipy.stats import f

if __name__ == '__main__':

    # 两样本的方差与样本数量
    var1, var2 = 0.04, 0.09
    df1, df2 = 16, 25

    alpha = 0.05
    rv = f(dfn= df1 - 1, dfd= df2 - 1)
    # 计算拒绝域
    fleft = rv.ppf(alpha / 2)
    fright = rv.ppf(1 - alpha / 2)
    # 计算统计量
    fstats = var1 / var2

    print("fleft: {}, fright: {}".format(round(fleft, 3), round(fright, 3)))
    print("fstats: ", round(fstats, 3))
    pval = rv.sf(fstats)

    if fstats > fright:
        pval = 2 * rv.sf(fstats)
    else:
        pval = 2 * rv.cdf(fstats)

    if fstats > fright or fstats < fleft:
        print("reject null hypothesis, pval is ", round(pval, 3))
    else:
        print("not reject null hypothesis, pval is ", round(pval, 3))

运行结果：

fleft: 0.37, fright: 2.437
fstats:  0.444
not reject null hypothesis, pval is  0.107

2. 单边检验

保持$\alpha =0.05$不变，修改原假设和备则假设$$H_0:{\delta }_A^2\ge {\delta }_B^2{H}_1:{\delta }_A^2<{\delta }_B^2$$，此时就变成了一个左侧检验。

代码实现：

from scipy.stats import f

if __name__ == '__main__':
    # 两样本的方差与样本数量
    var1, var2 = 0.04, 0.09
    df1, df2 = 16, 25

    alpha = 0.05
    rv = f(dfn=df1 - 1, dfd=df2 - 1)
    # 计算拒绝域
    fleft = rv.ppf(alpha)

    # 计算统计量
    fstats = var1 / var2

    print("fleft: {}".format(round(fleft, 3)))
    print("fstats: ", round(fstats, 3))
    pval = rv.cdf(fstats)


    if fstats < fleft:
        print("reject null hypothesis, pval is ", round(pval, 3))
    else:
        print("not reject null hypothesis, pval is ", round(pval, 3))

运行结果：

fleft: 0.437
fstats:  0.444
not reject null hypothesis, pval is  0.053