当上限为 Inf 时，R 中的集成失败

Question

我需要在数值上近似计算两个对数正态随机变量之和的对数方差。我想在 R 中执行此操作，这是一个示例：

################
## Setup data ##
################

sdX1 = 0.33
sdX2 = 0.70
muX1 = log(32765) - 0.5 * sdX1^2
muX2 = log(52650) - 0.5 * sdX2^2

####################################
## PDF for sum of 2 lognormal RVs ##
####################################

d2lnorm = function(z, muX1, muX2, sdX1, sdX2){

    #PDFs
    L1 = distr::Lnorm(meanlog = muX1, sdlog = sdX1)
    L2 = distr::Lnorm(meanlog = muX2, sdlog = sdX2)

    #Convlution integral
    L1plusL2 = distr::convpow(L1 + L2, 1)

    #Density function
    f.Z = distr::d(L1plusL2)

    #Evaluate
    return(f.Z(z))

}

############################################
## Expectation for sum of 2 lognormal RVs ##
############################################    

ex2lnorm = function(muX1, muX2, sdX1, sdX2){    

    #E(g(x)) = integral of g(x) * f(x) w.r.t x
    integrate(function(z) log(z) * d2lnorm(z, muX1 = muX1, muX2 = muX2, sdX1 = sdX1, sdX2 = sdX2), lower = 0, upper = +Inf)$value

}

##############
## Run code ##
##############

ex2lnorm(muX1, muX2, sdX1, sdX2)

但是，积分的计算结果为 0。如果我将 integrate 的上限更改为较大的数字，它就会起作用。但是，我是运行一堆模拟，我不能 fiddle 有上限让它每次都能工作。有没有其他方法可以使它始终如一地工作？

Answer 1

如果找到不为 0 的最大上限值，您可以做什么：

ex2lnorm = function(muX1, muX2, sdX1, sdX2){    

  f <- function(z) log(z) * d2lnorm(z, muX1 = muX1, muX2 = muX2, sdX1 = sdX1, sdX2 = sdX2)

  upper <- 2^10
  cond.lower <- FALSE
  repeat {
    new <- integrate(f, lower = 0, upper = upper)$value
    if (new > 0) cond.lower <- TRUE
    if (cond.lower && new == 0) break
    upper <- upper * 2
    prev <- new
  }

  prev
}

当上限为 Inf 时，R 中的集成失败

Integrate in R is failing when upper limit is Inf

r

integrate