计算多个数的几何平均值的有效方法

我需要计算大量数字的几何平均值,其值不受先验限制。天真的方法是

double geometric_mean(std::vector<double> const&data) // failure

{

auto product = 1.0;

for(auto x:data) product *= x;

return std::pow(product,1.0/data.size());

}

但是,这很可能由于累积的下溢或溢出而失败product(注意:long double并不能真正避免此问题)。因此,下一个选择是总结对数:

double geometric_mean(std::vector<double> const&data)

{

auto sumlog = 0.0;

for(auto x:data) sum_log += std::log(x);

return std::exp(sum_log/data.size());

}

这可行,但是需要std::log()每个元素,这可能很慢。

例如,通过分别跟踪累加的指数(等于)和尾数(等于)product

回答:

“拆分指数和尾数”解决方案:

double geometric_mean(std::vector<double> const & data)

{

double m = 1.0;

long long ex = 0;

double invN = 1.0 / data.size();

for (double x : data)

{

int i;

double f1 = std::frexp(x,&i);

m*=f1;

ex+=i;

}

return std::pow( std::numeric_limits<double>::radix,ex * invN) * std::pow(m,invN);

}

如果您担心ex可能会溢出,可以将其定义为double而不是a long

long,并invN在每一步乘以,但是使用这种方法可能会损失很多精度。

对于大型输入,我们可以将计算分为几个存储区:

double geometric_mean(std::vector<double> const & data)

{

long long ex = 0;

auto do_bucket = [&data,&ex](int first,int last) -> double

{

double ans = 1.0;

for ( ;first != last;++first)

{

int i;

ans *= std::frexp(data[first],&i);

ex+=i;

}

return ans;

};

const int bucket_size = -std::log2( std::numeric_limits<double>::min() );

std::size_t buckets = data.size() / bucket_size;

double invN = 1.0 / data.size();

double m = 1.0;

for (std::size_t i = 0;i < buckets;++i)

m *= std::pow( do_bucket(i * bucket_size,(i+1) * bucket_size),invN );

m*= std::pow( do_bucket( buckets * bucket_size, data.size() ),invN );

return std::pow( std::numeric_limits<double>::radix,ex * invN ) * m;

}

以上是 计算多个数的几何平均值的有效方法 的全部内容, 来源链接: utcz.com/qa/423240.html

回到顶部