Hi Ravi - Thank you. I'm wondering how the magic happens here with options(digits=10)? And the most important point, I do not quite understand why divergent error could occur. In this case, the function is analytically integrable. Any insight would be greatly appreciated so that I learn something out of this. Thank you.
Rob On Wed, Aug 11, 2010 at 11:46 PM, Ravi Varadhan <rvarad...@jhmi.edu> wrote: > Try this: > > options(digits=10) > >> integrate(f=powerLaw2, lower=0, upper=Inf, l1=1.8980185, l2=-0.0804259, >> c0=1, t0=259.78, rel.tol=1.e-10) > 0.01089019946 with absolute error < 3.7e-11 > > > Ravi. > ____________________________________________________________________ > > Ravi Varadhan, Ph.D. > Assistant Professor, > Division of Geriatric Medicine and Gerontology > School of Medicine > Johns Hopkins University > > Ph. (410) 502-2619 > email: rvarad...@jhmi.edu > > > ----- Original Message ----- > From: R_help Help <rhelp...@gmail.com> > Date: Wednesday, August 11, 2010 9:44 pm > Subject: [R] Need help to understand integrate function > To: r-help@r-project.org > > >> Hi, >> >> I'm running into a wall when trying to use the integrate function. I >> have the following setting: >> >> powerLaw2 <- function(x,l1,l2,c0,t0) { >> >> idx <- which(x <= 0); >> if (length(idx) > 0) { >> x[idx] <- 0; >> } >> >> xl <- (-l1+l2)*log(x/t0); >> L <- log(c0)-l1*log(x)-log(1+exp(xl)); >> L <- exp(L); >> return(L); >> } >> >> plCDF2 <- function(x,l1,l2,c0,t0) { >> >> print(c(l1,l2,c0,t0)); >> cdf <- >> integrate(f=powerLaw2,lower=x,upper=Inf,l1=l1,l2=l2,c0=c0,t0=t0,stop.on.error=T); >> return(cdf$value); >> >> } >> >> I know that as long as my l1 > 1 and l2 < 1 the integration above >> should converge. However, this doesn't seem to be the case when I call >> integrate function with the following parameter: >> >> >> Browse[1]> >> integrate(f=powerLaw2,lower=1e-09,upper=Inf,l1=1.8980185,l2=-0.0804259,c0=1,t0=259.78,stop.on.error=T) >> failed with message 'the integral is probably divergent' >> >> Would you please shed some light? And how could I prevent this >> situation? Thank you. >> >> Rgds, >> >> Rob >> >> ______________________________________________ >> r-h...@r-project.org mailing list >> >> PLEASE do read the posting guide >> and provide commented, minimal, self-contained, reproducible code. > ______________________________________________ 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.