производительность AdapteIntegrate против интеграции - PullRequest
Новая
       7

производительность AdapteIntegrate против интеграции

6 голосов
/ 28 октября 2011

Я хотел бы выполнить числовое интегрирование в одном измерении , где подынтегральное выражение имеет векторное значение .integrate() разрешает только скалярное интегрирование, поэтому мне нужно будет вызывать его несколько раз.Пакет cubature, кажется, хорошо подходит, но, похоже, работает плохо для 1D-интегралов.Рассмотрим следующий пример (скалярное подынтегральное выражение и 1D-интеграция): 

library(cubature)
integrand <- function(x, a=0.01) exp(-x^2/a^2)*cos(x)
Nmax <- 1e3
tolerance <- 1e-4

# using cubature's adaptIntegrate
time1 <- system.time(replicate(1e3, {
  a <<- adaptIntegrate(integrand, -1, 1, tol=tolerance, fDim=1, maxEval=Nmax)
}) )

# using integrate
time2 <- system.time(replicate(1e3, {
  b <<- integrate(integrand, -1, 1, rel.tol=tolerance, subdivisions=Nmax)
}) )

time1
user  system elapsed 
  2.398   0.004   2.403 
time2
user  system elapsed 
  0.204   0.004   0.208 

a$integral
> [1] 0.0177241
b$value
> [1] 0.0177241

a$functionEvaluations
> [1] 345
b$subdivisions
> [1] 10

Каким-то образом adaptIntegrate, похоже, использует гораздо больше оценок функций для аналогичной точности.Оба метода, очевидно, используют квадратуру Гаусса-Кронрода (1D случай: квадратурное правило Гаусса с 15 точками), хотя ?integrate добавляет «алгоритм Эпсилона Уинна».Это объяснит большую разницу во времени?

Я открыт для предложений об альтернативных способах работы с векторно-значными интегратами, таких как 

integrand <- function(x, a = 0.01) c(exp(-x^2/a^2), cos(x))
adaptIntegrate(integrand, -1, 1, tol=tolerance, fDim=2, maxEval=Nmax)
$integral
[1] 0.01772454 1.68294197

$error
[1] 2.034608e-08 1.868441e-14

$functionEvaluations
[1] 345

Спасибо.

1 Ответ

2 голосов
/ 23 марта 2013

В CRAN есть также пакет R2Cuba , который реализует несколько многомерных алгоритмов интеграции:

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

...