Интегрируйте функцию, возвращающую ошибку округления после предыдущей работы - PullRequest
Новая
       22

Интегрируйте функцию, возвращающую ошибку округления после предыдущей работы

1 голос
/ 30 мая 2019

При использовании интегрирования для интегрирования логнормальной функции плотности от 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

Мне должно быть возвращено значение, нет? Я изо всех сил пытаюсь понять, насколько незначительное изменение значений приведет к тому, что функция больше не будет интегрироваться.

1 Ответ

1 голос
/ 30 мая 2019

Прежде всего, я полагаю, что вы усложняете, когда определяете подынтегральную функцию, записывая все выражение, кажется, лучше использовать встроенную функцию 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)

enter image description here

Что касается проблемы интеграции, пакет 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)
...