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.