I have a random vector $x$ of length $N$, which is constructed as follows:
I know up to $k$ of the entries are non-zero.
The exact number of non-zero entries is based on a poisson distribution with parameter $\lambda$. i.e., Let $N_{nz}$ denote the number of non-zero entries, then
$\mathbb{P}(N_{nz}=i)=\frac{e^{-\lambda}\lambda^{i}/i!}{\sum_{j=0}^{k}{}e^{-\lambda}\lambda^{j}/j!}$
Once $N_{nz}=n_{nz}$ is realized, we uniformly pick $n_{nz}$ positions in $x$ for the non-zero entries . Moreover, each of the non-zero entry in $x$ follows an independent normal distribution with a different mean but same variance $N(h_i,\sigma^{2})$.
My aim is to compute the variance of $x$ defined as $\mathbb{E}(xx^{T})-(\mathbb{E}x)(\mathbb{E}x)^{T}.$
Before presenting my answer, let's establish some conventions.
In my answer, I will repeatedly use the Law of Total Probability for expected values. We are assuming that $\mathrm{Var}(x):=\mathbb{E}[xx^{T}]-(\mathbb{E}[x])(\mathbb{E}[x])^{T}$, So we just have to calculate $\mathbb{E}[xx^{T}]$ and $\mathbb{E}[x]$. Let's start with $\mathbb{E}[xx^{T}]$.
\begin{align} \mathbb{E}[xx^{T}]&=\mathbb{E}[xx^{T}\mathbf{1}_{\{N_{nz}>0\}}+xx^{T}\mathbf{1}_{\{N_{nz}=0\}}] \\ &=\mathbb{E}[\mathbb{E}[xx^{T}\mathbf{1}_{\{N_{nz}>0\}}+xx^{T}\mathbf{1}_{\{N_{nz}=0\}}|N_{nz}]] \\ &=\mathbb{E}\left[\mathbb{E}[xx^{T}\mathbf{1}_{\{N_{nz}>0\}}|N_{nz}]+\mathbb{E}[xx^{T}\mathbf{1}_{\{N_{nz}=0\}}|N_{nz}]\right] \\ &\qquad\qquad \mbox{Nore that if }N_{nz}=0\mbox{ then }x=\mathbf{0}\mbox{, then }\mathbb{E}[xx^{T}\mathbf{1}_{\{N_{nz}=0\}}|N_{nz}]=0.\\ &= \mathbb{E}\left[\mathbb{E}[xx^{T}\mathbf{1}_{\{N_{nz}>0\}}|N_{nz}]\right] \\ &=\mathbb{E}_{\color{red}{N_{nz}}}\left[\mathbf{1}_{\{N_{nz}>0\}}\mathbb{E}_{\color{red}{x,U}}\left[\sum_{i=1}^{N_{nz}}x_{u_{i}}^{2}\right]\right]\\ & \qquad\qquad\mbox{where }U:=(u_{1},\dots,u_{N_{nz}})\sim\mathrm{Uniform}(\mathcal{P}_{N_{nz}})\\ & \qquad\qquad\mbox{and } x_{u_{i}}\neq 0 \mbox{for each }i=1,\ldots,N_{nz}.\\ & \qquad\qquad \mbox{Also, note that for $u_{i}$ fixed } x_{ui}\sim N(h_{u^{i},\sigma^{2}}), \mbox{ then} \\ & \qquad\qquad x_{u_{i}}=\sigma Y_{u_{i}}+h_{i}\mbox{ where }\{Y_{i}\}_{i=1}^{N} \mbox{ are iid }N(0,1) \mbox{ random variables.}\\ & \qquad\qquad \mbox{Therefore, } x_{u_{i}}^{2}=\sigma^{2}Y_{u_{i}}^{2}+2h_{u_{i}}\sigma Y_{u_{i}}+h_{u_{i}}^{2}\\ &= \mathbb{E}_{\color{red}{N_{nz}}}\left[\mathbf{1}_{\{N_{nz}>0\}}\mathbb{E}_{\color{red}{Y,U}}\left[\sum_{i=1}^{N_{nz}}\sigma^{2}Y_{u_{i}}^{2}+2h_{u_{i}}\sigma Y_{u_{i}}+h_{u_{i}}^{2}\right]\right]\\ &= \mathbb{E}_{\color{red}{N_{nz}}}\left[\mathbf{1}_{\{N_{nz}>0\}}\mathbb{E}_{\color{red}{U}}\left[\mathbb{E}\left[\left.\sum_{i=1}^{N_{nz}}\sigma^{2}Y_{u_{i}}^{2}+2h_{u_{i}}\sigma Y_{u_{i}}+h_{u_{i}}^{2}\right|U\right]\right]\right]\\ &= \mathbb{E}_{\color{red}{N_{nz}}}\left[\mathbf{1}_{\{N_{nz}>0\}}\mathbb{E}_{\color{red}{U}}\left[\mathbb{E}_{\color{red}{Y}}\left[\sum_{i=1}^{N_{nz}}\sigma^{2}Y_{u_{i}}^{2}+2h_{u_{i}}\sigma Y_{u_{i}}+h_{u_{i}}^{2}\right]\right]\right]\\ &= \mathbb{E}_{\color{red}{N_{nz}}}\left[\mathbf{1}_{\{N_{nz}>0\}}\mathbb{E}_{\color{red}{U}}\left[\sum_{i=1}^{N_{nz}}\sigma^{2}\mathbb{E}_{\color{red}{Y}}\left[Y_{u_{i}}^{2}\right]+2h_{u_{i}}\sigma \mathbb{E}_{\color{red}{Y}}\left[Y_{u_{i}}\right]+h_{u_{i}}^{2}\right]\right]\\ & \qquad\qquad \mbox{Note that for $u_{i}$ fixed we have }Y_{u_{i}}^{2}\sim\mathcal{X}^{2}_{1} \mbox{ where}\\ &\qquad\qquad \mathcal{X}_{1}^{2} \mbox{ is Chi-squared distribution with 1 degree of freedom.}\\ &\qquad\qquad \mbox{Therefore, }\mathbb{E}_{\color{red}{Y}}\left[Y_{u_{i}}^{2}\right]=1. \mbox{In addition, we know that } \mathbb{E}_{\color{red}{Y}}\left[Y_{u_{i}}\right]=0.\\ &= \mathbb{E}_{\color{red}{N_{nz}}}\left[\mathbf{1}_{\{N_{nz}>0\}}\mathbb{E}_{\color{red}{U}}\left[\sum_{i=1}^{N_{nz}}(\sigma^{2}+h_{u_{i}}^{2})\right]\right]\\ &= \mathbb{E}_{\color{red}{N_{nz}}}\left[\mathbf{1}_{\{N_{nz}>0\}}\sum_{(u_{1},\ldots,u_{N_{nz}})\in \mathcal{P}_{N_{nz}}}\frac{1}{\binom{N}{N_{nz}}}\left(\sum_{i=1}^{N_{nz}}(\sigma^{2}+h_{u_{i}}^{2})\right) \right]\\ &= \mathbb{E}_{\color{red}{N_{nz}}}\left[\mathbf{1}_{\{N_{nz}>0\}}\left(N_{nz}\sigma^{2}+\frac{1}{\binom{N}{N_{nz}}}\sum_{(u_{1},\ldots,u_{N_{nz}})\in \mathcal{P}_{N_{nz}}}\left(\sum_{i=1}^{N_{nz}}h_{u_{i}}^{2}\right)\right) \right]\\ &= \mathbb{E}_{\color{red}{N_{nz}}}\left[ \mathbf{1}_{\{N_{nz}>0\}}\left( N_{nz}\sigma^{2}+\frac{1}{\binom{N}{N_{nz}}}\sum_{j=1}^{N}h_{j}\binom{N-1}{N_{nz}-1} \right)\right]\\ \end{align} \begin{align} &= \mathbb{E}_{\color{red}{N_{nz}}}\left[ \mathbf{1}_{\{N_{nz}>0\}}\left(N_{nz}\sigma^{2}+\frac{\binom{N-1}{N_{nz}-1} }{\binom{N}{N_{nz}}}\sum_{j=1}^{N}h_{j}^{2} \right)\right] \color{white}{--------.}\\ &= \mathbb{E}_{\color{red}{N_{nz}}}\left[\mathbf{1}_{\{N_{nz}>0\}}\left(N_{nz}\sigma^{2}+\frac{N_{nz} }{N}\sum_{j=1}^{N}h_{j}^{2} \right)\right]\\ &= \sum_{n_{nz}=0}^{k}\mathbf{1}_{\{n_{nz}>0\}}\left(n_{nz}\sigma^{2}+\frac{n_{nz} }{N}\sum_{j=1}^{N}h_{j}^{2} \right)\mathbb{P}[N_{nz}=n_{nz}]\\ &= \sum_{n_{nz}=1}^{k}\left(n_{nz}\sigma^{2}+\frac{n_{nz} }{N}\sum_{j=1}^{N}h_{j}^{2} \right)\mathbb{P}[N_{nz}=n_{nz}]\\ &= \sum_{n_{nz}=1}^{k}n_{nz}\left(\sigma^{2}+\frac{1 }{N}\sum_{j=1}^{N}h_{j}^{2} \right)\frac{e^{-\lambda}\lambda^{n_{nz}}/(n_{nz}!)}{\sum_{i=0}^{k}e^{-\lambda}\lambda^{i}/i!}\\ &= \left(\sigma^{2}+\frac{1 }{N}\sum_{j=1}^{N}h_{j}^{2} \right)\sum_{n_{nz}=1}^{k}n_{nz}\frac{e^{-\lambda}\lambda^{n_{nz}}/(n_{nz}!)}{\sum_{i=1}^{k}e^{-\lambda}\lambda^{i}/i!}\\ &=\left(\sigma^{2}+\frac{1 }{N}\sum_{j=1}^{N}h_{j}^{2} \right)\left(\frac{1}{\sum_{i=0}^{k}e^{-\lambda}\lambda^{i}/i!}\sum_{n_{nz}=1}^{k}n_{nz}\frac{e^{-\lambda}\lambda^{n_{nz}}}{n_{nz}!}\right)\\ &=\left(\sigma^{2}+\frac{1 }{N}\sum_{j=1}^{N}h_{j}^{2} \right)\left(\frac{1}{\sum_{i=0}^{k}e^{-\lambda}\lambda^{i}/i!}\sum_{n_{nz}=1}^{k}\frac{e^{-\lambda}\lambda^{n_{nz}}}{(n_{nz}-1)!}\right)\\ &=\left(\sigma^{2}+\frac{1 }{N}\sum_{j=1}^{N}h_{j}^{2} \right)\left(\frac{\lambda}{\sum_{i=0}^{k}e^{-\lambda}\lambda^{i}/i!}\sum_{\ell=0}^{k-1}\frac{e^{-\lambda}\lambda^{\ell}}{\ell!}\right)\\ &=\left(\sigma^{2}+\frac{1 }{N}\sum_{j=1}^{N}h_{j}^{2} \right)\left(\lambda -\lambda\frac{e^{-\lambda}\lambda^{k}/k!}{\sum_{i=0}^{k}e^{-\lambda}\lambda^{i}/i!}\right)\\ &=\left(\sigma^{2}+\frac{1 }{N}\sum_{j=1}^{N}h_{j}^{2} \right)\lambda\left(1-\mathbb{P}[N_{nz}=k]\right)\\ \end{align}
Now we are going to calculate $\mathbb{E}[x]:=[\mathbb{E}[x_{1}],\ldots,\mathbb{E}[x_{N}]]$, note that we only have to calculate $\mathbb{E}[x_{i}]$ for each $i=1,2,\ldots,N$.
\begin{align} \mathbb{E}[x_{i}] &= \mathbb{E}[x_{i}\mathbf{1}_{\{N_{nz}>0\}}+x_{i}\mathbf{1}_{\{N_{nz}=0\}}] \\ &= \mathbb{E}\left[\mathbb{E}[x_{i}\mathbf{1}_{\{N_{nz}>0\}}+x_{i}\mathbf{1}_{\{N_{nz}=0\}}|N_{nz}]\right] \\ &= \mathbb{E}\left[\mathbb{E}[x_{i}\mathbf{1}_{\{N_{nz}>0\}}|N_{nz}]+\mathbb{E}[x_{i}\mathbf{1}_{\{N_{nz}=0\}}|N_{nz}]\right] \\ &\qquad\qquad \mbox{Nore that if }N_{nz}=0\mbox{ then }x=\mathbf{0}\mbox{, then }\mathbb{E}[x_{i}\mathbf{1}_{\{N_{nz}=0\}}|N_{nz}]=0.\\ &=\mathbb{E}\left[\mathbb{E}\left[\left. x_{i}\mathbf{1}_{\{N_{nz}>0\}}\right|N_{nz}\right]\right]\\ &= \mathbb{E}_{\color{red}{N_{nz}}}\left[\mathbf{1}_{\{N_{nz}>0\}}\mathbb{E}_{\color{red}{X,U}}\left[X_{i}\mathbf{1}_{i\in U}+X_{i}\mathbf{1}_{i\notin U}\right]\right]\\ &\qquad\qquad \mbox{where } U:=(u_{1},\dots,u_{N_{nz}})\sim\mathrm{Uniform}(\mathcal{P}_{N_{nz}})\\ &= \mathbb{E}_{\color{red}{N_{nz}}}\left[\mathbf{1}_{\{N_{nz}>0\}}\mathbb{E}\left[\mathbb{E}\left[\left.X_{i}\mathbf{1}_{i\in U}+X_{i}\mathbf{1}_{i\notin U}\right|U\right]\right]\right]\\ &=\mathbb{E}_{\color{red}{N_{nz}}}\left[\mathbf{1}_{\{N_{nz}>0\}}\mathbb{E}_{\color{red}{U}}\left[\mathbb{E}_{\color{red}{X}}\left[X_{i}\mathbf{1}_{i\in U}\right]+\mathbb{E}_{\color{red}{X}}\left[X_{i}\mathbf{1}_{i\notin U}\right]\right]\right]\\ &\qquad\qquad \mbox{Note that if }i\in U \mbox{ then }x_{i}\sim N(h_{i},\sigma^{2}), \mbox{therefore}, \mathbb{E}_{\color{red}{X}}\left[X_{i}\mathbf{1}_{i\in U}\right]=h_{i}.\\ &\qquad\qquad \mbox{Note that if }i\notin U \mbox{ then }x_{i}=0, \mbox{therefore}, \mathbb{E}_{\color{red}{X}}\left[X_{i}\mathbf{1}_{i\notin U}\right]=0.\\ &=\mathbb{E}_{\color{red}{N_{nz}}}\left[\mathbf{1}_{\{N_{nz}>0\}}\mathbb{E}_{\color{red}{U}}\left[h_{i}\mathbf{1}_{i\in U}\right]\right]\\ &=\mathbb{E}_{\color{red}{N_{nz}}}\left[\mathbf{1}_{\{N_{nz}>0\}}h_{i}\mathbb{P}\left[i\in U\right]\right]\\ &=\mathbb{E}_{\color{red}{N_{nz}}}\left[\mathbf{1}_{\{N_{nz}>0\}}h_{i}\frac{\binom{N-1}{N_{nz}-1}}{\binom{N}{N_{nz}}}\right]\\ &=\frac{h_{i}}{N}\mathbb{E}_{\color{red}{N_{nz}}}\left[\mathbf{1}_{\{N_{nz}>0\}}N_{nz}\right]\\ &=\frac{h_{i}}{N}\sum_{n_{nz}=0}^{k}\mathbf{1}_{\{n_{nz}>0\}}n_{nz}\mathbb{P}[N_{nz}=n_{nz}]\\ &=\frac{h_{i}}{N}\sum_{n_{nz}=1}^{k}n_{nz}\frac{e^{-\lambda}\lambda^{n_{nz}}/(n_{nz}!)}{\sum_{j=0}^{k}e^{-\lambda}\lambda^{j}/j!}\\ &=\frac{h_{i}}{N\sum_{j=0}^{k}e^{-\lambda}\lambda^{j}/j!}\sum_{n_{nz}=1}^{k}\frac{n_{nz}e^{-\lambda}\lambda^{n_{nz}}}{n_{nz}!}\\ &=\frac{h_{i}}{N\sum_{j=0}^{k}e^{-\lambda}\lambda^{j}/j!}\sum_{n_{nz}=1}^{k}\frac{e^{-\lambda}\lambda^{n_{nz}}}{(n_{nz}-1)!}\\ &=\frac{h_{i}\lambda}{N\sum_{j=0}^{k}e^{-\lambda}\lambda^{j}/j!}\sum_{n_{nz}=1}^{k}\frac{e^{-\lambda}\lambda^{n_{nz}-1}}{(n_{nz}-1)!}\\ &=\frac{h_{i}\lambda}{N\sum_{j=0}^{k}e^{-\lambda}\lambda^{j}/j!}\sum_{\ell=0}^{k-1}\frac{e^{-\lambda}\lambda^{\ell}}{\ell!}\\ &=\frac{h_{i}\lambda}{N\sum_{j=0}^{k}e^{-\lambda}\lambda^{j}/j!}\left(\sum_{\ell=0}^{k}\frac{e^{-\lambda}\lambda^{\ell}}{\ell!}-\frac{e^{-\lambda}\lambda^{k}}{k!}\right)\\ &=\frac{h_{i}\lambda}{N}\left(1-\mathbb{P}[N_{nz}=k]\right) \end{align}
Therefore, we have \begin{align} (\mathbb{E}[x])(\mathbb{E}[x])^{T}&=\sum_{i=1}^{N} (\mathbb{E}[x_{i}])^{2}\\ &= \sum_{i=1}^{N} \frac{h_{i}^{2}\lambda^{2}}{N^2}\left(1-\mathbb{P}[N_{nz}=k]\right)^{2}\\ &= \lambda^{2}\left(\frac{1-\mathbb{P}[N_{nz}=k]}{N}\right)^{2}\sum_{i=1}^{N} h_{i}^{2}. \end{align}
So, we get the following
\begin{align} \mathrm{Var}(x)&:=\mathbb{E}[xx^{T}]-(\mathbb{E}[x])(\mathbb{E}[x])^{T}\\ &= \left(\sigma^{2}+\frac{1 }{N}\sum_{j=1}^{N}h_{j}^{2} \right)\lambda\left(1 -\mathbb{P}[N_{nz}=k]\right)-\lambda^{2}\left(\frac{1-\mathbb{P}[N_{nz}=k]}{N}\right)^{2}\sum_{i=1}^{N} h_{i}^{2}\\ &=\lambda(1-\mathbb{P}[N_{nz}=k])\left(\sigma^{2}+\frac{(N-\lambda(1-\mathbb{P}[N_{nz}=k]))}{N^{2}}\sum_{i=1}^{N} h_{i}^{2}\right). \end{align}