Expectation of $b^T \operatorname{sign}(Ab)$ I'm trying to compute the expectation of:
$$b^T \operatorname{sign}(Ab)$$
Where $b$ is a $n\times1$ vector of independent Bernoulli random variables:
$$P(b_i = 1) = 0.5,\quad P(b_i = -1) = 0.5$$
and $A$ is a fixed $n\times n$ matrix.
I think this might be hard to compute, since we need to find the probabilities of:
$$v = \operatorname{sign}(Ab)\quad P(v_i > 0).$$
And $v$ and $b$ are correlated...
However, computing $v_i$ corresponds to computing the sum:
$$v_i = \sum_j A_{ij} b_i.$$ 
I think there must be some method to find all values of $v_i$ that are interesting (near zero), since you are summing the values of $A$. This way it might be possible to find the probability density of $sign(v)$.
However we will also need to take the correlation into account.
If it's not possible or straightforward to compute, does anyone know how to bound the expectation? Or possibly, I will just sample the $b$'s to compute an empirical expectation, however I would like to know if we can make some gaurantee that the empirical expectation is near the real expectation.
 A: A closed-form solution looks nontrivial; in fact, an exact solution for arbitrary $A$ appears to be #P hard. But for fixed $A$ you can approximate the expectation with arbitrary precision using the law of large numbers. In addition, the central limit theorem (CLT) gives a rate of convergence. Let 
$$f(b) = b^T \text{sign}(Ab) = b^T v$$
Informally, the CLT states that for large $m$,
$$ \mathbb{E}[f(b)] \sim \frac{1}{m} \sum_{i=1}^m f(b^{(i)}) + \frac{\sigma \epsilon}{\sqrt{m}} $$
where $\sim$ denotes convergence in distribution as $m \rightarrow \infty$, the $b^{(i)}$ for $i=1,\ldots,m$ are a sequence of iid Bernoulli random vectors, $\sigma^2 = \text{Var}[f(b)]$ is a constant, and $\epsilon$ is a standard normal random variable. Beware that $\sigma^2$ depends on the entries of $A$ and could be very large. Also, since this is only an asymptotic result, arbitrarily large errors are still possible for finite $m$ and their probability is not guaranteed to follow the asymptotic distribution $\sigma \epsilon / \sqrt{m}$. It is probably not a good idea to rely on this method without some careful consideration of these issues.
Here are some detailed observations about special cases and general hardness:


*

*If $A$ is the identity matrix we have $f(b) = \sum_{i=1}^n b_i^2 = n$. This should be useful for approximating the solution when $A$ is a small perturbation of the identity. 

*Since you only need the expectation you should probably exploit the linearity of expectation:
\begin{align}
E[b^T v] 
&= \sum_{i=1}^n E[b_i v_i] \\
&= \sum_{i=1}^n E[b_i E[v_i|b_i]] \\
&= \tfrac{1}{2} \sum_{i=1}^n (-E[v_i|b_i=-1] + E[v_i|b_i=1]) 
\end{align}

*If the entries of $A$ are integers bounded uniformly by $|a_{ij}| \leq k$, the exact distribution of $(Ab)_i$ can be computed in $O(kn^2)$ time by iterative convolution. From this one can obtain $E[v_i|b_i]$, which in conjunction with observation (2) yields an $O(kn^3)$ algorithm for computing the expectation in the case of $A$ with bounded integer entries. 

*The algorithm from observation (3) can be extended to the case where the entries of $A$ have finite precision: for example, if the entries of $A$ are restricted to 0, +/- 0.1, +/- 0.2, etc. one can simply scale $A$ by a factor of 10 to obtain an integer matrix. This allows us to approximate the expectation with arbitrary precision, but the algorithm's complexity scales exponentially with the number of digits.

*From observation (2) we can see that $E[v_i|b_i]$ is needed to compute the expectation. This is trivial to compute if there is an entry $a_{ij}$ such that $$|a_{ij}| > \sum_{k \neq j} |a_{ik}|$$ But in general, computing this expectation requires us to count the number of ways that an arbitrary sequence of integers can be summed to exceed zero. This problem is known to be #P hard.

A: If all you want is a practical solution for any particular "value" of A, just employ stochastic simulation, and produce a confidence interval on your estimated expectation. 
I will not attempt to determine whether there is some nice variance reduction scheme available, which I suppose there is (such as control variables).  Because you can just very quickly blast the living daylights out of this by doing hundreds of millions of replications.
On the other hand, if you need to perform this to high accuracy in real-time for large values of n with new A matrices being provided multiple times per second ....
A: Take $n=2$. Then $\nu_1=a_{11}b_1+a_{12}b_2$. Let us examine $b_1sign(\nu_1)$. This is a function of $(b_1,b_2)$. Vector $(b_1,b_2)$ achieves values $\{(1,1),(1,-1),(-1,1),(-1,-1)\}$ with equal probability 1/4. Substituting these values into function $b_1sign(a_{11}b_1+a_{12}b_2)$ we can get the exact distribution of such function. 
We can replicate this process for $\nu_2$ and for general $n$, given the particular matrix $A$. It would not be hard to write a program which does that. Feeding it different matrices would give you an idea about theoretical behaviour of such function, however I am not sure that there is a exact analytic solution for arbitrary matrix $A$.
