Expectation of a quadratic form of Bernoulli random variables 
Let $X_i$ be independent Bernoulli random variables with success probabilities $p_i$. That is, $\mathbb{P}(X_i = 1) = p_i, \mathbb{P}(X_i = 0) = 1-p_i$. Is there a closed expression or an approximate formula for the following expectation?
  $$\mathbb{E}\left(\frac{\sum_i \sum_j X_i X_j B_{ij}}{(\sum_i X_i)^2}\right)$$
  where $B_{ij}$ are coefficients and $B_{ij} = B_{ji}, B_{ii} = 0, \forall i, j$.

In matrix form, this is a ratio of two quadratic forms (while the latter one has a power of 2)
$$\mathbb{E}\left(\frac{\mathbf{X}^T \mathbf{B} \mathbf{X}}{(\mathbf{X}^T \mathbf{X})^2}\right)$$
where $\mathbf{B}$ is a diagonal free symmetric matrix.
Intuitively, this value is the expected average value of a submatrix of $\mathbf{B}$ where the rows (and columns) are picked randomly.
It can be evaluated in exponential time by enumerating all the scenarios. But I am wondering if there are any methods that can evaluate this expectation in polynomial time.
Note that if we can evaluate 
$$\mathbb{E}\left(\sum_i \sum_j X_i X_j B_{ij} \Big| \sum_i X_i = m\right)$$
then we are done since the distribution of $\sum_i X_i$ is known (Poisson Binomial Disribution)
Note too that I have already figured out the following expectation
$$\mathbb{E}\left(\sum_i C_i X_i \Big| \sum_i X_i = m\right) = \frac{1}{P_k}\sum_{n=0}^N \sum_{j=1}^N \frac{A^{-nm+n}}{(N+1)}p_j C_j \prod_{\substack{k=1\\k\neq j}}^N (p_k A^n + (1-p_k))$$
$$P_k = \mathbb{P}\left(\sum_i X_i = m\right) = \frac{1}{N+1}\sum_{n=0}^N A^{-nm}\prod_{j=1}^N (p_j A^n + (1-p_j)),\quad A = e^{2i\pi/(N+1)}$$
using the method introduced in Fernandez (2010).
 A: I just found a closed-form expression for the expectation, which can be evaluated in polynomial time.
\begin{align}
  &\mathbb{E} \left[\left( \displaystyle\sum_{i=1}^N X_i \right)^{-2} \sum_{\substack{1\le i, j\le N \\ i\neq j}} X_iX_j B_{ij}   \right] \\
  &= \sum_{m=2}^N \mathbb{E} \left[\left( \displaystyle\sum_{i=1}^N X_i \right)^{-2} \sum_{\substack{1\le i, j\le N \\ i\neq j}} X_iX_j B_{ij}   \left|  \sum_{i=1}^N X_i = m \right. \right] \mathbb{P}\left(  \sum_{i=1}^N X_i = m  \right) \\
  &= \sum_{m=2}^N \mathbb{E} \left[\left( \displaystyle\sum_{i=1}^N X_i \right)^{-2} \sum_{\substack{1\le i, j\le N \\ i\neq j}} X_iX_j B_{ij}   \left|  \sum_{i=1}^N X_i = m \right. \right] \mathbb{P}\left(  \sum_{i=1}^N X_i = m  \right) \\
  &= \sum_{m=2}^N \frac{1}{m^2}  \sum_{\substack{1\le i, j\le N \\ i\neq j}} B_{ij} \mathbb{E} \left[ X_iX_j \left|  \sum_{i=1}^N X_i = m \right. \right] \mathbb{P}\left(  \sum_{i=1}^N X_i = m  \right) \\
  &= \sum_{m=2}^N \frac{1}{m^2}  \sum_{\substack{1\le i, j\le N \\ i\neq j}} B_{ij} \mathbb{P} \left(  X_i = 1, X_j = 1 \left|  \sum_{i=1}^N X_i = m \right.  \right) \mathbb{P}\left(  \sum_{i=1}^N X_i = m  \right) \\
  &= \sum_{m=2}^N \frac{1}{m^2}  \sum_{\substack{1\le i, j\le N \\ i\neq j}} B_{ij} \mathbb{P} \left( \left.  \sum_{i=1}^N X_i = m \right| X_i = 1, X_j = 1 \right) \mathbb{P}\left( X_i = 1, X_j = 1  \right) \\
  &= \sum_{m=2}^N \frac{1}{m^2}  \sum_{\substack{1\le i, j\le N \\ i\neq j}} B_{ij} p_i p_j \sum_{n=0}^N \frac{1}{(N+1)} A^{-nm+2n} \prod_{\substack{k=1\\k\neq i,j}}^N (p_k A^n + (1-p_k))  \\
\end{align}
Note that $m$ starts from 2 as when $m < 2$ the conditional expectation is 0.
In stead of $O(2^N)$, this formula can be evaluated in $O(N^5)$ time.
