Take the 2-minute tour ×
Mathematics Stack Exchange is a question and answer site for people studying math at any level and professionals in related fields. It's 100% free, no registration required.

I would like to know how to compute a probability function of a convolution of Negative Binomial distribution with Maple.

Here is an easy example of what I want to do : '

With(Statistics):

X[1]:=RandomVariable(NegativeBinomial(2,0.5)):
X[2]:=RandomVariable(NegativeBinomial(6,0.3)):

S:=X[1]+X[2]:

ProbabilityFunction(S,0);

If I ask for the Mean, it works fine with Mean(S) but when I ask for the Probability Function, Maple gives me "FAIL" as answer.

Thank you!

Jean-Philippe

share|improve this question
    
The PDF of the sum of random variables will be the convolution of the individual PDF's. Hence, take the Fourier Transform (FT) of the individual PDF's, multiply the FT's, and take the inverse FT of the product to obtain the convolution. –  Rod Carvalho Sep 16 '12 at 18:26
    
I want to generalize that situation... I want to be able to change de distribution easily... I'm sure there is a way to do it like I want but I don't know how... –  Jean-Philippe Le Cavalier Sep 16 '12 at 18:32

1 Answer 1

I am not really sure how to compute this with Maple, but here is how to get it in closed form.

Recall that the negative binomial distribution with parameters $n$ and $p$ has the following simple probability generating function: $$ \mathcal{P}_X\left(z\right) = \left( \frac{p}{1-(1-p)z} \right)^n $$ And that the probability generating function for the sum of two independent random variables $X_1+X_2$ is a product of individual generating functions: $$ \mathcal{P}_{X_1+X_2}(z) = \left( \frac{p_1}{1-(1-p_1)z} \right)^{n_1} \left( \frac{p_2}{1-(1-p_2)z} \right)^{n_2} $$ The probability mass function for $X_1+X_2$ is then read off as a series coefficient of PGF: $$ \mathbb{P}\left(X_1+X_2 = n\right) = [z^n] \mathcal{P}_{X_1+X_2}(z) $$ In particular: $$ \mathbb{P}\left(X_1+X_2=0\right) = p_1^{n_1} p_2^{n_2} = \mathbb{P}\left(X_1=0\right) \mathbb{P}\left(X_2=0\right) $$

By the way, here is what Mathematica gives for your problem:

In[14]:= Simplify[
 PDF[TransformedDistribution[
   x1 + x2, {Distributed[x1, NegativeBinomialDistribution[2, 1/2]], 
         Distributed[x2, NegativeBinomialDistribution[6, 3/10]]}], k]]

Out[14]= Piecewise[{{243*2^(-14 - k)*5^(-7 - k)*(66*5^(7 + k) - 
      15030*7^(3 + k) + 3*5^(7 + k)*k + 30673*7^(2 + k)*k - 
      4300*7^(2 + k)*k^2 + 60*7^(3 + k)*k^3 - 
            20*7^(2 + k)*k^4 + 2*7^(2 + k)*k^5), k >= 0}}, 0]

Compare with the series coefficient:

In[18]:= With[{p1 = 1/2, p2 = 3/10}, 
 Simplify[SeriesCoefficient[(p1^2*
      p2^6)/((1 + (-1 + p1)*z)^2*(1 + (-1 + p2)*z)^6), {z, 0, k}]]]

Out[18]= Piecewise[{{243*2^(-14 - k)*5^(-7 - k)*(66*5^(7 + k) - 
      15030*7^(3 + k) + 3*5^(7 + k)*k + 30673*7^(2 + k)*k - 
      4300*7^(2 + k)*k^2 + 60*7^(3 + k)*k^3 - 
            20*7^(2 + k)*k^4 + 2*7^(2 + k)*k^5), k >= 0}}, 0]
share|improve this answer
    
Thanks for your answer. In facts, I have a convolution of 500 negative binomial distributions... I don't use Mathematica, there is probably a way to do it with Maple... –  Jean-Philippe Le Cavalier Sep 16 '12 at 19:31

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

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