[R] Problem with Integrate for NEF-HS distribution

[xia wang]

I need to calcuate the cumulative probability for the Natural Exponential Family - Hyperbolic secant distribution with a parameter theta between -pi/2  and pi/2. The integration should be between 0 and 1 as it is a probability. 

The function "integrate" works fine when the absolute value of theta is not too large.  That is, the NEF-HS distribution is not too skewed.  However, once the theta gets large in absolute value, such as -1 as shown in the example below, "integrate" keeps give me error message for "non-finite function" when I put the lower bound as -Inf.  I suspect that it is caused by the very heavy tail of the distribution. 

Is there any way that I can get around of this and let "integrate" work for the skewed distribution?  Or is there any other function for integrating in R-package? Thanks a lot for your advice in advance!


> theta<--1
> sech<-function(X) 2/(exp(-X)+exp(X))
> integrand <- function(X) {0.5*cos(theta)*exp(X*theta)*sech(pi*X/2)}

> integrate(integrand, -3,1)
0.8134389 with absolute error < 7e-09
> integrate(integrand, -10,1)
0.9810894 with absolute error < 5.9e-06
> integrate(integrand, -15,1)
0.9840505 with absolute error < 7e-09
> integrate(integrand, -50,1)
0.9842315 with absolute error < 4.4e-05
> integrate(integrand, -100,1)
0.9842315 with absolute error < 3.2e-05
> integrate(integrand, -Inf,1)
Error in integrate(integrand, -Inf, 1) : non-finite function value


Xia
_________________________________________________________________
Be the filmmaker you always wanted to be—learn how to burn a DVD with Windows®.