PDF of $| \sum_{k=1}^{n}{|h_k||g_k|\exp\left( j \theta_k \right)} |^2$ for small values of $n$ and $Q$?

Given the following function of random variables

$$f = \left|\sum_{k=1}^{n}{|h_k||g_k|\exp\left( j \theta_k \right)} \right|^2,$$ where $$h_1, \cdots, h_n$$ and $$g_1, \cdots, h_n$$ are i.i.d. random variables following the complex Gaussian distribution $$\mathcal{CN}(0, \sigma_{h}^2)$$ and $$\mathcal{CN}(0, \sigma_{g}^2)$$ and $$\theta_1, \cdots, \theta_n$$ $$\in [\frac{-\pi}{Q}, \frac{\pi}{Q}]$$ are i.i.d. uniformly distributed random variables with probability density function (PDF) given by $$\frac{Q}{2\pi}$$, where $$Q$$ is a integer number greater than 0. Additionally, we assume that $$h_k$$, $$g_k$$ and $$\theta_k$$ are independent for all values of $$k$$.

UPDATED on 10/01/2020

Based on my simulations I know that the PDF of $$r$$, which is defined as $$r = \left|\sum_{k=1}^{n}{|h_k||g_k|\exp\left( j \theta_k \right)} \right|,$$ can be accurately approximated by a Gamma random variable even for small values of $$n$$, and $$Q$$, $$e.g.$$, $$n=1$$ and $$Q=4$$. Therefore, the PDF of $$f = r^2$$ can be approximated as $$P_{n}(f) = \frac{1}{2 \Gamma(\kappa)\Theta^{\kappa}} f^{\left(\frac{\kappa-2}{2} \right)} \exp\left(-\frac{\sqrt{f}}{\Theta}\right), f > 0.$$

However, I'm not being able to find the $$\kappa$$ and $$\Theta$$, $$i.e.$$, shape and scale, parameters of the Gamma random variable.

This problem arises from the study on wireless communications channels and is of great importance to the research community.

• in one of your earlier questions you mentioned drawing from a Nagakami distribution, rather than from a Gaussian distribution; is there a reason you switched? (I am asking, because that earlier question gave the same motivation from "wireless communication".) – Carlo Beenakker Jan 5 at 21:53
• @CarloBeenakker, thanks for your reply. Nakagami channels are only one of several possible channel models. I thought it would be easier to start from the simpler Rayleigh channel, which is described by the absolute values of $g$ and $h$. – Felipe Augusto de Figueiredo Jan 6 at 1:09
• Just to be sure for the definition of the law of $\theta_i$, Is it that $\theta_i\in [0,\frac{2\pi}{Q}]$? – RaphaelB4 Jan 6 at 11:13
• @RaphaelB4, $\theta_{k} \in [\frac{-\pi}{Q}, \frac{\pi}{Q}]$. I've just updated the quaestion with this information.Thanks! – Felipe Augusto de Figueiredo Jan 6 at 12:28

n $$\mathbf{\gg}$$ 1:

The variable $$f=r^2$$ is the square of the distance $$r$$ from the origin after $$n$$ steps of a random walk on the plane with random direction [$$\theta$$ uniformly in $$(0,2\pi)$$] and mean squared step size $$s^2$$ given by $$s^2=\mathbb{E}(|h|^2||g|^2)=\mathbb{E}(|h|^2|)\mathbb{E}(|g|^2)=\sigma_h^2\sigma_g^2.$$ For $$n\gg 1$$ the distribution of $$r$$ is a Maxwell distribution and hence the distribution of $$f$$ is exponential, $$P_n(r)=\frac{2r}{ns^2}e^{-\frac{r^2}{ns^2}}\Rightarrow P_n(f)=\frac{1}{ns^2}e^{-\frac{f}{ns^2}}.$$ So for $$n\gg 1$$ the distribution $$P_n(f)$$ is a Gamma distribution with shape $$k=1$$ and scale $$\theta=n\sigma_h^2\sigma_g^2$$.

n $$=$$ 1:

For $$n=1$$ one can calculate directly $$P_1(r)=\int_0^\infty d|g|\int_0^\infty d|h|\, \delta(r-|gh|)P(|h|)P(|g|)$$ $$=(\sigma_h\sigma_g)^{-2}\int_0^\infty dh\, \frac{4r}{h} e^{-h^2/\sigma_h^2} e^{-(r/h)^2/\sigma_g^2}=\frac{4r}{\sigma_g^2\sigma_h^2}K_0\left(\frac{2r}{\sigma_g\sigma_h}\right),$$ with $$K_0$$ a Bessel function.
(I have used that the absolute value $$|h|$$ when $$h$$ is $$\mathcal{CN}(0, 1)$$ has distribution $$p(|h|)=2|h|e^{-|h|^2}$$.)
So for $$n=1$$ the distribution of $$f=r^2$$ is $$P_1(f)=\frac{2}{\sigma_g^2\sigma_h^2}K_0\left(\frac{2\sqrt{f}}{\sigma_g\sigma_h}\right).$$

arbitrary n:

A closed form for $$P_n(r)$$ exists for arbitrary $$n$$, see The non-isotropic two-dimensional random walk, by B.C. Barber (1993):

$$P_n(f)=(\sigma_g\sigma_h)^{-n-1}\frac{2f^{(n-1)/2}}{\Gamma(n)}K_{n-1}\left(\frac{2\sqrt f}{\sigma_g\sigma_h}\right).$$

The large-$$f$$ decay is a stretched exponential $$\propto e^{-2\sqrt{f}}$$, so this cannot be precisely modeled by a Gamma distribution (which decays exponentially).

anisotropic random walk:

All of this was for an isotropic random walk, with $$\theta_k$$ at each step uniformly distributed in $$(-\pi,\pi)$$. The parameter $$Q>1$$ restricts the angle to $$(-\pi/Q,\pi/Q)$$, so the scattering is peaked in the forward direction. For $$n\gg 1$$ the anisotropy can be accounted for by a rescaling of the mean squared step size, $$s_Q^2=\frac{\sigma_h^2\sigma_g^2}{1-\langle \cos\theta\rangle},$$ where $$\langle\cdots\rangle$$ denotes the angular average. We have $$\langle\cos\theta\rangle=\frac{Q}{2\pi}\int_{-\pi/Q}^{\pi/Q}\cos\theta\,d\theta=\frac{Q}{\pi}\sin(\pi/Q),$$ so we arrive at the distribution $$P_n(f)=\frac{1}{ns_Q^2}e^{-\frac{f}{ns_Q^2}},\;\;s_Q^2=\frac{\sigma_h^2\sigma_g^2}{1-(Q/\pi)\sin(\pi/Q)}.$$ This is a Gamma distribution with shape $$k=1$$ and scale $$\theta=ns_Q^2$$.

For $$n=1$$ the distribution $$P_1(f)$$ is $$Q$$-independent, so the result $$\propto K_0$$ above still applies. For $$n>1$$ I presume the decay will still be a stretched exponential, but I do not have an exact result as for $$Q=1$$.

• for $n\gg 1$ the parameters of the Gamma distribution are shape $=1$ and scale $=n\sigma_h^2\sigma_g^2$. For $n=1$ I find a Bessel function distribution, which for small $f$ diverges logarithmically as $-\ln f$ and for large $f$ decays as a stretched exponential $e^{-\sqrt f}$. So this is definitely not a Gamma distribution, but if you insist you might take the shape $k=1/2$, which gives a reasonable fit (see plot). Interpolating between $k=1/2$ for $n=1$ and $k=1$ for $n\geq 10$, with $\theta=n\sigma_h^2\sigma_g^2$, should then cover the entire range of $n$. Is that satisfactory? – Carlo Beenakker Jan 6 at 11:17
• I've added the $Q>1$ case for $n\gg 1$, which is a Gamma distribution with shape $=1$ and scale $=n\sigma_h^2\sigma_g^2[1-(Q/\pi)\sin(\pi/Q)]^{-1}$. I have also added the exact result for all $n$ for $Q=1$ (in terms of a Bessel function $K_{n-1}$). – Carlo Beenakker Jan 6 at 18:49
• there is not really much more to say: in the limit $n\gg 1$ of many steps the random walk satisfies a diffusion equation, which contains only a single parameter, the diffusion constant $D$, or equivalently, the mean free path $l$. For isotropic scattering in $d$ dimensions $D$ is related to the mean-square-displacement $s^2$ by $s^2=2dD\tau$ (with $\tau$ the time between steps). For anisotropic scattering you should replace $D$ by $D/(1-\langle\cos\theta\rangle)$, or equivalently, replace $l$ by the transport mean free path. – Carlo Beenakker Jan 7 at 12:15
• certainly, once you are in the diffusion regime $n\gg 1$ all you need to know is the mean-square displacement, it does not matter how the step size is distributed. – Carlo Beenakker Jan 7 at 14:07
• if you have more and more forward scattering, for larger and larger $Q$, you will need larger and larger $n$ to reach the diffusion regime, simply because you need backscattering for diffusion and backscattering is suppressed for large $Q$; if $Q$ is large and $n$ is small a better approximation would be to ignore backscattering altogether, and just fix the angles $\theta_k$ at zero. – Carlo Beenakker Jan 9 at 15:28

I understand the answer I gave, which involves either exact or asymptotic results, is not quite what the OP was looking for, so let me try to suggest a different approach. The question is: please find parameters $$k$$ and $$\theta$$ such that $$p_{k,\theta}(r) = \frac{r^{k-1}}{ \Gamma(k)\theta^{k}} e^{-r/\theta}$$ approximates the distribution of the distance from the origin after $$n$$ steps in a certain anisotropic random walk on the two-dimensional plane. The random walk has step size distribution $$F(r)=\frac{4r}{\sigma_g^2\sigma_h^2}K_0\left(\frac{2r}{\sigma_g\sigma_h}\right)$$ and scattering angle uniformly distributed in $$(-\pi/Q,\pi/Q)$$.

The first step would be to eliminate the parameters $$\sigma_g$$ and $$\sigma_h$$ by rescaling $$r=R\sigma_g\sigma_h$$ The distribution $$P_{n,Q}(R)$$ of $$R$$ no longer depends on $$\sigma_g$$ or $$\sigma_h$$. It only contains the two parameters $$n$$ and $$Q$$, both $$\in\mathbb{N}$$, and we wish to fit this with the two-parameter distribution $$p_{k,\theta}(R)$$, with $$k,\theta\in\mathbb{R}$$.

The next thing to note is that in the context in which the problem appears, only a few integer values of $$Q$$ are of relevance. We also have exact asymptotic results for $$n\gg 1$$, so if we would take 100 different combinations of $$Q$$ and $$n$$ we would have covered much of the relevant parameter space. The fitting procedure could be simplified by fitting first for large $$R$$, to extract $$\theta$$, and then for small $$R$$, to extract $$k$$. I guess this could be automated to a large extent, or perhaps even done by hand.

As an example, I tried this for $$Q=1$$, when I know the exact $$P_{n,1}(r)$$. The Mathematica command

Manipulate[ LogPlot[{(4*r^n/Gamma[n])*BesselK[n-1, 2*r], (r^(k-1)/(Gamma[k]*t^k))*Exp[-r/t]},{r,0,20}], {n,2,10,1}, {k,0.1,5}, {t,0.1,5}]

provides for a convenient way to do the fit by sliding the parameters $$n$$, $$k$$, and $$t\equiv\theta$$. Here is a representative output for $$n=3,5,8$$:

I followed up on the suggestion by the OP to carry out the fit by adjusting the parameters $$k$$ and $$\theta$$ such that the first two moments of $$R$$ agree. For $$Q=1$$ we have the exact result $$\mathbb{E}(R)=\frac{\sqrt{\pi } \,\Gamma \left(n+\frac{1}{2}\right)}{2 \Gamma (n)},\;\;\mathbb{E}(R^2)=n,$$ to be compared with the Gamma-distribution result $$\mathbb{E}(R)=k\theta$$, $$\mathbb{E}(R^2)=k(1+k)\theta^2$$.
I show the comparison for $$n=5$$, when $$k=2.955$$ and $$\theta=0.654$$, see the plot below. It seems this fitting procedure is worse than the unconstrained fit (with $$k=3.6$$ and $$\theta=0.55$$) carried out above.

• Dear Carlo, thanks once more for the new answer. What I'm trying to do right now is to find $\kappa$ and $\Theta$ based on the Gamma's distribution definition of $\mu$ and $\sigma^2$, $i.e.$, $\mu = \kappa \Theta$ and $\sigma^2 = \kappa \Theta^2$. The only problem I've found is how to find the expectation of $r$. If you have any jint on that I would be very glad! – Felipe Augusto de Figueiredo Jan 13 at 19:20
• If we find two different moments for $Q \neq 1$, for example the mean and variance, we will be able to find $\kappa$ and $\Theta$. I will post an expansion for $r^2$. – Felipe Augusto de Figueiredo Jan 13 at 21:42
• By using the following identity $r^2 = | \sum_{i=1}^{N}{z_{i} e^{\theta_{1}}} |^2 = \sum_{i=1}^{N}{z_{i}^{2}} + 2\sum_{i=1}^{N}{\sum_{j=i+1}^{N}{z_{i}z_{j}\cos(\theta_{i}-\theta_{j})}}$ we can find the expectation of $r^2$. Then, if we find a different expectation we can use the two of them to figure out $\kappa$ and $\Theta$. – Felipe Augusto de Figueiredo Jan 13 at 23:59
• I tried this for $Q=1$. The result is not very encouraging, a global fit produces much better results than the fit of the first two moments. – Carlo Beenakker Jan 14 at 11:07
• for $Q=1$ I have an exact expression for $P(R)$ in terms of the $K_{n-1}$ Bessel function (see my other answer); from that I can immediately calculate all moments of $R$. – Carlo Beenakker Jan 14 at 14:38