统计学/数值方法/统计软件数值比较

介绍

统计计算需要额外的精度，并容易受到截断和抵消误差等错误的影响。这些错误是由于二进制表示和有限精度造成的，可能会导致不准确的结果。在这项工作中，我们将讨论统计软件的准确性，以及用于测量准确性的不同测试和方法，以及不同软件包的比较。

软件的准确性

准确性可以定义为结果的正确性。在使用统计软件包时，假设结果是正确的，以便对这些结果进行评论。但是，必须承认计算机有一些限制。主要问题是计算机系统提供的可用精度通常有限。很明显，统计软件无法提供超过这些限制的准确结果。但是，统计软件应该认识到它的限制，并清楚地表明这些限制已经达到。我们今天通常使用两种类型的精度

• 单精度
• 双精度

二进制表示和有限精度

二进制表示和有限精度会导致不准确 NIKHIL。计算机没有实数 - 我们用有限的近似值来表示它们。

示例：假设我们要在单精度浮点数中表示 0.1。结果将如下所示

0.1 = .00011001100110011001100110 = 0.99999964 (McCullough, 1998)

很明显，我们只能用二进制形式近似于 0.1。如果我们试图减去两个只有最低精度位不同的数字，这个问题就会加剧。例如，100000.1-100000 = .09375

使用单精度浮点数，我们只能表示 24 个有效二进制位，换句话说，6-7 个有效十进制位。在双精度中，可以表示 53 个有效二进制位，相当于 15-17 个有效十进制位。二进制表示的限制会创建五个不同的数值范围，会导致精度损失

• 负溢出
• 负下溢
• 零
• 正下溢
• 正溢出

溢出意味着值已变得太大而无法表示。下溢意味着操作数太小，太接近零，因此结果必须设置为零。单精度和双精度表示具有不同的范围。

二进制表示的结果

这些限制会导致不同情况下的不同错误

• 抵消误差 是由于减去两个几乎相等的数字而产生的。
• 累积误差 是在一系列计算中累积起来的连续舍入误差。在这种类型的误差中，有可能只有结果的最右边几位受到影响，或者结果没有一位是准确的。
• 二进制表示和有限精度的另一个结果是，两个在代数上等效的公式在数值上可能不等于。例如

$\sum _{n=1}^{10000}n^{-2}$ $\sum _{n=1}^{10000}(10001-n)^{-2}$

第一个公式按升序添加数字，而第二个公式按降序添加数字。在第一个公式中，最小的数字是在计算结束时才达到的，因此这些数字都丢失了舍入误差。误差是第二个的 650 倍。(McCullough, 1998)

• 截断误差 是由于二进制表示的限制而产生的近似误差。

示例

$\sin x=x-{\frac {x^{3}}{3!}}+{\frac {x^{5}}{5!}}-{\frac {x^{7}}{7!}}+\cdots$

sin(x) 的真值与通过对有限项求和而获得的结果之间的差异就是截断误差。(McCullough, 1998)

• 算法误差 是另一个不准确来源。计算一个量可能有多种方法，而这些不同的方法在准确性上可能会有所不同。例如，根据 Sawitzki (1994)，在单精度环境中，使用以下公式计算方差，

$S^{2}=(1/(1-n)(\sum x_{i}^{2}-n{\bar {x}}^{2}))$

测量准确度

由于计算机的限制，在计算统计值时会出现一些问题。我们需要一个指标来衡量计算值的准确程度。这个指标基于计算值 (q) 和真实值 (c) 之间的差异。一个常用的指标是对数相对误差或 LRE（正确有效数字的位数）（McCullough，1998）

$LRE=-\log _{10}{\left[|q-c|/|c|\right]}$

规则

• q 应该接近 c（小于 2 位数）。如果它们不接近，则将 LRE 设置为零
• 如果 LRE 大于 c 中的位数，则将 LRE 设置为 c 中的位数。
• 如果 LRE 小于 1，则将其设置为零。

测试统计软件

在本部分，我们将讨论两个旨在测量软件准确性的不同测试：Wilkinson 测试（Wilkinson，1985）和 NIST StRD 基准测试。

Wilkinson 的统计测验

在 Wilkinson 的统计测验中使用的 Wilkinson 数据集“NASTY”是由 Leland Wilkinson（1985）创建的。该数据集包含不同的变量，例如仅包含零的“Zero”，包含所有缺失值的“Miss”等等。NASTY 在其包含的值方面是一个合理的数据集。例如，“NASTY”中“Big”的值小于美国人口，而“Tiny”与许多工程值相当。另一方面，“统计测验”中的练习并非旨在合理。这些测试旨在检查统计计算中的一些特定问题。Wilkinson 的统计测验是一个入门级测试。

NIST StRD 基准测试

这些基准测试由美国国家标准与技术研究院设计，包含不同难度级别的数据集。其目的是测试统计软件在统计学不同主题和不同难度级别上的准确性。在“统计参考数据集”项目的网页中，有五组数据集

• 方差分析
• 线性回归
• 马尔可夫链蒙特卡罗
• 非线性回归
• 单变量汇总统计

在所有基准测试组中，有三种不同类型的数据集：低难度级别数据集、中等难度级别数据集和高难度级别数据集。通过使用这些数据集，我们将探讨统计软件是否能够为一些统计计算提供 15 位有效数字的准确结果。

NIST 提供了 11 个数据集，其中有 6 个数据集难度级别较低，2 个数据集难度级别中等，1 个数据集难度级别较高。为每个数据集提供了均值 (μ)、标准差 (σ) 和一阶自相关系数 (ρ) 的 15 位有效数字的认证值。

在 ANOVA 数据集组中，有 11 个数据集，难度级别分别为：4 个低、4 个中等和 3 个高。为每个数据集提供了处理间自由度、处理内自由度、平方和、均方差、F 统计量、$R^{2}$ 和残差标准差的 15 位有效数字的认证值。由于大多数认证值用于计算 F 统计量，因此，我们将仅比较其 LRE $\lambda _{F}$ 与统计软件包的结果。

为了测试统计软件的线性回归结果，NIST 提供了 11 个数据集，难度级别分别为 2 个低、2 个中等和 7 个高。对于每个数据集，我们都有系数估计、系数标准误差、残差标准差、$R^{2}$、线性回归方差分析表（包含残差平方和）的 15 位有效数字的认证值。我们将比较最小准确系数的 LRE $\lambda _{\beta }$、标准误差 $\lambda _{\sigma }$ 和残差平方和 $\lambda _{r}$。在非线性回归数据集组中，NIST 设计了 27 个数据集，难度级别分别为 8 个低、11 个中等和 8 个高。对于每个数据集，我们都有 NIST 提供的系数估计、系数标准误差、残差平方和、残差标准差和自由度的 11 位有效数字的认证值。

在非线性回归计算的情况下，我们应用曲线拟合方法。在这种方法中，我们需要起始值来初始化方程中的每个变量。然后我们生成曲线并计算收敛准则（例如，平方和）。然后我们调整变量以使曲线更接近数据点。有几种调整变量的算法

• Marquardt 和 Levenberg 方法
• 线性下降法
• 高斯-牛顿法

其中一种方法被反复应用，直到收敛准则的变化小于收敛容差。

NIST 还提供了两组起始值：起始 I（远离解的值），起始 II（接近解的值）。将起始 II 作为初始值更容易获得精确的解。因此，将优先考虑起始 I 解决方案。

其他重要设置如下

• 收敛容差（例如 1E-6）
• 解法（例如高斯-牛顿法或 Levenberg-Marquardt 法）
• 收敛准则（例如残差平方和 (RSS) 或参数差的最大值的平方）

我们也可以在数值导数和解析导数之间进行选择。

测试示例

测试软件包：SAS、SPSS 和 S-Plus

在本部分中，我们将讨论 M.D. McCullough 应用的三个统计软件包的测试结果。在 McCullough 的工作中，SAS 6.12、SPSS 7.5 和 S-Plus 4.0 被测试并比较了 NIST 提供的认证 LRE 值。比较将根据以下部分进行处理

• 单变量统计
• 方差分析
• 线性回归
• 非线性回归

单变量统计

表 1：SAS 单变量统计结果（McCullough，1998）

在 SAS 中计算的所有值似乎或多或少都是准确的。对于数据集 NumAcc1，由于观察次数不足，无法计算 p 值。对于数据集 NumAcc3（平均难度）和 NumAcc 4（高难度）计算标准差似乎会给 SAS 带来压力。

表 2：SPSS 单变量统计结果（McCullough，1998）

为均值和标准差计算的所有值似乎或多或少都是准确的。对于数据集 NumAcc1，由于观察次数不足，无法计算 p 值。计算数据集 NumAcc3 和 -4 的标准差似乎也会给 SPSS 带来压力。对于 p 值，SPSS 只用 3 位小数表示结果，这低估了第一个 p 值的准确性，而高估了最后一个 p 值的准确性。

表 3：S-Plus 单变量统计结果（McCullough，1998）

为均值和标准差计算的所有值似乎或多或少都是准确的。S-Plus 在计算数据集 NumAcc3 和 -4 的标准差时也遇到了问题。S-Plus 在计算 p 值方面表现不佳。

方差分析

表 4：SAS 方差分析结果（McCullough，1998）

结果

• SAS 只能解决低难度级别的方差分析问题。
• 对于中等或更高难度的数据集，F 统计量可以以非常差的性能和 0 位数字精度计算。
• SPSS 可以像 SAS 一样，为低难度级别的数据集显示准确的结果。
• SPSS 在计算方差分析方面的性能很差。
• 对于数据集“AtmWtAg”，SPSS 没有显示任何 F 统计量，这似乎比显示 0 精度结果更合理。
• S-Plus 在处理方差分析问题方面比其他软件包更出色。
• 即使对于更高难度的数据集，S-Plus 也能显示比其他软件包更准确的结果。但对于高难度数据集的结果仍然不够准确。
• S-Plus 可以以足够的精度解决中等难度的問題。

线性回归

表 5：SAS 线性回归结果（McCullough，1998）

SAS 无法为数据集 Filip 提供解决方案，Filip 是一个十次多项式。除了 Filip 之外，SAS 显示的结果或多或少都是准确的。但对于更高难度的数据集，性能似乎下降，尤其是在计算系数方面。

表 6：SPSS 线性回归结果（McCullough，1998）

SPSS 在处理“Filip”（一个 10 次多项式）时也遇到了问题。许多软件包无法计算其值。与 SAS 一样，SPSS 为高难度数据集提供了较低的准确性。

表 7：S-Plus 线性回归结果（McCullough，1998）

S-Plus 是唯一能够为数据集“Filip”提供结果的软件包。Filip 结果的准确性似乎并不差，而是处于中等水平。即使对于更高难度的数据集，S-Plus 也能计算出比其他软件包更准确的结果。只有数据集“Wrampler4”和“-5”的系数低于平均准确度。

非线性回归

表 8：SAS 非线性回归结果（McCullough，1998）

对于非线性回归，每个软件包都测试了两种设置组合，因为不同的设置会对结果产生影响。正如我们从表格中看到的那样，在 SAS 中，首选组合产生的结果比默认组合更好。在本表格中，使用默认组合产生的结果位于括号内。由于 NIST 为认证值提供了 11 位数字，因此我们寻找 LRE 值为 11 的值。

首选组合

• 方法：高斯-牛顿法
• 准则：PARAM
• 容差：1E-6

表 9：SPSS 非线性回归结果（McCullough，1998）

同样在 SPSS 中，首选组合显示出比默认选项更好的性能。所有问题都在初始值“起始 I”下得到解决，而在 SAS 中，更高难度级别的数据集是在起始 II 值下得到解决的。

首选组合

• 方法：Levenberg-Marquardt 法
• 准则：PARAM
• 容差：1E-12

表 10：S-Plus 非线性回归结果（McCullough，1998）

正如我们从表格中看到的那样，首选组合在 S-Plus 中也比默认组合更好。除“MGH10”以外的所有问题都在初始值“起始 I”下得到解决。我们可以说，S-Plus 在计算非线性回归方面比其他软件表现更好。

首选组合

• 方法：高斯-牛顿法
• 准则：RSS
• 容差：1E-6

比较结果

所有软件包都为单变量统计中的均值和标准差提供了准确的结果。测试的统计软件包之间没有太大区别。在方差分析计算中，SAS 和 SPSS 无法通过中等难度级别的問題，而 S-Plus 提供的结果比其他软件包更准确。但对于高难度数据集，它也产生了较差的结果。关于线性回归问题，所有软件包似乎都可靠。如果我们检查所有软件包的结果，我们可以说在计算非线性回归结果方面的成功很大程度上取决于所选择的选项。

其他重要结果如下

• S-Plus 从起始 II 开始解决了一次。
• SPSS 从未将起始 II 作为初始值使用，但有一次产生了 0 精度数字。
• SAS 使用了三次起始 II，每次都产生了 0 精度数字。

不同版本的 SPSS 比较

在本部分，我们将比较旧版和新版 SPSS，以查看旧版中的问题在新版中是否已解决。在本部分，我们比较了 SPSS 7.5 版和 SPSS 12.0 版。7.5 版的 LRE 值取自 B.D. McCullough 的文章（参见参考文献）。我们还将这些测试应用于 12.0 版，并计算了相应的 LRE 值。我们从每个难度组中选择了一个数据集，并在 12.0 版中应用了单变量统计、方差分析和线性回归。数据集的来源是 NIST 统计参考数据集档案。然后，我们使用 NIST 提供的认证值计算每个数据集的 LRE 值，以便比较 SPSS 的两个版本。

单变量统计

难度：低

我们的第一个数据集是 PiDigits，其难度级别较低，由 NIST 设计，用于检测计算单变量统计值时的缺陷。

PiDigits 的认证值如下：

• 样本均值：4.53480000000000
• 样本标准差：2.86733906028871

文件：Version12pd.jpg

表 11：SPSS 12.0 版对 PiDigits 数据集的分析结果

从表 13 中可以看出，SPSS 12.0 的结果与 NIST 提供的认证值匹配。因此，我们对均值和标准差的 LRE 为 $\lambda _{\mu }$: 15，$\lambda _{\delta }$: 15。在 7.5 版中，LRE 值为 $\lambda _{\mu }$: 14.7，$\lambda _{\delta }$: 15。（McCullough，1998）

难度：中等

第二个数据集是 NumAcc3，它来自 NIST 数据集，在单变量统计方面难度中等。NumAcc3 的认证值如下：

• 样本均值：1000000.2
• 样本标准差：0.1

文件：Version12na3.jpg

表 12：SPSS 12.0 版对 NumAcc3 数据集的分析结果

在表 14 中，我们可以看到计算出的均值与 NIST 的认证值相同。因此，我们对均值的 LRE 为 $\lambda _{\mu }$: 15。但是，标准差值与认证值不同。因此，计算标准差的 LRE 如下所示：

$\lambda _{\delta }$：-log10 |0,10000000003464-0,1|/|0,1| = 9.5

SPSS 7.5 版的 LRE 为 $\lambda _{\mu }$: 15，$\lambda _{\delta }$: 9.5。（McCullough，1998）

难度：高

单变量统计中的最后一个数据集是 NumAcc4，难度级别很高。NumAcc4 的认证值如下：

• 样本均值：10000000.2
• 样本标准差：0.1

文件：Version12na4.jpg

表 13：SPSS 12.0 版对 NumAcc4 数据集的分析结果

对于这个数据集，我们也没有遇到任何计算均值方面的问题。因此，LRE 为 $\lambda _{\mu }$: 15。但是，标准差值与认证值不匹配。因此，我们应该计算标准差的 LRE 如下所示：

$\lambda _{\delta }$：-log10 |0,10000000056078-0,1|/|0,1| = 8.3

SPSS 7.5 版的 LRE 为 $\lambda _{\mu }$: 15，$\lambda _{\delta }$：8.3（McCullough，1998）

对于我们测试的这一部分，可以说 SPSS 的两个版本之间没有差异。对于中等和高难度数据集，提供的标准差结果仍然具有平均准确性。

方差分析

难度：低

我们用于测试 SPSS 12.0 关于低难度级别问题的 dataset 是 SiRstv。SiRstv 的认证 F 统计量为 1.18046237440255E+00。

文件：Version12ss.jpg

表 14：SPSS 12.0 版对 SiRstv 数据集的分析结果

• LRE：$\lambda _{F}$：-log10 | 1,18046237440224- 1,18046237440255|/ |1,18046237440255| = 12,58
• SPSS 7.5 版的 LRE：$\lambda _{F}$：9,6（McCullough，1998）

难度：中等

我们用于测试中等难度级别问题的 dataset 是 AtmWtAg。AtmWtAg 的认证 F 统计量值为 1.59467335677930E+01。

文件：Version12awa.jpg

表 15：SPSS 12.0 版对 AtmWtAg 数据集的分析结果

• LREs : $\lambda _{F}$ : -log10 | 15,9467336134506- 15,9467335677930|/| 15,9467335677930| = 8,5
• SPSS v 7.5 的 LREs : $\lambda _{F}$: 缺失

难度：高

我们使用数据集 SmnLsg07 来测试高难度问题。SmnLsg07 的认证 F 值为 2.10000000000000E+01

文件:Version12sl07.jpg

表 16: SPSS 12.0 版本数据集 SmnLsg07 的结果

• LREs : $\lambda _{F}$ : -log10 | 21,0381922055595 - 21|/| 21| = 2,7
• SPSS v 7.5 的 LREs : $\lambda _{F}$: 0

在 12.0 版本中计算的 ANOVA 结果优于 7.5 版本中计算的结果。但是，准确度仍然太低。

线性回归

难度：低

我们的低难度数据集是用于线性回归的 Norris。Norris 的认证值如下：

• 样本残差平方和 : 26.6173985294224

• 

  表 17: Norris 的系数估计值 (www.itl.nist.gov)

• LREs : $\lambda _{r}$ : 9,9 $\lambda _{\beta }$ : 12,3 $\lambda _{\sigma }$ : 10,2
• SPSS v 7.5 的 LREs : $\lambda _{r}$: 9,9 , $\lambda _{\beta }$ : 12,3 , $\lambda _{\sigma }$ : 10,2 (McCullough, 1998)

文件:Version12nor.jpg

表 18: SPSS 12.0 版本数据集 Norris 的结果

难度：中等

我们使用数据集 NoInt1 来测试平均难度数据集的性能。回归模型如下

y = B1*x + e

NoInt1 的认证值 

• 样本残差平方和 : 127,272727272727
• 系数估计值 : 2.07438016528926，标准差 : 0.16528925619834E-0(www.itl.nist.gov)
• LREs: $\lambda _{r}$:12,8 $\lambda _{\beta }$: 15 $\lambda _{\sigma }$: 12,9
• SPSS v 7.5 的 LREs : $\lambda _{r}$: 12,8 , $\lambda _{\beta }$: 14,7 , $\lambda _{\sigma }$: 12,5 (McCullough, 1998)

文件:Version12nint.jpg

表 19: SPSS 12.0 版本数据集 NoInt1 的结果

难度：高

我们的高难度数据集是 NIST 设计的 Longley。

• 模型: y =B0+B1*x1 + B2*x2 + B3*x3 + B4*x4 + B5*x5 + B6*x6 +e
• LREs 
  • $\lambda _{r}$: -log10 |836424,055505842-836424,055505915|/ |836424,055505915| = 13,1
  • $\lambda _{\beta }$ : 15
  • $\lambda _{\sigma }$ : -log10 | 0,16528925619836E-01 – 0,16528925619834E-01|/ |0,16528925619834E-01| = 12,9
• SPSS v. 7.5 的 LREs : $\lambda _{r}$: 12,8 , $\lambda _{\beta }$ : 14,7 , $\lambda _{\sigma }$ : 12,5 (McCullough, 1998)

从计算得到的 LREs 可以得出结论，两个版本在进行线性回归方面的结果并没有太大差异。

结论

通过应用这些测试，我们试图确定软件是否可靠以及是否能提供准确的结果。然而，根据结果，我们可以说不同的软件包对于同一个问题会给出不同的结果，这可能导致我们对统计研究问题的错误解释。

具体来说，我们可以看到 SAS、SPSS 和 S-Plus 在解决线性回归问题方面比 ANOVA 问题做得更好。三者在 F 统计量的计算方面都表现不佳。

从两个不同版本的 SPSS 的比较结果中，我们可以得出结论，考虑到版本号之间的差异，SPSS v.12 和 v.7.5 所提供的结果的准确性之间的差异并不大。另一方面，SPSS v.12 在处理 ANOVA 问题方面比旧版本好得多。但是，它在更高难度的难题中仍然存在问题。

参考文献

• McCullough, B.D. 1998, 'Assessing The Reliability of Ststistical Software: Part I',The American Statistician, Vol.52, No.4, pp.358-366.
• McCullough, B.D. 1999, 'Assessing The Reliability of Ststistical Software: Part II', The American Statistician, Vol.53, No.2, pp.149-159
• Sawitzki, G. 1994, 'Testing Numerical Reliability of Data Analysis Systems', Computational Statistics & Data Analysis, Vol.18, No.2, pp.269-286
• Wilkinson, L. 1993, 'Practical Guidelines for Testing Statistical Software' in 25th Conference on Statistical Computing at Schloss Reisenburg, ed. P. Dirschedl& R. Ostermnann, Physica Verlag
• National Institute of Standards and Technology. (1 September 2000). The Statistical Reference Datasets: Archives, [Online], Available from: <http://www.itl.nist.gov/div898/strd/general/dataarchive.html> [10 November 2005].