Preprint. Under review.
K INKS VS . S MOOTHNESS : I DENTIFIABILITY OF R EAL A NALYTIC N ICA FOR L APLACE - LIKE S OURCES
arXiv:2609.21926v1 [cs.LG] 18 Sep 2026
Isaac Manring, Kejun Huang Department of Computer Science University of Florida Gainesville, FL 32611, USA {imanring, kejun.huang}@ufl.edu
A BSTRACT Many machine learning systems try to explain complex data - like images or financial time series - in terms of hidden, independent factors that generated them. Recovering the true underlying factors, rather than some scrambled version of them, is the central challenge of nonlinear Independent Component Analysis (nICA). We prove identifiability (exact recovery) up to trivial ambiguities for real analytic generating functions when source probability density functions have a finite number of discontinuities in the first derivative. The Laplace distribution is the most prominent example satisfying this assumption. Our proof relies on the contrast between kinks in the source distribution and the smoothness of real analytic functions. Real analytic functions comprise a broad class of generating mechanisms, and can be approximated with Normalizing Flows or Variational Autoencoders with standard activation functions (e.g., tanh, softplus, GELU), so our result applies with minimal changes to existing training pipelines. We perform experiments on real and synthetic data with both Normalizing Flows and Variational AutoEncoders demonstrating their identifiability properties. In experiments on CelebA data we recover several interpretable latent factors controlling unique attributes across the dataset.
1
I NTRODUCTION
Nonlinear Independent Component Analysis (nICA) has foundational theoretical implications for unsupervised and representation learning. For observations x(i) ∈ Rm , it seeks to recover statistically independent generative factors (also called sources) s(i) ∈ Rk , such that x(i) = f (s(i) ) + ϵ(i) ↔ g(x(i) − ϵ(i) ) = s(i) ,
(1)
where f is a nonlinear diffeomorphism, with its inverse g, and ϵ(i) is noise (see formal definitions in section 3). When m = k the noise component is generally dropped. This model is in general not able to recover the generative factors s(i) because of its non-identifiability (see Definition 4 for a formal definition). One of the central goals of nICA is to guarantee the recovery of the true latent factors given a set of assumptions. These guarantees are essential for interpretable latent representations. Additionally, they are useful for training stability and reproducibility (Hyvärinen et al., 2024; Kivva et al., 2022). Last but not least, D’Amour et al. (2022) point to non-identifiability as a cause of the drop in performance of deep learning systems at deployment. In this paper, we prove identifiability of nICA without auxiliary observations when the source distribution has a discontinuity in the first derivative (for example Laplacian) and f and g are real analytic. This is a substantial contribution to nICA theory as the assumption of analyticity is mild. As a caveat we note in subsection 4.2 that while there are several distributions on compact supports that are theoretically identifiable according to Theorem 1, they are not practically identifiable due to the universal approximation ability of real analytic diffeomorphisms on compact sets. Neural networks can easily be constrained to be real analytic by requiring the activation functions to be real analytic. Thus training for our method which we dub RAD (Real Analytic Decoders) requires 1
Preprint. Under review.
minimal changes to standard techniques, enabling easy adoption. We demonstrate the validity of the theoretical results through extensive simulation experiments. RAD also demonstrates good stability in training on real world data and recovers several interpretable latent variables for the CelebA dataset (Liu et al., 2015).
2
R ELATED W ORK
While nICA is still developing, approaches to guarantee identifiability generally fall into two categories. The first utilizes auxiliary observations such as class labels, time-series indices, or the preceding observation for auto-correlated data (Khemakhem et al., 2020; Hyvarinen et al., 2019; Sorrenson et al., 2020). The sources are assumed to be conditionally independent given the auxiliary variable. Hälvä & Hyvarinen (2020) extended this using a Hidden Markov Model so that the auxiliary variables can be unobserved too. The second approach does not assume any auxiliary information, but restricts the function class of f . Zheng et al. (2022); Zheng & Zhang (2023) proved identifiability assuming that the Jacobian of f had a certain sparsity pattern. Nguyen & Fu (2025) proposed the idea that in order to recover s(i) , each source should control a sufficiently distinct aspect of the observations. Mathematically, this was formulated as a sufficiently scattered condition on the rows of the Jacobian of f . This requires the observations to be in a much higher dimension than the sources. In a similar vein Buchholz et al. (2022) showed that if the Jacobian was orthonormal at every point (a conformal map), then identifiability could be achieved. Drawing from optimal transport theory, Huang et al. (2021) trained a unique neural network to recover the sources. The parameterization ensured that f was the gradient of a convex function (its Jacobian was symmetric positive semi-definite). Yang et al. (2022) shows identifiability for volume preserving flows when there are at least two distinct but overlapping classes in the data. In the most comparable work to ours, Kivva et al. (2022) considered real analytic mixture distributions while constraining the form of f to be piecewise affine. In several theorems they prove identifiability of increasing strength under assumptions of increasing strength. The result most comparable to Theorem 1 shows that with a prior of independent Gaussian Mixture Models (GMMs) and a piecewise affine f , si can be recovered up to permutation, scaling, and translation. While our proof does not resemble Kivva et al. (2022), the assumptions are interestingly duals of each other: we assume real analytic f and non-differentiable source distribution, while they assume piecewise affine f and real analytic source distribution. While at first glance RAD may seem similar to Sparse Autoencoders (SAE) because of the use of a Laplace source distribution, there are several important distinctions. First, SAEs generally only have three layers: linear encoder, ReLU, linear decoder (Huben et al., 2024). RAD in general has many more than this. Second, the latent dimension in SAEs is higher than the observation dimension, which is not allowed in RAD. Thus SAEs form a distinct line of research.
3
P RELIMINARIES
To begin, we introduce some preliminary definitions. We require f and g to be real analytic diffeomorphisms. Definition 1 (Diffeomorphism (Boumal, 2023)). A diffeomorphism is a bijective map f : U → V, where U ⊆ Rk , V ⊆ Rm are open sets such that both f and f −1 = g are smooth. Note that diffeomorphisms do not require m = k. If m > k (the data are higher dimensional than the generative factors), we simply require the data to live in a lower dimensional embedded submanifold (see Definition 3). Definition 2 (Real Analytic Function (Krantz & Parks, 2002)). A function f : U ⊆ Rm → Rn is real analytic if for every x0 ∈ U each component f i admits a convergent Taylor expansion about x0 whose sum is equal to f i in some neighborhood around x0 . The set of real analytic functions C ω is a subset of infinitely differentiable functions C ∞ . Many familiar functions, such as polynomials, exp(·), and cos(·), are real analytic. C ω is closed under 2
Preprint. Under review.
composition and inversion. If the functions u, v are real analytic at point x0 then u(x0 ) + v(x0 ) and u(x0 )v(y0 ) are real analytic at x0 . Additionally, u(x1 0 ) is real analytic if u(x0 ) ̸= 0. As mentioned earlier, in the case when m > k we adopt the assumption that the data approximately lie on a low-dimensional submanifold embedded in a higher dimension. This is the standard assumption of representation learning methods, such as VAEs, when the embedded dimension is lower than the observation dimension (Goodfellow et al., 2016). Because f is real analytic this submanifold must be real analytic also. Definition 3 (Real Analytic Submanifold (Krantz & Parks, 2002)). The set M ⊆ Rm is a kdimensional real analytic submanifold if for each p ∈ M there exist an open V with p ∈ V ⊆ Rm , a convex open U ⊆ Rk , and real analytic maps g : V → U, f : U → V such that M ∩ V = image f and g ◦ f is the identity on U. This assumption implies that there is redundant information in the observations. For example, in an image each pixel value is not free to take on any value as it must be similar to those around it for the image to convey meaningful information. This is why compression is possible. Because in practice x(i) almost never exactly lie on a submanifold, it is assumed that they are noisy. We define pf ,g (x) as the probability density of x given the model in Equation 1 parameterized by f and g with known source distribution ρ. We provide the standard definition of identifiability up to an ambiguity class below. Definition 4 (Identifiability). The model in Equation 1 is said to be identifiable up to equivalence class ≃ if ∀x ∈ M, pf̃ ,g̃ (x) = pf ,g (x) =⇒ ∀x ∈ M, g(x) ≃ g̃(x).
4
M AIN R ESULT
In this section we present the identifiability proof. First, using standard techniques, we transform the model in Equation 1 to be in terms of learning a self-map of the sources h = g̃ ◦ f . This step also eliminates the noise from ϵ in the case of m > k and requires Assumption 1, which comes directly from Khemakhem et al. (2020). Then in Theorem 1 we prove that h must be a signed permutation (scaling and translation are embedded in the knowledge of the source distribution). Assumption 1 ((Khemakhem et al., 2020)). The set {x ∈ M|φϵ (x) = 0} has measure zero, where φϵ is the characteristic function of the density of ϵ. Lemma 1. Suppose that ϵ satisfies Assumption 1 and we learn real analytic f̃ and g̃ such that g̃ ◦ f̃ is the identity and ∀x ∈ M, pf̃ ,g̃ (x) = pf ,g (x). Then h = g̃ ◦ f is a real analytic diffeomorphism and h(s) ∼ ρ. Proof. The first part of this proof is relatively trivial. From Definition 2.7.4 of (Krantz & Parks, 2002) it is clear that the composition g̃ ◦ f is real analytic in Rk . It is also clear from Definition 1 that the composition of diffeomorphisms is a diffeomorphism. The second part removes the noise ϵ. Using Assumption 1, Khemakhem et al. (2020) provided a proof for this in Appendix B.2.2 Step I using a deconvolution argument. We utilize the contrast between real analytic generating functions and non-differentiable source distributions to prove identifiability. The non-differentiability of the source distribution ρ must satisfy the following assumption. Assumption 2. ρ is a probability density function with connected support D ⊆ R. log ρ can be partitioned into a finite number of real analytic functions α0 , . . . , αd at distinct breakpoints β1 , . . . , βd for d ≥ 1. The following must hold, ′ 1. αi−1 (βi ) ̸= αi′ (βi ) ∀i = 1, . . . , d
2. βi − βi−1 > 0 ∀i = 1, . . . , d + 1 with the convention that β0 = inf D and βd+1 = sup D 3. For d odd, there are analytic continuations of α d−1 and α d+1 on an open neighborhood 2 2 around β d+1 . Similarly for d even, there are analytic continuations of α d −1 and α d +1 2 2 2 onto the interval [β d , β d +1 ]. 2
2
3
Preprint. Under review.
The last requirement in Assumption 2 is used in step 3 of the proof of Theorem 1 in order to use the Taylor expansion around a common point. Note that although distributions of the form ρ(x) ∝ exp(−|x|α ), α ∈ (0, 1) satisfy requirements 1 and 2 of Assumption 2, they fail on requirement 3. Several well known distributions satisfy this assumption including the Laplace (also called double exponential) distribution and the triangular distribution. The probability density function of the Laplace distribution is given by 1 1 ρL (s) = exp − |s − µ| , 2b b for location and scale parameters µ and b. The Laplace is well known because its maximum likelihood estimate uses the ℓ1 norm penalty which induces sparsity while allowing outliers. Derivatives of the Laplace distribution like Laplace Mixture Model, Log Laplace, and asymmetric Laplace also satisfy Assumption 2. A Spike-and-Slab distribution to induce sparsity as in Moran et al. (2022) also satisfies the assumption. The standard triangular distribution describes the probability density of the sum of two standard uniform random variables. However, in its general form it can also be asymmetric, 0, x<a 2(x−a) , a ≤x≤c ρ∆ (x) = (b−a)(c−a) 2(b−x) (b−a)(b−c) , c ≤ x ≤ b 0, b<x While this distribution satisfies Assumption 2, it may be difficult to use in practice (see subsection 4.2). Exponential of many piecewise defined polynomials and other non-standard distributions also satisfy Assumption 2. The essential characteristic is the discontinuity in the derivative of the probability density function which implies a sharp change in the probability law of the source. Larger discontinuities are harder to approximate with real analytic functions and thus may make recovery easier. Theorem 1. Suppose h is a real analytic diffeomorphism, s is iid according to ρ, and h(s) is also iid according to ρ. If ρ satisfies Assumption 2, then h must be the identity up to signed permutation and translation. P Proof. The proof proceeds in three steps. In the first step we show that i log ρ(si ) − P i log ρ(hi (s)) must be real analytic. In the second step, we show that because the kinks in log ρ(si ) must be canceled out by the kinks in log ρ(hπ1 (i) (s)) we can factorize hπ1 (i) as in Equation 2. Finally in step three, using Taylor expansion we show that ri (s) = ±1 in Equation 2. Because s ∼ ρ and h(s) ∼ ρ, the change of variables formula gives Y Y ρ(si ) = | det J h (s)| ρ(hi (s)) i
X i
log ρ(si ) −
i
X
log ρ(hi (s)) = log | det J h (s)|
i
Because h is invertible det J h (s) ̸= 0. Further due to the continuity of det J h (s) it always has the same sign. Therefore, the absolute value reduces to a constant multiple. Since h is real analytic, J h (s) is also real analytic. The determinant is a polynomial function, which is realPanalytic, as is the P logarithm. Therefore log det J h (s) is real analytic. This implies that γ(s) = i log ρ(si ) − i log ρ(hi (s)) must also be real analytic. Let Aij = {s|si = βj } and Blm = {s|hl (s) = βm }.PAt Aij there is a non-differentiable point in log ρ(s) as you vary siS . It must be canceled out by i log ρ(hi (s)) to maintain P the real analytic property. Thus Aij ⊂ l,m Blm because there can only be non-analytic points in i log ρ(hi (s)) S around l,m Blm because h is real analytic. Suppose that ∄ l, m, such that hl ≡ βm on an open neighborhood of s ∈ Aij . For real analytic hl − βm the measure of the zero-set is zero (Mityagin, 2015). Thus Blm ∩ Aij has a measure of zero. But if all Blm have measure zero in Aij then we reach a contradiction in the measure of Aij . 4
Preprint. Under review.
So ∃ l, m such that hl ≡ βm on an open neighborhood of Aij . By the real analytic identity theorem (Corollary 1.2.6 in (Krantz & Parks, 2002)) hl ≡ βm on all of Aij . Because hl is constant on Aij the directional derivatives v⊤∇hl (s) = 0, ∀s ∈ Aij , v i = 0 which implies ∇hl (s) ∈ span{ei }. Suppose that both hl − βj and hk − βj vanish on Aij , then ∇hl (s), ∇hk (s) ∈ span{ei }, ∀s ∈ Aij . However, this contradicts the inverse function theorem which says that the Jacobian of h is invertible. Therefore Aij = Bπ(i,j) for permutations π. π can be further simplified. For j ̸= j ′ , ∅ = Aij ∩Aij ′ = Bπ(i,j) ∩Bπ(i,j ′ ) . However, Blm ∩Bl′ m′ = ∅ only when l = l′ and m′ ̸= m. Therefore, we can factorize π(i, j) into π1 (i), π2 (j). This means that all kinks in ρ(si ) are canceled by a single dimension of h. For simplicity we decompose s into u = si and v be the rest of s. For a fixed v, h must pass through each kink exactly once or cancellation will not be achieved. Therefore because of continuity in hπ1 (i) (u, v) with respect to u π2 can be reduced to the identity or the reverse order permutation. Because hπ1 (i) − βπ2 (j) vanishes precisely when u − βj vanishes and ∇hπ1 (i) ̸= 0, we can use the Weierstrass Preparation Theorem (Theorem 6.1.3 (Krantz & Parks, 2002)) to factorize hπ1 (i) − βπ2 (j) hπ1 (i) (u, v) − βπ2 (j) = ri (u, v)(u − βj ) (2) where ri (u, v) ̸= 0 and is real analytic. Note here the Weierstrass polynomial is of degree 1 because ∇hi ̸= 0. Step 3 utilizes the uniqueness of the Taylor coefficients of γi (u) = log ρ(hπ1 (i) (u, v)) − log ρ(u) by expanding multiple segments of γi around a shared center point to prove that ri (s) ≡ ±1. Due to space constraints, we leave the full details for Appendix A. Using ri ≡ ±1 in Equation 2, we see that h is a signed permutation with a potential translation.
To learn real analytic diffeomorphisms we consider two well known neural network architectures: normalizing flows (Kobyzev et al., 2021) and variational autoencoders (Kingma & Welling, 2019). Since affine transformations are real analytic, constraining a neural network to be real analytic only requires the activation functions to also be real analytic. Normalizing Flows require invertible activation functions, which means that they must be monotonic. We classify several popular monotonic activation functions shown in Figure 1 by their membership in C ω . The activation functions below are monotonic and real analytic, x
−x
1. Tanh(x) = eex −e +e−x
2. Sigmoid(x) = 1+e1−x 3. Softplus(x) = log(1 + ex ) The following are monotonic but not real analytic,
Figure 1: Common monotonic activation functions x, x≥0 1. LeakyReLU(x) = αx, x < 0 αx − 1, x < −1 2. LeakyHardTanh(x) = x, |x| ≤ 1 βx + 1, x > 1
VAEs do not require invertible activation functions and so may also use more modern alternatives like Swish (Ramachandran et al., 2017) or GELU (Hendrycks & Gimpel, 2023), both of which are real analytic. 5
Preprint. Under review.
4.1
A SSUMPTION T IGHTNESS
In this section we discuss the tightness of the assumptions leading to identifiability in Theorem 1. Some of the most well known nICA identifiability counter-examples are the so called Measure Preserving Automorphisms (MPAs) (Hyvärinen & Pajunen, 1999). For any continuous source distribution ρ with mild regularity assumptions there exists a transformation q such that q(s) ∼ N (0, I). This q can be constructed by q = Φ−1 ◦ R where Φ is the Gaussian Cumulative Distribution Function (CDF) and R is the CDF of ρ. Since the normal distribution is invariant under orthogonal transformation, Qq(s) ∼ N (0, I). Finally the sources can be transformed back to their original distribution so that q −1 (Qq(s)) has the same distribution as s. This provides a clear counter-example when we relax the requirement of h to be real analytic. However, q cannot be real analytic if ρ is not real analytic because R would not be real analytic. MPAs also provide a counter-example when s is distributed according to a real analytic distribution. This includes obvious distributions like Gaussian, but also when there is a discontinuity in the PDF between the 0 and non-zero parts such as for the Gamma or Beta distributions. 4.2
L IMITATIONS
The Whitney Approximation Theorem (Krantz & Parks, 2002) states that real analytic functions can approximate any C r function for any r ≥ 0 to arbitrary precision on a compact set (closed and bounded in euclidean space). Further these approximation results can be attained by neural networks. Ishikawa et al. (2023) provides general results on the universal approximation ability of invertible neural networks (Normalizing Flows). With our architecture we impose the constraint that f ∈ C ω ⊂ C ∞ . So Ishikawa et al. (2023) demonstrates that our method is a universal function and universal density approximator. Puthawala et al. (2022) provides universal approximation results for deep latent neural networks with dimensionality reduction. Proving identifiability for universal function approximators is a two edged sword. On the one hand, they can approximate any generating function, but on the other hand they can also approximate identifiability counter-examples like Measure Preserving Automorphisms (MPAs). Therefore, while Theorem 1 indicates that RAD is identifiable for compactly supported source distributions, training a real analytic approximator is very difficult because it could approximate an MPA. Note that this does not apply to source distributions on non-compact supports like the Laplace. In particular, a model with triangular or truncated Laplace source distribution, where the tails are cut off, while theoretically identifiable, may be practically non-identifiable. In subsection 5.1 we tested these source distributions. Table 2 indicates that we were not able to recover the sources, confirming our hypothesis about the universal approximation problem. Incidentally, this is also a limitation of (Kivva et al., 2022) due to the universal approximation ability of piecewise affine functions. One cannot prove identifiability for a class of universal approximators and then rely on those functions to approximate any function not in the class. Rather the practitioner must know something substantive about the generating function in order to achieve identifiability.
5
E XPERIMENT
We run three experiments (Synthetic, Yahoo Stock, and CelebA) to verify the theory. The synthetic experiments have ground truth labels and thus are the primary means of identifiability verification. We utilize the standard Mean Correlation Coefficient (MCC) metric, which is the mean absolute correlation between learned and true sources after matching. Full details on experiments are in Appendix B. 5.1
S YNTHETIC DATA
In these experiments we generate synthetic data and train Normalizing Flows to recover the ground truth factors, in order to confirm Theorem 1. The synthetic data was generated with source distribution ρ and nonlinear functions fi as follows: 6
Preprint. Under review.
Table 1: Average MCC ± standard deviation for Laplace source model with different activation functions in the Normalizing Flow. Real analytic activation functions achieved substantially higher MCCs than non-real analytic activation functions. k=5 Activation Function Tanh Sigmoid Softplus LeakyHardTanh LeakyReLU
k = 10
Mixture A
Mixture B
Mixture A
Mixture B
0.944 ± 0.017 0.895 ± 0.030 0.833 ± 0.075 0.431 ± 0.198 0.143 ± 0.159
0.993 ± 0.003 0.987 ± 0.002 0.991 ± 0.002 0.660 ± 0.154 0.066 ± 0.041
0.930 ± 0.023 0.842 ± 0.085 0.765 ± 0.110 0.448 ± 0.057 0.257 ± 0.107
0.983 ± 0.006 0.987 ± 0.006 0.980 ± 0.007 0.540 ± 0.039 0.038 ± 0.011
iid
1. s(i) ∼ ρ 2. y (i) = As(i) h i⊤ (i) (i) 3. x = f1 (y 1 ), · · · , fk (y k ) For each trial a new random normal mixing matrix A was generated. To make training easier, we shrank the log singular values by 30%. This shrinks very high and very low singular values which were hard for the neural networks to learn. The√training would slow down when singular values were too high or low. Then we normalized by k so that the range of y would fall on the most nonlinear part of fi . Next we applied one of two randomly generated nonlinear fi , which we named Mixture A and Mixture B respectively. 1 y sin(y + 2πwi ) + 4 3 fi (y) = wi y + tanh(y)(1 − wi ) fi (y) =
where wi ∼ U (0, 1). Mixture A’s generating function includes the linear part in order to maintain invertibility. In the experiments below we train Normalizing flows with 6 dense invertible layers. Unless otherwise stated, the final layer has no activation applied. For each experiment we trained on 5 different randomly generated datasets using the same seeds across experiments. We calculate Mean Correlation Coefficient (MCC) and report mean and standard deviation across trials. In the first experiment, we compare identifiability of a Laplacian source model for different activation functions for both m = k = 5 and m = k = 10 with n = 10, 000 samples. Table 1 contains the results. Real analytic activation functions resulted in markedly higher MCCs for both Mixture A and Mixture B as Theorem 1 would indicate. In the second experiment, we compare identifiability of the model with different source distributions. We fixed the activation function as Tanh to keep g real analytic and because Tanh performed the best in the first experiment. Otherwise the architecture was the same. We considered higher dimensions m = k = 15, 25. Interestingly we found that a lower sample of n = 5, 000 would still achieve high MCCs. We tested six different source distributions: Laplace, Laplace Mixture Model (LMM), Gaussian, Gamma, Triangular, and Truncated Laplace. For Laplace and Gaussian source distributions, we set the mean to 0 and scale to 1. For the Gamma source distribution we used α = 5, β = 5 and added a Softplus activation to the final layer in g̃ to ensure positive values. The LMM was a mixture of two Laplace distributions with scales 1 and means 0 and 1 with a 0.75, 0.25 weighting respectively. Because the standard triangular distribution has support on [−1, 1] we applied a Tanh activation to the final layer. Finally the truncated Laplace distribution was a Laplace distribution truncated to the support [−2, 2]. So we applied a Tanh to the final layer and then multiplied by 2. According to Theorem 1 the Laplace, LMM, Triangular, and truncated Laplace models are identifiable while the Gaussian and Gamma are not. However, as mentioned in subsection 4.2, the Triangular, and truncated Laplace models have compactly supported source distributions and thus are not expected to be able to recover the sources. 7
Preprint. Under review.
Table 2: Average MCC ± standard deviation for different source distributions. The Laplace source distribution results in high MCCs even for higher dimensions. k = 15
k = 25
Distribution
Mixture A
Mixture B
Mixture A
Mixture B
Laplace LMM Gaussian Gamma Triangular Trunc. Laplace
0.933 ± 0.029 0.880 ± 0.045 0.514 ± 0.011 0.511 ± 0.029 0.453 ± 0.010 0.484 ± 0.033
0.989 ± 0.002 0.966 ± 0.027 0.520 ± 0.006 0.475 ± 0.025 0.433 ± 0.013 0.468 ± 0.025
0.857 ± 0.072 0.805 ± 0.048 0.417 ± 0.011 0.432 ± 0.021 0.385 ± 0.005 0.400 ± 0.006
0.986 ± 0.003 0.966 ± 0.019 0.437 ± 0.011 0.348 ± 0.031 0.385 ± 0.008 0.394 ± 0.012
Results reported in Table 2 show high MCCs for the Laplace and LMM models, and low MCCs for Gaussian, Gamma, Triangular, and Truncated Laplace models. The high MCCs of the LMM model demonstrate that even when there are multiple non-differentiable points the model is still identifiable. Finally we provide an example correlation plot between true and learned sources in Figure 2 comparing an identifiable Laplace model to a non-identifiable Gamma model.
Figure 2: Correlation plot between learned and ground truth sources for Laplace source (left) and Gamma source (right). The identifiable Laplace model is able to recover the true sources, while the Gamma model is not. 5.2
YAHOO S TOCK R ETURNS
Using yfinance1 , we downloaded daily closing stock prices for 15 companies for the period 199901-22 to 2026-08-31. See appendix for full data and training details. In order to obtain stationary data, we calculated the log returns for each stock, which is given by rs (t) = log cs (t) − log cs (t − 1) for closing price cs (t) of stock s at time t. We then trained Normalizing Flows with Tanh activations to generate the returns. Because there is no ground truth to validate our results, following Kivva et al. (2022) we trained 10 models starting from distinct random seeds and compared latent sources between models. For all 45 distinct pairs of models we calculated the MCC between the learned sources. This was done for several source distributions as reported in Table 3. The results demonstrate that a Laplace prior made training much more stable across seeds, in accordance with Theorem 1. Additionally, to compare against the results from Kivva et al. (2022) we trained a piecewise affine Normalizing Flow with independent Gaussian Mixture Model source distributions, which we term GMM. Each component had a mixture of two Gaussians that were jointly learned with the Normalizing Flow. We tried LeakyReLU and LeakyHardTanh activations with a 6-layer network. We report 1
https://ranaroussi.github.io/yfinance/
8
Preprint. Under review.
Table 3: Average ± standard error (SE) of MCC between model pairs trained on the Yahoo stock returns dataset for different source distributions. The Laplace distribution is the only one that satisfies Assumption 2 and thus has the highest average MCC. Source Distribution
Avg MCC ± SE
Laplace Gaussian Exponential Gamma GMM
0.871 ± 0.012 0.516 ± 0.003 0.447 ± 0.003 0.492 ± 0.003 0.480 ± 0.044
the LeakyReLU because it performed slightly better. The poor performance reported in Table 3 does not disprove the results from (Kivva et al., 2022), but it does demonstrate that this dataset does not satisfy the assumptions of (Kivva et al., 2022). 5.3
C ELEBA
Lastly, we performed a larger scale experiment on the CelebA dataset (Liu et al., 2015) using the framework provided by Subramanian (2020). We trained a VAE with Laplacian prior and a latent space of 128. The encoder and decoder were CNNs with softplus and tanh activations to ensure the network is real analytic. See Appendix C for details. The posterior distribution parameters are provided by g̃. We modeled the posterior with a Laplace distribution since the prior distribution is also a Laplace. For sample i index j ! (i) (i) |sj − aj | (i) (i) (i) (i) p(sj |x ) = Lap(aj , bj ) ∝ exp − . (i) bj Nawa & Nadarajah (2024) provide the exact form of the KL Divergence between two Laplace distributions, b1 b0 |a0 − a1 | |a0 − a1 | DKL (Lap(a0 , b0 )∥ Lap(a1 , b1 )) = log + exp − + − 1. b0 b1 b0 b1 Following (Bojanowski et al., 2017) we augmented the standard mean squared error (MSE) reconstruction loss with Laplacian Pyramid loss, which helps avoid the excessive blurriness resulting from MSE reconstruction loss. Bojanowski et al. (2017) gives the Laplacian Pyramid loss as, X Lap1 (x, x′ ) = 22j |Lj (x) − Lj (x′ )|1 , j
where Lj is the j-th level of the Laplacian pyramid (Ling & Okada, 2006). The results shown in Appendix C demonstrate that many of the resulting latent variables have a clear interpretation across different images.
6
C ONCLUSION AND F UTURE W ORK
Real analytic functions comprise a large class of practically useful functions. Proving identifiability for real analytic transformations of sources with non-differentiable probability density functions is a significant step forward in nICA research. Additionally, unlike many nICA methods (Nguyen & Fu, 2025; Zheng et al., 2022; Gresele et al., 2021), RAD does not require expensive Jacobian regularization terms in the objective function. In fact, it can be trained with existing techniques such as Normalizing Flows and Variational Autoencoders. However, there are still many open problems to be investigated, including the sample complexity of RAD and other nICA methods, and the trade-off between universal approximation and source recovery ability in nICA. Additionally, many nICA methods use VAEs, but assume that the likelihood can be perfectly learned. This is only possible when the KL-divergence between the true posterior and the variational distribution from the encoder is zero. If the encoder has infinite capacity this can happen, but in practice there is a gap between the ELBO and the likelihood. Future work would study how a non-zero gap affects identifiability results. 9
Preprint. Under review.
7
LLM U SAGE D ISCLOSURE
While all proofs were written by the authors, an LLM was employed to ideate on proof techniques for some parts of steps 2 and 3 of the proof of Theorem 1. The authors take responsibility for the correctness and comprehensibility of the proofs. LLMs were also used to provide feedback and suggest improvements for the work. LLMs were not used to implement any significant portion of the experiments or propose novel ideas. Further LLMs were not used to write the paper apart from the feedback already mentioned.
8
R EPRODUCIBILITY S TATEMENT
Theorem 1, one of the main contributions, rests on the assumptions explicitly stated in section 4. The proof of the theorem is contained in section 4 and Appendix A. To reproduce experiments, we will release the code upon acceptance. Training details are provided in Appendix B and Appendix C.
R EFERENCES Piotr Bojanowski, Armand Joulin, David Lopez-Paz, and Arthur Szlam. Optimizing the latent space of generative networks. In International Conference on Machine Learning, 2017. URL https: //api.semanticscholar.org/CorpusID:2019311. Nicolas Boumal. An Introduction to Optimization on Smooth Manifolds. Cambridge University Press, 2023. Simon Buchholz, Michel Besserve, and Bernhard Schölkopf. Function classes for identifiable nonlinear independent component analysis. Advances in Neural Information Processing Systems, 35: 16946–16961, 2022. Alexander D’Amour, Katherine Heller, Dan Moldovan, Ben Adlam, Babak Alipanahi, Alex Beutel, Christina Chen, Jonathan Deaton, Jacob Eisenstein, Matthew D. Hoffman, Farhad Hormozdiari, Neil Houlsby, Shaobo Hou, Ghassen Jerfel, Alan Karthikesalingam, Mario Lucic, Yian Ma, Cory McLean, Diana Mincu, Akinori Mitani, Andrea Montanari, Zachary Nado, Vivek Natarajan, Christopher Nielson, Thomas F. Osborne, Rajiv Raman, Kim Ramasamy, Rory Sayres, Jessica Schrouff, Martin Seneviratne, Shannon Sequeira, Harini Suresh, Victor Veitch, Max Vladymyrov, Xuezhi Wang, Kellie Webster, Steve Yadlowsky, Taedong Yun, Xiaohua Zhai, and D. Sculley. Underspecification presents challenges for credibility in modern machine learning. J. Mach. Learn. Res., 23(1), January 2022. ISSN 1532-4435. Ian Goodfellow, Yoshua Bengio, and Aaron Courville. Deep learning. MIT press, 2016. Luigi Gresele, Julius Von Kügelgen, Vincent Stimper, Bernhard Schölkopf, and Michel Besserve. Independent mechanism analysis, a new concept? Advances in neural information processing systems, 34:28233–28248, 2021. Hermanni Hälvä and Aapo Hyvarinen. Hidden markov nonlinear ica: Unsupervised learning from nonstationary time series. In Conference on Uncertainty in Artificial Intelligence, pp. 939–948. PMLR, 2020. Dan Hendrycks and Kevin Gimpel. Gaussian error linear units (gelus), 2023. URL https:// arxiv.org/abs/1606.08415. Chin-Wei Huang, Ricky T. Q. Chen, Christos Tsirigotis, and Aaron Courville. Convex potential flows: Universal probability distributions with optimal transport and convex optimization. In International Conference on Learning Representations, 2021. URL https://openreview. net/forum?id=te7PVH1sPxJ. Robert Huben, Hoagy Cunningham, Logan Riggs Smith, Aidan Ewart, and Lee Sharkey. Sparse autoencoders find highly interpretable features in language models. In The Twelfth International Conference on Learning Representations, 2024. URL https://openreview.net/forum? id=F76bwRSLeK. 10
Preprint. Under review.
Aapo Hyvarinen, Hiroaki Sasaki, and Richard Turner. Nonlinear ica using auxiliary variables and generalized contrastive learning. In The 22nd international conference on artificial intelligence and statistics, pp. 859–868. Pmlr, 2019. Aapo Hyvärinen, Ilyes Khemakhem, and Ricardo Monti. Identifiability of latent-variable and structural-equation models: from linear to nonlinear. Annals of the Institute of Statistical Mathematics, 76(1):1–33, 2024. Aapo Hyvärinen and Petteri Pajunen. Nonlinear independent component analysis: Existence and uniqueness results. Neural Networks, 12(3):429–439, 1999. ISSN 0893-6080. doi: https: //doi.org/10.1016/S0893-6080(98)00140-3. URL https://www.sciencedirect.com/ science/article/pii/S0893608098001403. Isao Ishikawa, Takeshi Teshima, Koichi Tojo, Kenta Oono, Masahiro Ikeda, and Masashi Sugiyama. Universal approximation property of invertible neural networks. J. Mach. Learn. Res., 24(1), January 2023. ISSN 1532-4435. Ilyes Khemakhem, Diederik Kingma, Ricardo Monti, and Aapo Hyvarinen. Variational autoencoders and nonlinear ica: A unifying framework. In International conference on artificial intelligence and statistics, pp. 2207–2217. PMLR, 2020. Diederik P Kingma and Max Welling. An introduction to variational autoencoders. Foundations and Trends® in Machine Learning, 12(4):307–392, 2019. Bohdan Kivva, Goutham Rajendran, Pradeep Ravikumar, and Bryon Aragam. Identifiability of deep generative models without auxiliary information. In Proceedings of the 36th International Conference on Neural Information Processing Systems, NIPS ’22, Red Hook, NY, USA, 2022. Curran Associates Inc. ISBN 9781713871088. Ivan Kobyzev, Simon J.D. Prince, and Marcus A. Brubaker. Normalizing flows: An introduction and review of current methods. IEEE Transactions on Pattern Analysis and Machine Intelligence, 43(11):3964–3979, November 2021. ISSN 1939-3539. doi: 10.1109/tpami.2020.2992934. URL http://dx.doi.org/10.1109/TPAMI.2020.2992934. Steven G. Krantz and Harold R. Parks. A Primer of Real Analytic Functions. Birkhäuser Boston, MA, 2nd edition, 2002. URL https://doi.org/10.1007/978-0-8176-8134-0. Haibin Ling and K. Okada. Diffusion distance for histogram comparison. In 2006 IEEE Computer Society Conference on Computer Vision and Pattern Recognition (CVPR’06), volume 1, pp. 246– 253, 2006. doi: 10.1109/CVPR.2006.99. Ziwei Liu, Ping Luo, Xiaogang Wang, and Xiaoou Tang. Deep learning face attributes in the wild. In Proceedings of International Conference on Computer Vision (ICCV), December 2015. Boris Mityagin. The zero set of a real analytic function, 2015. URL https://arxiv.org/ abs/1512.07276. Gemma Elyse Moran, Dhanya Sridhar, Yixin Wang, and David Blei. Identifiable deep generative models via sparse decoding. Transactions on Machine Learning Research, 2022. ISSN 28358856. URL https://openreview.net/forum?id=vd0onGWZbE. Victor Nawa and Saralees Nadarajah. Exact expressions for kullback–leibler divergence for univariate distributions. Entropy, 26(11), 2024. ISSN 1099-4300. doi: 10.3390/e26110959. URL https://www.mdpi.com/1099-4300/26/11/959. Hoang-Son Nguyen and Xiao Fu. Diverse influence component analysis: A geometric approach to nonlinear mixture identifiability, 2025. URL https://arxiv.org/abs/2510.17040. Michael Puthawala, Matti Lassas, Ivan Dokmanic, and Maarten De Hoop. Universal joint approximation of manifolds and densities by simple injective flows. In Kamalika Chaudhuri, Stefanie Jegelka, Le Song, Csaba Szepesvari, Gang Niu, and Sivan Sabato (eds.), Proceedings of the 39th International Conference on Machine Learning, volume 162 of Proceedings of Machine Learning Research, pp. 17959–17983. PMLR, 17–23 Jul 2022. URL https: //proceedings.mlr.press/v162/puthawala22a.html. 11
Preprint. Under review.
Prajit Ramachandran, Barret Zoph, and Quoc V. Le. Searching for activation functions, 2017. URL https://arxiv.org/abs/1710.05941. Peter Sorrenson, Carsten Rother, and Ullrich Köthe. Disentanglement by nonlinear ica with general incompressible-flow networks (gin). In International Conference on Learning Representations, 2020. URL https://openreview.net/forum?id=rygeHgSFDH. A.K Subramanian. Pytorch-vae. https://github.com/AntixK/PyTorch-VAE, 2020. Xiaojiang Yang, Yi Wang, Jiacheng Sun, Xing Zhang, Shifeng Zhang, Zhenguo Li, and Junchi Yan. Nonlinear ICA using volume-preserving transformations. In International Conference on Learning Representations, 2022. URL https://openreview.net/forum?id= AMpki9kp8Cn. Yujia Zheng and Kun Zhang. Generalizing nonlinear ica beyond structural sparsity. Advances in Neural Information Processing Systems, 36:13326–13355, 2023. Yujia Zheng, Ignavier Ng, and Kun Zhang. On the identifiability of nonlinear ica: Sparsity and beyond. Advances in neural information processing systems, 35:16411–16422, 2022.
A
M AIN T HEOREM C ONTINUED P ROOF S TEP 3
We can also decompose γ as X
γ(u, v) = γi (u, v) +
log ρ(hj (u, v)) −
X
log ρ(v l )
l
j̸=π1 (i)
γi (u, v) = log ρ(ri (u, v)(u − βj ) + βπ2 (j) ) − log ρ(u). P P Near u = βj , s ∈ Blm , ∀(l, m) ̸= (i, j), l̸=π1 (i) log ρ(hl (s)) − l log ρ(v l ) is real analytic, so γi must also be real analytic. At this point for notational simplicity, we focus on a fixed v and therefore drop v from the notation. We consider two cases which require slightly different techniques. The cases correspond to when hπ1 (i) is the identity of u and when hπ1 (i) reverses u. Case 1: ri (u) > 0 Because γi (u) is real analytic, γi′ (u) must be continuous across u = βj . ′ lim+ γi′ (u) = lim+ απ2 (j) (ri (u)(u − βj ) + βπ2 (j) ) − αj (u) u→βj
u→βj
= lim+ απ′ 2 (j) (ri (u)(u − βj ) + βπ2 (j) )(ri′ (u)(u − βj ) + ri (u)) − αj′ (u) u→βj = απ′ 2 (j) (βπ2 (j) )ri (βj ) − αj′ (βj ) The same reasoning can be done for u → βj− . Continuity at u = βj gives, ′ απ′ 2 (j) (βπ2 (j) )ri (βj ) − αj′ (βj ) = απ′ 2 (j)−1 (βπ2 (j) )ri (βj ) − αj−1 (βj )
ri (βj ) =
′ αj′ (βj ) − αj−1 (βj ) ′ ′ απ2 (j) (βπ2 (j) ) − απ2 (j)−1 (βπ2 (j) )
(3)
∂h (β )
j l = ri′ (βj )(βj − βj ) + ri (βj ) > 0, π2 is the identity from the continuity argument But since ∂u in step 2. Therefore ri (βj ) = 1. P∞ l Let the power series of ri around u = βj be ri (u) = 1 + l=n cl (u − βj ) , where coefficients c1 , . . . , cn−1 = 0 for n ≥ 1, where n = 1 means no coefficients are 0. Since αj is real analytic it P∞ (j) can be described by its power series αj (x) = m=1 am (x − βj )m . Therefore !m ∞ ∞ X X αj (βj + (u − βj )ri (u)) = a(j) (u − βj ) + cl (u − βj )l+1 m
m=1
u − βj +
∞ X
l=n
!m cl (u − βj )l+1
= (u − βj )m + mcn (u − βj )m+n + O((u − βj )m+n+1 )
l=n (j)
αj (βj + (u − βj )ri (u)) − αj (u) = a1 cn (u − βj )n+1 + O((u − βj )n+2 ) 12
Preprint. Under review.
We can do the same thing for αj−1 to get (j−1)
αj−1 (βj + (u − βj )ri (u)) − αj−1 (u) = a1
cn (u − βj )n+1 + O((u − βj )n+2 )
Because γi (u) is real analytic, as u → βj the coefficients of (u − βj )n+1 must be equal. But, (j−1) (j) a1 ̸= a1 , so cn = 0. We can do this inductively to show that all Taylor coefficients of ri are 0 except for the constant term. Since ri is real analytic, this means that ri ≡ 1. Case 2: ri (u) < 0 When d is odd, for j = d+1 2 π2 (j) = j due to the fact that π2 is a reverse ordering. The rest of the arguments mirror Case 1, resulting in ri ≡ −1. However, for d even, the argument does not mirror exactly. Because ri ̸= 0 and ri is continuous, it must have the same sign on its domain. There exists δ ∈ [β d , β d +1 ], such that hπ1 (i) (δ, v) = δ. So 2 2 using the Weierstrass Preparation Theorem as before we factorize hπ1 (i) hπ1 (i) (s) − δ = ri (s)(si − δ) P∞ P∞ (j) l The power series of ri and αj are ri (u) = c0 + l=n cl (u − δ) and αj (x) = m=1 am (x − δ)m . !m ∞ ∞ X X l+1 (π2 (j)) απ2 (j) (ri (u)(u − δ) + δ) = am c0 (u − δ) + cl (u − δ) m=1
c0 (u − δ) +
∞ X
l=n
!m cl (u − δ)
l+1
m−1 m = cm cn (u − δ) 0 (u − δ) + mc0
n+m
n+m+1 + O (u − δ)
l=n
=
∞ X
απ2 (j) (ri (u)(u − δ) + δ) − αj−1 (u) (π2 (j)) n+1 n+2 (j−1) m 2 (j)) m (a(π c − a )(u − δ) + a c (u − δ) + O (u − δ) n m 0 m 1
m=1
This is precisely γi for u ∈ (βj−1 , βj ). For each interval we can do this same expansion. Since γi is real analytic, any power series expansion around δ must have the same coefficients. From the coefficients of the first order terms we have, (π (j))
a1 2
(j−1)
(π (j)−1)
= a1 2
c0 − a1
c0 =
(j)
c0 − a1
(j)
(j−1)
(π (j)−1)
− a1 2
a1 − a1 a1 2
(π (j))
With doing this same calculation for the intervals (βπ2 (j)−1 , βπ2 (j) ) and (βπ2 (j) , βπ2 (j)+1 ) we arrive at
c0 =
(π (j))
(π (j)−1)
− a1 2
(j)
(j−1)
a1 2
a1 − a1 1 c0 = =⇒ c0 = −1 c0
c0 cannot be positive because we assumed ri < 0. The n + 1 order term must also have equal coefficients for all intervals. Thus for any j the (βj , βj+1 ) interval has n + 1 order coefficient, (π (j)−1)
2 an+1
(j)
(π (j)−1)
(−1)n+1 − an+1 + a1 2
cn
(4)
If n is odd we equate the n + 1 coefficient for intervals (βj−1 , βj ) and (βj , βj+1 ) and then for intervals (βπ2 (j)−1 , βπ2 (j) ) and (βπ2 (j) , βπ2 (j)+1 ) to get the following two expressions for cn , 13
Preprint. Under review.
(π (j)−1)
cn =
2 (an+1
(π (j))
a1 2 (j)
cn =
(π (j))
2 − an+1
(j)
(j−1)
) − (an+1 − an+1 ) (π (j)−1)
− a1 2
(j−1)
(π (j)−1)
2 (an+1 − an+1 ) − (an+1
(π (j))
2 − an+1
)
(j) (j−1) a1 − a1
cn = −cn =⇒ cn = 0 If n is even we equate the coefficients in Equation 4 for intervals (βj−1 , βj ) and (βπ2 (j) , βπ2 (j)+1 ) to get, (j−1)
cn =
(π (j))
2 (an+1 − an+1
B
E XPERIMENT D ETAILS
B.1
S YNTHETIC
(π (j))
2 )(−1) − (an+1
(π (j))
a1 2
(j−1)
− a1
(j−1)
− an+1 )
=0
For the first experiment with m = k = 5, 10 we generated 10,000 data points to train on for 500 epochs. For m = k = 15, 25 we generated 5,000 data points to train on for 1000 epochs. The optimizer was Adam with learning rate 1e-3 and batches of 256 samples. B.2
S TOCK NF
The stock data was from the period 1999-01-22 to 2026-08-31. We chose large companies from a wide variety of market sectors (Table 4). For each model we trained for 750 epochs with Adam optimizer, learning rate set at 1e-2, with a decline on plateau scheduler, and batches of 256 samples. We arbitrarily chose the seeds 9-18 to initialize the 10 models with. Table 4: Company names of stock data in the Yahoo Stock experiment. For easier replication we provide the Yahoo Finance Symbol over the company names. AAPL Apple Inc. ALL The Allstate Corp. CVX Chevron Corp. XOM ExxonMobil Corp.
C
MSFT Microsoft Corp. JPM JPMorgan Chase & Co. WMT Walmart Inc. AMD Advanced Micro Devices Inc.
NVDA NVIDIA Corp. PFE Pfizer Inc. CAT Caterpillar Inc. BRK-A Berkshire Hathaway Inc.
AMZN Amazon Inc. JNJ Johnson & Johnson DIS The Walt Disney Co.
C ELEBA
The Variational Autoencoder architecture used convolutional layers for the encoder and convolutional transpose layers for the decoder. The kernel size was 3, with 32, 64, 128, 256, 512 filters for the 5 convolutional layers. We used the Softplus activation for hidden layers and the Tanh for the final layer of the decoder. We used a latent dimension of k = 128 and trained for 100 epochs with a learning rate of 0.005 on one B200 GPU.
14
Preprint. Under review.
Figure 3: Head orientation
Figure 4: Bangs
Figure 5: Skin Color
15
Preprint. Under review.
Figure 6: Smile/Broad face
16