以前工作后集成函数返回舍入错误答案

【问题标题】：Integrate function returning roundoff error after working previously以前工作后集成函数返回舍入错误
【发布时间】：2019-05-31 04:12:59
【问题描述】：

当使用积分来积分 2000 -> Inf 的对数正态密度函数时，我返回一个错误。我之前使用过一个非常相似的方程，没有任何问题。

我尝试禁用错误停止，并将 rel.tol 设置得更低。我对 r 还很陌生，所以如果他们都没有做任何事情，我深表歉意。

> integrand = function(x) {(x-2000)*(1/x)*(1/(.99066*((2*pi)^.5)))*exp(-((log(x)-7.641)^2)/((2*(.99066)^2)))}
> integrate(integrand,lower=2000,upper=Inf)
1854.002 with absolute error < 0.018
#returns value fine

> integrand = function(x) {(x-2000)*(1/x)*(1/(1.6247*((2*pi)^.5)))*exp(-((log(x)-9.0167)^2)/((2*(1.6247)^2)))}
> integrate(integrand,lower=2000,upper=Inf)
Error in integrate(integrand, lower = 2000, upper = Inf) : 
  roundoff error is detected in the extrapolation table
#small change in the mu and sigma in the lognormal density function results in roundoff error

> integrate(integrand,lower=1293,upper=Inf)
29005.08 with absolute error < 2
#integrating on lower bound works fine, but having lower=1294 returns a roundoff error again
> integrate(integrand,lower=1294,upper=Inf)
Error in integrate(integrand, lower = 1294, upper = Inf) : 
  roundoff error is detected in the extrapolation table

我应该得到一个值，不是吗？我很难看出稍微改变值会导致函数不再集成。

【问题讨论】：

为什么是第一个因素(x - 2000)？
我正在寻找带有可扣除的预期值，将数据拟合到对数正态分布。期望值公式是从密度函数的免赔额到无穷大的积分乘以 x 减去免赔额。将免赔额视为平均值。类似于中心时刻，但我的值都是正数，非零。
我可能弄错了，但是当我设置所有出现的 x 减去 2000 年时，我得到的值与免赔额为 0 相同。我也尝试使用 6000 并将下限也设置为 6000结果相同。我不精通统计或保险，但我认为这不是预期的结果。
看看答案是不是你想要的。
效果很好，解决了我的问题。谢谢。

标签： r integrate

【解决方案1】：

首先，我相信您在通过写下整个表达式来定义被积函数时会变得复杂，使用内置的dlnorm 函数似乎更好。

g <- function(x, deduce, meanlog, sdlog){
  (x - deduce) * dlnorm(x, meanlog = meanlog, sdlog = sdlog)
}

curve(g(x, deduce = 2000, meanlog = 9.0167, sdlog = 1.6247), 
      from = 1294, to = 1e4)

至于集成问题，包cubature 通常在integrate 失败时会做得更好。以下所有内容均产生结果，没有错误。

library(cubature)

cubintegrate(g, lower = 1293, upper = Inf, method = "pcubature", 
             deduce = 2000, meanlog = 9.0167, sdlog = 1.6247)

cubintegrate(g, lower = 1294, upper = Inf, method = "pcubature", 
             deduce = 2000, meanlog = 9.0167, sdlog = 1.6247)

cubintegrate(g, lower = 2000, upper = Inf, method = "pcubature", 
             deduce = 2000, meanlog = 9.0167, sdlog = 1.6247)

【讨论】：