Hi, It is quite easy to integrate ((q+x)^n)*exp(-2x) over x from 0 to Inf, if you are familiar with the following basic Laplace transform results: (a) L_s[x^n] = \Gamma(n+1)/s^(n+1) (b) L_s[f(x+q)] = exp(-qs) F(s)
Using (a) and (b) and substituting s=2, you get your desired integral, I: I = exp(-2q) \Gamma(n+1) / 2^(n+1). Hope this helps, 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] -------------------------------------------------------------------------- > -----Original Message----- > From: [EMAIL PROTECTED] [mailto:r-help- > [EMAIL PROTECTED] On Behalf Of Matthias Kohl > Sent: Friday, December 23, 2005 10:09 AM > To: Bickel, David > Cc: [email protected]; Duncan Murdoch > Subject: Re: [R] convolution of the double exponential distribution > > Duncan Murdoch schrieb: > > >On 12/22/2005 7:56 PM, Bickel, David wrote: > > > > > >>Is there any R function that computes the convolution of the double > >>exponential distribution? > >> > >>If not, is there a good way to integrate ((q+x)^n)*exp(-2x) over x from > >>0 to Inf for any value of q and for any positive integer n? I need to > >>perform the integration within a function with q and n as arguments. The > >>function integrate() is giving me this message: > >> > >>"evaluation of function gave a result of wrong length" > >> > >> > > > >Under the substitution of y = q+x, that looks like a gamma integral. > >The x = 0 to Inf range translates into y = q to Inf, so you'll need an > >incomplete gamma function, such as pgamma. Be careful to get the > >constant multiplier right. > > > >Duncan Murdoch > > > >______________________________________________ > >[email protected] mailing list > >https://stat.ethz.ch/mailman/listinfo/r-help > >PLEASE do read the posting guide! http://www.R-project.org/posting- > guide.html > > > > > > Hi, > > you can use our package "distr". > > require(distr) > ## define double exponential distribution > loc <- 0 # location parameter > sca <- 1 # scale parameter > > rfun <- function(n){ loc + scale * ifelse(runif(n) > 0.5, 1, -1) * rexp(n) > } > body(rfun) <- substitute({ loc + scale * ifelse(runif(n) > 0.5, 1, -1) * > rexp(n) }, > list(loc = loc, scale = sca)) > > dfun <- function(x){ exp(-abs(x-loc)/scale)/(2*scale) } > body(dfun) <- substitute({ exp(-abs(x-loc)/scale)/(2*scale) }, list(loc > = loc, scale = sca)) > > pfun <- function(x){ 0.5*(1 + sign(x-loc)*(1-exp(-abs(x-loc)/scale))) } > body(pfun) <- substitute({ 0.5*(1 + > sign(x-loc)*(1-exp(-abs(x-loc)/scale))) }, > list(loc = loc, scale = sca)) > > qfun <- function(x){ loc - scale*sign(x-0.5)*log(1 - 2*abs(x-0.5)) } > body(qfun) <- substitute({ loc - scale*sign(x-0.5)*log(1 - 2*abs(x-0.5)) > }, > list(loc = loc, scale = sca)) > > D1 <- new("AbscontDistribution", r = rfun, d = dfun, p = pfun, q = qfun) > plot(D1) > > D2 <- D1 + D1 # convolution based on FFT > plot(D2) > > hth, > Matthias > > -- > StaMatS - Statistik + Mathematik Service > Dipl.Math.(Univ.) Matthias Kohl > www.stamats.de > > ______________________________________________ > [email protected] mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide! http://www.R-project.org/posting- > guide.html ______________________________________________ [email protected] mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
