Hello.
I am trying to obtain probability density function (PDF) H(x) resulting from the the convolution between two probability density functions E(x) and F(x) - the objective is to obtain the PDF of the sum between two random variables. I want to do it using the Symbolic Toolbox, calculating the convolution as the product of the Laplace transforms of E(x) and F(x). Here is what I am doing.
>> syms x
>> E(x) = 2.*(1 - 0.01)*dirac(x) + rectangularPulse(0, 300, x).*0.01./300;
>> F(x) = 2.*(1 - 0.01)*dirac(x) + rectangularPulse(0, 400, x).*0.01./400;
I am checking that the integral of these functions is 1, as should be for a PDF (integrating up to 1000 is enough):
>> int(E,x,0,1000)
ans = 1
>> int(F,x,0,1000)
ans = 1
If I calculate the Laplace transform of one of these PDFs and then apply the inverse Laplace transform, I get the same function, as expected:
>> ilaplace(laplace(E,x,s),s,x)
ans = (99*dirac(x))/50 - heaviside(x - 300)/30000 + 1/30000
>> int(ilaplace(laplace(E,x,s),s,x), 0, 1000)
ans = 1
OK, now I generate the Laplace transform of the convolution between E(x) and F(x), as G(s):
>> syms s
>> LE=laplace(E,x,s)
LE = 1/(30000*s) - exp(-300*s)/(30000*s) + 99/50
>> LF=laplace(F,x,s)
LF = 1/(40000*s) - exp(-400*s)/(40000*s) + 99/50
>> G=LE*LF
G = (1/(30000*s) - exp(-300*s)/(30000*s) + 99/50)*(1/(40000*s) - exp(-400*s)/(40000*s) + 99/50)
Now, I generate the resulting PDF, by inverting the Laplace transform:
>> H=ilaplace(G,s,x)
H = x/1200000000 - (33*heaviside(x - 300))/500000 - (99*heaviside(x - 400))/2000000 + (9801*dirac(x))/2500 - (heaviside(x - 300)*(x - 300))/1200000000 - (heaviside(x - 400)*(x - 400))/1200000000 + (heaviside(x - 700)*(x - 700))/1200000000 + 231/2000000
And this is where I find the problem. H(x) is not a PDF, since its integral is greater than 1:
>> int(H, x, 0, 1000)
ans = 19999/10000
Is there someone able to tell me what I am doing wrong and how I should correct it? Thank you in advance.
