0
$\begingroup$

Hi I was trying to compute the convolution of Weibull distribution as follows with parameter t, lambda(scale parameter) and k:

f.x <- function(x,lambda,k){
  dweibull(x,k,lambda)
}

F1_t = function(t,lambda,k){
  cum(t,lambda,k)
}
F2_t = function(t,lambda,k){
  sapply(t, function(y) integrate(function(x) F1_t(t-x,lambda,k)*f.x(x,lambda,k), lower = 0, upper = y)$value)
}

F3_t = function(t,lambda,k){
  sapply(t, function(y) integrate(function(x) F2_t(t-x,lambda,k)*f.x(x,lambda,k), lower = 0, upper = y)$value)
}
F4_t = function(t,lambda,k){
  integrate(function(x) F3_t(t-x,lambda,k)*f.x(x,lambda,k),0,t,tol=.Machine$double.eps^.1,subdivisions = 200000)$value

}

Where Fn_t are the convolution definition of the n-th convolution,The first three convolutions work fine but as I evaluated F4_t with

F4_t(1,1,2)

I got the following warning from R:

Error in integrate(function(x) F2_t(t - x, lambda, k) * f.x(x, lambda,  : 
  maximum number of subdivisions reached 

I know that this is because of non-convergence of the integration function so I have tried to adjust both tolerance(tried:rel.tol=.Machine$double.eps^.1) and subdivision(tried: 20000000000) value of the integrate command but the same problem still appears, so I wonder if there is more efficient way to calculate convolution.

$\endgroup$
  • 1
    $\begingroup$ Yes: use a Fast Fourier Transform (FFT). Various posts on CV illustrate how this can be done in R, such as my answer at stats.stackexchange.com/a/49444/919. $\endgroup$ – whuber Nov 24 '19 at 17:08

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.