假设检验

假设检验概览

在总体分布函数完全未知，或者只知其形式、但不知道其参数的情况下，为了推断总体的某些未知的特性，提出的关于总体的假设。

算法流程

假设检验的算法流程如下：

1. 确定显著性水平
2. 确定原假设
3. 确定备选假设
4. 确定检验统计量
5. 使用样本数据计算检验统计量的值
6. 计算检验统计量发生的概率
7. 对比发生的概率和显著性水平，确定是否拒绝原假设。

在实际使用的过程中，因为步骤6是需要根据统计量分布计算的，如使用计算机计算，或者查表比较麻烦，而步骤5可以直接根据样本数据计算。所以实际中步骤5、6、7会做如下修改：

1. 根据统计量分布和显著性水平计算拒绝原假设的临界值，确定拒绝域和接受域
2. 使用样本数据计算检验检验统计量的值
3. 根据检测统计量的值在拒绝域或接受域内来决定是否拒绝原假设。

假设检验的难点主要是指在步骤2和步骤4。在学习假设检验之前，参照附录一先理解一下p_value。

定义
第一类错误：原假设为真，结果拒绝了原假设。
第二类错误：原假设为假，结果接受了原假设。

如何确定原假设

假设检验中，一般情况下原假设是为了实现目标，在基于先验知识的情况下做出的假设，所以假设检验算法本身是倾向于接受原假设的，只有数据给的理由足够充分，才会拒绝原假设。基于这一点，有以下结论：

结论一：如果原假设没有设定好，原假设本身是错的，数据量又不够的情况下，检测结果会比较容易犯第二类错误(附录二)。在数据量已定的情况下，为了使得整个假设检验更有意义，需要根据原理来设置合理的原假设。

结论二：如果数据量很少的话，并没有太大的意义做假设检验。假设检验的结果很容易犯错（附录二）。

结论三：当数据量比较大的时候，需要设定更大的显著水平。否则，只要真实结果与原假设有一点点差，就很容易拒绝原假设（附录三）。

结论四：扩大数据量能够降低犯第二类错误的概率。换个角度，为了控制第二类错误的概率，需要更多的数据量。（多大的数据量是合适的，就引申出功效函数）

确定检验统计量

假设检验结果是根据p_value值来判断的。p_value值，也是一个概率值。概率值常见两种计算方式，一种是使用频率估计概率，一种是基于数据分布函数计算得到。假设检验是检验判断对错的概率，如果使用频率估计概率的思想，就需要相同事件出现很多次，每次判断后得到反馈结果，然后基于判断正确的频率来估计概率，而这在实际中是不现实的。所以假设检验中的p_value是通过数据分布函数计算得到的。

假设检验需要分布函数，那么就需要构造满足某种分布的统计量。而该统计量的值基于样本数据和假设条件计算得到。

不同的数据分布前提，不同的实际情况下，有不同的数据分布。可以从不同的角度构造统计量。统计量有好坏之分，可以基于统计量的无偏性、有效性、相合性等来比较统计量的优劣。

构造统计量一般情况下还有一个前提：尽可能的利用已知的先验知识。在数据量已定的前提下，利用了更多的先验知识，更容易拒绝原假设。能在不提高犯第一列错误的前提下，降低犯第二类错误的概率（附录四）。

在假设检验中，正态总体的假设检验是比较多用的（附录五）。

假设检验案例

基础知识（p142）
定理一： 独立同正态分布的n个随机变量的和还是正态分布。

定理二： 独立同正态分布的n个随机变量的均值和方差相互独立，且方差服从自由度为n-1的卡方分布

定理三： 独立同正态分布的n个随机变量的均值与标准差可以构造t分布。（直接由t分布的定义和定理一、定理二得到）

定理四： 一组独立同正态分布的m个随机变量的方差与另一组独立同正态分布的n个随机变量的方差可以构造F分布。（由定理二和F分布的定义直接得到）；一组独立同正态分布的m个随机变量的均值、方差与另一组独立同正态分布的n个随机变量的均值、方差，当方差相等时，可以构造t分布。由定理一、定理二和t分布的定义可以得到。

• 基于定理一，可以做方差已知的关于均值的显著性检验；
• 基于定理二，可以做方差未知的关于均值的显著性检验；
• 基于定理四，可以做方差相同的两正态总体的均值差的显著性检验；
• 基于定理二，可以做正态总体的方差的显著性检验；
• 基于定理四，可以做两正态总体的方差是否有差异的显著性检验；

附录

函数准备

import numpy as np
from scipy import stats

def z_test(X, t_mean=0, std=1, test_type='both', is_print=True):
    """Z检验
    
    总体方差已知情况下，正态分布总体的均值的检验
    
    Parameters
    ----------
    X: list or np.array, 样本数据
    t_mean: float, 检测均值
    std: float, 总体方差
    test_type: 检测类型，双边检验/左侧单边检验/右侧单边检验
    
    Return
    ------
    float, p_value值
    """
    n = len(X)
    T = np.sqrt(n) * (np.mean(X) - t_mean) / std
    if test_type == 'both':
        p_value = 2 * (1 - stats.norm.cdf(np.abs(T))) # 正态分布就是在这里用到的
    elif test_type == 'right':
        p_value = stats.norm.cdf(T)
    else:
        p_value = 1 - stats.norm.cdf(T)
    if is_print:
        print('t_mean:{}, std:{}, test_type:{} | p_value:{:.5}'.format(
            t_mean, std, test_type, p_value))
    return p_value

def t_test(X, t_mean=0, test_type='both', is_print=True):
    """t检验
    
    总体方差未知情况下，正态分布总体的均值的检验
    
    Parameters
    ----------
    X: list or np.array, 样本数据
    t_mean: float, 检测均值
    test_type: 检测类型，双边检验/左侧单边检验/右侧单边检验
    
    Return
    ------
    float, p_value值
    """
    n = len(X)
    T = (np.mean(X) - t_mean) / (np.std(X) / np.sqrt(n-1))
    if test_type == 'both':
        p_value = 2 * (1 - stats.t.cdf(np.abs(T), n-1)) # t分布就是在这里用到的 
    elif test_type == 'right':
        p_value = stats.t.cdf(T, n-1)
    else:
        p_value = 1 - stats.t.cdf(T, n-1)
    if is_print:
        print('t_mean:{}, test_type:{} | p_value:{:.5}'.format(
            t_mean, test_type, p_value))
    return p_value

附录一

假设数据取自同一个分布，使用相同的显著水平，相同的检测方法做检验。

for i in range(5):
    _ = z_test(np.random.normal(0, 1, 100), t_mean=0, std=1, is_print=True)

结果如下：
t_mean:0, std:1, test_type:both | p_value:0.83079
t_mean:0, std:1, test_type:both | p_value:0.042919
t_mean:0, std:1, test_type:both | p_value:0.25172
t_mean:0, std:1, test_type:both | p_value:0.80495
t_mean:0, std:1, test_type:both | p_value:0.91295

每次的检验结果是可能不同的，对比理解p_value的含义。对比不同的次数下p_value的表现。

n = 5
p_value = [np.mean([z_test(np.random.normal(0, 1, 100), t_mean=0, std=1, is_print=False) for i in range(n)]) for j in range(10)]
print('n={}, p_value_mean={:.3}, p_value_var={:.3}'.format(n, np.mean(p_value), np.var(p_value)))
n = 10
p_value = [np.mean([z_test(np.random.normal(0, 1, 100), t_mean=0, std=1, is_print=False) for i in range(n)]) for j in range(10)]
print('n={}, p_value_mean={:.3}, p_value_var={:.3}'.format(n, np.mean(p_value), np.var(p_value)))
n = 20
p_value = [np.mean([z_test(np.random.normal(0, 1, 100), t_mean=0, std=1, is_print=False) for i in range(n)]) for j in range(10)]
print('n={}, p_value_mean={:.3}, p_value_var={:.3}'.format(n, np.mean(p_value), np.var(p_value)))

结果如下：
n=5, p_value_mean=0.48, p_value_var=0.0157
n=10, p_value_mean=0.504, p_value_var=0.00723
n=20, p_value_mean=0.499, p_value_var=0.00193

当相同情况做的次数多了以后，检测的p_value的均值是越来越稳定的。之后的对比不同检验情况下p_value的差异的时候设置n为100。

附录二

样本取自相同的总体分布，假设总体分布为标准正态分布$N(0,1)$。对比相同的数据量下，Z检验对关于均值为0和均值为0.1的检测结果。

sample_num = 10
t_mean = 0
p_value = np.mean([z_test(np.random.normal(0, 1, sample_num), t_mean=t_mean, std=1, is_print=False) for i in range(100)])
print('sample num:{}, t_mean:{}, p_value={}'.format(sample_num, t_mean, p_value))
t_mean = 0.1
p_value = np.mean([z_test(np.random.normal(0, 1, sample_num), t_mean=t_mean, std=1, is_print=False) for i in range(100)])
print('sample num:{}, t_mean:{}, p_value={}'.format(sample_num, t_mean, p_value))

sample_num = 100
t_mean = 0
p_value = np.mean([z_test(np.random.normal(0, 1, sample_num), t_mean=t_mean, std=1, is_print=False) for i in range(100)])
print('sample num:{}, t_mean:{}, p_value={}'.format(sample_num, t_mean, p_value))
t_mean = 0.1
p_value = np.mean([z_test(np.random.normal(0, 1, sample_num), t_mean=t_mean, std=1, is_print=False) for i in range(100)])
print('sample num:{}, t_mean:{}, p_value={}'.format(sample_num, t_mean, p_value))

sample_num = 1000
t_mean = 0
p_value = np.mean([z_test(np.random.normal(0, 1, sample_num), t_mean=t_mean, std=1, is_print=False) for i in range(100)])
print('sample num:{}, t_mean:{}, p_value={}'.format(sample_num, t_mean, p_value))
t_mean = 0.1
p_value = np.mean([z_test(np.random.normal(0, 1, sample_num), t_mean=t_mean, std=1, is_print=False) for i in range(100)])
print('sample num:{}, t_mean:{}, p_value={}'.format(sample_num, t_mean, p_value))

结果如下：
sample num:10, t_mean:0, p_value=0.5072029252323097
sample num:10, t_mean:0.1, p_value=0.46138045248769893
sample num:100, t_mean:0, p_value=0.4875769624054251
sample num:100, t_mean:0.1, p_value=0.4245074008824946
sample num:1000, t_mean:0, p_value=0.4773629769262692
sample num:1000, t_mean:0.1, p_value=0.023915295736835533

对比不同的样本数量下，不同的均值的检测结果

sample_num = 10
t_mean = 0
p_value = np.mean([z_test(np.random.normal(0, 1, sample_num), t_mean=t_mean, std=1, is_print=False) for i in range(100)])
print('sample num:{}, t_mean:{}, p_value={}'.format(sample_num, t_mean, p_value))
t_mean = 0.5
p_value = np.mean([z_test(np.random.normal(0, 1, sample_num), t_mean=t_mean, std=1, is_print=False) for i in range(100)])
print('sample num:{}, t_mean:{}, p_value={}'.format(sample_num, t_mean, p_value))

sample_num = 5
t_mean = 0
p_value = np.mean([z_test(np.random.normal(0, 1, sample_num), t_mean=t_mean, std=1, is_print=False) for i in range(100)])
print('sample num:{}, t_mean:{}, p_value={}'.format(sample_num, t_mean, p_value))
t_mean = 1
p_value = np.mean([z_test(np.random.normal(0, 1, sample_num), t_mean=t_mean, std=1, is_print=False) for i in range(100)])
print('sample num:{}, t_mean:{}, p_value={}'.format(sample_num, t_mean, p_value))

结果如下：
sample num:10, t_mean:0, p_value=0.47722880090175634
sample num:10, t_mean:0.5, p_value=0.23667633388385417
sample num:5, t_mean:0, p_value=0.5310308791352258
sample num:5, t_mean:1, p_value=0.11377033525596417

附录三

样本取自相同的总体分布，假设总体分布为标准正态分布$N(0,1)$。对比不同的数据量的情况下，z检验对关于均值为0.1的检测结果差异。

sample_num = 10
p_value = np.mean([z_test(np.random.normal(0, 1, sample_num), t_mean=0.1, std=1, is_print=False) for i in range(100)])
print('sample num:{}, p_value={}'.format(sample_num, p_value))
sample_num = 100
p_value = np.mean([z_test(np.random.normal(0, 1, sample_num), t_mean=0.1, std=1, is_print=False) for i in range(100)])
print('sample num:{}, p_value={}'.format(sample_num, p_value))
sample_num = 1000
p_value = np.mean([z_test(np.random.normal(0, 1, sample_num), t_mean=0.1, std=1, is_print=False) for i in range(100)])
print('sample num:{}, p_value={}'.format(sample_num, p_value))

结果如下：
sample num:10, p_value=0.4848865002720213
sample num:100, p_value=0.3826647304964832
sample num:1000, p_value=0.026368291730919216

样本量越多，检测结果越容易拒绝原假设。

附录四

在样本量已定的前提下。当总体方差已知的前提下，对比t_test和z_test的检测结果。

sample_num = 10
t_p_value = []
z_p_value = []
for j in range(1000):
    X = np.random.normal(0, 0.5, sample_num)
    t_p_value.append(t_test(X, t_mean=0.2, is_print=False))
    z_p_value.append(z_test(X, t_mean=0.2, std=0.5, is_print=False))
t_p_value = np.mean(t_p_value)
z_p_value = np.mean(z_p_value)
print('sample num:{}, t_p_value={}, z_p_value={}'.format(sample_num, t_p_value, z_p_value))

结果如下：
sample num:10, t_p_value=0.31059238958063545, z_p_value=0.30291339323045574

当总体方差已知的前提下，t_test和z_test都可以做正态总体均值的假设检验，但是因为z_test使用了总体方差的信息，一般情况下，当原假设出错的时候，z_test的p值会比较小。

附录五

中心极限定理：假设随机变量$X_1,X_2,...,X_n$相互独立，服从同一分布。且具有数学期望和方差：$E(X_k)=\mu ,D(X_k)={\sigma }^2>0(k=1,2,3...)$，则随机变量之和${\sum }_{k=1}^n{X}_k$的标准化变量
$$Y_n=\frac{{\sum }_{k=1}^n{X}_k-E({\sum }_{k=1}^n{X}_k)}{\sqrt{D({\sum }_{k=1}^n{X}_k)}}=\frac{{\sum }_{k=1}^n{X}_k-n\mu }{\sqrt{n}\sigma }$$
的分布函数$F_n(x)$对于任意x满足
lim ⁡ n → + ∞ F n ( x ) = lim ⁡ n → + ∞ P { ∑ k = 1 n X k − n μ n σ ≤ x } = ∫ − ∞ x 1 2 π e − t 2 / 2 d t = Φ ( x ) \lim\limits_{n \rightarrow +\infty}F_n(x)=\lim\limits_{n \rightarrow +\infty}P\{ {\sum_{k=1}^nX_k-n\mu \over \sqrt n\sigma}\leq x\}=\int_{-\infty}^x \frac{1}{\sqrt {2\pi}}e^{-t^2/2}{\rm d}t=\Phi(x) n→+∞lim​Fn​(x)=n→+∞lim​P{n ​σ∑k=1n​Xk​−nμ​≤x}=∫−∞x​2π ​1​e−t2/2dt=Φ(x)

# 总体数据为y
y = np.random.choice(list(range(100)), 1000)
# 样本数据为从总体y中独立抽取的50个数据并取均值,重复抽取100次
X = [np.mean(np.random.choice(y, 50)) for i in range(100)]
# 使用shapiro做X的正态性检验
# 为了得到更加可靠的数据，做50次检验取均值
p_value = np.mean([stats.shapiro([np.mean(np.random.choice(y, 50)) for i in range(100)])[1] for i in range(50)])
print('{:.5}'.format(p_value))

结果如下：
0.43548

检测结果发现，这样取出来的值服从正态分布，跟中心极限定理一致。因为每个X的均值的期望与等于总体y的均值。所以可以根据X的值来做总体y的均值的显著性检验。

很容易犯的一个错误是直接对从总体y中随机选择的数据做t检验。