GENERAL INCOMPRESSIBLE-FLOW NETWORKS (GIN)
Year: 2020

Abstract:
A central question of representation learning asks under which conditions it is possible to reconstruct the true latent variables of an arbitrarily complex generative process. Recent breakthrough work by Khemakhem et al. (2019) on nonlinear ICA has answered this question for a broad class of conditional generative processes. We extend this important result in a direction relevant for application to real-world data. First, we generalize the theory to the case of unknown intrinsic problem dimension and prove that in some special (but not very restrictive) cases, informative latent variables will be automatically separated from noise by an estimating model. Furthermore, the recovered informative latent variables will be in one-to-one correspondence with the true latent variables of the generating process, up to a trivial component-wise transformation. Second, we introduce a modification of the RealNVP invertible neural network architecture (Dinh et al., 2016) which is particularly suitable for this type of problem: the General Incompressibleflow Network (GIN). Experiments on artificial data and EMNIST demonstrate that theoretical predictions are indeed verified in practice. In particular, we provide a detailed set of exactly 22 informative latent variables extracted from EMNIST.

heading: 1 INTRODUCTION
text: Deep latent-variable models promise to unlock the key factors of variation within a dataset, opening a window to interpretation and granting the power to manipulate data in an intuitive fashion. The theory of identifiability in linear independent component analysis (ICA) (Comon, 1994) tells us when this is possible, if we restrict the model to a linear transformation, but until recently there was no corresponding theory for the highly nonlinear models needed to manipulate complex data. This changed with the recent breakthrough work by Khemakhem et al. (2019), which showed that under relatively mild conditions, it is possible to recover the joint data and latent space distribution, up to a simple transformation in the latent space. The key requirement is that the generating process is conditioned on a variable which is observed along with the data. This condition could be a class label, time index of a time series, or any other piece of information additional to the data. They interpret their theory as a nonlinear version of ICA.
This work extends this theory in a direction relevant for application to real-world data. The existing theory assumes knowledge of the intrinsic problem dimension, but this is unrealistic for anything but artificially generated datasets. Here, we show that in the special case of Gaussian latent space distributions, the intrinsic problem dimension can be discovered. The important latent variables are organically separated from noise variables by the estimating model. Furthermore, the variables discovered correspond to the true generating latent variables, up to a trivial component-wise translation and scaling. Very similar results exist for other members of the exponential family with two parameters, such as the beta and gamma distributions.
We introduce a variant of the RealNVP (Dinh et al., 2016) invertible neural network: the General Incompressible-flow Network (GIN). The flow is called incompressible in reference to fluid dynamics, since it preserves volumes: the Jacobian determinant is simply unity. We emphasise its generality and increased expressive power in comparison to previous volume-preserving flows, such as NICE (Dinh et al., 2014). As already noted in Khemakhem et al. (2019), flow-based generative models are a natural fit for the theory of nonlinear ICA, as are the variational autoencoders (VAEs) (Kingma &
Welling, 2013) used in that work. For us, major advantages of invertible architectures over VAEs are the ability to specify volume preservation and directly optimize the likelihood, and freedom from the requirement to specify the dimension of the model’s latent space. An INN always has a latent space of the same dimension as the data. In addition, the forward and backward models share parameters, saving the effort of learning separate models for each direction.
In summary, our work makes the following contributions:
• We extend the theory of nonlinear ICA to allow for unknown intrinsic problem dimension. Doing so, we find that this dimension can be discovered and a one-to-one correspondence between generating and estimated latent variables established.
• We propose as an implementation an invertible neural network obtained by modifying the RealNVP architecture. We call our new architecture GIN: the General Incompressible-flow Network.
• We demonstrate the viability of the model on artificial data and the EMNIST dataset. We extract 22 meaningful variables from EMNIST, encoding both global and local features.

heading: 2 RELATED WORK
text: The basic goals of nonlinear ICA stem from the original work on linear ICA. An influential formulation, as well as the first identifiability results, were given in Comon (1994). These stated the conditions which allow the generating latent variables to be discovered, when the mixing function is a linear transformation. However, it was shown in Hyvärinen & Pajunen (1999) that this approach to identifiability does not extend to general nonlinear functions.
The first identifiability results in nonlinear ICA came in Hyvärinen & Morioka (2016) and Hyvärinen & Morioka (2017), applied to time series, and implemented via a discriminative model and semisupervised learning. A more general formulation, valid for other forms of data, was given in Hyvärinen et al. (2018) and the theory was extended to generative models in Khemakhem et al. (2019), where experiments were implemented by a VAE.
Many authors have addressed the general problem of disentanglement, and proposed models to learn disentangled features. Prominent among these is β-VAE (Higgins et al., 2017) and its variations (e.g. Chen et al., 2018) which augment the standard ELBO loss with tunable hyperparameters to encourage disentanglement. There are also attempts to modify the GAN framework (Goodfellow et al., 2014) such as InfoGAN (Chen et al., 2016), which tries to maximize the mutual information between some dimensions of the latent space and the observed data. Many of these approaches are unsupervised. However, as pointed out and empirically demonstrated in Locatello et al. (2018), unsupervised models without conditioning in the latent space are in general unidentifiable.
Several unsupervised VAE models implement conditioning in the latent space by means of Gaussian mixtures (Johnson et al., 2016; Dilokthanakul et al., 2016; Zhao et al., 2019). Our work differs mainly by (i) only considering supervised tasks, therefore being safely covered by the theory of Khemakhem et al. (2019), and (ii) enforcing volume-preservation, not possible in a VAE.
The invertible neural networks in this work build upon the NICE framework (Dinh et al., 2014) and its extension in RealNVP (Dinh et al., 2016). A similar network design to ours is Ardizzone et al. (2019). This is a conditioned INN based on RealNVP, however the conditioning information is applied as a parameter of the network. The authors find in experiments with MNIST that meaningful variables are present in their latent space, but are rotated such that they are not aligned with the axes of the space. In this work, the conditioning is only present as a parameter of the latent space distributions. As a result, it is covered by the theory of Khemakhem et al. (2019) and its extension here, which results in non-rotated, meaningful latent variables.

heading: 3 THEORY

heading: 3.1 EXISTING THEORY OF NONLINEAR ICA
text: This section adapts theoretical results from Khemakhem et al. (2019) to the context of invertible neural networks. Suppose the existence of the following three random variables: a latent generating
variable z ∈ Rn, a condition u ∈ Rm and a data point x ∈ Rd, where n ≤ d. If n < d, also suppose the existence of a noise variable ε ∈ Rd−n. The variables z and ε make up the generating latent space, and are not necessarily known, whereas u and x are observable and therefore always known.
The distribution of z is a factorial member of the exponential family with k sufficient statistics, conditioned on u. In its most general form the distribution can be written as
pz(z|u) = n∏
i=1
Qi(zi)
Zi(u) exp  k∑ j=1 Ti,j(zi)λi,j(u)  (1) where the Ti,j are the sufficient statistics, the λi,j their coefficients and Zi the normalizing constant. Qi is called the base measure, which in many cases is simply 1.
The distribution of εmust not depend on u and its probability density function must be always finite.
The variable x is the result of an arbitrarily complex, invertible, deterministic transformation from the generating latent space to the data space: x = f(z, ε). This can alternatively be formulated as an injective, non-deterministic transformation from the lower-dimensional z-space to the higherdimensional x-space: x = fε(z).
In general, an observed dataset D will consist only of instances of x and u. The task of nonlinear ICA is to disentangle the data to recover the generating latent variables z, as well as the form of the function f and its inverse. We can try to achieve this with a sufficiently general, invertible function approximator g, which maps from a latent variable w ∈ Rd to the data space: x = g(w; θ). Here θ denotes the parameters of g. Note that the dimensions of the latent space and data space are the same due to the invertibility of g. We assume that w follows a conditionally independent exponential probability distribution, conditioned on the same condition u as z:
pw(w|u) = d∏
i=1
Q′i(wi)
Z ′i(u) exp  k∑ j=1 T ′i,j(wi)λ ′ i,j(u)  . (2) The coefficients of the sufficient statistics λ′i,j(u) are not restricted, which means that a given variable i of the estimated latent space is allowed to lose its dependence on u. If this occurs, we will consider this variable to encode noise, playing a role equivalent to ε in the generating latent space.
In addition to this specification of the generative and estimating models, some conditions are necessary to ensure the latent variables can be recovered. The most important of these concerns the variability of λi,j(u) under u. However, as long as the λi,j are randomly and independently generated, and there are at least nk + 1 distinct conditions u, this condition is almost surely fulfilled. See Appendix A.1 for further details.
In the limit of infinite data and perfect convergence, the estimating model will give the same conditional likelihood to all data points as the true generating model:
pT,λ,f,ε(x|u) = pT′,λ′,g(x|u). (3) If this is the case, the vector of sufficient statistics T from the generating latent space will be related to that of the estimated latent space by an affine transformation:
T(z) = AT′(w) + c, (4)
where A is some constant, full-rank, nk × dk matrix and c ∈ Rnk some constant vector. The relationship holds for all values of z and w. The proof of this result can be found in Appendix A.
This relationship is a generalization of that derived in Khemakhem et al. (2019). In that version, they assume knowledge of n, the dimension of the generating latent space, and give their estimating model a latent space of the same dimension. In this context the matrix A defines a relationship between latent spaces of the same dimension, making it square and invertible.

heading: 3.2 NONLINEAR ICA WITH A GAUSSIAN LATENT SPACE
text: Given the relationship in equation (4), we can ask whether any stronger results hold for particular special cases. We hope in particular to induce sparsity in the matrix A, so that the estimated latent space is related to the true one by as simple a transformation as possible. We show (proof in Appendix B) that when both the generating and estimated latent spaces follow a Gaussian distribution, the generating latent space variables are recovered up to a trivial translation and scaling. Furthermore, the dimension of the generating latent space is recovered.
More precisely, each generating latent variable zi is related to exactly one estimated latent variable wj , for some j, as zi = aiwj + bi, for some constant ai and bi. Furthermore, each estimated latent variable wj is related to at most one zi. If the estimating latent space has higher dimension than the dimension of the generating latent space, some estimating latent variables are not related to any generating latent variables and so must encode only noise. This is the dimension discovery mechanism, since the estimated latent space organically splits into informative and non-informative parts, with the dimension of the informative part equal to the unknown intrinsic dimension of the generating latent space. Very similar results can be derived for all common continuous two-parameter members of the exponential family, including the gamma and beta distributions. See Appendix C.

heading: 3.3 VOLUME PRESERVATION VIA INCOMPRESSIBLE FLOW
text: In a multivariate Gaussian distribution with diagonal covariance, the standard deviations are directly proportional to the principal axes of the ellipsoid defining a surface of constant probability. In classical PCA, this property is used to assign importance to each principal component based on its standard deviation. In PCA we can think of the rotation of the data into the basis defined by the principal components as the latent space. In a deep latent variable model the transformation between the latent space and the data is significantly more complex, but if this transformation preserves volumes like the rotation in PCA, it will retain the desirable correspondence between the standard deviation of a latent space variable and its importance in explaining the data. Because invertible neural networks (normalizing flows) have a tractable Jacobian determinant, they open the possibility to constrain it to unity. This is equivalent to volume preservation and is the rationale behind the General Incompressible-flow Network.

heading: 4 EXPERIMENTS
text: Experiments on artificial datasets confirm the theory for a normally distributed latent space, as well as identifying potential causes of failure. Experiments on the EMNIST dataset (Cohen et al., 2017) demonstrate the ability of GIN to estimate independent and interpretable latent variables from realworld data.

heading: 4.1 MODEL DESCRIPTION
text: The GIN model is similar in form to RealNVP (Dinh et al., 2016) and shares its flexibility, but retains the volume-preserving properties of the NICE framework (Dinh et al., 2014).
RealNVP coupling layers splitD-dimensional input x into two parts, x1:d and xd+1:D where d < D. The output of the layer is the concatentation of y1:d and yd+1:D with
y1:d = x1:d (5) yd+1:D = xd+1:D exp(s(x1:d)) + t(x1:d) (6)
where addition, multiplication and exponentiation are applied component-wise. The logarithm of the Jacobian determinant of a coupling layer is simply the sum of the scaling function s(x1:d). Volume preservation is achieved by setting the final component s(xd) to the negative sum of the previous components, so the total sum of s(x1:d) is zero (in contrast, NICE sets s(x1:d) ≡ 0). Hence the Jacobian determinant is unity and volume is preserved. The network is free to allow the volume contribution of some dimensions of the output of any coupling layer to grow, but only by shrinking the other dimensions of the output in direct proportion. As well as enforcing a strong correspondence between the importance of a latent variable and its standard deviation, we believe volume-preservation has a regularizing effect, since it is a very strong constraint, comparable to orthonormality in linear transformations.

heading: 4.2 OPTIMIZATION
text: The experiments in this section deal with labeled data, where each data point belongs to one of M different classes. This class label is used as the condition u associated with each data point x. In the estimated latent space, all data instances with the same label should belong to the same Gaussian distribution. Hence we are learning a Gaussian mixture in the estimated latent space, with M mixture components. Since the distribution for each class in the estimated latent space is required to be factorial, the variance of each mixture component is diagonal, and we can write σ2i (u) for the variance in the i-th dimension.
Given a set of data, condition pairs D = {(x(1), u(1)), . . . , (x(N), u(N))} and model g, parameterized by θ (where g maps from the latent space to the data space) we can construct a loss from the log-likelihood. Using the change of variables formula, we have log p(x|u) = log p(w|u), with no Jacobian term since the transformation w = g−1(x;θ) is volume-preserving. To maximize the likelihood of D, we minimize the negative log-likelihood of w in the estimated latent space:
L(θ) = E(x,u)∈D
[ 1
n n∑ i=1 ( (g−1i (x;θ)− µi(u;θ))2 2σ2i (u;θ) + log(σi(u;θ)) )] . (7)

heading: 4.3 ARTIFICIAL DATA

heading: 4.3.1 EXPERIMENT 1: CONDITIONS OF THEORY FULFILLED
text: Samples are generated in two dimensions, conditioned on five different cluster labels, see Fig. 2. The means of the clusters are chosen independently from a uniform distribution on [−5, 5] and variances from a uniform distribution on [.5, 3].1 This data is then concatenated with independent Gaussian noise in eight dimensions to make a ten-dimensional generating latent space, where only the first two variables are informative. The noise is scaled by 0.01 to be small in comparison to the informative dimensions. The latent space samples are then passed through a RealNVP network
1These values are the same as in an artificial data experiment in Khemakhem et al. (2019). The difference is that the data will be projected into a higher-dimensional space and that we will not assume that the estimating latent space is two-dimensional, this will instead be inferred by the model.
with 8 fully connected coupling blocks with randomly initialized weights to produce the observed data. This acts as a highly nonlinear mixing which can only be successfully treated with nonlinear methods.
GIN is used as the estimating model, with 8 fully connected coupling blocks (full details in Appendix D). Training converges quickly and stably using the Adam optimizer (Kingma & Ba, 2014) with initial learning rate 10−2 and other values set to the usual recommendations. Batch size is 1,000 and the data is augmented with Gaussian noise (σ = 0.01) at each iteration. After convergence of the loss, the learning rate is reduced by a factor of 10 and trained again until convergence.
Over a number of experiments we made the following observations:
• The model converges stably and gives importance (quantified by standard deviation) to only two variables in the estimated latent space, provided there is sufficient overlap between the mixture components in the generating latent space. • Where there is not enough overlap in the generating latent space, the model cannot rec-
ognize common variables across all the different classes, and tends to split one genuine dimension of variation into two or more in its estimated latent space. This appears to be a problem of finite data. We have observed that when this behaviour occurs, and if the gap between mixture components is not too large, it can be prevented by increasing the number of samples so that the space between the mixture components is better filled (see Fig. 6 and 7 in Appendix E). This is consistent with the theory, where equation (4) is true asymptotically. Since the latent space distributions are members of the exponential family, they have support across the entire domain of the latent space, hence gaps can never remain in the limit of infinite samples. • Choice of learning rate is important. If the initial learning rate is too low, training gets stuck
in bad local optima, where single true variables are split into several latent dimensions.

heading: 4.3.2 EXPERIMENT 2: CONDITIONS OF THEORY NOT FULFILLED
text: Samples are generated as in Experiment 1, but with only three mixture components. Since there are two sufficient statistics per dimension, and two dimensions of variation, according to the theory we need at least nk + 1 = 5 distinct conditions u for equation (4) to hold (see section 3.1). Therefore, we might not expect successful experiments. Nonetheless, we observe essentially the same results as for the previous experiments (see Fig. 8 in Appendix E), with the same caveats regarding gaps between the mixture components in the generating latent space. This suggests that the conditions derived in Khemakhem et al. (2019), although sufficient for disentanglement, are not necessary.

heading: 4.4.1 EXPERIMENT

heading: 4.4 EMNIST
text: The data comes from the EMNIST Digits training set of 240,000 images of handwritten digits with labels (Cohen et al., 2017). EMNIST is a larger version of the well-known MNIST dataset, and also includes handwritten letters. Here we use only the digits. The digit label is used as the condition u, hence we construct a Gaussian mixture in the estimated latent space with 10 mixture components. According to the theory (see 3.1), these are only enough conditions to guarantee identifiability if there are only four informative latent variables in the generating latent space. We expect that the true number of informative variables is somewhat higher, so as in Experiment 2 above, we are operating outside of the guarantees of the theory. In addition, the true generative process of human handwriting may not exactly fulfill our method’s assumptions, and we might still lack data on rare or subtle variations, despite the large size of the dataset in comparison to MNIST.
The estimating model is a GIN which uses convolutional coupling blocks and fully connected coupling blocks to transform the data to the latent space (full details in Appendix D). Optimization is with the Adam optimizer, with initial learning rate 3e-4. Batch size is 240 and the data is augmented with Gaussian noise (σ = 0.01) at each iteration. The model is trained for 45 epochs, then for a further 50 epochs with the learning rate reduced by a factor of 10.

heading: 4.4.2 RESULTS
text: The model encodes information into 22 broadly interpretable latent space variables. The eight variables with the highest standard deviation encode global information, such as slant, line thickness and width which affect all digits more or less equally. The remaining 14 meaningful variables encode more local information, which does not affect all digits equally. We interpret all other variables as encoding noise. This is motivated mainly by our observation of an abrupt end to any observable change in the reconstructed digits when any variable after the first 22 is changed. This can be seen clearly in Fig. 15 in Appendix F. Some variables, particularly the global ones, are not entirely disentangled from another. Nevertheless, the results are compelling and suggest that the assumptions required by the theory are approximately met.
We observed some global features which are not usually seen in disentanglement experiments on MNIST (e.g. Chen et al., 2016; Dupont, 2018). These experiments usually obtain digit slant, width and line thickness as the major global independent variables. We too observe slant and line thickness as independent variables (variables 2 and 5), but find width to be split into two variables, one governing the width of the upper half of a digit (variable 1) and the other the lower half (variable 8). This makes sense, since these can in fact vary independently. We also observe the height of a digit (variable 3) (see Fig. 4). This does not usually appear in disentanglement experiments, possibly because it is too subtle a variation for those experiments, but possibly because it is not present or not discoverable in the smaller MNIST dataset but is in EMNIST.
Local features are also usually not observed in such experiments, so the variables which control these are particularly interesting. These variables modify only a region of the digit, leaving the rest of the digit untouched. In addition, digits which do not have the feature which is being modified in that region are left alone. Examples include variable 13, which changes the orientation of the top-left stroke in 2, 3 and 7, and variable 16, which modifies only the lower-right stroke of a 2 (see Fig. 5). The full set of local and global variables can be seen in Fig. 10 to Fig. 15 in Appendix F.

heading: 4.4.3 NOTE ON THE CHOICE OF CONDITIONING VARIABLE
text: The concept of a “true generative process” rests on the assumption that conditioning variables u act as causes for the latent variables z, such that the elements of the latter become independent given the former. In the case of handwritten digits, u represents a person’s decision to write a particular digit, and z determines the hand motions to produce this digit. Modelling u in terms of digit labels is a plausible proxy for the unknown causal variables in the human brain, and the success of the resulting model suggests that the assumptions are approximately fulfilled. Identifying promising conditioning variables u in less obvious situations is an important open problem for future work.

heading: 4.4.4 POTENTIAL IMPROVEMENTS
text: Increasing the number of conditioning variables would bring this work closer in line with the existing theory. At the current value of ten, the theory only applies if there are at most four informative generating latent variables, which is almost certainly too low a number. One option to increase the number of conditions is to relax the labels from hard to soft, i.e. compute posterior probabilities of class membership given a data example. This could be achieved by information distillation as in Hinton et al. (2015). Another possibility is to split existing labels into sub-labels, for example making a distinction between sevens with and without crossbars. This may also aid the generation of more realistic samples, since we observed no sevens with crossbars in our generated samples, even though they make up approximately 5% of the sevens in the dataset.

heading: 5 CONCLUSION AND OUTLOOK
text: We have expanded the theory of nonlinear ICA to cover problems with unknown intrinsic dimension and demonstrated an implementation with GIN, a new volume-preserving modification of the RealNVP invertible network architecture. The variables discovered by GIN in EMNIST are interpretable and detailed enough to suggest that the assumptions made about the generating process are approximately true. Furthermore, our experiments with EMNIST demonstrate the viability of applying models inspired by the new theory of nonlinear ICA to real-world data, even when not all of the conditions of the theory are met.
It is not clear if the methods from this work will scale to larger problems. However, given the recent advances of similar flow-based generative models in density estimation on larger datasets such as CelebA and ImageNet (e.g. Kingma & Dhariwal, 2018; Ho et al., 2019), it is a plausible prospect. In addition, it is not clear whether the method can successfully be extended to the context of semisupervised learning, or ultimately, unsupervised learning.

heading: A PROOF OF IDENTIFIABILITY
text: This section reproduces a proof from Khemakhem et al. (2019) with some modifications to adapt it to the context of invertible neural networks (normalizing flows).
Define the domain of fε as Z = Z1 × · · · × Zn. Define the vector of sufficient statistics T(z) = (T1,1(z1), . . . , Tn,k(zn)) and the vector of their coefficients λ(u) = (λ1,1(u), . . . , λn,k(u)).
Theorem 1 Assume we observe data distributed according to the generative model defined in Sec. 3.1. Further suppose the following:
(i) The sufficient statistics Ti,j(z) are differentiable almost everywhere and their derivatives dTi,j(z)/dz are nonzero almost surely for all z ∈ Zi and all 1 ≤ i ≤ n and 1 ≤ j ≤ k.
(ii) There exist nk + 1 distinct conditions u0, . . . , unk such that the matrix
L = (λ(u1)− λ(u0), . . . ,λ(unk)− λ(u0)) (8) of size nk × nk is invertible.
Then the sufficient statistics of the generating latent space are related to those of the estimated latent space by the following relationship:
T(z) = AT′(w) + c (9)
whereA is a constant, full-rank nk × dk matrix and c ∈ Rnk a constant vector.
Proof: The conditional probabilities pT,λ,f,ε(x|u) and pT′,λ′,g(x|u) are assumed to be the same in the limit of infinite data. By expanding these expressions via the change of variables formula and taking the logarithm we find
log pT,λ(z|u) + log p(ε) + log ∣∣det Jf−1(x)∣∣ = log pT′,λ′(w|u) + log ∣∣det Jg−1(x)∣∣ (10)
where J is the Jacobian matrix. Let u0, . . . , unk be the conditions from (ii) above. We can subtract this expression for u0 from the expression for some condition ul. The Jacobian terms and the term involving ε will vanish, since they do not depend on u:
log pz(z|ul)− log pz(z|u0) = log pw(w|ul)− log pw(w|u0). (11) Since the conditional probabilities of z and w are exponential family members, we can write this expression as
n∑ i=1 log Zi(u0) Zi(ul) + k∑ j=1 Ti,j(z)(λi,j(ul)− λi,j(u0))  =
d∑ i=1 log Z ′i(u0) Z ′i(ul) + k∑ j=1 T ′i,j(w)(λ ′ i,j(ul)− λ′i,j(u0))  (12) where the base measures Qi have cancelled out, since they do not depend on u. Defining λ̄(u) = λ(u)− λ(u0) and writing the above in terms of inner products, we find〈
T(z), λ̄ 〉 + ∑ i log Zi(u0) Zi(ul) = 〈 T′(w), λ̄′ 〉 + ∑ i log Z ′i(u0) Z ′i(ul)
(13)
for 1 ≤ l ≤ nk. Combining the nk expressions into a single matrix equation we can write this in terms of L from condition (ii) above and an analogously defined L′ ∈ Rdk×nk:
LTT(z) = L′TT′(w) + b (14)
where bl = ∑
i log Z′i(u0)Zi(ul) Zi(u0)Z′i(ul)
. SinceLT is invertible, we can multiply this expression by its inverse from the left to get
T(z) = AT′(w) + c (15)
whereA = L−TL′T and c = L−T b. It now remains to show thatA has full rank.
According to a lemma from Khemakhem et al. (2019), there exist k distinct values z1i to z k i such that (
dTi dzi (z1i ), . . . , dTi dzi (zki ) )
are linearly independent in Rk, for all 1 ≤ i ≤ n. Define k vectors zl = (zl1, . . . , z l n) from the points given by this lemma. Take the derivative of equation (15) evaluated at each of these k vectors and concatentate the resulting Jacobians as Q = [JT(z1), . . . , JT(zk)]. Each Jacobian has size nk × n, hence Q has size nk × nk and is invertible by the lemma and the fact that each component of T is univariate. We can construct a corresponding matrix Q′ made up of the Jacobians of T′(g−1 ◦ fε(z)) evaluated at the same points and write
Q = AQ′. (16)
SinceQ is full rank, bothA andQ′ must also be full rank.
A.1 REMARK REGARDING CONDITION (II)
If the coefficients of the sufficient statistics λ are generated randomly and independently, then condition (ii) is almost surely fulfilled. In this case, we can ignore the dependence on u to consider the λ(ul) as independent random variables. Then condition (ii) states that instances of these random variables are in general position in Rnk, which is true almost surely.
A.2 REMARK REGARDING NOISE VARIABLES IN THE ESTIMATED LATENT SPACE
Any variables in the estimated latent space whose distribution does not depend on the conditioning variable u is considered as encoding only noise. Such variables will be cancelled out in equation (12) and the corresponding column of L′T will contain only zeros. This means the corresponding column of A will also contain only zeros, so any variation in such a noise variable has no effect on z. The reverse is also true: if a column of A contains only zeros, the corresponding variable in the estimated latent space does not depend on u and must encode only noise.

heading: B SPARSITY IN UNMIXING MATRIX: GAUSSIAN DISTRIBUTION
text: Suppose samples z from the generating latent space follow a conditional Gaussian distribution. Suppose the estimating model g faithfully reproduces the observed conditional density p(x|u) and samples w from its latent space also follow a conditional Gaussian distribution. Then we can apply equation (4) to relate the generating and estimating latent spaces.
The sufficient statistics of a normal distribution with free mean and variance are z and z2. Hence the relationship between the latent spaces becomes(
z z2
) = A ( w w2 ) + c (17)
where the squaring is applied element-wise. We can writeA in block matrix form as
A =
( A(1) A(2)
A(3) A(4)
) (18)
and c as
c =
( c(1)
c(2)
) (19)
Then:
z = A(1)w +A(2)w2 + c(1) (20)
z2 = A(3)w +A(4)w2 + c(2) (21)
so we can write for each dimension i of z zi = ∑ j A (1) ij wj + ∑ j A (2) ij w 2 j + c (1) i (22)
z2i = ∑ j A (3) ij wj + ∑ j A (4) ij w 2 j + c (2) i (23)
In order to compare the equations, we need to square (22). To do so, we will have to square the second term on the right hand side, involving w2j . There is no matching term in (23), so we have to set all entries ofA(2) to zero. In more detail:∑
j
A (2) ij w 2 j 2 = ∑ j ∑ j′ A (2) ij A (2) ij′w 2 jw 2 j′ (24)
= ∑ j (A (2) ij ) 2w4j + ∑ j 6=j′ A (2) ij A (2) ij′w 2 jw 2 j′ (25)
The first term with w4j matches no term in (23), so we have to set A (2) ij = 0 for all i and j. This simplifies the earlier equation: zi = ∑ j A (1) ij wj + c (1) i (26)
The square of the first term on the right hand side involves terms with wjwj′ cross terms:∑ j A (1) ij wj 2 = ∑ j A (1) ij w 2 j + ∑ j 6=j′ A (1) ij A (1) ij′wjwj′ (27)
so we have to setA(1)ij A (1) ij′ = 0 for all j 6= j′. This means that the i-th row ofA(1) can have at most one nonzero entry. It must also have at least one nonzero entry, since if the row were all zero, a row ofA would be all zero (sinceA(2) = 0), butA has full rank. Since there are as many or fewer rows than columns (n ≤ d), each row ofA is linearly independent, so it is not possible for one to be zero. Hence each row of A(1) has exactly one nonzero entry. Moreover, no two rows of A(1) have their nonzero entries in the same column. If they did, the two rows would not be linearly independent, but they must be sinceA has full rank. Therefore we can write
zi = aiwj + bi (28)
where ai = A (1) ij and bi = c (1) i . That is, the generating latent variable zi is linearly related to some latent variable of the estimating model wj . This estimated latent variable is uniquely associated with zi and any estimated latent variables not associated with a generating latent variable zi (in the case d > n) encode no information about the generating latent space. So the model has decoded the original latent variables z up to an affine transformation and permutation as a subset of variables in its estimated latent space and has encoded no information (only noise) into the remaining latent variables.

heading: C SPARSITY IN UNMIXING MATRIX: TWO-PARAMETER EXPONENTIAL FAMILY MEMBERS
text: The results of Appendix B can be extended to other members of the exponential family with 2 parameters. In the general case, writing Ti,1 = T1 and Ti,2 = T2 for all i, equations (22) and (23) become
T1(zi) = ∑ j A (1) ij T1(wj) + ∑ j A (2) ij T2(wj) + c (1) i (29)
T2(zi) = ∑ j A (3) ij T1(wj) + ∑ j A (4) ij T2(wj) + c (2) i (30)
which, with the definition t = T2 ◦ T−11 becomes T1(zi) = ∑ j A (1) ij T1(wj) + ∑ j A (2) ij t(T1(wj)) + c (1) i (31)
t(T1(wj)) = ∑ j A (3) ij T1(wj) + ∑ j A (4) ij t(T1(wj)) + c (2) i . (32)
We can combine these equations to get
t ∑ j A (1) ij T1(wj) + ∑ j A (2) ij t(T1(wj)) + c (1) i  = ∑ j A (3) ij T1(wj) + ∑ j A (4) ij t(T1(wj)) + c (2) i .
(33) This equation has many summed terms on the right. Suppose that t has a convergent Taylor expansion in some region of its domain. We can take this expansion of the term on the left and compare coefficients. Since t cannot be linear, there will be terms in the expansion of order two or higher. As in the Gaussian case, these polynomial terms create cross terms which are impossible to reconcile with those on the right hand side of the equation. The only consistent solutions can be found by setting all coefficients except those of functions of wj (for some j) to zero:
t ( A
(1) ij T1(wj) +A (2) ij t(T1(wj)) + c (1) i
) = A
(3) ij T1(wj) +A (4) ij t(T1(wj)) + c (2) i (34)
We can further simplify this equation by examining higher-order terms, but need to know the form of t to do so. In any case, we can now write equation (29) as
T1(zi) = A (1) ij T1(wj) +A (2) ij t(T1(wj)) + c (1) i (35)
showing that each generating latent variable zi is related to exactly one estimated latent variable wj . As in the Gaussian case, we can use the full rank property of A to see that each estimated latent variable is associated with at most one generating latent variable, and any estimated latent variables not associated with any generating latent variables must encode only noise. The task now is to check the form of t for each two-parameter member of the exponential family, to see what further constraints we can derive from equation (34). The results are stated in Table 1.

heading: D NETWORK ARCHITECTURE
text: The estimating model g is built in the reverse direction for practical purposes, so the models described here are g−1 which maps from the data space to the latent space. The type and number of coupling blocks for the different experiments are shown below. The affine coupling function is the concatentation of the scale function s and the translation function t, computed together for efficiency, as in Kingma & Dhariwal (2018). It is implemented as either a fully connected network (MLP) or convolutional network, with the specified layer widths and a ReLU activation after all but the final layers. For the convolutional coupling blocks, the splits are along the channel dimension. The scale function s is passed through a clamping function 2 tanh(s), which limits the output to the range (-2, 2), as in Ardizzone et al. (2019). Two affine coupling functions are applied per block, as described in Dinh et al. (2016). Downsampling increases the number of channels by a factor of 4 and decreases the image width and height by a factor of 2, done in a checkerboard-like manner, as described in Jacobsen et al. (2018). The dimensions are permuted in a random but fixed way before application of each fully connected coupling block and likewise the channels for the convolutional coupling blocks. The network for the artificial data experiments has 4,480 learnable parameters and the network for the EMNIST experiments has 2,620,192 learnable parameters.
D.1 NOTE ON OPTIMIZATION METHOD
In the experiments described in this paper, the mean and variance of a mixture component was updated at each iteration as the mean and variance of the transformations to latent space of all data points belonging to that mixture component in a minibatch of data D:
µi(u′;θ)← ED:u=u′(g−1i (x;θ)) (36) σ2i (u ′;θ)← VarD:u=u′(g−1i (x;θ)). (37)
Hence the parameters of the mixture components would change in each batch, according to the data present. The notation µi(u;θ) does not indicate that µ is directly parameterized by θ and learned, instead it indicates that µ is a function of g−1 which is parameterized by θ. A change in θ will also change the output of µ, given the same mini-batch of data D. The same holds for σi(u;θ). We expect that specifying the means and variances as learnable parameters updated in tandem with the model weights would have worked equally well.

heading: E FIGURES FROM THE ARTIFICIAL DATA EXPERIMENTS

heading: F EMNIST FIGURES

Review summaries:
This paper builds upon the recent theoretical framework on nonlinear ICA, put forward in recent work Khemakhem et al. (2019) that draw a lot of attention. The latter work provides an extension of the basic nonlinear ICA that is closely related to a VAE with a conditional factorized prior, essentially introducing side information (with assumed extra observables u) to resolve non-identifiability issues.