在R中处理非常小的数字 [英] Dealing with very small numbers in R

问题描述

我需要计算一个非常小的数字列表，例如

(0.1)^1000，0.2^(1200)，

然后将其规格化，使其合计为1 即

A1=0.1^1000， A2=0.2^1200

我想计算一下 A1‘=a1/(a1+a2)， A2‘=a2(a1+a2)。

我遇到了下溢问题，因为我得到A1=0。我怎么才能绕过这个问题呢？ 从理论上讲，我可以处理日志，然后log(A1)=1000*log(0.l)将是一种表示a1而不会出现下溢问题的方法--但为了标准化，我需要 Log(a1+a2)-我无法计算，因为我不能直接表示a1。

我正在用R编程--据我所知，在c#中没有Decimal这样的数据类型 让你得到比双精度更好的值。

如有任何建议，我们将不胜感激，谢谢

推荐答案

从数学角度讲，其中一个数字将是APPX。零，另一个。你们之间的数字差异很大，所以我甚至想知道这是否有意义。

一般来说，要做到这一点，您可以使用位于R底层的logspace_addC函数中的思想。您可以将logxpy ( =log(x+y) )Whenlx = log(x)和ly = log(y)定义为：

logxpy <- function(lx,ly) max(lx,ly) + log1p(exp(-abs(lx-ly)))

这意味着我们可以使用：

> la1 <- 1000*log(0.1)
> la2 <- 1200*log(0.2)

> exp(la1 - logxpy(la1,la2))
[1] 5.807714e-162

> exp(la2 - logxpy(la1,la2))
[1] 1

如果您有更多的数字，也可以递归调用此函数。请注意，1仍然是1，而不是1减5.807...e-162。如果你真的需要更高的精度，并且你的平台支持LONG DOUBLE类型，你可以用C或C++编写所有代码，然后在稍后返回结果。但如果我是对的，R目前只能处理正常的双打，所以最终当结果显示时，你会再次失去精度。

---

编辑：

为您计算：

log(x+y) = log(exp(lx)+exp(ly))
         = log( exp(lx) * (1 + exp(ly-lx) )
         = lx + log ( 1 + exp(ly - lx)  )

现在您只需将最大的作为lx，然后您就会得到logxpy()中的表达式。

编辑2：为什么要取最大值？很容易，以确保您在EXP(lx-ly)中使用负数。如果lx-ly变得太大，则exp(lx-ly)将返回inf。这不是一个正确的结果。Exp(ly-lx)将返回0，这将产生更好的结果：

假设lx=1且ly=1000，则：

> 1+log1p(exp(1000-1))
[1] Inf
> 1000+log1p(exp(1-1000))
[1] 1000

这篇关于在R中处理非常小的数字的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！