Representation of arbitrary probability distributions by latent variable models with gaussian posteriors and priors A probability distribution of a random variable $\mathbf x$ may be expressed as
$ P(x)=\int P(x|z) P(z) dz$
where $\mathbf z$ is a latent random variable.
In models such as Variational Autoencoders (VAEs) it is common to let $\mathbf z \sim \mathcal{N}(0,I)$ and $P(x|z) = \mathcal{N}\left(x;\mu(z),\Sigma(z)\right)$. However, it is not obvious to me whether those kind of models are able to represent any arbitrary probability distribution $P(x)$ by varying $\mu(z)$ and $\Sigma(z)$, and if that is the case, why it is possible.
Do you have any intuition or references where this is proven?
Thank you very much
 A: Background: we have a prior $p(z) = \mathcal{N}(z|0,I)$, a decoder 
$ p_\theta(x|z) = \mathcal{N}(x|\mu_\theta(z),\Sigma_\theta(z))$, and an encoder $q_\phi(z|x)=\mathcal{N}(z|\mu_\phi(x),\Sigma_\phi(x))$.
Often, $\Sigma_\theta(z)$ is very small, so the decoder is deterministic. 
This means the decoder is essentially just a neural network mapping from an $n$-D ball to the data space. Clearly, this is a problem if the true distribution $q(x)$ has a different topology, e.g. forms a torus, or is strongly separated multimodal distribution. People are aware of this problem: e.g.,  Falorsi et al, Explorations in homeomorphic variational autoencoding.
There are other problems, however. Since $\widehat{x} = \mu_\theta(z)$ is computed with an artificial neural network, we have the guarantees of the universal approximation theorem. For instance, given some complex image of a face, we should be able to reconstruct it. But in practice, VAE outputs are blurry (compared to GANs). Recall that the (Beta) VAE loss is 
$$ \mathcal{L}(\theta,\phi) = 
-\mathbb{E}_{q_\phi(z|x)}[ \log p_\theta(x|z) ] + 
\beta\, \mathcal{D}_\text{KL}[q_\phi(z|x)\mid\mid p(z)] $$
This blurriness is due to (1) the normal noise assumption of an $L_2$ likelihood and (2) the regularization effect of the KL penalty. In other words, the form of the loss in standard VAEs forms a practical, though not theoretical I think, bound on what the network can represent.
There is another issue with representation capability, but with the encoder rather than the decoder. As you likely know, the standard VAE uses a variational approximation to the inferred posterior $q_\phi(z|x)$, via a parameterizing a mean and variance. 
This is a very severe limitation. Consider an ambiguous $x$, where the ideal posterior is a complex multimodal distribution over $z$. The encoder is obviously highly limited in its ability to match it. This is a general problem with using a weak posterior approximation in variational inference.
This problem is in addition to the Bayesian regularization induced by the second term in the loss function, which further limits the expressivity of the posterior via penalizing its distance from the unexpressive prior.
The former issue has hard theoretical limits, whereas the latter is more practical in nature. 
