Second order approximations of expected exponentiated sums of dependent Bernoulli random variables

Question

I would like to approximate $\mathbb{E}[e^y]$ for $y=\sum_i^n c_i x_i$, and $x_i\sim \text{Bernoulli}(p_i)$ ($p_i$'s distinct) using $p_i$ and $\text{cov}(x_i,x_j)$. The $c_i\in \mathbb{R}$ are fixed, and the $x_i$ are not independent.

---

Currently, I am estimating $\mathbb{E}[e^y]$ by assuming the $x_i$ are independent. That is, I take $\mathbb{E}[e^y]\approx \prod_i[ (1-p_i)+p_i e^{c_i}]$, which only requires me to know $p_i=\mathbb{E}[x_i]$.

However, in addition to knowing $p_i$, I also know $\text{cov}(x_i,x_j)$, but I am not using this information. I am searching for a way to improve my estimation of $\mathbb{E}[e^y]$ by making use of the covariances.

I have tried a few things. Each has problems which I will attempt to illustrate with a toy example. First a bit more notation: Let $\mu=\mathbb{E}[y]$, $\sigma^2=\text{var}[y]$, and take $\hat{y}=\sum_i^n c_i\hat{x}_i$ with $\hat{x}_i\sim \text{Bernoulli}(p_i)$ all independent.

Now a toy example. Suppose $y=7x_1+5x_2-8x_3$, $P(x_1=1,x_2=1,x_3=0)=P(x_1=0,x_2=1,x_3=1)=0.2$, and all other assignments of $x_1,x_2,x_3$ occur with probability $0.1$. Then, $\text{cov}(x_1,x_3)=-0.05$, $\text{cov}(x_1,x_2)=\text{cov}(x_2,x_3)=0$, and $\mathbb{E}[e^y]=32,681.07$.

Methods of approximation:

1) Ignore the dependencies: $\mathbb{E}[e^y]\approx \mathbb{E}[\text{exp}(\hat{y})]$. This approach only uses first order information ($p_i$). On the toy example we get $\mathbb{E}[\text{exp}(\hat{y})]=24,553.48$.

2) Taylor series approximation: $\mathbb{E}[e^y]\approx e^\mu(1+\tfrac12 \sigma^2)$. This approach gives terrible approximations for large $c$'s. On the toy example we get $e^\mu(1+\tfrac12 \sigma^2)=7.697877$.

3) Log-normal approximation: $\mathbb{E}[e^y]\approx e^{\mu+\tfrac12 \sigma}$. This approach overshoots very badly which is problematic for my application. On the toy example $e^{\mu+\tfrac12 \sigma}=165,585,792$. But note that $\mathbb{E}[e^y]$ is bounded between $\text{exp}({\sum_i \mathbb{1}_{\{c_i<0\}}c_i}) \le \mathbb{E}[e^y] \le \text{exp}({\sum_i \mathbb{1}_{\{c_i>0\}}c_i})=162,754.8$.

4) Treat $x_i$ as independent, but correct the second moment: Note that $\mathbb{E}[e^y]=\sum_{k=0}^\infty \frac{1}{k!}\mathbb{E}[y^k]$, then we might consider $\mathbb{E}[e^y]\approx \mathbb{E}[\text{exp}(\hat{y})]+\tfrac12 \mathbb{E}[y^2]-\tfrac12 \mathbb{E}[\hat{y}^2]$. I'm not sure this is well principled, and I empirically observe very small corrections when the coefficients $c_i$ are large (as in the toy example). In this case, this approximation gives a value of $24,556.28$.

---

More details:

The above problem arises from the following. Imagine an $n\times m$ matrix $R$ of indicator vectors, i.e. each row is an independent random vector in $\{0,1\}^m$ that takes at most a single 1 (all other positions being zero), the probabilities of which are known. Then $y=\tfrac12 \text{tr}(R^T C R)$ where $C$ is a symmetric matrix of constants with zeros on the main diagonal. In the above formulation $x_i$ is an off diagonal element of the matrix $RR^T$. I presented the problem in that way because I thought it makes the problem more clear. I've given these additional details in case they are important in some way.

For what it's worth, in my problems $n$ is usually $O(10^1)-O(10^2)$ and some of the $c_i$ may be large in magnitude (e.g. $\ge O(10^1)$).

Answer 1

Without further restrictions on the dependence structure, knowing the covariance matrix is not generally very useful in computing or approximating $\mathbb{E}[e^y]$.

If the joint distribution of $(x_1,x_2,\cdots,x_n)$ was multivariate normal then knowledge of the covariance matrix $\text{cov}(x_i,x_j)$, combined with knowledge of the individual $\mathbb{E}[x_i]$, would completely determine the value of $\mathbb{E}[e^y]$. However, in this case the joint distribution of $(x_1,x_2,\cdots,x_n)$ is not determined fully by knowledge of the covariances, except of course in the special case where $n=2$.

Here's an example of what can happen in the case $n=3$. For $n>3$ there are more 'degrees of freedom' so similar things can happen but the notation becomes more complicated. Let $p_{ijk}$ denote $P(x_1=i,x_2=j,x_3=k)$. Suppose that $p_{ijk}=f(i+j+k)$ for some function $f$ i.e. the joint probability of a sequence of $x_i$ values depends only on the number of non-zero entries. For simplicity suppose that $p_1=p_2=p_3=\frac{1}{2}$, which means that $f(1)+2f(2)+f(3)=\frac{1}{2}$. Covariances are determined by the value of $\mathbb{E}[x_1 x_2]=P(x_1=1,x_2=1)=f(2)+f(3)$; let's fix this value at $\frac{1}{4}$. We also require that the probabilities sum to $1$, so $f(0)+3f(1)+3f(2)+f(3)=1$.

You can check that $(f(0),f(1),f(2),f(3))=(\frac{1}{4},0,\frac{1}{4},0)$ and $(f(0),f(1),f(2),f(3))=(0,\frac{1}{4},0,\frac{1}{4})$ both satisfy all the above constraints. In particular they both lead to the same covariances. However, in the case where all the $c_i$ are large and positive, it is clear that $\mathbb{E}[e^y]$ is dominated by the term $f(3)e^{c_1+c_2+c_3}$, which will be very different for these two possible functions.

So, to get any further using the covariances you'd need to have further restrictions on the dependence structure in the $x_i$ variables, such that the covariances restrict higher order dependencies.

Alternative approach if you can calculate the joint probabilities (can be used regardless of the signs of the $c_i$): look for the values of $(x_1,x_2,x_3)$ which maximise $\sum c_i x_i$; suppose that these values have joint probability $p$. Then approximate $\mathbb{E}[e^y]\approx p\exp(\sum c_i x_i)$. In your example I get $0.2\times \exp(12)\approx 32551$, which doesn't seem too bad.