Neural Autoregressive Flows

Introduction

Invertible transformations with a tractable Jacobian, also known as normalizing flows, are useful tools in many machine learning problems, for example: (1) In the context of deep generative models, training necessitates evaluating data samples under the model’s inverse transformation (Dinh et al., 2016). Tractable density is an appealing property for these models, since it allows the objective of interest to be directly optimized; whereas other mainstream methods rely on alternative losses, in the case of intractable density models (Kingma & Welling, 2013; Rezende et al., 2014), or implicit losses, in the case of adversarial models (Goodfellow et al., 2014). (2) In the context of variational inference (Rezende & Mohamed, 2015), they can be used to improve the variational approximation to the posterior by parameterizing more complex distributions. This is important since a poor variational approximation to the posterior can fail to reflect the right amount of uncertainty, and/or be biased (Turner & Sahani, 2011), resulting in inaccurate and unreliable predictions. We are thus interested in improving techniques for normalizing flows.

Recent work by Kingma et al. (2016) reinterprets autoregressive models as invertible transformations suitable for constructing normalizing flows. The inverse transformation process, unlike sampling from the autoregressive model, is not sequential and thus can be accelerated via parallel computation. This allows multiple layers of transformations to be stacked, increasing expressiveness for better variational inference (Kingma et al., 2016) or better density estimation for generative models (Papamakarios et al., 2017). Stacking also makes it possible to improve on the sequential conditional factorization assumed by autoregressive models such as PixelRNN or PixelCNN (Oord et al., 2016), and thus define a more flexible joint probability.

We note that the normalizing flow introduced by Kingma et al. (2016) only applies an affine transformation of each scalar random variable. Although this transformation is conditioned on preceding variables, the resulting flow can still be susceptible to bad local minima, and thus failure to capture the multimodal shape of a target density; see Figure 1 and 2.

We propose replacing the conditional affine transformation of Kingma et al. (2016) with a more rich family of transformations, and note the requirements for doing so. We determine that very general transformations, for instance parametrized by deep neural networks, are possible. We then propose and evaluate several specific monotonic neural network architectures which are more suited for learning multimodal distributions. Concretely, our method amounts to using an autoregressive model to output the weights of multiple independent transformer networks, each of which operates on a single random variable, replacing the affine transformations of previous works.

Empirically, we show that our method works better than the state-of-the-art affine autoregressive flows of Kingma et al. (2016) and Papamakarios et al. (2017), both as a sample generator which captures multimodal target densities with higher fidelity, and as a density model which more accurately evaluates the likelihood of data samples drawn from an unknown distribution.

We also demonstrate that our method is a universal approximator on proper distributions in real space, which guarantees the expressiveness of the chosen parameterization and supports our empirical findings.

Background

A (finite) normalizing flow (NF), or flow, is an invertible function $f_{\theta}:\mathcal{X}\rightarrow\mathcal{Y}$ used to express a transformation between random variables We use $\operatorname{\mathbf{x}}$ and $\operatorname{\mathbf{y}}$ to denote inputs and outputs of a function, not the inputs and targets of a supervised learning problem. . Since $f$ is invertible, the change of variables formula can be used to translate between densities $p_{Y}(\operatorname{\mathbf{y}})$ and $p_{X}(\operatorname{\mathbf{x}})$:

The determinant of $f$’s Jacobian appears on the right hand side to account for the way in which $f$ can (locally) expand or contract regions of $X$, thereby lowering or raising the resulting density in those regions’ images in $Y$. Since the composition of invertible functions is itself invertible, complex NFs are often formed via function composition (or “stacking”) of simpler NFs.

Applying the change of variables formula from Equation 1 to the right hand side of Equation 2 yields:

Thus for efficient training, the following operations must be tractable and cheap:

Sampling $\operatorname{\mathbf{x}}\sim p_{X}(\operatorname{\mathbf{x}})$

Computing $\operatorname{\mathbf{y}}=f(\operatorname{\mathbf{x}})$

Computing the gradient of the log-determinant of the Jacobian of $f$

Research on constructing NFs, such as our work, focuses on finding ways to parametrize flows which meet the above requirements while being maximally flexible in terms of the transformations which they can represent. Note that some of the terms of of Equation 3 may be constant with respect to $\theta$ There might be some other parameters other than $\theta$ that are learnable, such as parameters of $p_{X}$ and $p_{target}$ in the variational inference and maximum likelihood settings, respectively. and thus trivial to differentiate, such as $p_{X}(\operatorname{\mathbf{x}})$ in the maximum likelihood setting.

Affine autoregressive flows (AAFs) Our terminology differs from previous works, and hence holds the potential for confusion, but we believe it is apt. Under our unifying perspective, NAF, IAF, AF, and MAF all make use of the same principle, which is an invertible transformer conditioned on the outputs of an autoregressive (and emphatically not an inverse autoregressive) conditioner. , such as inverse autoregressive flows (IAF) (Kingma et al., 2016), are one particularly successful pre-existing approach. Affine autoregressive flows yield a triangular Jacobian matrix, so that the log-determinant can be computed in linear time, as the sum of the diagonal entries on log scale. In AAFs, the components of $\operatorname{\mathbf{x}}$ and $\operatorname{\mathbf{y}}$ are given an order (which may be chosen arbitrarily), and $y_{t}$ is computed as a function of $x_{1:t}$. Specifically, this function can be decomposed via an autoregressive conditioner, $c$, and an invertible transformer, $\tau$, as Dinh et al. (2014) use $m$ and $g^{-1}$ to denote $c$ and $\tau$, and refer to them as the “coupling function” and “coupling law”, respectively.:

It is possible to efficiently compute the output of $c$ for all $t$ in a single forward pass using a model such as MADE (Germain et al., 2015), as pointed out by Kingma et al. (2016).

with $\sigma$ produced by an exponential nonlinearity. Kingma et al. (2016) use:

with $\sigma$ produced by a sigmoid nonlinearity. Such transformers are trivially invertible, but their relative simplicity also means that the expressivity of $f$ comes entirely from the complexity of $c$ and from stacking multiple AAFs (potentially using different orderings of the variables) Permuting the order of variables is itself a normalizing flow that does not expand or contract space and can be inverted by another permutation. . However, the only requirements on $\tau$ are:

The transformer $\tau$ must be invertible as a function of $x_{t}$.

$\frac{dy_{t}}{dx_{t}}$ must be cheap to compute.

This raises the possibility of using a more powerful transformer in order to increase the expressivity of the flow.

Neural Autoregressive Flows

We propose replacing the affine transformer used in previous works with a neural network, yielding a more rich family of distributions with only a minor increase in computation and memory requirements. Specifically,

is a deep neural network which takes the scalar $x_{t}$ as input and produces $y_{t}$ as output, and its weights and biases are given by the outputs of $c(x_{1:t-1})$ We’ll sometimes write $\tau_{c}$ for $\tau(c(x_{1:t-1}),\cdot)$. (see Figure 4(a)). We refer to these values $\phi$ as pseudo-parameters, in order to distinguish them from the statistical parameters of the model.

We now state the condition for NAF to be strictly monotonic, and thus invertible (as per requirement 1):

Using strictly positive weights and strictly monotonic activation functions for $\tau_{c}$ is sufficient for the entire network to be strictly monotonic.

Meanwhile, $\frac{dy_{t}}{dx_{t}}$ and gradients wrt the pseudo-parameters Gradients for pseudo-parameters are backpropagated through the conditioner, $c$, in order to train its parameters. can all be computed efficiently via backpropagation (as per requirement 2).

Whereas affine transformers require information about multimodality in $y_{t}$ to flow through $x_{1:t-1}$, our neural autoregressive flows (NAFs) are able to induce multimodality more naturally, via inflection points in $\tau_{c}$, as shown in Figure 5. Intuitively, $\tau_{c}$ can be viewed as analogous to a cumulative distribution function (CDF), so that its derivative corresponds to a PDF, where its inflection points yield local maxima or minima.

In this work, we use two specific architectures for $\tau_{c}$, which we refer to as deep sigmoidal flows (DSF) and deep dense sigmoidal flows (DDSF) (see Figure 4(b), 4(c) for an illustration). We find that small neural network transformers of 1 or 2 hidden layers with 8 or 16 sigmoid units perform well across our experiments, although there are other possibilities worth exploring (see Section 3.3). Sigmoids contain inflection points, and so can easily induce inflection points in $\tau_{c}$, and thus multimodality in $p(y_{t})$. We begin by describing the DSF transformation, which is already sufficiently expressive to form a universal approximator for probability distributions, as we prove in section 4.

The DSF transformation resembles an MLP with a single hidden layer of sigmoid units. Naive use of sigmoid activation functions would restrict the range of $\tau_{c}$, however, and result in a model that assigns 0 density to sufficiently large or small $y_{t}$, which is problematic when $y_{t}$ can take on arbitrary real values. We address this issue by applying the inverse sigmoid (or “logit”) function at the output layer. To ensure that the output’s preactivation is in the domain of the logit (that is, $(0,1)$), we combine the output of the sigmoid units via an attention-like (Bahdanau et al., 2014) softmax-weighted sums:

Since all of the sigmoid activations are bounded between 0 and 1, the final preactivation (which is their convex combination) is as well. The complete DSF transformation can be seen as mapping the original random variable to a different space through an activation function, where doing affine/linear operations is non-linear with respect to the variable in the original space, and then mapping it back to the original space through the inverse activation.

When stacking multiple sigmoidal transformation, we realize it resembles an MLP with bottleneck as shown by the bottom left of Figure 4. A more general alternative is the deep dense sigmoidal flow (DDSF), which takes the form of a fully connected MLP:

for $1\leq l\leq L$ where $h_{0}=x$ and $y=h_{L}$; $d_{0}=d_{L}=1$. We also require $\sum_{j}w_{ij}=1$, $\sum_{j}u_{kj}=1$ for all $i,k$, and all parameters except $b$ to be positive.

We use either DSF (Equation 8) or DDSF (Equation 19) to define the transformer function $\tau$ in Equation 4. To compute the log-determinant of Jacobian in a numerically stable way, we need to apply log-sum-exp to the chain rule

We elaborate more on the numerical stability in parameterization and computation of logarithmic operations in the supplementary materials.

2 Efficient Parametrization of Larger Transformers

Multi-layer NAFs, such as DDSF, require $c$ to output $\mathcal{O}(d^{2})$ pseudo-parameters, where $d$ is the number of hidden units in each layer of $\tau$. As this is impractical for large $d$, we propose parametrizing $\tau$ with $\mathcal{O}(d^{2})$ statistical parameters, but only $\mathcal{O}(d)$ pseudo-parameters which modulate the computation on a per-unit basis, using a technique such as conditional batch-normalization (CBN) (Dumoulin et al., 2016). Such an approach also makes it possible to use minibatch-style matrix-matrix products for the forward and backwards propagation through the graph of $\tau_{c}$. In particular, we use a technique similar to conditional weight normalization (CWN) (Krueger et al., 2017) in our experiments with DDSF; see appendix for details.

3 Possibilities for Alternative Architectures

Finally, we emphasize that in general, $\tau$ need not be expressed as a neural architecture; it only needs to satisfy the requirements of invertibility and differentiability given at the end of section 2.

NAFs are Universal Density Approximators

In this section, we prove that NAFs (specifically DSF) can be used to approximate any probability distribution over real vectors arbitrarily well, given that $\tau_{c}$ has enough hidden units output by generic neural networks with autoregressive conditioning. Ours is the first such result we are aware of for finite normalizing flows.

Our result builds on the work of Huang et al. (2017), who demonstrate the general universal representational capability of inverse autoregressive transformations parameterized by an autoregressive neural network (that transform uniform random variables into any random variables in reals). However, we note that their proposition is weaker than we require, as there are no constraints on the parameterization of the transformer $\tau$, whereas we’ve constrained $\tau$ to have strictly positive weights and monotonic activation functions, to ensure it is invertible throughout training.

The idea of proving the universal approximation theorem for DSF (1) in the IAF direction (which transforms unstructured random variables into structured random variables) resembles the concept of the inverse transform sampling: we first draw a sample from a simple distribution, such as uniform distribution, and then pass the sample though DSF. If DSF converges to any inverse conditional CDF, the resulting random variable then converges in distribution to any target random variable as long as the latter has positive continuous probability density everywhere in the reals. (2) For the MAF direction, DSF serves as a solution to the non-linear independent component analysis problem (Hyvärinen & Pajunen, 1999), which disentangles structured random variables into uniformly and independently distributed random variables. (3) Combining the two, we further show that DSF can transform any structured noise variable into a random variable with any desired distribution.

We define the following notation for the pre-logit of the DSF transformation (compare equation 8):

where $\mathcal{C}=(w_{j},b_{j},\tau_{j})_{j=1}^{n}$ are functions of $x_{1:1-t}$ parameterized by neural networks. Let $b_{j}$ be in $(r_{0},r_{1})$; $\tau_{j}$ be bounded and positive; $\sum_{j=1}^{n}w_{j}=1$ and $w_{j}>0$. See Appendix F and G for the proof.

where $\operatorname{\mathcal{C}}_{t}=(a_{tj},b_{tj},\tau_{tj})_{j=1}^{n}$ are functions of $x_{1:t-1}$, such that $Y_{n}\doteq G_{n}(X)$ converges in distribution to $Y$.

where $\operatorname{\mathcal{C}}_{t}=(a_{tj},b_{tj},\tau_{tj})_{j=1}^{n}$ are functions of $x_{1:t-1}$, such that $Y_{n}\doteq H_{n}(X)$ converges in distribution to $Y$.

where $\operatorname{\mathcal{C}}_{t}=(a_{tj},b_{tj},\tau_{tj})_{j=1}^{n}$ are functions of $x_{1:t-1}$, such that $Y_{n}\doteq K_{n}(X)$ converges in distribution to $Y$.

Related work

Neural autoregressive flows are a generalization of the affine autoregressive flows introduced by Kingma et al. (2016) as inverse autoregressive flows (IAF) and further developed by Chen et al. (2016) and Papamakarios et al. (2017) as autoregressive flows (AF) and masked autoregressive flows (MAF), respectively; for details on their relationship to our work see Sections 2 and 3. While Dinh et al. (2014) draw a particular connection between their NICE model and the Neural Autoregressive Density Estimator (NADE) (Larochelle & Murray, 2011), (Kingma et al., 2016) were the first to highlight the general approach of using autoregressive models to construct normalizing flows. Chen et al. (2016) and then Papamakarios et al. (2017) subsequently noticed that this same approach could be used efficiently in reverse when the key operation is evaluating, as opposed to sampling from, the flow’s learned output density. Our method increases the expressivity of these previous approaches by using a neural net to output pseudo-parameters of another network, thus falling into the hypernetwork framework (Ha et al., 2016; Bertinetto et al., 2016; Brabandere et al., 2016).

There has been a growing interest in normalizing flows (NFs) in the deep learning community, driven by successful applications and structural advantages they have over alternatives. Rippel & Adams (2013), Rezende & Mohamed (2015) and Dinh et al. (2014) first introduced normalizing flows to the deep learning community as density models, variational posteriors and generative models, respectively. In contrast to traditional variational posteriors, NFs can represent a richer family of distributions without requiring approximations (beyond Monte Carlo estimation of the KL-divergence). The NF-based RealNVP-style generative models (Dinh et al., 2016, 2014) also have qualitative advantages over alternative approaches. Unlike generative adversarial networks (GANs) (Goodfellow et al., 2014) and varational autoencoders (VAEs) (Kingma & Welling, 2013; Rezende et al., 2014), computing likelihood is cheap. Unlike autoregressive generative models, such as pixelCNNs (Oord et al., 2016), sampling is also cheap. Unfortunately, in practice RealNVP-style models are not currently competitive with autoregressive models in terms of likelihood, perhaps due to the more restricted nature of the transformations they employ.

Several promising recent works expand the capabilities of NFs for generative modeling and density estimation, however. Perhaps the most exciting example is Oord et al. (2017), who propose the probability density distillation technique to train an IAF (Kingma et al., 2016) based on the autoregressive WaveNet (van den Oord et al., 2016) as a generative model using another pretrained WaveNet model to express the target density, thus overcoming the slow sequential sampling procedure required by the original WaveNet (and characteristic of autoregressive models in general), and reaching super-real-time speeds suitable for production. The previously mentioned MAF technique (Papamakarios et al., 2017) further demonstrates the potential of NFs to improve on state-of-the-art autoregressive density estimation models; such highly performant MAF models could also be “distilled” for rapid sampling using the same procedure as in Oord et al. (2017).

Other recent works also find novel applications of NFs, demonstrating their broad utility. Loaiza-Ganem et al. (2017) use NFs to solve maximum entropy problems, rather than match a target distribution. Louizos & Welling (2017) and Krueger et al. (2017) apply NFs to express approximate posteriors over parameters of neural networks. Song et al. (2017) use NFs as a proposal distribution in a novel Metropolis-Hastings MCMC algorithm.

Finally, there are also several works which develop new techniques for constructing NFs that are orthogonal to ours (Tomczak & Welling, 2017, 2016; Gemici et al., 2016; Duvenaud et al., 2016; Berg et al., 2018).

Experiments

Our experiments evaluate NAFs on the classic applications of variational inference and density estimation, where we outperform IAF and MAF baselines. We first demonstrate the qualitative advantage NAFs have over AAFs in energy function fitting and density estimation (Section 6.1). We then demonstrate the capability of NAFs to capture a multimodal Bayesian posterior in a limited data setting (Section 6.2). For larger-scale experiments, we show that using NAF instead of IAF to approximate the posterior distribution of latent variables in a variational autoencoder (Kingma & Welling, 2013; Rezende et al., 2014) yields better likelihood results on binarized MNIST (Larochelle & Murray, 2011) (Section 6.3). Finally, we report our experimental results on density estimation of a suite of UCI datasets (Section 6.4).

First, we demonstrate that, in the case of marginally independent distributions, affine transformation can fail to fit the true distribution. We consider a mixture of Gaussian density estimation task. We define the modes of the Gaussians to be laid out on a 2D meshgrid within the range , and consider 2, 5 and 10 modes on each dimension. While the affine flow only produces a single mode, the neural flow matches the target distribution quite well even up to a 10x10 grid with 100 modes (see Figure 6).

1.2 Convergence

We then repeat the experiment that produces Figure 1 and 2 16 times, smooth out the learning curve and present average convergence result of each model with its corresponding standard deviation. For affine flow, we stack 6 layers of transformation with reversed ordering. For DSF and DDSF we used one transformation. We set $d=16$ for both, $L=2$ for DDSF.

2 Sine Wave experiment

Here we demonstrate the ability of DSF to capture multimodal posterior distributions. To do so, we create a toy experiment where the goal is to infer the posterior over the frequency of a sine wave, given only 3 datapoints. We fix the form of the function as $y(t)=\sin(2\pi f\cdot t)$ and specify a Uniform prior over the frequency: $p(f)=U()$. The task is to infer the posterior distribution $p(f|T,Y)$ given the dataset (T,Y)=((0,\ThisStyle{5\mathord{\stretchto{\raisebox{2.3pt}{\SavedStyle/}}{}}6},\ThisStyle{10\mathord{\stretchto{\raisebox{2.3pt}{\SavedStyle/}}{}}6}),(0,0,0)), as represented by the red crosses of Figure 8 (left). We assume the data likelihood given the frequency parameter to be $p(y_{i}|t_{i},f)=\mathcal{N}(y_{i};y_{f}(t_{i}),0.125)$, where the variance $\sigma^{2}=0.125$ represents the inherent uncertainty of the data. Figure 8 (right) shows that DSF learns a good posterior in this task.

3 Amortized Approximate Posterior

We evaluate NAF’s ability to improve variational inference, in the context of the binarized MNIST (Larochelle & Murray, 2011) benchmark using the well-known variational autoencoder (Kingma & Welling, 2013; Rezende et al., 2014) (Table 1). Here again the DSF architecture outperforms both standard IAF and the traditional independent Gaussian posterior by a statistically significant margin.

4 Density Estimation with Masked Autoregressive Flows

We replicate the density estimation experiments of Papamakarios et al. (2017), which compare MADE (Germain et al., 2015) and RealNVP (Dinh et al., 2016) to their proposed MAF model (using either 5 or 10 layers of MAF) on BSDS300 (Martin et al., 2001) as well as 4 UCI datasets (Lichman, 2013) processed as in Uria et al. (2013). Simply replacing the affine transformer with our DDSF architecture in their best performing architecture for each task (keeping all other settings fixed) results in substantial performance gains, and also outperforms the more recent Transformation Autoregressive Networks (TAN) Oliva et al. (2018), setting a new state-of-the-art for these tasks. Results are presented in Table 2.

Conclusion

In this work we introduce the neural autoregressive flow (NAF), a flexible method of tractably approximating rich families of distributions. In particular, our experiments show that NAF is able to model multimodal distributions and outperform related methods such as inverse autoregressive flow in density estimation and variational inference. Our work emphasizes the difficulty and importance of capturing multimodality, as previous methods fail even on simple toy tasks, whereas our method yields significant improvements in performance.

Acknowledgements

We would like to thank Tegan Maharaj, Ahmed Touati, Shawn Tan and Giancarlo Kerg for helpful comments and advice. We also thank George Papamakarios for providing details on density estimation task’s setup.

References

Appendix A Exclusive KL View of the MLE

The maximum likelihood principle requires us to minimize the following metric:

This means we want to transform the empirical distribution $p_{data}$, or $p_{X}$, to fit the target density (the usually unstructured, base distribution $p_{0}$), as explained in Section 2.

Appendix B Monotonicity of NAF

Here we show that using strictly positive weights and strictly monotonically increasing activation functions is sufficient to ensure strict monotonicity of NAF. For brevity, we write monotonic, or monotonicity to represent the strictly monotonically increasing behavior of a function. Also, note that a continuous function is strictly monotonically increasing exactly when its derivative is greater than 0.

Suppose we have an MLP with $L+1$ layers: $h_{0},h_{1},h_{2},...,h_{L}$, $x=h_{0}$ and $y=h_{L}$, where $x$ and $y$ are scalar, and $h_{l,j}$ denotes the $j$-th node of the $l$-th layer.For $1\leq l\leq L$, we have

for some monotonic activation function $A_{l}$, positive weight vector $w_{l,j}$, and bias $b_{l,j}$. Differentiating this yields

which is greater than 0 since $A$ is monotonic (and hence has positive derivative), and $w_{l,j,k}$ is positive by construction. Thus any unit of layer $l$ is monotonic with respect to any unit of layer $l-1$, for $1\leq l\leq L$.

Now, suppose we have $J$ monotonic functions $f_{j}$ of the input $x$. Then the weighted sum of these functions $\sum_{j=1}^{J}u_{j}f_{j}$ with $u_{j}>0$ is also monotonic with respect to $x$, since

Finally, we use induction to show that all $h_{l,j}$ (including $y$) are monotonic with respect to $x$.

The base case ($l=1$) is given by Equation 27.

Suppose the inductive hypothesis holds, which means $h_{l,j}$ is monotonic with respect to $x$ for all $j$ of layer $l$. Then by Equation 28, $h_{l+1,k}$ is also monotonic with respect to $x$ for all k.

Thus by mathematical induction, monotonicity of $h_{l,j}$ holds for all $l$ and $j$. ∎

Appendix C Log Determinant of Jacobian

As we mention at the end of Section 3.1, to compute the log-determinant of the Jacobian as part of the objective function, we need to handle the numerical stability. We first derive the Jacobian of the DDSF (note that DSF is a special case of DDSF), and then summarize the numerically stable operations that were utilized in this work.

Again defining $x=h_{0}$ and $y=h_{L}$, the Jacobian of each DDSF transformation can be written as a sequence of dot products due to the chain rule:

For notational convenience, we define a few more intermediate variables. For each layer of DDSF, we have

C.2 Numerically Stable Operations

Since the Jacobian of each DDSF transformation is chain of dot products (Equation 29), with some nested multiplicative operations (Equation 53), we calculate everything in the log-scale (where multiplication is replaced with addition) to avoid unstable gradient signal.

To ensure the summing-to-one and positivity constraints of $u$ and $w$, we let the autoregressive conditioner output pre-activation $u\_$ and $w\_$, and apply softmax to them. We do the same for $a$ by having the conditioner output $a\_$ and apply softplus to ensure positivity. In Equation 53, we have

where $x^{*}=\max_{i}\{x_{i}\}$ and $\delta$ is a small value such as $10^{-6}$.

C.2.2 Logarithmic dot product

Appendix D Scalability and Parameter Sharing

As discussed in Section 3.3, a multi-layer NAF such as DDSF requires the autoregressive conditioner $c$ to output many pseudo-parameters, on the order of $\mathcal{O}(Ld^{2})$, where $L$ is the number of layers of the transformer network ($\tau$), and $d$ is the average number of hidden units per layer. In practice, we reduce the number of outputs (and thus the computation and memory requirements) of DDSF by instead endowing $\tau$ with some learned (non-conditional) statistical parameters. Specifically, we decompose $w\_$ and $u\_$ (the preactivations of $\tau$’s weights, see section C.2.1) into pseudo-parameters and statistical parameters. Take $u\_$ for example:

where $v^{(l+1)}$ is a $d_{l+1}\times d_{l}$ matrix of statistical parameters, and $\eta$ is output by $c$. See figure 9 for a depiction.

The linear transformation before applying sigmoid resembles conditional weight normalization (CWN) (Krueger et al., 2017). While CWN rescales the weight vectors normalized to have unit L2 norm, here we rescale the weight vector normalized by softmax such that it sums to one and is positive. We call this conditional normalized weight exponentiation. This reduces the number of pseudo-parameters to $\mathcal{O}(Ld)$.

Appendix E Identity Flow Initialization

In many cases, initializing the transformation to have a minimal effect is believed to help with training, as it can be thought of as a warm start with a simpler distribution. For instance, for variational inference, when we initialize the normalizing flow to be an identity flow, the approximate posterior is at least as good as the input distribution (usually a fully factorized Gaussian distribution) before the transformation. To this end, for DSF and DDSF, we initialize the pseudo-weights $a$ to be close to $1$, the pseudo-biases $b$ to be close to .

This is achieved by initializing the conditioner (whose outputs are the pseudo-parameters) to have small weights and the appropriate output biases. Specifically, we initialize the output biases of the last layer of our MADE (Germain et al., 2015) conditioner to be zero, and add $\operatorname{softplus}^{-1}(1)\approx 0.5413$ to the outputs of which correspond to $a$ before applying the softplus activation function. We initialize all conditioner’s weights by sampling from from $\textnormal{Unif}(-0.001,0.001)$. We note that there might be better ways to initialize the weights to account for the different numbers of active incoming units.

Appendix F Lemmas: Uniform Convergence of DSF

We want to show the convergence result of Equation 21. To this end, we first show that DSF can be used to universally approximate any strictly monotonic function. The is the case where $x_{1:t-1}$ are fixed, which means $\operatorname{\mathcal{C}}(x_{1:t-1})$ are simply constants. We demonstrate it using the following two lemmas.

(Step functions universally approximate monotonic functions) Define:

For brevity, we write $s_{j}(x)=s(x-b_{j})$. For any $\epsilon>0$, we choose $n=\lceil\frac{1}{\epsilon}\rceil$, and divide the range $(0,1)$ into $n+1$ evenly spaced intervals: $(0,y_{1})$, $(y_{1},y_{2})$, $...$, $(y_{n},1)$. For each $y_{j}$, there is a corresponding inverse value since $S$ is strictly monotonic, $x_{j}=S^{-1}(y_{j})$. We want to set $S_{n}^{*}(x_{j})=y_{j}$ for $1\leq j\leq n-1$ and $S_{n}^{*}(x_{n})=1$. To do so, we set the bias terms $b_{j}$ to be $x_{j}$. Then we just need to solve a system of $n$ linear equations $\sum_{j^{\prime}=1}^{n}w_{j^{\prime}}\cdot s_{j^{\prime}}(x_{j})=t_{j}$, where $t_{j}=y_{j}$ for $1\leq j<n$, $t_{0}=0$ and $t_{n}=1$. We can express this system of equations in the matrix form as $\mathbf{S}\mathbf{w}=\mathbf{t}$, where:

where $\delta_{\omega\geq\eta}=1$ whenever $\omega\geq\eta$ and $\delta_{\omega\geq\eta}=0$ otherwise. Then we have $\mathbf{w}=\mathbf{S}^{-1}\mathbf{t}$. Note that ${\bf S}$ is a lower triangular matrix, and its inverse takes the form of a Jordan matrix: $(\mathbf{S}^{-1})_{i,i}=1$ and $(\mathbf{S}^{-1})_{i+1,i}=-1$. Additionally, $t_{j}-t_{j-1}=\frac{1}{n+1}$ for $j=1,...,n-1$ and is equal to $\frac{2}{n+2}$ for $j=n$. We then have $S_{n}^{*}(x)=\mathbf{t}^{T}\mathbf{S}^{-T}\mathbf{s}(x)$, where $\mathbf{s}(x)_{j}=s_{j}(x)$; thus

where $C_{v}(z)=\sum_{k}\delta_{z\geq v_{k}}$ is the count of elements in a vector that $z$ is no smaller than. ∎

Note that the additional constraint that $\mathbf{w}$ lies on an $n-1$ dimensional simplex is always satisfied, because

See Figure 10 for a visual illustration of the proof. Using this result, we now turn to the case of using sigmoid functions instead of step functions.

(Superimposed sigmoids universally approximate monotonic functions) Define:

With the same constraints and definition in Lemma 1, given any $\epsilon>0$, there exists a positive integer $n$, real constants $w_{j}$, $\tau_{j}$ and $b_{j}$ for $j=1,...,n$, where additionally $\tau_{j}$ are bounded and positive, such that $|S_{n}(x)-S(x)|<\epsilon\quad\forall x\in(r_{0},r_{1})$.

Let $\epsilon_{1}=\frac{1}{3}\epsilon$ and $\epsilon_{2}=\frac{2}{3}\epsilon$. We know that for this $\epsilon_{1}$, there exists an $n$ such that $\left|S_{n}^{*}-S\right|<\epsilon_{1}$.

We chose the same $w_{j}$, $b_{j}$ for $j=1,...,n$ as the ones used in the proof of Lemma 1, and let $\tau_{1},...,\tau_{n}$ all be the same value denoted by $\tau$.

Take $\kappa=\min_{j\neq j^{\prime}}|b_{j}-b_{j^{\prime}}|$ and $\tau=\frac{\kappa}{\sigma^{-1}(1-\epsilon_{0})}$ for some $\epsilon_{0}>0$. Take $\Gamma$ to be a lower triangular matrix with values of 0.5 on the diagonal and 1 below the diagonal.

Since the product $\Gamma\cdot w$ represents the half step points of $S_{n}^{*}$ at $x=b_{j}$’s, the result above entails $\left|S_{n}(x)-S_{n}^{*}(x)\right|<\epsilon_{2}=2\epsilon_{1}$ for all $x$. To see this, we choose $\epsilon_{0}=\frac{1}{2(n+1)}$. Then $S_{n}$ intercepts with all segments of $S_{n}^{*}$ except for the ends. We choose $\epsilon_{2}=2\epsilon_{1}$ since the last step of $S_{n}^{*}$ is of size $\frac{2}{n+1}$, and thus the bound also holds true in the vicinity where $S_{n}$ intercepts with the last step of $S_{n}^{*}$.

Now we show that (the pre-logit) DSF (Equation 21) can universally approximate monotonic functions. We do this by showing that the well-known universal function approximation properties of neural networks (Cybenko, 1989) allow us to produce parameters which are sufficiently close to those required by Lemma 2.

where $t\in[1,m]$, and $\mathcal{C}_{t}=(w_{tj},b_{tj},\tau_{tj})_{j=1}^{n}$ are functions of $x_{1:1-t}$ parameterized by neural networks, with $\tau_{tj}$ bounded and positive, $b_{tj}\in[r_{0},r_{1}]$, $\sum_{j=1}^{n}w_{tj}=1$, and $w_{tj}>0$

The idea is first to show that we can find a sequence of parameters, $\operatorname{\mathcal{C}}_{n}(x_{1:t-1})$, that yield a good approximation of the target function, $S(x_{t},x_{1:t-1})$. We then show that these parameters can be arbitrarily well approximated by the outputs of a neural network, $\operatorname{\mathcal{C}}_{k}(x_{1:t-1})$, which in turn yield a good approximation of $S$.

From Lemma 2, we know that such a sequence $\operatorname{\mathcal{C}}_{n}(x_{1:t-1})$ exists, and furthermore that we can, for any $\epsilon$, and independently of $S$ and $x_{1:t-1}$ choose an $N$ large enough so that:

Now, $\operatorname{\mathcal{S}}$ has bounded derivative wrt $\operatorname{\mathcal{C}}$, and is thus uniformly continuous, so long as $\tau$ is greater than some positive constant, which is always the case for any fixed $\operatorname{\mathcal{C}}_{n}$, and thus can be guaranteed for $\operatorname{\mathcal{C}}_{k}$ as well (for large enough $k$). Uniform continuity allows us into translate the convergence of $\operatorname{\mathcal{C}}_{k}\rightarrow\operatorname{\mathcal{C}}_{n}$ to convergence of $\operatorname{\mathcal{S}}_{n}(x_{t},\operatorname{\mathcal{C}}_{k}(x_{1:t-1}))\rightarrow\operatorname{\mathcal{S}}_{n}(x_{t},\operatorname{\mathcal{C}}_{n}(x_{1:t-1}))$, since for any $\epsilon$, there exists a $\delta>0$ such that

Combining this with Equation 54, we have for all $x_{t}$ and $x_{1:t-1}$, and for all $n\geq N$ and $k\geq K$

Appendix G Proof of Universal Approximation of DSF

Given an arbitrary ordering, let $F$ be the CDFs of $Y$, defined as $F_{t}(y_{t},x_{1:t-1})=\Pr(Y_{t}\leq y_{t}|x_{1:t-1})$. According to Theorem 1 of Hyvärinen & Pajunen (1999), $F(Y)$ is uniformly distributed in the cube $^{m}$. $F$ has an upper triangular Jacobian matrix, whose diagonal entries are conditional densities which are positive by assumption. Let $G$ be a multivariate and multivariable function where $G_{t}$ is the inverse of the CDF of $Y_{t}$: $G_{t}(F_{t}(y_{t},x_{1:t-1}),x_{1:t-1})=y_{t}$.

According to Lemma 3, there exists a sequence of functions in the given form $(\operatorname{\mathcal{S}}_{n})_{n\geq 1}$ that converge uniformly to $\sigma\circ G$. Since uniform convergence implies pointwise convergence, $G_{n}=\sigma^{-1}\circ S_{n}$ converges pointwise to $G$, by continuity of $\sigma^{-1}$. Since $G_{n}$ converges pointwise to $G$ and $G(X)=Y$, by Lemma 4, we have $Y_{n}\xrightarrow{d}Y$

Given an arbitrary ordering, let $H$ be the CDFs of $X$:

Due to Hyvärinen & Pajunen (1999), $y_{1},...y_{m}$ are independently and uniformly distributed in $(0,1)^{m}$.

According to Lemma 3, there exists a sequence of functions in the given form $(\operatorname{\mathcal{S}}_{n})_{n\geq 1}$ that converge uniformly to $H$. Since $H_{n}=\operatorname{\mathcal{S}}_{n}$ converges pointwise to $H$ and $H(X)=Y$, by Lemma 4, we have $Y_{n}\xrightarrow{d}Y$

Given an arbitrary ordering, let $H$ be the CDFs of $X$ defined the same way in the proof for Proposition 3, and let $G$ be the inverse of the CDFs of $Y$ defined the same way in the proof for Proposition 2. Due to Hyvärinen & Pajunen (1999), $H(X)$ is uniformly distributed in $(0,1)^{m}$, so $G(H(X))=Y$. Since $H_{t}(x_{t},x_{1:t-1})$ is monotonic wrt $x_{t}$ given $x_{1:t-1}$, and $G_{t}(H_{t},H_{1:t-1})$ is monotonic wrt $H_{t}$ given $H_{1:t-1}$, $G_{t}$ is also monotonic wrt $x_{t}$ given $x_{1:t-1}$, as

According to Lemma 3, there exists a sequence of functions in the given form $(\operatorname{\mathcal{S}}_{n})_{n\geq 1}$ that converge uniformly to $\sigma\circ G\circ H$. Since uniform convergence implies pointwise convergence, $K_{n}=\sigma^{-1}\circ S_{n}$ converges pointwise to $G\circ H$, by continuity of $\sigma^{-1}$. Since $K_{n}$ converges pointwise to $G\circ H$ and $G(H(X))=Y$, by Lemma 4, we have $Y_{n}\xrightarrow{d}Y$

Appendix H Experimental Details

For the experiment of amortized variational inference, we implement the Variational Autoencoder (Kingma & Welling, 2013). Specifically, we follow the architecture used in Kingma et al. (2016): the encoder has three layers with $$feature maps. We use resnet blocks (He et al., 2016) with$3\times 3$convolution filters and a stride of$2$to downsize the feature maps. The convolution layers are followed by a fully connected layer of size$450$as a context for the flow layers that transform the noise sampled from a standard normal distribution of dimension$32. The decoder is symmetrical with the encoder, with the strided convolution replaced by a combination of bilinear upsampling and regular resnet convolution to double the feature map size. We used the ELUs activation function (Clevert et al., 2015) and weight normalization (Salimans & Kingma, 2016) in the encoder and decoder. In terms of optimization, Adam (Kingma et al., 2015) is used with learning rate fined tuned for each inference setting, and Polyak averaging (Polyak & Juditsky, 1992) was used for evaluation with\alpha=0.998$which stands for the proportion of the past at each time step. We also consider a variant of Adam known as Amsgrad (Reddi et al., 2018) as a hyperparameter. For vanilla VAE, we simply apply a resnet dot product with the context vector to output the mean and the pre-softplus standard deviation, and transform each dimension of the noise vector independently. We call this linear flow. For IAF-affine and IAF-DSF, we employ MADE (Germain et al., 2015) as the conditioner$c(x_{1:t-1})$, and we apply dot product on the context vector to output a scale vector and a bias vector to conditionally rescale and shift the preactivation of each layer of the MADE. Each MADE has one hidden layer with$1920$hidden units. The IAF experiments all start with a linear flow layer followed by IAF-affine or IAF-DSF transformations. For DSF, we choose$d=16$.

For the experiment of density estimation with MAF, we followed the implementation of Papamakarios et al. (2017). Specifically for each dataset, we experimented with both $5$ and $10$ flow layers, followed by one linear flow layer. The following table specifies the number of hidden layers and the number of hidden units per hidden layer for MADE: