Consider the sum squared deviation $\text{SS}$ of $N$ observations $\hat{x}_1,\ldots,\hat{x}_N$ of a poisson random variable with mean $\mu$, i.e. $$ \text{SS} = \sum_{m=0}^{\infty} \left( \sum_{n=1}^N \frac{\delta_{\hat{x}_n m}}{N} - P_m \right)^2 $$ where $P_m$ is the probability mass function of the poisson distribution.
To compute the variance of this quantity it is necessary to evaluate expectation values involving up to four products of kronecker deltas i.e. $$ \langle \delta_{\hat{x}_{n} m} \delta_{\hat{x}_{n'} m} \delta_{\hat{x}_{\tilde{n}} \tilde{m}} \delta_{\hat{x}_{\tilde{n}'} \tilde{m}}\rangle $$ I can see that $\langle \delta_{\hat{x}_{n} m}\rangle = P_m$, and that $\langle \delta_{\hat{x}_{n} m} \delta_{\hat{x}_{n'} m'} \rangle = \begin{cases} P_m \delta_{mm'} & n'=n \\ P_m P_{m'} & n' \neq n \end{cases}$, but I am getting lost in the combinatorics beyond this point.
How do we calculate the three and four term products?