Convolution of uniform PDF and normal PDF in Matlab

1.2k Views Asked by At

I have a problem with Matlab and probability density functions:

$x$ is a random variable, uniformly distributed between $[-0.5,0.5]$, so its probability density function is $p_x(t)=\text{rect}(t)$.

$w$ is also a random variable, randomly distributed with zero-mean and variance of 1. $p_w(t)=\frac{1}{\sqrt{2\pi}}\exp(\frac{-t^2}{2})$.

$x$ and $w$ are uncorrelated.

I have to find $p_{x+w}(t)$ and I know that this is equal to the convolution between $p_x(t)$ and $p_w(t)$: $$ \int\frac{1}{\sqrt{2\pi}}\exp\Big(\frac{-(t-\tau)^2}{2}\Big)\text{rect}(\tau)d\tau = \int_{-0.5}^{0.5} \frac{1}{\sqrt{2\pi}}\exp\Big(\frac{-(t-\tau)^2}{2}\Big)d\tau \\= -\int_{t-0.5}^{t+0.5} \frac{1}{\sqrt{2\pi}}\exp\Big(\frac{-z^2}{2}\Big)dz = Q(t-0.5)-Q(t+0.5) $$

My problem occurs when I plot this result with matlab:

plot(t,qfunc(t-0.5)-qfunc(t+0.5));

It's different from plotting the convolution of the two PDF's:

plot(t,decimate(conv(rectpuls(t),normpdf(t,0,1)),2));

Is there a way to understand which one of the two plots is the correct one?

Thank you.

2

There are 2 best solutions below

0
On

Comment:

Let $U \sim \mathsf{Unif}(-.5, .5)$ and $Z \sim \mathsf{Norm}(0,1).$ Then $X = U + Z.$ I suppose that in your final result you use $Q$ for the standard normal CDF, often written $\Phi.$ However, what you have is a non-positive function. I think it should be $\Phi(x + .5)-\Phi(x - .5).$ [Thus, your analysis is correct except for a simple mistake in algebra.]

The following demonstration in R statistical software with a million observations $X$ shows a histogram consistent with that density function:

 u = runif(10^6, -.5, .5);  z = rnorm(10^6);  x = u + z
 hist(x, prob=T, br=50, col="skyblue2")
    curve(pnorm(x+.5)-pnorm(x-.5), col="red", lwd=2, add=T)

enter image description here

Note: In R runif and rnorm sample from uniform and normal distributions and pnorm is a normal CDF. I'm sorry not to use Matlab, but I do not have access to it.

0
On

There is a perfect agreement in fact if, instead of decimating, you simply add the 'same' parameter that forces the convolution result to have the same size as the original, instead of twice its size, the reason for which you were using a decimation by 2), as written in the program below:

clear all;close all;hold on
t=-5:0.01:5;
a=1/sqrt(2*pi);
normpdf=@(t)(a*exp(-t.^2/2));
plot(t,qfunc(t-0.5)-qfunc(t+0.5),'b');
plot(t,0.01*conv(rectpuls(t),normpdf(t),'same'),'r');