[R] R numerical integration

Hans W Borchers hwborchers at googlemail.com
Sat Mar 24 08:53:00 CET 2012


casperyc <casperyc <at> hotmail.co.uk> writes:

> Is there any other packages to do numerical integration other than the
> default 'integrate'?
> Basically, I am integrating:
>
>  integrate(function(x) dnorm(x,mu,sigma)/(1+exp(-x)),-Inf,Inf)$value
>
> The integration is ok provided sigma is >0.
> However, when mu=-1.645074 and sigma=17535.26 It stopped working.
> On the other hand, Maple gives me a value of 0.5005299403.

Using `integrate()` to integrate from -1e-8 to 1e-8 will give quite a correct
result, while integrating from -1e-10 to 1e-10 will return 0.
It seems more appropriate to transform the infinite into a finite interval.
Function `quadinf` in package pracma does this automatically, applying the
`integrate` routine to the finite interval [-1, 1].

    library(pracma)
    quadinf(fun, -Inf, Inf, tol=1e-10)
    # [1] 0.4999626

I am astonished to see the result from Maple as this does not appear to be
correct --- Mathematica, for instance, returns 0.499963.

> It is an important line of the coding that I am doing and I am looking for
> some package that is able to do numerical integration efficiently (fast and
> accurate to a tol=1e-4). I have tried 'cubature', which does not give me
> anything even after 10 minutes.
>
> Thanks. casper



More information about the R-help mailing list