Создание более эффективного моделирования monte carlo

Итак, я написал этот код, который должен эффективно оценивать область под кривой функции, определенной как h (x). Моя проблема заключается в том, что мне нужно оценить область с точностью до 6 знаков после запятой, но алгоритм, который я определил в оценкеN, кажется, слишком тяжел для моей машины. По сути, вопрос в том, как я могу сделать следующий код более эффективным? Есть ли способ избавиться от этого цикла?Создание более эффективного моделирования monte carlo

h = function(x) { 
    return(1+(x^9)+(x^3)) 
} 
estimateN = function(n) { 
    count = 0 
    k = 1 
    xpoints = runif(n, 0, 1) 
    ypoints = runif(n, 0, 3) 
    while(k <= n){ 
    if(ypoints[k]<=h(xpoints[k])) 
     count = count+1 
    k = k+1 
    } 
    #because of the range that im using for y 
    return(3*(count/n)) 
} 
#uses the fact that err<=1/sqrt(n) to determine size of dataset 
estimate_to = function(i) { 
    n = (10^i)^2 
    print(paste(n, " repetitions: ", estimateN(n))) 
} 

estimate_to(6)

источник

2012-11-19 Khodeir

'estimate_to (6)' может быть немного слишком жадными: ваш текущий алгоритм, вероятнее всего, R запустить из памяти пытается выделить числовые векторы длины 1e12. – flodel

Если вам действительно нужны симуляции 1e12, вам придется переписать свой алгоритм, чтобы найти компромисс между временем вычислений и использованием памяти. – flodel

Лучший способ сделать интеграцию в Монте-Карло (много) более эффективной, т. Е. Получить ту же точность в меньших итерациях, - использовать выборку важности, например Metropolis Monte Carlo. В вашем случае точки с 'x' ближе к' 1.0' вносят больший вклад в значение интеграла, чем те, которые ближе к '0.0'. –

Заменить этот код:

count = 0 
k = 1 
while(k <= n){ 
if(ypoints[k]<=h(xpoints[k])) 
    count = count+1 
k = k+1 
}

С этой линией:

count <- sum(ypoints <= h(xpoints))

источник

2012-11-19 02:01:44

Я думаю, что красноречиво суммирует силу векторизации. –

Если это действительно эффективность вы стремитесь, integrate это на несколько порядков быстрее (не говоря уже о более эффективной памяти) для этой проблемы.

integrate(h, 0, 1) 

# 1.35 with absolute error < 1.5e-14 

microbenchmark(integrate(h, 0, 1), estimate_to(3), times=10) 

# Unit: microseconds 
#    expr  min   lq  median   uq  max neval 
# integrate(h, 0, 1)  14.456  17.769  42.918  54.514  83.125 10 
#  estimate_to(3) 151980.781 159830.956 162290.668 167197.742 174881.066 10

источник

2014-06-13 23:51:29 jbaums

Создание более эффективного моделирования monte carlo

ответ

Смежные вопросы