Интеграция Монте-Карло - как найти ошибку

Мне нужно применить интегрирование Монте-Карло, чтобы оценить следующую функцию в [0,1]:

f(x) = exp(-ax) cos(bx)

где a=0.3060734 и b=0.11221230. Но мне нужно использовать четыре разных варианта Монте-Карло: удачный или промах, грубый, выборка по важности и контроль дисперсии. Все они находятся на странице 392 из этой книги: https://edisciplinas.usp.br/pluginfile.php/5168099/mod_resource/content/1/Julio%20Stern.pdf. Однако после оценки интеграла f (x) в каждом варианте мне нужно вычислить относительную ошибку ( | g* - g | / g ) < 1%, где g - действительное значение интеграла (неизвестно), а g* - оценочное значение. Как рассчитать погрешность? Я думал об использовании дисперсии каждого варианта, но я не уверен, как это сделать и сделать это менее 1%. У меня уже есть код, рассчитывающий каждую оценку и дисперсию. OBS: Я должен использовать Python или R.


python r statistics montecarlo
person carolina131    schedule 24.03.2020    source источник

Ответы (1)


arrow_upward
0
arrow_downward

Вот пример (с R) для моделирования Монте-Карло для интеграции

a <- 0.3060734
b <- 0.11221230
f <- function(x) exp(-a*x)*cos(b*x)

MCsim <- function(n) {
  gs <- sapply(n,function(k) mean(runif(k) <= f(runif(k))))
  g <- integrate(f,0,1)$value
  rel_err <- abs(gs-g)/g
  data.frame(n,gs,g,rel_err)
}

n <- 10**(1:7)
res <- MCsim(n)

такой, что

> res
      n        gs         g      rel_err
1 1e+01 0.8000000 0.8597815 6.953106e-02
2 1e+02 0.8700000 0.8597815 1.188497e-02
3 1e+03 0.8590000 0.8597815 9.089768e-04
4 1e+04 0.8570000 0.8597815 3.235149e-03
5 1e+05 0.8606900 0.8597815 1.056639e-03
6 1e+06 0.8596650 0.8597815 1.355245e-04
7 1e+07 0.8597802 0.8597815 1.536980e-06
person ThomasIsCoding    schedule 24.03.2020
comment
Спасибо за ваш ответ! Я понял, что вы сделали, но в моем случае я не могу использовать реальное значение og, которое вы получили в rel_err ‹- abs (gs-g) / g, мне нужно использовать дисперсию четырех разных вариантов метода чтобы вычислить ошибку, но я не знаю как. - person carolina131; 25.03.2020