I'm trying to compute using matlab the mutual information for an $ \infty $-PAM input (the limit of a very dense PAM constellation) for a range of snr and I got stuck.
I'm working with a real-valued scalar Gaussian channel with input power constrained to unity and standard Gaussian noise $ N \sim (0, 1) $. The output of the channel is given by, $$ y = \sqrt{ \gamma } x + n $$
I know that under the power input constraint the input is uniformly distributed on $ \left[ -\sqrt{ 3 }, \sqrt{ 3 } \right] $.
I found that in HIGH-SNR regime, there exists a gap between the Shannon capacity and the capacity achieved with $ \infty $-PAM around 0.254 bits. I know this is called the shaping loss and represents the loss of using an uniform input rather than a Gaussian input.
I know how to compute that. But, I got stuck trying to compute the capacity for a range of snr or $ \gamma $. I tried using the formula, $$ I(x;y) = h(y) - h(n) $$ with $ f_{y}(y) = \int f_{x}(x)f_{y|x}(y|x) dx $ and I got $$ f_{y} (y) = \frac{ 1 }{ 2 \sqrt{ 3 \gamma } } \frac{ 1 }{ 2 } \int_{ \frac{ y - \sqrt{3} }{ \sqrt{2} } }^{ \frac{ y + \sqrt{3} }{ \sqrt{2} } } \frac{ 2 }{ \sqrt{ \pi } } e^{ -t^{2} } dt \\ = \frac{ 1 }{ 2 \sqrt{ 3 \gamma } } \frac{ 1 }{ 2 } \left[ \text{erf} \left( \frac{ y + \sqrt{ 3 \gamma } }{ \sqrt{ 2 } } \right) - \text{erf} \left( \frac{ y - \sqrt{ 3 \gamma } }{ \sqrt{ 2 } } \right) \right] $$
I tried to compute $ h(y) = \int_{-\infty}^{\infty} f_{y} (y) \log_{2} f_{y} (y) $ and the integral does not converge. I tried with int$(\cdot)$ and integral$(\cdot)$ methods in matlab.
I also tried approximating the $ \text{erf}(\cdot) $ function with Taylor series expansion $$ erf(x) = \frac{ 2 }{ \sqrt{ \pi } } \sum_{ n = 0 }^{ \infty } \frac{ (-1)^{n} x^{ 2n + 1 } }{ n!(2n + 1) } $$ but it didn't work.
I think there would be another approach to solve this problem, but I don't know it. So,
Could you help me? How can I obtain a more manageable integral? Is my approach correct?
I appreciate any help!! Thanks.
Note that the integrand has the general form of $$g(z,a)\log_2 g(z,a)$$
where $g(z,a)=\mathrm{erf}(z+a) - \mathrm{erf}(z-a)$. Try to prove that the integrand $\sim0$ after, say $z>a+\delta$, $2\leq\delta\leq 3$ ($\delta =3$ could be good enough). After this you should be integrate it as a proper integral, i..e $\int_{-{a+\delta}}^{a+\delta}...$ without stability problem.
Note, that there is
erffunction in matlab, so you don't need to use it's Taylor expansion.If you want to try another integration function in matlab you can try
trapz?