6

我在 PHP 中实现毕达哥拉斯平均值,算术和几何平均值是小菜一碟,但我很难想出一个可靠的调和平均值实现。

这是WolframAlpha 的定义

WolframAlpha 的谐波平均定义


这是 PHP 中的等效实现:

function harmonicMeanV1()
{
    $result = 0;
    $arguments = func_get_args();

    foreach ($arguments as $argument)
    {
        $result += 1 / $argument;
    }

    return func_num_args() / $result;
}

现在,如果任何参数是0这将引发除以 0 警告,但由于与n -11 / n相同并且优雅地返回常量而不引发任何错误,我可以将其重写为以下内容(如果它仍然会引发错误没有论据,但现在让我们忽略它):pow(0, -1)INF

function harmonicMeanV2()
{
    $arguments = func_get_args();
    $arguments = array_map('pow', $arguments, array_fill(0, count($arguments), -1));

    return count($arguments) / array_sum($arguments);
}

两种实现在大多数情况下都可以正常工作(例如v1v2WolframAlpha),但是如果1 / n i系列的总和 为 0 ,它们会失败,我应该得到另一个除以 0 的警告,但我不...

考虑以下集合:(-2, 3, 6WolframAlpha说它是一个复数无穷大):

  1 / -2    // -0.5
+ 1 / 3     // 0.33333333333333333333333333333333
+ 1 / 6     // 0.16666666666666666666666666666667

= 0

但是,我-2.7755575615629E-17的两个实现都返回总和v1v2)而不是0.

虽然 CodePad 上的返回结果是-108086391056890000我的开发机器(32 位)说它是-1.0808639105689E+17,但它仍然不像我所期望的0那样。INF我什至尝试调用is_infinite()返回值,但它false按预期返回。

我还找到了作为 PECL 扩展stats_harmonic_mean()的一部分的函数stats,但令我惊讶的是,我得到了完全相同的错误结果:-1.0808639105689E+17如果有任何参数是00则返回但不检查系列的总和如您所见在第 3585 行

3557    /* {{{ proto float stats_harmonic_mean(array a)
3558       Returns the harmonic mean of an array of values */
3559    PHP_FUNCTION(stats_harmonic_mean)
3560    {
3561        zval *arr;
3562        double sum = 0.0;
3563        zval **entry;
3564        HashPosition pos;
3565        int elements_num;
3566    
3567        if (zend_parse_parameters(ZEND_NUM_ARGS() TSRMLS_CC, "a",  &arr) == FAILURE) {
3568            return;
3569        }
3570        if ((elements_num = zend_hash_num_elements(Z_ARRVAL_P(arr))) == 0) {
3571            php_error_docref(NULL TSRMLS_CC, E_WARNING, "The array has zero elements");
3572            RETURN_FALSE;
3573        }
3574    
3575        zend_hash_internal_pointer_reset_ex(Z_ARRVAL_P(arr), &pos);
3576        while (zend_hash_get_current_data_ex(Z_ARRVAL_P(arr), (void **)&entry, &pos) == SUCCESS) {
3577            convert_to_double_ex(entry);
3578            if (Z_DVAL_PP(entry) == 0) {
3579                RETURN_LONG(0);
3580            }
3581            sum += 1 / Z_DVAL_PP(entry);
3582            zend_hash_move_forward_ex(Z_ARRVAL_P(arr), &pos);   
3583        }
3584    
3585        RETURN_DOUBLE(elements_num / sum);
3586    }
3587    /* }}} */

这看起来像一个典型的浮动精度错误,但我无法真正理解原因,因为单个计算非常精确:

Array
(
    [0] => -0.5
    [1] => 0.33333333333333
    [2] => 0.16666666666667
)

gmp是否可以在不恢复到/bcmath扩展名的情况下解决此问题?

4

2 回答 2

4

你是对的。您找到的数字是浮点运算特性的产物。

增加更多的精度不会帮助你。你所做的只是移动球门柱。

底线是计算是以有限的精度完成的。这意味着在某些时候,中间结果将被四舍五入。然后,该中间结果不再准确。错误通过计算传播,并最终成为您的最终结果。当精确结果为零时,您通常会得到 1e-16 左右的双精度数字结果。

每当您的计算涉及分母不是 2 的幂的分数时,都会发生这种情况。

解决它的唯一方法是用整数或有理数(如果可以的话)表示计算,并使用任意精度整数包进行计算。这就是 Wolfram|Alpha 所做的。

请注意,几何平均值的计算也不是微不足道的。尝试 20 次 1e20 的序列。由于数字都相同,结果应该是 1e20。但你会发现结果是无限的。原因是这 20 个数字的乘积 (10e400) 超出了双精度浮点数的范围,因此设置为无穷大。无穷大的第 20 个根仍然是无穷大。

最后,一个元观察:毕达哥拉斯的意思是真的只对正数有意义。3 和 -3 的几何平均值是多少?是想象的吗??您链接到的维基百科页面上的不等式链只有在所有值都是正数时才有效。

于 2012-05-04T05:30:37.083 回答
3

是的,这是浮点精度的问题。-1/2 可以精确表示,但 1/3 和 1/6 不能。因此,当您将它们相加时,您不会得到零。

您可以使用您提到的使用公分母方法(您发布的 H2 和 H3 公式),但这只是将罐子踢了一点,一旦产品总和,您仍然会得到不准确的结果学期开始四舍五入。

无论如何,您为什么要采用可能为负数的数字的调和平均值?这是一个固有的不稳定计算(H(-2,3,6+epsilon) 对于非常小的 epsilon 变化很大)。

于 2012-05-04T04:22:18.140 回答