Численное интегрирование в R. Unirroot Ошибка

Question

fun = function(lambda){
  G = function(s){exp(lambda*(s-1))-s}

  roots = uniroot.all(G,c(0,1))
  probExtinction = roots[2]
  integrand = probExtinction*exp((-9)*lambda)*lambda^(14)

  return(integrand)
}

ExpectedProb = integrate(fun,0,Inf)

Ошибка, которую я получаю, когда запускаю этот код:

Error in uniroot(f, lower = xseq[i], upper = xseq[i + 1], ...) : 
  f() values at end points not of opposite sign
In addition: Warning messages:
1: In lambda * (s - 1) :
  longer object length is not a multiple of shorter object length
2: In if (is.na(f.lower)) stop("f.lower = f(lower) is NA") :
  the condition has length > 1 and only the first element will be used
3: In if (is.na(f.upper)) stop("f.upper = f(upper) is NA") :
  the condition has length > 1 and only the first element will be used

Что странно, когда я запускаю только этот кодкод Я не получаю uniroot-error:

lambda = 2 
G = function(s){exp(lambda*(s-1))-s}
roots = uniroot.all(G,c(0,1))
probExtinction = roots[2])

Кто-нибудь знает, как я могу это исправить?

Answer 1

У меня постоянно возникала такая же проблема, как и у вас.Проблема где-то с функцией integrate (которой у меня мало опыта).Вместо этого я просто отобрал кучу гамма-переменных и оценил функцию s (лямбда) для каждой.Возьмем среднее значение и до:

library(rootSolve)

a = 15
b = 9 
n = 10000
lambda = rgamma(n, a, b)    # lambda ~ Gamma(15, 9)

# 's' is s(lambda)
s = sapply(lambda, function(l){
    G = function(s) exp(l*(s-1))-s  # solve G(s) = s
    min(uniroot.all(G, c(0, 1)))    # smallest root
    })  

mean(s)
# plot(cumsum(s) / 1:n)

Среднее будет вашей оценкой интеграла (в отличие от выполнения mean(exp(rnorm(1000))) для оценки среднего экспоненциальных нормалей).Вы можете построить совокупное среднее значение, чтобы увидеть сходимость.