How should I derive a PDF from the difference of a couple of given two parameter Weibull distributions ? The convolution integrals become very difficult to evaluate. PS : The parameters are not necessarily the same since the idea is to get a distribution for a general case.

  • $\begingroup$ please add more information to clarify your question $\endgroup$
    – Antoine
    Commented Jun 15, 2016 at 11:49
  • $\begingroup$ Do they have the same values of their parameters or not? $\endgroup$
    – whuber
    Commented Jun 15, 2016 at 12:13
  • $\begingroup$ Since we're considering a general case, it'll be better to assume that they do not have the same parameters. However, I don' think it'll make much difference since the convolution integral is equally difficult to compute. $\endgroup$ Commented Jun 15, 2016 at 12:20
  • $\begingroup$ It makes a profound difference! When either both shape parameters are equal and equal to certain special values (such as $1$ or $2$), or both scale parameters are the same, a closed form answer is possible. Otherwise it does look difficult. $\endgroup$
    – whuber
    Commented Jun 15, 2016 at 12:37
  • $\begingroup$ I am sorry. Your argument sounds correct. It looks tough to get a closed form solution in other cases. $\endgroup$ Commented Jun 15, 2016 at 18:28

1 Answer 1


This seems to be a difficult question. I tried to find a closed form answer with the help of maple, but maple gives up! So I will use the same numerical approach as used in Difference of two i.i.d. lognormal random variables. The convolution formulas used are developed there, so will not be redone here. If $X, Y$ are independent weibull distributed (possibly with different parameters) then the density of $D=X-Y$ is $$ f_D(t) = \int_0^\infty f(t+y) g(y) \; dy $$ with $f$ the density of $X$, $g$ the density of $Y$.

The density of the difference $D$ can then be calculated via numerical integration:

dDIFF <-  function(x, shape1, shape2=shape1, scale1=1, scale2=scale1) {
    res  <- x
    for (xx in seq(along=x)) {
        res[xx]  <-  integrate(function(y) dweibull(y+x[xx], shape1, scale1) * dweibull(y, shape2, scale2),  lower=0.0,  upper=+Inf)$value

and we can make a plot of that, first in a symmetrical case:

plot of density by numerical integration, symmetric case

(note the cusp at zero), and then an asymmetrical case:

plot of density by numerical integration, asymmetrical case

This numerical solution is very fast and should be enough for most practical purposes! If you have to evaluate it very often, maybe make a spline approximation?


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.