经验功率估计 (skbio.stats.power

本模块的目的是提供正态和非正态分布数据的经验、事后功率估计。它还为子样本数据提供支持,以便于分析。

其基本原理是基于子抽样和蒙特卡罗模拟。假设有一些种群, \(K_{{1}}, K_{{2}}, ... K_{{n}}\) 有一些财产, \(\mu\) 这样的话 \(\mu_{{1}} \neq \mu_{{2}} \neq ... \neq \mu_{{n}}\) . 对于每个群体,一个样本, \(S\) 可以用参数绘制, \(x\) 在哪里? \(x \approx \mu\) 对于样品,我们可以用一个测试, \(f\) ,来证明这一点 \(x_{{1}} \neq x_{{2}} \neq ... \neq x_{{n}}\) .

既然我们知道了 \(\mu_{{1}} \neq \mu_{{2}} \neq ... \neq \mu_{{n}}\) ,我们知道我们应该拒绝无效假设。如果我们不能拒绝无效假设,我们就犯了第二类错误,我们的结果是假阴性。我们可以通过反复对总体进行二次抽样并观察我们看到假阴性的频率,来估计不同采样深度下II型错误的频率。如果我们对每个子取样深度重复几次,并且改变我们使用的深度,我们就可以开始近似计算我们使用的样本数量和假阴性率之间的关系,也就是测试的统计能力。

为了从功率不足的数据中生成完整的功率曲线 statsmodels.stats.power 包可用于求解效果大小。效果大小可用于外推数据的幂曲线。

此模块中的大多数函数都接受一个统计测试函数,该函数获取样本列表并返回p值。然后在一系列子样本上对测试进行评估。

取样可采用两种方式。对于任何一组样品,我们都可以选择画画 \(n\) 对每个样本进行随机观察。或者,如果元数据可用,则可以基于一组控制类别匹配样本,以便从可用匹配集随机抽取成对样本。

功能

subsample_power(test, samples[, draw_mode, ...])

子样本数据迭代计算功率

subsample_paired_power(test, meta, cat, ...)

使用具有匹配元数据的样本迭代地估计幂

confidence_bound(vec[, alpha, df, axis])

假设正态分布,计算置信区间

paired_subsamples(meta, cat, control_cats[, ...])

绘制样本列表 cat 与之匹配 control_cats

示例

假设我们想测试两个随机变量之间的关系, inddep . 让我们使用随机子抽样来估计α为0.1、0.01和0.001的测试的统计能力。

为了控制伪随机数的生成,我们将使用一个种子。将这些函数用于自己的数据时,不需要包含步骤。

>>> import numpy as np
>>> np.random.seed(20)
>>> ind = np.random.randint(0, 20, 15)
>>> ind
array([ 3, 15,  9, 11,  7,  2,  0,  8, 19, 16,  6,  6, 16,  9,  5])
>>> dep = (3 * ind + 5 + np.random.randn(15) * 5).round(3)
>>> dep
array([ 15.617,  47.533,  28.04 ,  33.788,  19.602,  12.229,   4.779,
        36.838,  67.256,  55.032,  22.157,   7.051,  58.601,  38.664,
        18.783])

让我们定义一个测试,它将绘制一个样本对列表,并确定它们是否相关。我们将使用 scipy.stats.pearsonr 它接受两个数组并返回一个相关系数和一个表示两个分布相关的概率的p值。

>>> from scipy.stats import pearsonr
>>> f = lambda x: pearsonr(x[0], x[1])[1]

现在,让我们使用随机抽样来估计我们对第一个分布的测试的能力。

>>> samples = [ind, dep]
>>> print("%.3e" % f(samples))
3.646e-08

subsample_power ,我们可以通过设置 draw_mode “匹配”。我们也可以设置临界值,这样我们就可以估计临界值的功率 \(\alpha = 0.05\) ,临界值为0.01,临界值为0.001。

>>> from skbio.stats.power import subsample_power
>>> pwr_100, counts_100 = subsample_power(test=f,
...                                       samples=samples,
...                                       max_counts=10,
...                                       min_counts=3,
...                                       counts_interval=1,
...                                       draw_mode="matched",
...                                       alpha_pwr=0.1,
...                                       num_iter=25)
>>> pwr_010, counts_010 = subsample_power(test=f,
...                                       samples=samples,
...                                       max_counts=10,
...                                       min_counts=3,
...                                       counts_interval=1,
...                                       draw_mode="matched",
...                                       alpha_pwr=0.01,
...                                       num_iter=25)
>>> pwr_001, counts_001 = subsample_power(test=f,
...                                       samples=samples,
...                                       max_counts=10,
...                                       min_counts=3,
...                                       counts_interval=1,
...                                       draw_mode="matched",
...                                       alpha_pwr=0.001,
...                                       num_iter=25)
>>> counts_100
array([3, 4, 5, 6, 7, 8, 9])
>>> pwr_100.mean(0)
array([ 0.484,  0.844,  0.932,  0.984,  1.   ,  1.   ,  1.   ])
>>> pwr_010.mean(0)
array([ 0.044,  0.224,  0.572,  0.836,  0.928,  0.996,  1.   ])
>>> pwr_001.mean(0)
array([ 0.   ,  0.016,  0.108,  0.332,  0.572,  0.848,  0.956])

基于这种功率估计,随着我们增强我们对没有犯I型错误和识别出假阳性的信心,我们需要确信我们没有犯II型错误的样本数量会增加。