<div dir="ltr">Apenas uma correção... a velocidade da função de integração na minha máquina é 0.02 segundos e combinado com a de otimização é 0.2 segundos.</div><div class="gmail_extra"><br><br><div class="gmail_quote">Em 13 de agosto de 2014 22:27, Alexandre Castagna <span dir="ltr"><<a href="mailto:alexandre.castagna@gmail.com" target="_blank">alexandre.castagna@gmail.com</a>></span> escreveu:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr">Prezados,<div><br></div><div>Preciso aplicar repetidas vezes uma integral em um distribuição de probabilidade e por isso preciso de um método muito rápido (são milhares de repetições). Porém não existe uma função matemática que descreva a distribuição, os valores estão tabulados.</div>

<div><br></div><div>A função de integração que construí é muito lenta (0.2 segundos na minha máquina); queria reduzir essa velocidade em 10 ou 100 vezes. Apresento a função com os valores tabulados:</div><div><br></div><div>

<div>    angles <- log(c(0.001, 0.01, 0.1, 0.12589, 0.15849, 0.19953, 0.25119, 0.31623, 0.39811, </div><div>    0.50119, 0.63096, 0.79433, 1, 1.2589, 1.5849, 1.9953, 2.5119, 3.1623, </div><div>    3.9811, 5.0119, 6.3096, 7.9433, 10, 15, 20, 25, 30, 35, 40, 45, 50, 55, </div>

<div>    60, 65, 70, 75, 80, 85, 90, 95, 100, 105, 110, 115, 120, 125, 130, 135, </div><div>    140, 145, 150, 155, 160, 165, 170, 175, 180)*pi/180)</div><div><br></div><div>    mc.petold <- log(c(869105, 39180, 1.76661e3, 1.29564e3, 9.50172e2, 6.99092e2, 5.13687e2, </div>

<div>    3.76373e2, 2.76318e2, 2.18839e2, 1.44369e2, 1.02241e2, 7.16082e1, </div><div>    4.95803e1, 3.39511e1, 2.28129e1, 1.51622e1, 1.00154e1, 6.57957, 4.29530, </div><div>    2.80690, 1.81927, 1.15257, 4.89344e-1, 2.44424e-1, 1.47151e-1, 8.60848e-2, </div>

<div>    5.93075e-2, 4.20985e-2, 3.06722e-2, 2.27533e-2, 1.69904e-2, 1.31254e-2, </div><div>    1.04625e-2, 8.48826e-3, 6.97601e-3, 5.84232e-3, 4.95306e-3, 4.29232e-3, </div><div>    3.78161e-3, 3.40405e-3, 3.11591e-3, 2.91222e-3, 2.79696e-3, 2.68568e-3, </div>

<div>    2.57142e-3, 2.47603e-3, 2.37667e-3, 2.32898e-3, 2.31308e-3, 2.36475e-3, </div><div>    2.50584e-3, 2.66183e-3, 2.83472e-3, 3.03046e-3, 3.09206e-3, 3.15366e-3))</div><div><br></div><div>cum.beta.ph.p.petzold <- function(psi)</div>

<div>   {</div><div>    Integral <- integrate(function(x) { exp(approx(angles, </div><div>                 mc.petold, xout = log(x))$y)*sin(x) }, </div><div>                 lower = 1.745329e-05, upper = psi)$value*2*pi</div>

<div>    return(Integral)</div><div>   }</div></div><div><br></div><div>cum.beta.ph.p.petzold(psi = pi)<br></div><div><br></div><div>Onde psi deve ser entrado em radianos. A integral de "0" a pi é apenas aproximadamente 1, por conta de uma simplificação dos valores e pelo limite inferior de integração em 0.001º (1.745329e-05 rad), já que o valor é infinito em 0. A normalização para 1 é fácil de ser realizada.</div>

<div><br></div><div>A função é aplicada milhares de vezes em um programa de simulação com números aleatórios, onde para tirar uma amostra aleatória dessa distribuição, eu uso:</div><div><br></div><div>randon.n <- runif(1)<br>

</div><div>optimize( function(x) abs(cum.beta.ph.p.petzold(x)-randon.n), interval = c(0, pi), tol = 1e-3)$minimum</div><div> </div><div>Agradeço se alguém puder oferecer uma sugestão.</div><div><br></div><div>Obrigado,</div>
<span class="HOEnZb"><font color="#888888">
<div>-- <br>Alexandre Castagna Mourão e Lima<br>Laboratório de Radioecologia e Mudanças Globais - UERJ<br>Rua São Francisco Xavier, 524, Maracanã<br>Pav. Haroldo Lisboa da Cunha, Sub-solo
</div></font></span></div>
</blockquote></div><br><br clear="all"><div><br></div>-- <br>Alexandre Castagna Mourão e Lima<br>Laboratório de Radioecologia e Mudanças Globais - UERJ<br>Rua São Francisco Xavier, 524, Maracanã<br>Pav. Haroldo Lisboa da Cunha, Sub-solo
</div>