警惕 SMOOTH 函数的行为
原文链接：https://www.nv5geospatialsoftware.com/Learn/Blogs/Blog-Details/beware-the-behavior-of-smooth
16554 为本文评分：
1.5
警惕 SMOOTH 函数的行为
匿名 2013年5月30日，星期四
假设我有一个包含异常值的数组。一个非常大的异常值：
IDL> a = [1.0, 1.0, 2.0, 3.0, 4.0, 1.0d18, 4.0, 3.0, 2.0, 1.0, 1.0]
我想使用滑动平均（或称框状滤波、顶帽滤波，具体叫法取决于您学习的背景）滤波器对这个双精度数组进行平滑处理，例如使用 IDL 的 SMOOTH 函数。为了效率，SMOOTH 函数将核宽度除以相邻值之间差值的滑动总和。以下是将 SMOOTH 应用于数组 a、滤波器宽度为 3 的结果：
IDL> print, smooth(a, 3)
1.0000000 1.3333333 2.0000000 3.0000000 3.3333333e+017
3.3333333e+017 3.3333333e+017 0.00000000 -1.0000000 -1.6666667
1.0000000
哇。输入数组是对称的，为什么输出不对称？此外——这令人担忧——输入数组由正数组成，那么这些数字的任何子集的平均值怎么会是负数呢？答案在于计算机上浮点数的表示方式。为了理解原因，我们可以使用 MACHAR 函数返回的信息：
IDL> m = machar(/double)
IDL> help, m
** 结构体 DMACHAR，13 个标签，长度 72，数据长度 68:
IBETA LONG 2
IT LONG 53
IRND LONG 5
NGRD LONG 0
MACHEP LONG -52
NEGEP LONG -53
IEXP LONG 11
MINEXP LONG -1022
MAXEXP LONG 1024
EPS DOUBLE 2.2204460e-016
EPSNEG DOUBLE 1.1102230e-016
XMIN DOUBLE 2.2250739e-308
XMAX DOUBLE 1.7976931e+308
使用 IBETA 字段（给出构造数字所用的基数，2，意料之中）和 IT 字段（给出数字尾数中使用的基数为 2 的位数），两个双精度数字之间可解析的最大距离必然由下式给出：
IDL> mrd = double(m.ibeta)^m.it
IDL> print, mrd
9.0071993e+015
因此，当 SMOOTH 尝试对距离大于 MRD 的两个数字求差时，由于精度损失，可能会发生糟糕的事情。让我们将这些信息应用到上面的例子中。定义略高于和略低于 MRD 设定的阈值的值：
IDL> below = double(m.ibeta)^(m.it-1)
IDL> above = double(m.ibeta)^(m.it+1)
IDL> print, below, above
4.5035996e+015 1.8014399e+016
并将它们代入上面使用的数组中：
IDL> b = [1.0, 1.0, 2.0, 3.0, 4.0, below, 4.0, 3.0, 2.0, 1.0, 1.0]
IDL> c = [1.0, 1.0, 2.0, 3.0, 4.0, above, 4.0, 3.0, 2.0, 1.0, 1.0]
现在对这两个数组应用 SMOOTH 并评估结果：
IDL> print, smooth(b, 3)
1.0000000 1.3333333 2.0000000 3.0000000 1.5011999e+015
1.5011999e+015 1.5011999e+015 3.0000000 2.0000000 1.3333333
1.0000000
IDL> print, smooth(c, 3)
1.0000000 1.3333333 2.0000000 3.0000000 6.0047995e+015
6.0047995e+015 6.0047995e+015 3.3333333 2.3333333 1.6666667
1.0000000
注意，SMOOTH 应用于数组 b 时效果良好，但应用于 c 时则不然——结果不对称。
IDL 帮助页面中 SMOOTH 函数的 关于平滑处理大范围数据的注意事项 部分对此行为进行了更好的讨论，并提供了针对这种情况的解决方案。
更新（2013-07-01）：我忘了提，这篇博文的灵感以及帮助中的说明，都源于几年前与 DRDC Atlantic 的 Carmen Lucas 的讨论。感谢 Carmen 向我指出这个意外行为。