Variance of summation of products of all pairs of Bernoulli random variables

Question

I have $n$ independent Bernoulli random variables with parameter $p$: $\{X_1, \dots, X_n\}$ and I want to find the variance of the sum of the products of all pairs of these random variables. Specifically, I want to calculate

$Var\left[\sum_{ij \in [n]^2} (c_i X_i \cdot X_j)\right]$ where $c_{ij} \geq 0$ are constants.

Even a somewhat tight upper bound (tighter than $n^2$) is also helpful.

Answer 1

For simplicity, we symmetrize $(c_{ij})$ by letting $s_{ij} = \frac{1}{2}(c_{ij} + c_{ji})$. Also, define $Y$ by

$$ Y = \sum_{(i,j) \in [n]^2} c_{ij} X_i X_j = \sum_{(i,j) \in [n]^2} s_{ij} X_i X_j . $$

Then

\begin{align*} \mathbf{Var}(Y) = \mathbf{Cov}(Y, Y) = \sum_{(i,j), (k,l) \in [n]^2} s_{ij}s_{kl} \mathbf{Cov}(X_i X_j, X_k X_l). \end{align*}

Now by noting that $X_i^2 = X_i$, we have

\begin{align*} \mathbf{Cov}(X_i X_j, X_k X_l) &= \mathbf{E}[X_i X_j X_k X_l] - \mathbf{E}[X_i X_j] \mathbf{E}[X_k X_l] \\ &= p^{|\{i, j, k, l\}|} - p^{|\{i, j\}|}p^{|\{k, l\}|}. \end{align*}

Using the idea as in this answer, we can classify all the different scenarios for the values of $|\{i, j, k, l\}|$, $|\{i, j\}|$, and $|\{k, l\}|$. Consequently,

\begin{align*} \mathbf{Var}(Y) &= \Biggl( \sum_{i} s_{ii}^2 \Biggr) p(1-p) + 4 \Biggl( \sum_{i, j \text{ distinct}} s_{ii}s_{ij} \Biggr) p^2(1-p) \\ &\quad + 2 \Biggl( \sum_{i, j \text{ distinct}} s_{ij}^2 \Biggr) p^2(1-p^2) + 4 \Biggl( \sum_{i, j, k \text{ distinct}} s_{ij}s_{jk} \Biggr) p^3(1 - p). \end{align*}

Let $M = \max_{(i,j) \in [n]^2} s_{ij}$. Then an upper bound can be obtained:

\begin{align*} \mathbf{Var}(Y) &\leq M^2 \biggl[ \begin{gathered} np(1-p) + 4n(n-1) p^2(1-p) \\ + 2n(n-1)p^2(1-p^2) + 4n(n-1)(n-2) p^3(1-p) \end{gathered} \biggr] \\ &= M^2 np(1-p) [(2n-2)(2n-3)p^2 + 6(n-1)p + 1] \end{align*}

Note that this bound is also tight, in the sense that if $(s_{ij})$ is a constant array, then the equally holds for the above inequality. In particular, $\mathbf{Var}(Y)$ is in general of order $n^3$.

Below is a comparison of the values of $\mathbf{Var}(Y)$ by a direct computation and using the formula. Mathematica 13 is used to produce this.

n = 4;
(* Direct computation *)
dist = Table[x[k] \[Distributed] BernoulliDistribution[p], {k, n}];
Sum[Expectation[x[i]*x[j]*x[k]*x[l], dist] - Expectation[x[i]*x[j], dist] Expectation[x[k]*x[l], dist], {i, n}, {j, n}, {k, n}, {l, n}]
(* Using formula *)
n p (1 - p) + 4 n (n - 1) p^2 (1 - p) + 2 n (n - 1) p^2 (1 - p^2) + 4 n (n - 1) (n - 2) p^3 (1 - p) // Expand
Clear[n, dist];