Let $\varepsilon\sim N(0, I_D)$ be a $D$-dimensional random vector, distributed normally with mean $0$ and covariance given by the identity matrix of size $D$. In some computations I'm doing in my research, the following expectation arises: $$\mathbb{E}_{\varepsilon\sim N(0, I_D)}\left[\frac{\exp(\varepsilon_i b_i)}{\sum_{k=1}^D \exp(a_k + \varepsilon_kb_k) )}\right],$$ where $a, b\in\mathbb{R}^D$ are some $D$-dimensional vectors.
I wonder if this expectation has a closed form. It looks gnarly, though. The fact that passing a Gaussian through a sigmoid also has to be approximated gives me even less hope.