# How to estimate normalization constant during optimization of complex parameterized PDF using MLE?

If I have some data in $\mathbb{R}^N$ space, with $N$ beiling large and I want to estimate the density function of this data, $\mathbf{P}_{data}(x)$, using the fantastic property of neural networks that allows approximation of very complex functions...I will then start with some parameterized formulation $\mathbf{P}_{model}(x|{\theta})$, with ${\theta}$ being the weights in the neural network.

It seems of course I could use the traditional MLE formulation to estimate model parameters, ${\theta}^* = \mathrm{argmax}_{\theta}\prod_{i=1}^m \left[\mathbf{P}_{model}(x|{\theta})\right]$.

The problem I see though is that MLE traditionally relies on a known PDF being used....satisfying the property that $\int_x\mathbf{P}_{model}=1$. Thus I'm not sure how to do an optimization where the normalizing constant changes throughout the optimization, and/or the integral of $\int_x\mathbf{P}_{model}$ is very hard. I also can't quite see how the normalizing constant could be pulled in to the backprop equation for weight updates when it is seems like a separate constraint (unless this can be formulated as a constrained optimization problem)?

Doing some reading, I've seen some hints that importance sampling (or some kind of MCMC method) could be used here, but I've not been able to figure out how to implement it, or if its even tractable on a typical deep learning machine (>= 8 cores, 1 high end GPU, >100 GB RAM, etc).

How can I solve this problem...is it even possible?

There are many algorithms for training "unnormalized" models (e.g., score matching, minimum probability flow, noise contrastive estimation, contrastive divergence). In the following, I explain the basic idea behind persistent contrastive divergence (PCD).

Assuming our distribution is positive, we express it in terms of an energy function $E$ and normalization constant $Z$:

$$\textbf{P}_\text{model}(x \mid \theta) = \frac{e^{-E_\theta(x)}}{Z_\theta}$$

The maximum likelihood gradient is \begin{align} \nabla \mathbb{E}_\text{data}[\log \textbf{P}_\text{model}(x \mid \theta)] &= \mathbb{E}_\text{data}[-\nabla E_\theta(x)] - \nabla \log Z_\theta \\ &= \mathbb{E}_\text{model}[\nabla E_\theta(x)] - \mathbb{E}_\text{data}[\nabla E_\theta(x)] \end{align} since $$\nabla \log Z_\theta = \frac{1}{Z_\theta} \nabla \int e^{-E_\theta(x)} \, dx = - \int \frac{e^{-E_\theta(x)}}{Z_\theta} \nabla E_\theta(x) \, dx = -\mathbb{E}_\text{model}[\nabla E_\theta(x)].$$ The problem with this gradient is that we cannot easily estimate an expectation over the model distribution, since it is difficult to sample from an unnormalized model. We can try to run an MCMC sampler to sample from the model, but this can take a long time and we have to do it for every gradient step. So in PCD we initialize our sampler with the samples from the last gradient step, which if the model hasn't changed much should already be close to the desired distribution. We then only apply one or a few MCMC updates before using the samples. A commonly used sampler is Hamiltonian Monte Carlo.

One example paper which used deeper nets to represent the energy function and which used PCD is Learning Deep Energy Models (Ngiam et al., ICML, 2011).

• I'm still carefully reading the references and will likely have follow up questions, but this is most definitely the answer I was looking for.
– JPJ
Commented Nov 19, 2017 at 19:48
• Found out Ian Goodfellow has about 25 pages specifically devoted to this in his e-book here
– JPJ
Commented Nov 22, 2017 at 6:46
• The way I understand this after examining resources.....the negative phase, $\nabla \log Z_\theta$, can't be solved by importance sampling since we don't know $Z_\theta$, but we leverage the fact that MCMC works on unnormalized models and can still draw samples from $e^{-E_\theta(x)}$, correct?
– JPJ
Commented Dec 6, 2017 at 5:05
• There is a biased form of importance sampling you could use to estimate the second term in principle, but it won't work well for any interesting, high-dimensional model. There's attempts to continuously track $Z$ during training as well (Desjardins et al., 2011), but the question of a good proposal distribution remains. Commented Dec 6, 2017 at 17:02

This is just supplementary information (that won't fit in a comment) after studying the topic much further based on the accepted answer.

I've seen this particular question come up extensively under the area of Restricted Boltzmann Machines (RBMs) with physics analogies like here.

There is still an issue with tracking progress in such a scenario because the partition function $$Z_\theta$$ still cannot be estimated within reason. Directly from the link:

Tracking Progress

RBMs are particularly tricky to train. Because of the partition function Z of Eq. (1), we cannot estimate the log-likelihood $$\log(P(x))$$ during training. We therefore have no direct useful metric for choosing the optimal hyperparameters.

Several options are available to the user.

Inspection of Negative Samples

Negative samples obtained during training can be visualized. As training progresses, we know that the model defined by the RBM becomes closer to the true underlying distribution, $$p_{train}(x)$$. Negative samples should thus look like samples from the training set. Obviously bad hyperparameters can be discarded in this fashion.

Visual Inspection of Filters

The filters learnt by the model can be visualized. This amounts to plotting the weights of each unit as a gray-scale image (after reshaping to a square matrix). Filters should pick out strong features in the data. While it is not clear for an arbitrary dataset, what these features should look like, training on MNIST usually results in filters which act as stroke detectors, while training on natural images lead to Gabor like filters if trained in conjunction with a sparsity criteria.

Proxies to Likelihood

Other, more tractable functions can be used as a proxy to the likelihood. When training an RBM with PCD, one can use pseudo-likelihood as the proxy. Pseudo-likelihood (PL) is much less expensive to compute, as it assumes that all bits are independent.

Others have used annealed importance sampling (AIS) to estimate the partition function directly, although I'm not exactly sure how well it scales to high dimensions it is claimed to be invented for such cases.

DBNs

AIS