#数据分析 #统计学 #统计数据 #分析 #抽样

stats-ci

一个用于在样本数据上计算和操作置信区间的纯Rust库

13个版本

0.1.1 2024年4月7日
0.1.0 2023年6月29日
0.0.13 2023年6月28日
0.0.11 2023年5月27日
0.0.5 2023年4月27日

#99数学

Download history 2/week @ 2024-05-17 1/week @ 2024-05-24 2/week @ 2024-06-14

811 每月下载次数

MIT/Apache

220KB
2.5K SLoC

MIT license Apache 2.0 license Docs Tests Dependencies Downloads Latest crates.io

stats-ci

一个用于在样本数据上计算和操作置信区间的纯Rust库。

此包提供方便的手段,以帮助在多种情况下计算置信区间。它旨在用于分析实验数据和测量。实验数据受实验误差的影响,分析时必须考虑这种误差。在几种统计方法中,置信区间提供的信息易于解释。

动机源于个人需求,似乎没有任何crate能够提供计算此类区间既简单又全面的解决方案。一个例外是crate criterion,它可以计算其测量的置信区间,但并未导出此类功能。

此crate提供以下情况下轻松有效地计算样本数据置信区间的手段

  • mean 对均值(算术、调和、几何)的置信区间(数值数据),
  • quantile 对分位数(例如,中位数)的置信区间(任意有序数据),
  • proportion 比例的置信区间。
  • comparison 比较的置信区间(成对或不成对观察值)。

这是通过类型 Confidence 来表达置信水平,以及类型 Interval 来表示置信区间来实现的。

此crate目前不支持以下内容

  • 回归参数的置信区间。
  • 其他统计量的置信区间(例如,方差等)
  • 卡方检验

本库的文档提供了几个简单示例,展示了如何使用每个功能。

用法

本库位于 crates.io 上,可以通过简单地将 stats-ci 添加到项目 Cargo.toml 的依赖项中来使用(在 crates.io 上检查最新版本号,并替换下面的 {最新版本}

[dependencies]
stats-ci = "{ latest version }"

该库仍相对不稳定,从一个小版本到下一个版本可能会出现破坏性更改。

置信区间(概述)

统计学基于这样一个观点:实验数据是对某个理论(但尚未知的)概率分布的随机抽样。假设,如果抽样数据可以由无限数量的测量组成,则数据将精确描述该理论分布。然而,在现实情况下,可用的测量数量非常有限。这引发了如何根据有限的样本数量估计理论分布的参数(例如,平均值)以及准确度的问题。

在此背景下,置信区间为这两个问题提供了一个统计学的答案。大致来说,想法是选择一个期望的 置信水平(或者结论不太可能是虚假的几率有多高),估计的参数不再表示为一个单一值,而是表示为一个区间。该区间的宽度取决于置信水平、数据的变异性和样本数量。解释它的方法是,参数的理论值(例如,平均值,中位数)可以是该区间的任何值。95% 的置信水平意味着,如果在完全相同的条件下重复整个实验,大约 95% 的每次实验中获得的置信区间将包含参数的理论值。

考虑一个简单的实验,以下是从某个未知分布中采样的随机数据:

let data = [
    10.6, 6.6, 26.7, 0.4, 5.7, 0.3, 1.1, 5.0, 8.4, 1.4, 15.1, 0.3,
    20.4, 1.2, 28.4, 10.7, 0.4, 10.1, 4.5, 7.1, 4.3, 37.4, 0.9, 10.1,
    12.6, 21.7, 21.9, 2.0, 8.4, 9.3
];
  • 原始分布的平均值是多少?
    样本均值为 9.76667,但它与理论平均值有多接近?对这组数据的平均值计算 95% 的置信区间得到 [6.18467, 13.34866],这意味着理论平均值可以是这个区间内的任何数字(以 95% 的置信度)。所得区间很宽,因此估计不太精确。这是实验误差相当大的良好证据,并且基于均值为 9.76667 的确切值得出结论和推断是非常危险的。保持相同的置信水平,唯一减少此区间的方法是增加样本量,即进行额外的实验。

  • 原始分布的中位数是多少?
    对中位数计算 95% 的置信区间得到 [4.3, 10.6]。在这种情况下,区间的两个界限都是观察到的值。

  • 理论分布
    在这个例子中,数据实际上是从参数 λ = 0.1(平均值 = 1/λ = 10,中位数 = ln(2)/λ = 6.93147…)的指数分布中获取的。在这种情况下,我们可以验证,理论平均值和中位数确实分别包含在其各自的置信区间内。

示例

本库使得基于样本数据计算置信区间变得容易,适用于各种情况,包括平均值、分位数、比例和比较。

本库围绕一个类型 Confidence 来表达置信水平,以及一个类型 Interval 来表示置信区间。区间是通用的,可以为各种类型实例化,而不仅仅是通常的浮点或整数类型。

平均值和中位数

根据上面讨论的示例,可以对平均值或分位数(例如,中位数)进行区间计算,如下所述

// 1. import the crate
use stats_ci::*;
// 2. collect the data
let data = [
    10.6, 6.6, 26.7, 0.4, 5.7, 0.3, 1.1, 5.0, 8.4, 1.4, 15.1, 0.3,
    20.4, 1.2, 28.4, 10.7, 0.4, 10.1, 4.5, 7.1, 4.3, 37.4, 0.9,
    10.1, 12.6, 21.7, 21.9, 2.0, 8.4, 9.3
];
// 3. define the confidence level (for 95% confidence)
let confidence = Confidence::new(0.95);

// 4a. compute the interval for the arithmetic mean
if let Ok(ci) = mean::Arithmetic::ci(confidence, &data) {
    // display the interval
    println!("{}% c.i. for the mean = {}", confidence.percent(), ci);
    if ! ci.contains(&10.) {
        println!("Does NOT contains the theoretical mean!");
    }
}
// 4b. compute the interval for the median (i.e., 0.5-quantile)
if let Ok(ci) = quantile::ci(confidence, &data, 0.5) {
    // display the interval
    println!("{}% c.i. for the median = {}", confidence.percent(), ci);
    if ! ci.contains(&6.93147) {
        println!("Does NOT contains the theoretical median!");
    }
}

同样,可以按以下方式计算几何平均数和调和平均数的置信区间。

# use stats_ci::*;
let data = [ 10.6, 6.6, /* ... */ ];
let confidence = Confidence::new(0.95);
let ci = mean::Geometric::ci(confidence, &data);
let ci = mean::Harmonic::ci(confidence, &data);

比例

置信区间也可以用于计算比例。例如,当试图估计某个通信通道的损失率(即,消息丢失数与发送消息总数的比例)时,这种情况就会发生。如果知道损失的数量(在crate中标记为"success")和总消息数(标记为"population"),则可以简单地计算。

use stats_ci::*;
let confidence = Confidence::new(0.95);
let population = 10_000; // total number of sent messages
let successes = 89;      // number of lost messages
let ci = proportion::ci(confidence, population, successes).unwrap();
println!("Loss rate: {}", ci);
// > Loss rate: [0.007238518896069928, 0.010938644303608623]
//
// which means that the loss rate is estimated (95% confidence) to
// be between 0.7238% and 1.0939%.

// One-sided confidence
let confidence = Confidence::new_lower(0.95);
let ci = proportion::ci(confidence, population, successes).unwrap();
println!("Loss rate less than: {}", ci);
// > Loss rate less than: [0, 0.010583156571857643]
//
// which means that the loss rate is likely (95% confidence) to be
// less than 1.05832%.

一些便利函数使得处理数据迭代器或数组以及“成功”条件变得更加容易。

use stats_ci::*;
let data = [
    8, 11, 4, 18, 17, 9, 20, 3, 10, 14, 12, 7, 13, 16, 1, 6, 5, 2,
    15, 19
];
let confidence = Confidence::new(0.95);
let ci = proportion::ci_if(confidence, &data, |&x| x <= 10).unwrap();
println!("ci: {}", ci);
// > ci: [0.2992980081982124, 0.7007019918017876]
//
// yields the estimated proportion of numbers that are less or
// equal to 10, based on data obtained from random sampling.

比较

置信区间的一个常见用途是用于比较数据组。例如,当比较两个系统(比如系统A和系统B)时,例如确定新系统B是否比现有系统A带来改进。在这种情况下,存在两种情况

配对观察

只有当系统B的每个测量值都可以与系统A的相应测量值配对时,才会出现配对观察。这发生在实验被精心设计,使得两个系统在完全相同的条件下,针对完全相同的输入进行测量时。

use stats_ci::*;
// Zinc concentration in water samples from a river
// (from <https://online.stat.psu.edu/stat500/lesson/7/7.3/7.3.2>)
let data_bottom_water = [
    0.430, 0.266, 0.567, 0.531, 0.707, 0.716, 0.651, 0.589, 0.469,
    0.723,
];
let data_surface_water = [
    0.415, 0.238, 0.390, 0.410, 0.605, 0.609, 0.632, 0.523, 0.411,
    0.612,
];
let confidence = Confidence::new(0.95);
let ci = comparison::Paired::ci(
    confidence,
    &data_bottom_water,
    &data_surface_water
).unwrap();

置信区间是两个平均数的差异,可以按以下方式解释

  • 如果两个界限都是严格正的,则系统A的平均值高于系统B。此外,平均值的增加至少是区间下限,最多是区间上限。
  • 如果两个界限都是严格负的,则系统A的平均值低于系统B。
  • 否则,如果零包含在区间内,则系统A和系统B的平均值没有显著差异。无论值如何,在此情况下声称改进(或退化)是无效的。

非配对观察

在所有其他情况下,当观察值无法一一配对时,会出现非配对观察。例如,在上面的例子中,两个独立用户组被要求使用系统A或系统B,收集的数据用于比较。在这种情况下,测量值不能一对一匹配,并且一个系统的测量值数量可能与另一个系统不同。在这种情况下计算平均数差异的区间更为复杂,并且需要更多的测量来实现等效的精度。

使用这个crate,非配对观察与配对观察非常相似

use stats_ci::*;
// Gain in weight of 19 female rats between 28 and 84 days after birth.
// 12 were fed on a high protein diet and 7 on a low protein diet.
// (from <https://www.statsdirect.co.uk/help/parametric_methods/utt.htm>)
let data_high_protein = [
    134., 146., 104., 119., 124., 161., 107., 83., 113., 129., 97., 123.,
];
let data_low_protein = [70., 118., 101., 85., 107., 132., 94.];
let confidence = Confidence::new(0.95);
let ci = comparison::Unpaired::ci(
    confidence,
    &data_high_protein,
    &data_low_protein
).unwrap();

置信区间的解释与上述配对观察相同

  • 如果两个界限都是严格正的,则系统A的平均值比系统B高。
  • 如果两个界限都是严格负的,则系统A的平均值比系统B低。
  • 如果零包含在区间内,则系统A和系统B的平均值没有显著差异。

高级:增量统计

虽然到目前为止描述的简单接口对于大多数用例来说是足够的,但crate提供了一种稍微灵活的接口,使得在获得区间之后可以扩展样本数量。

以下示例显示了在计算算术平均数时如何进行此操作,但该方法也适用于其他类型的区间;请参阅API文档以获取更多信息和其他示例。点1到3与原始示例相同

// 1. import the crate
use stats_ci::*;
// 2. collect the data
let data = [
    10.6, 6.6, 26.7, 0.4, 5.7, 0.3, 1.1, 5.0, 8.4, 1.4, 15.1, 0.3,
    20.4, 1.2, 28.4, 10.7, 0.4, 10.1, 4.5, 7.1, 4.3, 37.4, 0.9,
    10.1, 12.6, 21.7, 21.9, 2.0, 8.4, 9.3
];
// 3. define the confidence level (for 95% confidence)
let confidence = Confidence::new(0.95);

// 4. create a statistics object
let mut stats = mean::Arithmetic::new();
// 5. add data
stats.extend(&data).unwrap();
// shortcut: combines 4. and 5.
let mut stats = mean::Arithmetic::from_iter(&data).unwrap();

// 6. compute the confidence interval over the mean for some
//    confidence level
let ci = stats.ci_mean(confidence).unwrap();

// 7. add more data
stats.extend(&[ 10.7, 9.8, /**/ ]).unwrap();
// 8. compute the new confidence interval
let ci = stats.ci_mean(confidence).unwrap();
// 9. or maybe with a different confidence level
let ci = stats.ci_mean(Confidence::new(0.99)).unwrap();
//10. and get other statistics
let mean = stats.sample_mean();
let variance = stats.sample_variance();
let std_dev = stats.sample_std_dev();
let std_err = stats.sample_sem();

请注意,只有点5和7在数据非常大时可能是成本昂贵的操作。

此接口在以下情况下很有用

  • 当你有一串数据,不想保留所有值时。
  • 当你想继续收集数据,直到你有足够的统计显著性(例如,区间宽度相对于平均值的宽度小于某个值)时。
  • 当你想继续收集数据,并在截止日期前提高准确性时。
  • 当您希望在单次数据遍历中计算多个置信水平或不同类型的置信区间时。
  • 当您需要通过并发执行来提高性能时(参见下一节)。

高级:并行执行

如果性能是问题,该库允许通过库 rayon 实现简单的并行执行。

// 1. import the crates
use stats_ci::*;
use rayon::prelude::*;
// 2. (as before:) collect the data
let data = [
    10.6, 6.6, 26.7, 0.4, 5.7, 0.3, 1.1, 5.0, 8.4, 1.4, 15.1, 0.3,
    20.4, 1.2, 28.4, 10.7, 0.4, 10.1, 4.5, 7.1, 4.3, 37.4, 0.9,
    10.1, 12.6, 21.7, 21.9, 2.0, 8.4, 9.3, /* … many more … */
];
// 3. (as before:) define the confidence level
let confidence = Confidence::new(0.95);
// 4. create and compute the statistics object
let stats = data
    .clone()
    .par_iter()
    .map(|&x| mean::Arithmetic::from_iter(&[x]).unwrap())
    .reduce(|| mean::Arithmetic::new(), |s1, s2| s1 + s2);
// 5. (as before:) compute the confidence interval
let ci = stats.ci_mean(confidence).unwrap();

注意,处理只有几千个样本时并行化几乎没有意义。

更多示例

您可以从该库的 API 文档 中找到更多信息和其他示例。

统计学 / 计算

  • 均值区间使用 Student t 分布,适用于约 100_000 个值,超过此值计算将切换到正态分布。
  • 比例和分位数区间依赖于 Wilson 分数方法,该方法比教科书通常介绍的 Wald 分数方法在统计上更稳定。
  • 该库使用补偿求和(Kahan 求和)以避免在求和非常大的数据时累积舍入误差。

库功能

该库有两个功能

  • approx (默认) 启用区间之间的近似比较。添加了对库 approx 的依赖。
  • serde 功能将库 serde 作为依赖项,并为 ConfidenceInterval 以及均值的区间增量状态提供序列化和反序列化功能。
stats-ci = { version = "{ latest version }", features = ["serde"] }

参考文献

许可

许可方式为以下之一

任选其一。

贡献

除非您明确声明,否则根据Apache-2.0许可定义,您有意提交以包含在作品中的任何贡献都应作为上述双重许可,不附加任何额外条款或条件。

贡献

我将乐意并仔细考虑您提供的任何建设性评论。特别是,我将考虑以下方面的建设性反馈:界面和计算的正确性、代码可读性、通用性和效率。

目前,以下是我待办事项列表中的内容

  • [功能] 回归参数的置信区间。
  • [统计] 审查/修复统计测试
  • [API] 减少恐慌代码
  • [重构] 重构错误结果

依赖项

~6.5MB
~124K SLoC