I'm trying to calculate the convolution of two probability density functions (PDFs) defined on the real axis:
a normal distribution (with \[Sigma] > 0 and \[Mu] a real number):
f[x_, \[Mu]_, \[Sigma]_] := 1/(Sqrt[2 \[Pi]] \[Sigma]) Exp[-(1/2) ((x - \[Mu])/\[Sigma])^2] and a Breit-Wigner distribution:
g[x_, m_, \[CapitalGamma]_] := 1 / ((x - m)^2 + 1/4 * \[CapitalGamma]^2) * \[CapitalGamma] / (2 * Pi) with \[CapitalGamma] > 0 and m > 0.
I checked that for both of them the integral from minus infinity to infinity is one.
I'm interested in calculating the convolution of them which I do with Fourier and inverse Fourier transforms:
conv = Assuming[{\[CapitalGamma] > 0 && m > 0 && \[Sigma] > 0 && \[Mu] \[Element] Reals}, InverseFourierTransform[ FourierTransform[f[x, \[Mu], \[Sigma]], x, w] * FourierTransform[g[x, m, \[CapitalGamma]], x, w], w, x] ] This gives me a closed form solution:
(E^(-((2 m - 2 x - I \[CapitalGamma] + 2 \[Mu])^2/( 8 \[Sigma]^2))) (1 - I Erfi[(2 m - 2 x - I \[CapitalGamma] + 2 \[Mu])/( 2 Sqrt[2] \[Sigma])]) + E^(-((2 m - 2 x + I \[CapitalGamma] + 2 \[Mu])^2/( 8 \[Sigma]^2))) (1 + I Erfi[(2 m - 2 x + I \[CapitalGamma] + 2 \[Mu])/( 2 Sqrt[2] \[Sigma])]))/(4 \[Pi] \[Sigma]) I expect the convolution product to have unit area again. I managed to manipulate the analytic expression for the indefinite integral of the convolution product but I find that after taking the limits I get 1/Sqrt[2*Pi] ($\approx 0.398942$) rather than one. A numerical example with bounds chosen quite large:
N[ NIntegrate[ Re[conv /. { \[CapitalGamma] -> 2.4952, m -> 91.1876, \[Mu] -> 0, \[Sigma] -> 2}] , {x, -1000, 1000}] ] also gives a numerical value (0.398623) close to 1/Sqrt[2*Pi].
Is this expected or did I do a mistake ?
Checks:
I checked that when Fourier and inverse Fourier transforming the original PDFs f and g I get the original PDFs f and g without any additional factor.
Direct convolution: Using the Convolve[] function I need to do some variable substitution to get a closed form result:
convStd = Assuming[{\[CapitalGamma] > 0 && m > 0 && \[Sigma] > 0 && \[Mu] \[Element] Reals && x \[Element] Reals && y \[Element] Reals}, Convolve[f[x, \[Mu], \[Sigma]] /. x -> u + \[Mu], g[x, m, \[CapitalGamma]] /. x -> u + \[Mu], u, y] ] /. y -> x While I did not manage to get a closed form expression for the indefinite integral, a numerical example:
N[ NIntegrate[ Re[convStd /. { \[CapitalGamma] -> 2.4952, m -> 91.1876, \[Mu] -> 0, \[Sigma] -> 2}] , {x, -1000, 1000}] ] shows that the area of the convolution product is close to one (0.999199)