正确的方法来获得置信区间
我有一个一维数据数组:
a = np.array([1,2,3,4,4,4,5,5,5,5,4,4,4,6,7,8])
为此,我想获得68%的置信区间(即1
sigma)。
此答案中的第一条评论指出,可以使用scipy.stats.norm函数通过以下方式实现scipy.stats.norm.interval
此目的:
from scipy import statsimport numpy as np
mean, sigma = np.mean(a), np.std(a)
conf_int = stats.norm.interval(0.68, loc=mean,
scale=sigma)
但是这篇文章中的评论指出,获得置信区间的 实际 正确方法是:
conf_int = stats.norm.interval(0.68, loc=mean, scale=sigma / np.sqrt(len(a)))
也就是说,sigma除以样本大小的平方根:np.sqrt(len(a))
。
问题是:哪个版本正确?
回答:
从具有正mu和std偏差sigma的正态分布 的68%置信区间为
stats.norm.interval(0.68, loc=mu, scale=sigma)
的68%置信区间 从均值mu和std偏差sigma为正态分布
stats.norm.interval(0.68, loc=mu, scale=sigma/sqrt(N))
从直觉上讲,这些公式很有意义,因为如果您举起一罐果冻豆并要求大量人猜测果冻豆的数量,则每个人的偏差可能都很大-
相同的标准偏差sigma
-但猜测的平均值将在估计实际数量方面做得非常好,这可以通过平均收缩量的标准偏差来反映1/sqrt(N)
。
如果单个抽签具有方差sigma**2
,则根据Bienaymé公式,N
不相关 抽签的总和具有方差N*sigma**2
。
平均值等于总和除以N。将随机变量(如总和)乘以常数后,方差乘以常数的平方。那是
Var(cX) = c**2 * Var(X)
所以均值的方差等于
(variance of the sum)/N**2 = N * sigma**2 / N**2 = sigma**2 / N
因此平均值的标准偏差(即方差的平方根)等于
sigma/sqrt(N).
这是sqrt(N)
分母中的原点。
这是一些基于Tom的代码的示例代码,该代码演示了上述声明:
import numpy as npfrom 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_a67.78% of the means are in conf_int_b
请注意,如果您conf_int_b
为 样本定义mean
并sigma
基于样本进行估计a
,则平均值可能不conf_int_b
符合所需的频率。
如果您从分布中获取 样本 并计算样本均值和标准差,
mean, sigma = a.mean(), a.std()
请注意,不能保证这些将等于 总体 均值和标准差,并且我们 假设 总体是正态分布的-这些不是自动给出的!
如果您抽样并想要 估计 总体平均值和标准偏差,则应使用
mean, sigma = a.mean(), a.std(ddof=1)
因为此sigma值是总体标准偏差的无偏估计量。
以上是 正确的方法来获得置信区间 的全部内容, 来源链接: utcz.com/qa/425358.html