用scipy获得置信区间的正确方法 [英] Correct way to obtain confidence interval with scipy

问题描述

我有一个一维数据数组:

I have a 1-dimensional array of data:

a = np.array([1,2,3,4,4,4,5,5,5,5,4,4,4,6,7,8])

我想获得 68% 的置信区间(即:1 西格玛).

for which I want to obtain the 68% confidence interval (ie: the 1 sigma).

这个答案中的第一条评论指出，这可以使用 scipy.stats.norm.interval 来实现 来自 scipy.stats.norm 函数，通过:

The first comment in this answer states that this can be achieved using scipy.stats.norm.interval from the scipy.stats.norm function, via:

from scipy import stats
import numpy as np
mean, sigma = np.mean(a), np.std(a)

conf_int = stats.norm.interval(0.68, loc=mean, 
    scale=sigma)

但是这篇文章中的评论指出获得置信区间的实际正确方法是:

But a comment in this post states that the actual correct way of obtaining the confidence interval is:

conf_int = stats.norm.interval(0.68, loc=mean, 
    scale=sigma / np.sqrt(len(a)))

即sigma除以样本大小的平方根:np.sqrt(len(a)).

that is, sigma is divided by the square-root of the sample size: np.sqrt(len(a)).

问题是:哪个版本是正确的?

The question is: which version is the correct one?

推荐答案

单次抽取的 68% 置信区间来自正态分布均值 mu 和标准差 sigma 为

The 68% confidence interval for a single draw from a normal distribution with mean mu and std deviation sigma is

stats.norm.interval(0.68, loc=mu, scale=sigma)

N 的均值 的 68% 置信区间来自正态分布平均 mu 和标准偏差 sigma 是

The 68% confidence interval for the mean of N draws from a normal distribution with mean mu and std deviation sigma is

stats.norm.interval(0.68, loc=mu, scale=sigma/sqrt(N))

直觉上，这些公式是有道理的，因为如果你拿着一罐软糖，让很多人猜软糖的数量，每个人可能会偏离很多——同样的标准偏差sigma -- 但是猜测的平均值可以很好地估计实际数字，这反映在均值的标准偏差缩小了 1/sqrt(N).

Intuitively, these formulas make sense, since if you hold up a jar of jelly beans and ask a large number of people to guess the number of jelly beans, each individual may be off by a lot -- the same std deviation sigma -- but the average of the guesses will do a remarkably fine job of estimating the actual number and this is reflected by the standard deviation of the mean shrinking by a factor of 1/sqrt(N).

如果单次抽签有方差 sigma**2，则通过 Bienaymé 公式，N uncorrelated 抽奖的总和有方差N*sigma**2.

If a single draw has variance sigma**2, then by the Bienaymé formula, the sum of N uncorrelated draws has variance N*sigma**2.

均值等于总和除以 N.当您将随机变量(如总和)乘以常数时，方差乘以常数平方.那是

The mean is equal to the sum divided by N. When you multiply a random variable (like the sum) by a constant, the variance is multiplied by the constant squared. That is

Var(cX) = c**2 * Var(X)

所以均值的方差等于

(variance of the sum)/N**2 = N * sigma**2 / N**2 = sigma**2 / N

因此均值的标准差(即方差的平方根)等于

and so the standard deviation of the mean (which is the square root of the variance) equals

sigma/sqrt(N).

这就是分母中sqrt(N)的由来.

以下是一些示例代码，基于 Tom 的代码，用于演示上述声明:

Here is some example code, based on Tom's code, which demonstrates the claims made above:

import numpy as np
from scipy import stats

N = 10000
a = np.random.normal(0, 1, N)
mean, sigma = a.mean(), a.std(ddof=1)
conf_int_a = stats.norm.interval(0.68, loc=mean, scale=sigma)

print('{:0.2%} of the single draws are in conf_int_a'
      .format(((a >= conf_int_a[0]) & (a < conf_int_a[1])).sum() / float(N)))

M = 1000
b = np.random.normal(0, 1, (N, M)).mean(axis=1)
conf_int_b = stats.norm.interval(0.68, loc=0, scale=1 / np.sqrt(M))
print('{:0.2%} of the means are in conf_int_b'
      .format(((b >= conf_int_b[0]) & (b < conf_int_b[1])).sum() / float(N)))

印刷品

68.03% of the single draws are in conf_int_a
67.78% of the means are in conf_int_b

请注意，如果您使用 mean 和 sigma 的估计值定义 conf_int_b基于样本 a，平均值可能不会落入 conf_int_b 与所需的频率.

Beware that if you define conf_int_b with the estimates for mean and sigma based on the sample a, the mean may not fall in conf_int_b with the desired frequency.

如果你从分布中抽取一个样本并计算样本均值和标准差，

If you take a sample from a distribution and compute the sample mean and std deviation,

mean, sigma = a.mean(), a.std()

请注意，不能保证这些会等于人口均值和标准差，我们假设人口是正态分布的——这些不是自动给定的！

be careful to note that there is no guarantee that these will equal the population mean and standard deviation and that we are assuming the population is normally distributed -- those are not automatic givens!

如果您抽取样本并希望估计总体均值和标准偏差，你应该使用

If you take a sample and want to estimate the population mean and standard deviation, you should use

mean, sigma = a.mean(), a.std(ddof=1)

因为 sigma 的这个值是总体标准差的无偏估计.

since this value for sigma is the unbiased estimator for the population standard deviation.

这篇关于用scipy获得置信区间的正确方法的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！