Try the following:
sech<-function(X) 2/(exp(-X)+exp(X))
your.integrand <- function(X) {0.5*cos(theta)*exp(X*theta)*sech(pi*X/2)}
my.integrand <- function(X) {0.5 * cos(theta) * 2 / (exp(-pi*X/2 - X*theta)
+ exp(pi*X/2 - X*theta)) }
theta <- -1
my.integrand(-800)
your.integrand(-800)
Do you see the problem?
Ravi.
----------------------------------------------------------------------------
-------
Ravi Varadhan, Ph.D.
Assistant Professor, The Center on Aging and Health
Division of Geriatric Medicine and Gerontology
Johns Hopkins University
Ph: (410) 502-2619
Fax: (410) 614-9625
Email: [EMAIL PROTECTED]
Webpage: http://www.jhsph.edu/agingandhealth/People/Faculty/Varadhan.html
----------------------------------------------------------------------------
--------
_____
From: xia wang [mailto:[EMAIL PROTECTED]
Sent: Tuesday, August 26, 2008 12:59 PM
To: Ravi Varadhan
Cc: r-help@r-project.org
Subject: RE: [R] Problem with Integrate for NEF-HS distribution
Thanks so much. It works well in my MCMC sampling. May I know why
re-writing the integrand can solve the problem? I have been thinking it was
the skewness of the distribution brought the error. Thanks!
Xia
> From: [EMAIL PROTECTED]
> To: [EMAIL PROTECTED]; r-help@r-project.org
> Subject: RE: [R] Problem with Integrate for NEF-HS distribution
> Date: Tue, 26 Aug 2008 11:59:44 -0400
>
> Hi,
>
> Simply re-write your integrand as follows:
>
>
> integrand <- function(X) {0.5 * cos(theta) * 2 / (exp(-pi*X/2 - X*theta) +
> exp(pi*X/2 - X*theta)) }
>
> Now the problem goes away.
>
> > theta <- -1
> > integrate(integrand, -Inf, 1)
> 0.9842315 with absolute error < 1.2e-05
>
> This would work for all theta such that abs(theta) < -pi/2.
>
> Ravi.
>
>
>
----------------------------------------------------------------------------
> -------
>
> Ravi Varadhan, Ph.D.
>
> Assistant Professor, The Center on Aging and Health
>
> Division of Geriatric Medicine and Gerontology
>
> Johns Hopkins University
>
> Ph: (410) 502-2619
>
> Fax: (410) 614-9625
>
> Email: [EMAIL PROTECTED]
>
> Webpage: http://www.jhsph.edu/agingandhealth/People/Faculty/Varadhan.html
>
>
>
>
----------------------------------------------------------------------------
> --------
>
>
> -----Original Message-----
> From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED]
On
> Behalf Of xia wang
> Sent: Tuesday, August 26, 2008 1:01 AM
> To: r-help@r-project.org
> Subject: [R] Problem with Integrate for NEF-HS distribution
>
>
>
> 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
> WindowsR.
>
> ______________________________________________
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
_____
See what people are saying about Windows Live. Check out featured posts.
Check It Out!
<http://www.windowslive.com/connect?ocid=TXT_TAGLM_WL_connect2_082008>
[[alternative HTML version deleted]]
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
