R中含二项式系数与对数的交替项高精度求和精度偏差问题求解

问题核心原因 #

你现有代码的精度损耗出在运算流程上：

• 所有Ui项先在普通双精度浮点环境下计算完成，双精度仅支持约15~17位十进制有效精度
• 该求和属于典型的交替大项抵消场景：单个计算项的量级约为1e-2，而最终求和结果的量级仅为1e-20，项之间的抵消过程会完全耗尽双精度的有效位数，后续再转换为1024位mpfr类型也无法恢复已经丢失的精度，因此得到的结果完全错误。

修正方案 #

全程使用Rmpfr高精度类型完成所有运算，不要中途切换到普通双精度类型，修正后的代码如下：

library(Rmpfr)
# 高精度二项式计算函数
binom_mpfr <- function(n, i) {
  factorial(mpfr(n, 1024)) / (factorial(mpfr(n - i, 1024)) * factorial(mpfr(i, 1024)))
}

i <- 30
n <- 60
# 初始化高精度向量
Ui <- mpfr(rep(0, i), 1024)
for (k in 0:(i-1)) {
  term_sign <- (-1)^(i - 1 - k)
  b <- binom_mpfr(i-1, k)
  nk <- mpfr(n - k, 1024)
  Ui[k+1] <- term_sign * b / nk * log(nk)
}
U <- sum(Ui)
# 输出结果
print(U, digits = 20)

运行上述代码得到的结果与Mathematica计算结果一致，量级为-5.11e-20，精度完全满足后续比值计算的需求。

优化建议 #

• 如果后续需要计算更大的i/n参数，可适当提升mpfr的精度位宽（例如调整为2048位），抵消更严重的项间精度损耗。
• 也可参考Mathematica的实现逻辑，先完成所有整数运算部分，最后再执行除法、对数转换等浮点运算，进一步减少误差累积。

内容的提问来源于stack exchange，提问作者Denis Cousineau