Sample Complexity of Interventional Causal Representation Learning

Authors: emre  acartürk
burak  varıcı
karthikeyan  shanmugam
ali  tajer
Venue: NeurIPS_cc_2024_Conference
Published: 2024-05-15 05:18:15

Abstract: Consider a data-generation process that transforms low-dimensional _latent_ causally-related variables to high-dimensional _observed_ variables. Causal representation learning (CRL) is the process of using the observed data to recover the latent causal variables and the causal structure among them. Despite the multitude of identifiability results under various interventional CRL settings, the existing guarantees apply exclusively to the _infinite-sample_ regime (i.e., infinite observed samples). This paper establishes the first sample-complexity analysis for the finite-sample regime, in which the interactions between the number of observed samples and probabilistic guarantees on recovering the latent variables and structure are established. This paper focuses on _general_ latent causal models, stochastic _soft_ interventions, and a linear transformation from the latent to the observation space. The identifiability results ensure graph recovery up to ancestors and latent variables recovery up to mixing with parent variables. Specifically, ${\cal O}((\log \frac{1}{\delta})^{4})$ samples suffice for latent graph recovery up to ancestors with probability $1 - \delta$, and ${\cal O}((\frac{1}{\epsilon}\log \frac{1}{\delta})^{4})$ samples suffice for latent causal variables recovery that is $\epsilon$ close to the identifiability class with probability $1 - \delta$.

Main paper: 
 introduction: 
 The observed data generated in a wide range of technological, social, and biological domains has high-dimensional, complex, and often unexplainable structures. Nevertheless, sometimes, such complex structures can be explained by latent generating factors that often form lower-dimensional structures. The field of causal representation learning (CRL) is motivated by this premise and focuses on disentangling the causally-related latent generating factors that underlie a given high-dimensional observable dataset. In particular, given a high-dimensional dataset, the objective of CRL is to learn (i) the underlying latent causal variables and (ii) the causal relationships among the latent variables. The learned representations then provide an explainable structure for the observed data and facilitate informed reasoning for the downstream tasks [1] .Formally, CRL consists of a data generation and data transformation pipeline as follows. There exists a set of high-level latent variables Z ∈ R n that are causally related. The causal interactions are captured by a directed acyclic graph (DAG) G. The latent variables Z go through an unknown transformation g and generate the observed variables X ∈ R d , i.e., X = g(Z). The objective of CRL is to use the observed variables X and learn the latent graph G and latent variables Z. The primary question facing CRL is identifiability, which refers to establishing the conditions under which recovering Z and G are possible.Identifiability is known to be impossible without proper inductive biases or additional supervision [2, 3] . One approach to address this issue is performing interventions -leading to interventional CRL -which has gained significant recent attention [4, 5, 6, 7, 8, 9, 10] . Specifically, interventions on latent causal variables create additional statistical diversity in the observed data. Accessing to enough interventional environments enables identifiability [1] . Despite the recent advances in establishing identifiability guarantees for interventional CRL, all the existing guarantees focus on the asymptotic infinite-sample regime, where one assumes access to an infinite number of observable samples X.In this paper, we establish conditions for finite-sample non-asymptotic identifiability and performance guarantees for interventional CRL. We focus on CRL under the general (non-parametric) latent causal models, soft interventions, and linear transformations from the latent to the observation space. This setting is well-studied in the infinite-sample regime [4, 5, 7, 11] , where the existing identifiability results show that the causal graph can be recovered up to ancestral nodes, and latent variables can be learned up to mixing with their ancestors. In this paper, we provide the probabilistic finite-sample counterparts of these identifiability guarantees and establish the first sample complexity analysis for interventional CRL in the finite-sample regime.Our CRL approach falls in the category of score-based CRL [5] , based on which our sample complexity analysis consists of two steps. In the first step, we delineate a general sample complexity for any desired consistent score estimator. Subsequently, we specialize these general results by adopting the reproducing kernel Hilbert space (RKHS)-based score estimator [12] . We establish the following identifiability guarantees for any desired pair of constants ϵ, δ ∈ R + .• Latent graph recovery: Using RHKS-based score estimator, O((log 1 δ ) 4 ) samples suffice to recover the transitive closure of the latent graph with probability at least 1 -δ.• Latent variables recovery: Using the same score estimator, O(( 1 ϵ log 1 δ ) 4 ) samples suffice to ensure that the mean squared error in the estimated latent causal variables is at most ϵ 2 with probability at least 1 -δ.• Dependence on model dimension: We further establish the explicit and implicit dependence of sample complexity on the dimensions of the latent and observable spaces, n and d, respectively. The precise characterizations of these expressions involve additional model-dependent constants specified in Section 5.• Improved guarantees: Finally, we note that our latent variables recovery result, where we show recovery up to mixing with parents, which improves upon the results in existing literature that guarantee recovering up to mixing with ancestors.Methodology. We offer novel finite-sample CRL algorithms and characterize the sample complexity guarantees achievable by these algorithms. We design our algorithms using the properties of score functions (i.e., the gradient of the log density). This approach is inspired by the score-based CRL framework [5, 8] . However, the algorithms are significantly different, led by the need to integrate finite-sample score estimation routines and their imperfections. In our analysis, first, we provide sample complexity upper bounds for a score-based framework that uses any desired consistent score difference estimator. Then, we adopt a specific score estimator [12] and provide explicit sample complexity guarantees.Related work. All the existing studies on interventional CRL focus on the infinite-sample regime.We briefly review the studies that adopt a similar CRL model, i.e., linear transformations and singlenode soft interventions. In this setting, complete identifiability results for linear latent models are established in [7] . Similar results are shown for nonlinear latent causal models in [4, 11] and for linear non-Gaussian models in [13] . In the most closely related setting to ours, identifiability up to ancestors is shown to be possible using soft interventions on general latent causal models [8] . Other studies on interventional CRL include using do interventions for polynomial transformations [6, 14] and hard interventions for general transformations [10, 8] . On a partially related problem, error rates for using score functions on observed causal discovery are provided in [15] .Notations. For n ∈ N, we define [n] ≜ {1, . . . , n}. Vectors are represented by lower case bold letters, and element i of vector v is denoted by v i . Matrices are represented by upper case bold letters, and we denote row i and column j of matrix A by A i,: and A :,j , respectively. The row permutation matrix for permutation π of [n] is denoted by P π . Random variables and their realizations are presented by upper and lower case letters, respectively. We denote the Moore-Penrose pseudoinverse of a matrix A by A † . Given a symmetric matrix A ∈ R d×d , we denote the vector of eigenvalues of A ordered in ascending order by λ(A) ∈ R d and the matrix of eigenvectors byQ(A) ∈ R d×d such that A = Q(A) • diag(λ(A)) • Q(A) ⊤ .For any matrix A, we denote the rank, column, and null spaces of A by rank(A), col(A), and null(A), respectively.2 Finite-sample data-generation processThe existing studies on analyzing the identifiability guarantees of CRL inevitably necessitate the availability of an infinite number of samples as otherwise, perfect identifiability is rendered impossible. In this paper, our objective is identifiability analysis assuming only a finite number of samples are available -a significant departure from the infinite-sample regime. In this section, we specify a general data-generating process consisting of a latent causal space and an unknown linear transformation that maps the latent variables to observed variables.Latent causal model. We have a latent space of n causally related random variables. We denote the latent causal variables by Z ≜ [Z 1 , . . . , Z n ] ⊤ . The causal relationships among the variables of Z are represented by a directed acyclic graph G with n nodes where the i-th node represents Z i . We denote the parents, children, and ancestors of a node i ∈ [n] in G by pa(i), ch(i), and an(i), respectively. We denote the probability density function (pdf) of Z by p Z , which is assumed to be well-defined without any zeros over lower-dimensional manifolds and to have full support on R n . Given the DAG G, p Z factorizes according toEQUATIONWe assume that the conditional pdfs {p i : i ∈ [n]} are continuously differentiable with respect to z.We call a permutation π = (π 1 , . . . , π n ) of [n] a valid causal order if for all i, j ∈ [n], i ∈ pa(j) implies π i < π j . Without loss of generality, we assume that (1, . . . , n) is a valid causal order.Transformation model. The latent variables Z are mapped to the observed variables X ≜ [X 1 , . . . , X d ] ⊤ through an unknown linear transformation G ∈ R d×n , where d ≥ n. Specifically,EQUATIONwhere G has full column rank. We denote the pdf of X by p X . Note that owing to the generation process of X, the pdf p X is supported on an n-dimensional subspace embedded in R d . We denote the L 2 norm of any function f with finite variance under p X by ∥f ∥ 2 p X ≜ E p X ∥f (x)∥ 2 2 . Intervention model. We consider CRL under interventions, and in particular focus on soft interventions, as the most general form of interventions. Hence, in addition to the observational model specified by (2) , we consider n single-node interventional environments {E m : m ∈ [n]}. We assume that the node intervened in the environment E m is unknown and denote it by I m . Leaving one node unintervened renders identifiability impossible [7] . Hence, we inevitably have{I m : m ∈ [n]} = [n].Applying a soft intervention on node i changes its causal mechanism specified by the observational conditional pdfp i (z i | z pa(i) ) to a distinct interventional conditional pdf q i (z i | z pa(i) ).The conditional pdfs {q i : i ∈ [n]} are assumed to be continuously differentiable with respect to z. Subsequently, the pdf of Z in the interventional environment E m , denoted by p m , factorizes according toEQUATIONTo distinguish the observational and interventional data, we denote the latent and observed random variables in environment E m by Z m and X m , respectively. Interventions do not alter the transformation matrix G, indicating that similarly to (2), we have X m = G • Z m . Throughout the paper, we adopt the convention that E 0 refers to the observational environment.Score functions. The score function associated with a pdf is defined as the gradient of its logarithm. We denote the score functions associated with the observational distributions of Z and X by s Z and s X , respectively, i.e.,EQUATIONSimilarly, we denote the score functions associated with the distributions of Z and X in environmentEQUATIONIn our algorithm and analysis, we will analyze the score variations across different interventions. For this purpose, we define the score difference functions of Z and X between environments E m and E 0 asEQUATIONTo ensure that interventions affect the target node and its parents distinctly, we adopt the following assumption. This assumption holds for a large class of models, e.g., additive noise models under stochastic hard interventions [5, Lemma 2] . Assumption 1 ([5, Assumption 1] ). For any m ∈ [n], and for all k ∈ pa(I m ), the term[d m Z ] k /[d m Z ] I m is not a constant function of z.Finite sample data. We assume that we have only N data samples under each of the environments. We denote the N samples of Z m and X m byEQUATIONAccordingly, for each sample j, we concatenate all latent and observational samples under different environments to create sample matricesEQUATIONFinally, we define the set of samplesX N ≜ {X j : j ∈ [N ]} and Z N ≜ {Z j : j ∈ [N ]}. 
finite-sample identifiability objectives: 
 The CRL framework specified in Section 2 has two objectives: use the N observational samples X N = {X j : j ∈ [N ]} to recover the causal graph G and the latent causal variables Z N = {Z j : j ∈ [N ]}. In this section, we formalize the inference rules and their associated fidelity metrics for recovering G and Z N .Latent graph recovery. Define G as the set of all DAGs on n nodes. We define Ĝ as a generic estimator of the latent graph G, i.e., Ĝ : R d×(n+1) N → G .To quantify the accuracy of the estimate Ĝ(X N ), we provide the following probably approximately correct (PAC) guarantee, which is the probabilistic counterpart of the standard identifiability up to ancestor definition in the interventional CRL literature [5, 7, 11] .Definition 1 (δ-PAC graph recovery). Graph estimate Ĝ(X N ) achieves δ-PAC latent graph recovery if, with probability at least 1 -δ, the transitive closures of Ĝ(X N ) and G are isomorphic.Latent variables recovery. We investigate a two-step estimator for the latent variables Z N . First, we define a generic estimator for the pseudoinverse of G asEQUATIONThen, given an estimate H(X N ) for G † , we estimate {Z j : j ∈ [N ]} according toEQUATIONWe provide the following definition to quantify the fidelity of these estimates with respect to the ground truth variables. Definition 2 ((ϵ, δ)-PAC variables recovery). The estimate H(X N ) achieves (ϵ, δ)-PAC latent variables recovery if the estimated causal variables { Ẑj : j ∈ [N ]} satisfyEQUATIONwhere for all i ̸ ∈ pa(j), C pa satisfies (C pa ) i,j = 0, and, with probability at least 1 -δ, we have ∥C err ∥ 2 ≤ ϵ.We note that Definition 2 is the probabilistic counterpart of the standard identifiability up to ancestor definition in interventional CRL literature.Noisy score model. The score-based CRL framework uses properties of score function differences to build the estimators for the latent graph and variables. Therefore, when we only have access to finite data, we need to estimate the score functions. We denote a generic score function estimator for environment E m by ŝm X (x; X N ) :EQUATION) When the dependence is clear from the context, we will drop the explicit dependence of ŝm X on X N . Similarly to (6) , we define the estimated score difference functionsEQUATIONIn this paper, we focus on consistent score difference estimators { dm X : m ∈ [n]}, i.e., they satisfy convergence in probability. Specifically, for any ϵ, δ > 0, there exists N (ϵ, δ) ∈ N such thatEQUATIONNotably, we note that many score estimators are known to be consistent [12, 16, 17, 18] . 
finite-sample crl algorithms: 
 In this section, we design finite-sample interventional CRL algorithms through which we establish finite-sample identifiability guarantees as well as bounds on the associated sample complexities.Our algorithms fall within the score-based category of CRL algorithms [5] . The main intuition of this framework is that score functions in the observed space contain all the information needed for recovering the latent graph G and the inverse transform G † . Specifically, two metrics are pivotal for retrieving the latent space:1. Score differences: As shown in Lemma 1, the nonzero entries of the latent score differences d m Z encode the graph structure and the observed score differences d m X are generated using the inverse transform G † . 2. Score difference correlations: As shown in Lemma 2, the column space of the correlation matrix of the score differences contains crucial information about the latent graph and the inverse transform, which we will leverage to form estimates for the latent graph G and the inverse transform G † .To proceed, for m ∈ [n] we denote the correlation matrices of d m X and dmEQUATIONNext, we formalize two critical properties of score differences and their correlation matrix. Specifically, the following lemma states that the sparsity pattern of the latent score differences d m Z exposes the structure of the latent graph, and the observed score differences d m X preserve this information through the pseudoinverse of G. Lemma 1 ( [5] ). Score function differences in the latent space satisfy, for all m ∈ [n],EQUATIONFurthermore, the score differences in the observed domain are given byEQUATIONThe next lemma specifies that the structure of the column space of correlation matrix R m X is heavily constrained by the graph structure and the inverse transformationG † . Lemma 2. For any m ∈ [n], we have col(R m X ) ⊆ span (G † ) ⊤ :,i : i ∈ pa(I m ) .While these two properties enable recovering the latent graph and variables, the ground truth score functions d m Z (z) and correlation matrix R m X are unknown. We only have access to their noisy counterparts, as estimated from the observed data. In our algorithms, we use methods with soft decision rules that can counter the effects of the errors introduced by the score estimation procedure. These are facilitated by the following three key definitions. for j ∈ (t + 1, n) do 6: M t,j ← V j \ ({π t } ∪ ĉh(π t ))▷ set for determining whether π t → π j 7:if col( m∈Vt Rm ; η) ⊥ γ null( m∈Mt,j Rm ; η) then 8:Add π t → π j and π t → l to Ĝ for all l ∈ ĉh(π j ) in Ĝ 9: Return Ĝ Definition 3 (Approximate column space). We define the η-approximate column space of a positive semidefinite matrix A ∈ R d×d , denoted by col(A; η), as the span of the eigenvectors of A associated with the eigenvalues that are strictly greater than η.Definition 4 (Approximate null space). We define the η-approximate null space of A, denoted by null(A; η), as the orthogonal complement of col(A; η).Definition 5 (Approximate subspace orthogonality). We say two subspaces A and B with orthonormal bases A and B are γ-approximately orthogonal, denoted byA ⊥ γ B, if ∥B ⊤ A∥ 2 ≤ γ.In these definitions, setting η = γ = 0 yields the standard definitions of column and null spaces and orthogonality. Next, we describe our algorithms for latent graph recovery (Algorithms 1 and 2) and latent variables recovery (Algorithm 3). We note that the algorithms for latent graph and latent variable recovery algorithms are fully decoupled, enabling the independent recovery of the graph and the latent variables.Algorithm 1 -Causal order estimation. In this algorithm, we estimate a permutation π of [n] such that I • π is a valid causal order. This serves as an intermediate step for estimating the latent graph.For this purpose, we use the noisy correlation matrices { Rm X : m ∈ [n]} to identify the leaf nodes with no children. In particular, the key property is that when I k is a leaf node, by carefully selecting the threshold η, the approximate column space of the term( m∈[n] Rm X -Rk X ) has dimension n -1.Precisely, with high probability, we haveEQUATION)After finding a leaf node, we iteratively identify the youngest node among the remaining set of nodes.Leveraging this, we construct the permutation π, which consists of nodes ordered from the eldest to the youngest. In Lemma 3, we show that I • π is a valid causal order with high probability.Algorithm 2 -Latent graph estimation. In Algorithm 2, we use the causal order found in Algorithm 1 and correlation matrices { Rm X : m ∈ [n]} to form a graph estimate Ĝ. We build this graph iteratively by considering a candidate edge π t → π j for all possible (t, j) pairs starting from a leaf node t. The key property that we leverage is that we can form two subsets V t and M t,j of [n] using π,Algorithm 3 Inverse transform estimation 1: Input: Rm for all m ∈ [n], η ≥ 0 2: H ← 0 n×d 3: for t ∈ (1, . . . , n) do 4: C ← col( Rπt ; η) 5:H πt,: ← v for v ∈ C and ∥v∥ 2 = 1 6: Return H t, and j such that, for sufficiently small γ, with high probability the following approximation holds.EQUATION)We check each edge π t → π j and add the detected edges to the estimated graph Ĝ. In Theorem 1, we show that this procedure guarantees a PAC latent graph recovery. Algorithm 3 -Inverse transform estimation. In Algorithm 3, we build our inverse transform estimate H one row at a time. The key property we use is that for any m ∈ [n] and unit vector v ∈ col( Rm X ), the following error term is small with high probabilityEQUATION)In other words v ⊤ • X is approximately equal to Z I m up to mixing with {Z i : i ∈ pa(I m )}, which conforms to our latent variables recovery objective. Note that in the noise-free setting, an exact equality holds due to Lemma 2 and property G † •G = I n . Based on this observation, in our algorithm, we construct our inverse transform estimate by setting row m of H to a unit vector v ∈ col( Rm X ) for all m ∈ [n]. In Theorem 2, we show that Algorithm 3 achieves a PAC latent variables recovery. 
sample complexity analysis: 
 In this section, we analyze the sample complexity of our CRL algorithm to establish an achievable sample complexity for CRL.Threshold selection. A critical step in our algorithm designs is choosing thresholds η and γ specified in Definitions 3-5. These thresholds determine the approximate ranks of noisy matrices and approximate orthogonality between subspaces, respectively. Specifically, we select η such that the approximate rank of Rm X tracks that of R m X with high probability. To proceed, we denote the minimum nonzero eigenvalue among arbitrary R m X sums byEQUATIONwhere λ(A) ∈ R d denotes the vector of eigenvalues of A. In our algorithms, we let η ∈ (0, η * ) and show that this choice ensures that eigenvalues of the null and column spaces can be separated via η.For the approximate orthogonality test in Definition 5, we show that γ * defined below serves as a bound on how close the non-orthogonal column spaces of R m X sums become orthogonal.EQUATION)Similarly, we let γ ∈ (0, γ * ) in our algorithms. For any choice of η and γ, we define η * ≜ min{η, η * -η} and γ * ≜ min{γ, γ * -γ}. We note that we do not have to know η * and γ * a priori, and can include routines to estimate them in practice. These estimates do not have to be highly accurate since even rough estimates suffice to choose reliable thresholds η and γ. Specifically, based on estimates for η * and γ * , we can choose "safe" thresholds for η and γ, e.g., one-fourth of the estimates for η * and γ * , so that (i) with high probability we satisfy the requirement η ∈ (0, η * ) and γ ∈ (0, γ * ), and (ii) we can avoid collecting excessive samples and compromising the sample complexity bounds. We provide the details of constructing such estimates of η * and γ * in Appendix G.Sample complexity results. We provide two sets of sample complexity analyses for the latent graph and the latent variables. We estimate the latent graph in two steps: (i) causal order estimation (Algorithm 1) and (ii) latent graph estimation (Algorithm 2). Algorithm 1 only employs approximate rank tests. Therefore, we first show that the approximate rank of the sum of correlation matrices{ Rm X : m ∈ [n]} tracks that of {R m X : m ∈ [n]} with high probability. To proceed, we define the following constants that appear repeatedly in our analysis:EQUATION)Lemma 3. Let η ∈ (0, η * ). For any δ > 0, N rank (δ) samples suffice to ensure that with probability at least 1 -δ,EQUATIONwhereEQUATIONNote that the function N (•, •) is an inherent property of the score difference estimator specified in (15) , and it is monotonically decreasing in its second argument. Using this lemma, we can show that Algorithm 1 returns a permutation π such that I • π is a valid causal order with high probability. Lemma 4. Let η ∈ (0, η * ). Under Assumption 1, for any δ > 0, N rank (δ) samples suffice to ensure that with probability at least 1 -δ, I • π is a valid causal order, where π is the output of Algorithm 1.Next, we show that Algorithm 2 achieves latent graph recovery with high probability.Theorem 1 (Sample complexity -Graph). Let η ∈ (0, η * ) and γ ∈ (0, γ * ). Under Assumption 1, for any δ > 0, N G (δ) samples suffice to ensure that the output Ĝ(X N ) collectively generated by Algorithms 1 and 2 satisfies δ-PAC graph recovery, whereEQUATIONWe note that constants β and β min are sample-independent. Then, (27) specifies the sample complexity as a function of the hyperparameters η, γ, target reliability δ ∈ (0, 1), as well as the latent dimension n. Next, we state the complementary result for recovering the latent variables.Theorem 2 (Sample complexity -Variables). Let η ∈ (0, η * ). For any ϵ > 0 and δ > 0, N Z (ϵ, δ) samples suffice to ensure that the output H(X N ) of Algorithm 3 satisfies (ϵ, δ)-PAC causal variables recovery, whereEQUATIONTheorems 1 and 2 collectively specify an extent of identifiability achievable in the finite-sample regime. Importantly, we note that since we can establish (ϵ, δ)-PAC identifiability guarantees for any vanishing ϵ, δ > 0 using finite samples, Algorithms 2 and 3 are consistent estimators of the latent graph and the inverse transformation up to the corresponding equivalence classes.Finally, we note that the latent variables are recovered up to mixing with parents, as specified in Theorem 2. This result is a refinement of the existing latent variable recovery results under soft interventions in existing literature [5, 7] , which recover the latent variables up to mixing with ancestors.RKHS-based score estimator. So far, we have characterized the sample complexity for any consistent score difference estimator. Next, we specialize the results by adopting the RKHS-based score estimator of [12] . We adopt this particular choice since it has known non-asymptotic sample complexity properties. To our knowledge, it is the only score estimator equipped with such a guarantee. To use the sample complexity of this score estimator in our paper, we first show that it is consistent in the sense of (15) . For the formal statement of the following property and its attendant assumptions, we refer to Appendix E. Lemma 5 (Informal). Assume that s X and s m X satisfy the conditions of the RKHS-based score estimator in [12] . Then, for any given δ and ϵ, the convergence specified in (15) holds whenEQUATIONwhere κ and C are sample independent constants that depends only on p X , p m X for m ∈ [n] and the structure of the RKHS.We use Lemma 5 to customize the general sample complexity bounds in Theorems 1 and 2 to the RKHS-based score estimator. For this purpose, we first provide our sample complexity result for recovering the latent graph. Theorem 3 (RKHS-based sample complexity -Graph). Let η ∈ (0, η * ) and γ ∈ (0, γ * ). Under Assumption 1 and the conditions of Lemma 5, N G (δ) samples suffice to ensure that Ĝ(X N ) collectively generated by Algorithms 1 and 2 satisfies δ-PAC graph recovery, whereN G (δ) = max C ϵ G , 2 √ 2κ 2 4• log 8n δ, and ϵ G ≜ min βη * n , βη * γ * 2n , β min .(30) Remark 1 (RKHS-based error bound -Graph). Theorem 3 implies that using N samples, the output Ĝ(X N ) of Algorithms 1 and 2 satisfies δ G (N )-PAC graph recovery, whereEQUATION)Similarly, we leverage Lemma 5 to specialize the general sample complexity in Theorem 2 to the RKHS-based score estimator. Theorem 4 (RKHS-based sample complexity -Variables). Let η ∈ (0, η * ). Under the conditions of Lemma 5, N Z (ϵ, δ) samples suffice to ensure that the output H(X N ) of Algorithm 3 satisfies (ϵ, δ)-PAC causal variables recovery, whereEQUATIONEQUATION)We note that the results in Theorems 3 and 4 show a more explicit dependence of the sample complexity on the latent space dimension n. We observe that the sample complexity of latent graph recovery explicitly depends on δ and the latent dimension n according to O((n log n δ ) 4 ). Similarly, the sample complexity of latent variables recovery explicitly depends on ϵ, δ, and the latent dimension n according to O((√ n ϵ log n δ ) 4). We note that constants β, β min , η * , and γ * are model parameters, and their scaling behavior in terms of n and d are investigated numerically in the next section. 
experiments: 
 In this section, we perform numerical assessments of our analyses to provide complementary insight into the sample complexity results of Section 5. Specifically, we evaluate the variations of the model constants with respect to problem dimensions n and d. Furthermore, we also evaluate performance variations of the finite-sample algorithm in terms of the variations in the score estimation error. For this purpose, we focus on the mean squared error (MSE) of the score estimator, specified byE∥ dm X -d m X ∥ 2 p X .These evaluations facilitate a better understanding of the properties of the CRL problem that are not explicit in our theoretical queries. 2Experimental details. We consider problem dimensions n ∈ {3, 5, 10} and d ∈ {n, 15} and generate G using Erdős-Rényi model with density 0.5 on n nodes. We adopt linear Gaussian models as the latent causal model. We consider N ∈ {10 2.5 , 10 3 , 10 3.5 , 10 4 , 10 4.5 , 10 5 } samples, and generate 100 latent models for each triplet (N, n, d). We generate N samples of Z from each environment for each latent model. Transformation G ∈ R d×n is randomly sampled under full-rank and bounded condition number constraints, and the observed variables are generated as X = G • Z. Due to the Gaussian structure of Z, the observed variables X also have multivariate Gaussian distributions with score function s X (x) = -Θ • x, where Θ is the precision matrix of X. Therefore, we estimate s X using the sample precision matrix estimated from N samples, ΘN , as ŝX (x;X N ) = -ΘN • x.Results. The sample complexity results depend on dimension n and d through their explicit presence as well as their implicit effect on model-dependent constants β, β min , η * , and γ * , as established in Theorems 1 and 2. First, we observe that all these parameters are independent of d. This is also expected theoretically since we can project down the d dimensional observed space to the n dimensional col(G). In Figure 1a , we illustrate the dependence of these constants on data dimensions n and d. We observe that γ * remains constant, while β and β min are decreasing as n increases, but the decay is not steep. On the other hand, we observe that η * decays exponentially with n. γ * is implicitly controlled by the condition number, which is kept bounded. We recall that β min is the minimum among {∥d m X ∥ p X : m ∈ [n]}, and similarly, β is the inverse of the maximum of {∥d m X ∥ p X : m ∈ [n]}, so they are expected to become smaller as n increases. Finally, for η * , we note that the definition in ( 22) is an order statistic among exponentially many, 2 n , subsets.Next, we note that all the imperfections in the graph and variable recovery are due to the imperfections in score estimates. In fact, this is the bottleneck in our decisions: the key impact of finite samples versus infinite samples is the degradation in the quality of score estimates. To have a direct insight into how the score estimation noise power translates into decision imperfection, we assess the quality of graph recovery success rate versus varying degrees of mean score error of score difference estimatesE[∥ dm X -d m X ∥ 2 p X ] for different triplets (N, n, d).In Figure 1b , each data point corresponds to a different sample size N , and the legend entries {3, 5, 10} denote the latent dimension n. This figure illustrates the scatter plot of the error probability δ versus MSE for different parameter combinations and the corresponding local regression curves. It is observed that in the low MSE regime, 1 -δ decays linearly with respect to log E[∥ dm X -d m X ∥ 2 p X ] and it plateaus as the MSE increases. This trend is due to our algorithms' high sensitivity to errors in estimating approximate matrix ranks. 
conclusion: 
 In this paper, we have established the first sample complexity results for interventional causal representation learning. Specifically, we have characterized upper bounds on the sample complexity of latent graph and variables recovery in terms of the finite sample performance of any consistent score difference estimator. We have, furthermore, adopted a particular score estimator [12] to derive explicit sample complexity statements.Our sample complexity results are given for partial identifiability (up to ancestors) in the soft intervention and linear transformation setting. Establishing sample complexity results for perfect identifiability via hard interventions or considering general transformations are important future directions. 
broader impacts: 
 Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?Answer: [NA] Justification: This paper is mostly theoretical and does not pose potential negative societal impacts.Guidelines:• The answer NA means that there is no societal impact of the work performed.• If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact.• Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations. • The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is a direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster. • The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional or unintentional) misuse of the technology. • If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of ML). 
safeguards: 
 Question: Does the paper describe safeguards that have been put in place for responsible release of data or models that have a high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?Answer: [NA] Justification: This paper is mostly theoretical and does not pose such risks.Guidelines:• The answer NA means that the paper poses no such risks.• Released models that have a high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters. • Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images. • We recognize that providing effective safeguards is challenging, and many papers do not require this, but we encourage authors to take this into account and make a best faith effort.12. Licenses for existing assets Question: Are the creators or original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?Answer: [Yes] Justification: Portions of the publicly available code of [5] , available under Apache 2.0 license, is adopted in the code of our experiments.Guidelines:• The answer NA means that the paper does not use existing assets.• The authors should cite the original paper that produced the code package or dataset.• The authors should state which version of the asset is used and, if possible, include a URL. • The name of the license (e.g., CC-BY 4.0) should be included for each asset.• For scraped data from a particular source (e.g., website), the copyright and terms of service of that source should be provided.• If assets are released, the license, copyright information, and terms of use in the package should be provided. For popular datasets, paperswithcode.com/datasets has curated licenses for some datasets. Their licensing guide can help determine the license of a dataset. • For existing datasets that are re-packaged, both the original license and the license of the derived asset (if it has changed) should be provided. • If this information is not available online, the authors are encouraged to reach out to the asset's creators. • We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors to adhere to the NeurIPS Code of Ethics and the guidelines for their institution. • For initial submissions, do not include any information that would break anonymity (if applicable), such as the institution conducting the review. 
 
 Appendix: Learning: Appendices A Proof structureWe structure the proofs of our main results as follows.• Appendix B -Infinite sample guarantees: First, we prove identifiability results in the noise-free setting, that is, using the ground truth R m X matrices instead of using the estimated Rm X . These intermediate results lay the groundwork for the proofs for the finite sample setting.• Appendix C -Bounded noise guarantees: Next, we show that if estimation error for the score difference functions, ∥ dm X -d m X ∥ p X , is bounded, then the noise-free identifiability guarantees transfer to the noisy setting.• Appendix D -Sample complexity guarantees: We use the noisy identifiability guarantees to prove that for any consistent score difference estimator, we can guarantee (ϵ, δ)-PAC guarantees for CRL identifiability objectives.• Appendix E -RKHS-based score estimator: Finally, we use the finite sample guarantees provided by the RKHS-based score estimator of [12] to adapt the sample complexity statements for a generic estimator to explicit guarantees in terms of model-dependent constants.To be used throughout the proof, we define the correlation matrices ofEQUATION)We recall the correlation matrices R m X and Rm X specified in (16) . For a subset M ⊆ [n], we use the following shorthand notations for sums of the matrices in {R m X : m ∈ M} and { RmEQUATIONIn this section, we provide the identifiability guarantees of Algorithms 1-3 with no error and probability 1 when using {R m X : m ∈ [n]} as input. Denote the image of a function f : R n → R k for any k ∈ Z + by im(f ), that is,EQUATIONFirst, we note a property of the column space of correlation matrices. Lemma 6. For a continuous function f : R n → R k , span of the image of f equals toEQUATIONProof: Note that two subspaces are equal if and only if their orthogonal complements are equal. We denote the orthogonal complement of span im(f ) byEQUATION)The orthogonal complement of the column space of a symmetric matrix is the null space, which we denote byEQUATION)Since any correlation matrix is positive semidefinite, using Proposition 2, (39) can be written asEQUATION)Since f is a continuous function, so is y ⊤ • f for any y ∈ R k . We can then use [5, Proposition 2], which is stated below for the sake of conciseness. Proposition 1 ([5, Proposition 2] ). Consider two continuous functions f, g : R n → R. Then, for anyEQUATION)Under this proposition, (41) can be written asEQUATIONwhich is the same as the definition of S 1 in (38). This concludes the proof. Corollary 1. For a continuous function f : R d → R k , the span of the image of f • G equals to the column space ofEQUATIONWe can specialize this lemma forEQUATIONandEQUATIONWe can state (17) from Lemma 1 equivalently asEQUATIONRecall from (18) in Lemma 1 thatEQUATIONThen, due to (47), for any given z ∈ R n and corresponding x = G • z, we haveEQUATIONSince d m X (x) lies in the subspace span{(G † ) i,: : i ∈ pa(I m )} for any x ∈ col(G), we haveEQUATIONDue to (46), we haveEQUATIONwhich concludes the proof.We can use Lemma 2 to show that in the infinite sample regime, Algorithm 3 works correctly. Specifically, let us show the following. Lemma 7 (Infinite samples -Encoder). For any v ∈ col(R m X ), we haveEQUATIONEQUATIONIn other words, any v ∈ col(R m X ) can be expressed asEQUATIONfor some constants {c i : i ∈ pa(I m )}. Then, we haveEQUATION)Since G has full column rank, the pseudoinverseEQUATIONwhere e i denotes the i-th standard basis vector. Substituting this, (55) becomesEQUATIONNote that i-th entry of a row vector w ∈ R n is given by w • e i . Therefore, we haveEQUATION)This concludes the proof.In this section, we adapt the workhorse lemmas of [5] to our setting. We first note that due to (18) and the definitions of R m X and R m Z in ( 16) and (34), we haveEQUATIONwhere G † is a full row rank matrix. Therefore, for any B ⊆EQUATIONNext, we note that for a function f : R n → R k , the notation R(f ) in [5] is used to denote the "rank" of a function, and it corresponds to dim span im(f ) in the notation of this paper. Using Proposition 3 and Lemma 6, we can see that the notation R(F) for a set of functions F can be given asEQUATION= dim col(EQUATIONUsing ( 66) and (62), we can adopt [5, Lemmas 5 and 6] using {R m X : m ∈ [n]} matrices. We note that (59) and (60) together with Lemma 2 imply that, for any ancestrally closed A ⊆ [n] and k ∈ A for which there exists j ∈ A \ {k} such that I k ∈ pa(I j ), we haveEQUATIONWe note that due to the propertyEQUATIONEQUATIONWe will show that Algorithm 1 gives a valid causal order in the infinite sample regime using Lemma 8. Lemma 9 (Infinite samples -Causal order). Setting η = 0, I • π is a valid causal order, where π is the output of Algorithm 1 using {R m X : m ∈ [n]}.Proof: The proof of this lemma follows closely to the proof of [5, Lemma 7] . Note that setting η = 0 means that approximate rank and column/null spaces are actually equal to the standard rank and column/null space definitions. We note that I • π is a valid causal order if and only if, for any t ∈ [n], we have de(I πt ) ⊆ I πi : i > t .(70) We prove (70) by induction as follows.In the base case, for t = n, we have V t = [n]. Therefore, using Lemma 8, we see thatEQUATIONTherefore, for t = n, rank(REQUATIONSince { I πi : i > t } = ∅ for t = n, the base case is proved by choosing any such k as π t .For the induction step at time t ∈ {n -1, . . . , 1}, we assume that, for any u ∈ {t + 1, . . . , n}, we have de(I πu ) ⊆ { I πi : i > u }. By construction, V t = [n] \ { π u : u > t }. Therefore, I Vt is ancestrally closed: For any k ∈ I Vt , k cannot be a descendant of I [n]\Vt = { I πu : u > t } due to the induction hypothesis. In other words, pa(I Vt ) = I Vt . Therefore, we can use Lemma 8 asEQUATIONTherefore, for any t ∈ { n -1, . . . , 1 }, rank(REQUATION)This is equivalent toEQUATIONTherefore, choosing any such k as π t satisfies the induction hypothesis for t as well, and the proof by induction is concluded.Lemma 10 (Infinite samples -Graph). Setting η = γ = 0, the collective output of Algorithms 1 and 2 using {R m X : m ∈ [n]} is equal to the graph isomorphism of the transitive closure of G by permutation I.Proof: The proof of this lemma follows closely to the proof of [5, Lemma 8] . Note that setting η = 0 means that approximate ranks and column/null spaces are equal to the standard definitions. Similarly, γ = 0 means that approximate orthogonality is equivalent to the standard orthogonality definition. We first note that the transitive closure of two DAGs are equal if and only if their transitive reductions are equal. Therefore, the lemma statement is equivalent toEQUATIONWe will prove this result by induction over the transitive reduction edges. In the rest of the proof, using Lemma 9, we will leverage that I • π is a valid causal order.In the base case, we have t = n -1, and thus necessarily j = n. The only possible path between I πt to I πj is therefore a transitive reduction edge. Since we did not add any children to node π t yet, the set M t,j is equal to V j \ {π t }. Noting that I Vj is ancestrally closed and the only possible child of I πt is I πj when I • π is a valid causal order, we can use Lemma 8 to obtainEQUATIONSpecifically, using (67), we can specify the exact column spaces we can obtain from M t,j asEQUATION)Then, using (68), the null space of R Mt,j X is given byEQUATION)The test in Algorithm 2 involves checking whether the column space of R Vt X and the null space of R Mt,j X are orthogonal. Note that I Vt is ancestrally closed, therefore, using (67), we getEQUATION)Next, we note that the orthogonality test which uses orthonormal bases of subspaces (Definition 5) is equal to the following. Let A ∈ R d×k1 and B ∈ R d×k2 be full column rank matrices. Then,EQUATION)This holds since A • A † is the orthogonal projection operator onto col(A). Therefore, given orthonormal bases A 0 and B 0 for A and B, we haveEQUATION) which is precisely the metric used in Definition 5. To use this property, let us first consider the matrix G :,A which is formed by stacking the independent column vectors in set {G :,i : i ∈ A} for some subset A ⊆ [n] . Since G has full column rank, we see that the pseudoinverse of G :,A is exactly given by GEQUATIONwhere the matrix G † A,: is formed by stacking the independent row vectors in set {(G † ) i,: : i ∈ A}. Using this notation, if we haveEQUATIONthen the orthogonality metric is equal toEQUATION)Due to the construction of sets V t and M t,j , we have I Vt \ I Mt,j = I πt . Therefore, the only entry in I Vt that may be missing from pa(I Mt,j ) is I πt . Therefore, previous equation is equal toEQUATION)In summary, if (84) holds, the orthogonality inspected in Algorithm 2 happens if and only ifEQUATION)In the base case, since t = n -1, this holds if and only if I πt ∈ pa tr (I πj ), and the base case is done.For the induction hypothesis, we assume that for any u ∈ { n -1, . . . , t + 1 }, all the transitive reductive edges has been correctly identified, that is, for any i ∈ [n], π u ∈ pa tr (π i ) if and only if I πu ∈ pa tr (I πi ). We show that the induction hypothesis is also satisfied for t and j ∈ [n] using induction on j.In the base case, we have j = t + 1. This case works exactly the same way as the t = n -1 and j = n case. For the induction hypothesis, assume that for any u ∈ { n + 1, . . . , j -1 }, all the transitive reductive edges between I πt and I πu have been correctly identified, that is, π t ∈ pa tr (π u ) if and only if I πt ∈ pa tr (I πu ). We prove that the hypothesis also holds for j as follows. First, we claim that I Mt,j ∪{πt} is an ancestrally closed set. Note that I Vj is an ancestrally closed set by construction. Since the induction hypothesis ensures that all transitive reduction relations between I πu and I πi for any u ∈ { t + 1, . . . , n -1 } and i ∈ [n] are recovered, this also means that transitive closure relations, that is, descendants, are also recovered. Namely, the induction hypothesis impliesEQUATIONNext, the induction hypothesis for j implies that for any k ∈ {t + 1, . . . , j -1 }, we have also identified all the transitive reductive edges I πt → I π k . Therefore, for set M t,j ∪ {π t } = V j \ ĉh(π t ), all the elements in ĉh(π t ) has the property that their descendants have been fully identified, and therefore are excluded from the set ( Ĝ keeps the transitive closure), which means the set I Mt,j ∪{πt} is ancestrally closed. Since this set is ancestrally closed, we haveEQUATION+ 2 • E p X d m X (x) • ( dm X (x) -d m X (x)) ⊤ 2 .(100)Note that the first term is the spectral norm of a covariance matrix, and is upper bounded byEQUATION(101)EQUATION)For the second term, we can upper bound the spectral norm by Frobenius norm asEQUATION)Using Cauchy-Schwarz inequality,EQUATION(106)EQUATION)Then, the rank statement about RMEQUATION)This statement can be equivalently written asEQUATIONNext, we find an upper bound on χ in terms of ∥R m X -R m X ∥ 2 using Weyl's inequality (Lemma 12). χEQUATION)Therefore, for any η ∈ (0, η * ), leveraging the condition (120), we have (i) If |M| = 1 and M = {m}, we can use (123) to show that the conditionEQUATION)EQUATION)(ii) For any generic M, we can use (124) to show that the conditionEQUATIONEQUATION)This concludes the proof.Since we have a sufficient condition to ensure that approximate ranks of { RM X : M ⊆ [n]} are equal to the ranks of {R M X : M ⊆ [n]}, now we can state a sufficient condition for finding a valid causal order. Lemma 14 (Bounded noise -Causal order). Let η ∈ (0, η * ). IfEQUATIONthen I • π is a valid causal order, where π is the output of Algorithm 1.Proof: We note that Algorithm 1 only uses approximate rank statements for { RM X : M ⊆ [n]}. Therefore, if the approximate rank statements are equal to the ranks of {R M X : M ⊆ [n]}, then, the infinite sample guarantees of Algorithm 1 provided in Lemma 9 immediately transfer to the noisy R m X regime. Next, we note that Algorithms 2 and 3 use the direction information of the approximate column and null spaces of { RM X : M ⊆ [n]} as well. Therefore, we need to bound the "distance" of the estimated subspaces from their noise-free counterparts. The main results we will use for this purpose are the Davis-Kahan sin θ theorems [19] , specialized for approximate column and null spaces.Lemma 15 ([19, sin θ theorem]). Consider two k × k real symmetric matrices A and Â. Denote the orthonormal bases of the column and null spaces of A by A 1 and A 0 , respectively, such thatEQUATIONwhere C 1 is a symmetric matrix with the same eigenvalues as A. Similarly, denote the orthonormal bases for the approximate column and null space of Â by Â1 and Â0 , such thatEQUATIONwhere eigenvalues of Ĉ0 are that of the approximate null space of Â, and the eigenvalues of Ĉ1 are that of the approximate column space of Â. DefineEQUATIONwhere λ max and λ min denote the maximum and minimum eigenvalue of a symmetric matrix. Then, given that the approximate rank of Â is equal to the rank of A and χ > 0, we haveEQUATIONandEQUATION)Lemma 16 ([19, Symmetric sin θ theorem]). Under the setting described in the sin θ theorem above, we haveEQUATIONandEQUATION)We note that χ in these sin θ theorems can be bounded using Weyl's inequality (Lemma 12) and the minimum nonzero eigenvalue of matrix A as follows.EQUATIONTherefore,EQUATIONUsing R M X for A and RM X for Â, and using the definition of η * in (110), we obtainEQUATIONIn order to ensure that the approximate rank of RM X is equal to the rank of R M X , it suffices thatEQUATIONdue to (121). If this holds, we getEQUATIONTherefore, the left-hand sides of (133)-( 136) are upper bounded byEQUATIONWe can use Davis-Kahan symmetric sin θ theorem (Lemma 16) to prove that encoder estimation error is upper bounded by the error in Rm X . Lemma 17 (Bounded noise -Encoder). Let η ∈ (0, η * ). If (125) holds, then the output H of Algorithm 3 satisfiesEQUATIONwhere for all i ̸ ∈ pa(j), C pa satisfies (C pa ) i,j = 0, andEQUATIONProof: We prove this result by considering each row m of H separately. In Algorithm 3, we select H m,: from col( Rm X ; η) for each m ∈ [n]. We start by noting that in Lemma 7, we have shown that using the noise-free score difference correlation matrices, we getEQUATIONWe can decompose any v ∈ col( Rm X ; η) asEQUATIONwhere A 0 and A 1 are the orthonormal bases of the null and column spaces of R m X in line with the notation of the Davis-Kahan sin θ theorem. Note that A 1 • A ⊤ 1 is the orthogonal projection operator onto the column space of R m X , therefore, using (145), we haveEQUATIONthe "erroneous" entries in v ⊤ • G can be written asEQUATIONNote that if we choose v ∈ col( Rm X ; η) as row m of H, we haveEQUATION) Since v ∈ col( Rm X ; η), we have Â1 • Â⊤ 1 • v = v,where Â1 is the orthonormal basis of the approximate column space of Rm X . Using this and ∥v∥ 2 = 1, we can write (149) asEQUATIONUsing Davis-Kahan sin θ theorem (Lemma 15), when (125) holds, we can further upper bound the m-th row of the error term byEQUATION)Note that the spectral norm is upper bounded by the Frobenius norm. Therefore,EQUATIONNext, we investigate the approximate orthogonality test used for Algorithm 2 under noisy score difference correlation matrices. Lemma 18. Let A, B, Â, B be k × k symmetric matrices, such thatEQUATIONEQUATIONwhere[A 0 , A 1 ], [B 0 , B 1 ], [ Â0 , Â1 ], [ B0 , B1] are all k × k orthogonal matrices, and shapes of A 0 -Â0 and B 0 -B0 match. Given that the spectra of (i) C 0 and Ĉ1 , (ii) C 1 and Ĉ0 , (iii) D 0 and D1 , and (iv) D 1 and D0 are all separated with gap χ, we haveEQUATIONEQUATIONProof: We first prove a perturbation bound on the product of orthogonal projectors.EQUATIONUsing Davis-Kahan symmetric sin θ theorem (Lemma 16) yieldsEQUATION)We use this result in two separate cases: col(A 0 ) ⊥ col(B 0 ) and col(A 0 ) ̸ ⊥ col(B 0 ).(i) Case 1: col(A 0 ) ⊥ col(B 0 ). In this case, (164) immediately yields (158).(ii) Case 2:col(A 0 ) ̸ ⊥ col(B 0 ). Let us define ∥A ⊤ 0 • B 0 ∥ 2 = µ ̸ = 0.We can use (164) and triangle inequality to getÂ0 • 0 • B0 • B⊤ 0 2 + 1 χ Â -A 2 + B -B 2 (165) ≥ Â0 • Â⊤ 0 • B0 • B⊤ 0 2 + Â0 • Â⊤ 0 • B0 • B⊤ 0 -A 0 • A ⊤ 0 • B 0 • B ⊤ 0 2(166)EQUATION)Rearranging the terms yields (159).Finally, we show that the graph estimation is correct under bounded noise. Lemma 19 (Bounded noise -Graph). Let η ∈ (0, η * ) and γ ∈ (0, γ * ), whereγ * ≜ min t∈[n]G :,t 2 • (G † ) t,: 2 , and γ * ≜ min{γ, γ * -γ} .(169)Under Assumption 1, the output Ĝ of Algorithm 2 is the graph isomorphism of the transitive closure of the true latent graph G ifEQUATIONProof: Algorithm 2 uses approximate column and null spaces of R M X for different M ⊆ [n]. Therefore, the generic sufficient condition (113), that is,EQUATIONensures that the approximate ranks of all RM X matrices investigated in the algorithm will be equal to the rank of R M X . To ensure that the estimated graph will be correct, we must guarantee that the approximate orthogonality tests will also give correct results. Using Lemma 18, we can show that the test metric is upper bounded in terms of max m∈[n] ∥ Rm X -R m X ∥ 2 .We do so by settingA = R Vt X , Â = RVt X , B = R Mt,j X , B = RMt,j X; and A 0 as the orthonormal basis of the column space of A, and B 0 as the null space of B. As shown in (141), (113) ensures that the spectral gap χ in Lemma 18 statement is at least χ > η * . In Lemma 10, we show that in the infinite sample, noise-free regime, when investigating (t, j), the test metric is given byEQUATION)Since the approximate orthogonality test needs to distinguish even the weakest instance of nonorthogonality from orthogonality, we defineEQUATION)Using Lemma 18, we see that for (t, j) for which I t → I j exists in G, if (113) holds, thenEQUATIONSimilarly, for (t, j) for I t → I j does not exist in G, if (113) holds, thenEQUATIONTherefore, for a choice of γ ∈ (0, γ * ) to yield correct orthogonality decisions, it suffices thatEQUATIONor, equivalently,EQUATIONIn other words, when using γ ∈ (0, γ * ) as the approximate orthogonality threshold, we can correctly identify orthogonal subspaces from non-orthogonal ones whenEQUATION)Due to Lemma 10, if all the approximate rank and orthogonality tests are correct, then the output Ĝ of Algorithm 2 is the graph isomorphism of the transitive closure of the true latent graph G.In this section, we will transform the bounded noise identifiability guarantees in Appendix C to sample complexity statements about a generic consistent score difference estimator. Specifically, recall the consistent score difference estimator model from (15) .EQUATION)Recall the definitions of the constants in (24),β ≜ 4 max m∈[n] d m X p X -1andEQUATIONFirst, we derive the sufficient samples for (ϵ, δ)-approximating Rm X using Lemma 11. Lemma 20. For any ϵ > 0 and δ > 0, using N R (ϵ, δ) samples suffice to ensure thatEQUATIONwhereEQUATIONProof: In Lemma 11, in order to bound the error in Rm X for a specific m ∈ [n], we requireEQUATION)To ensure this for all m ∈ [n], it suffices thatEQUATIONUnder this condition, Lemma 11 shows that the error in correlation matrices R m X is bounded by the error in d m X as RmEQUATIONEQUATIONand col(A + B) = col(A) + col(B) , (216) where + denotes the Minkowski sum between sets, which is defined, for any two subsets A and B of a space, asEQUATIONProof: Using Proposition 2, we can write null(A + B) asEQUATION)Since A and B are both PSD, for any v ∈ R k , we have v ⊤ Av ≥ 0 and v ⊤ Bv ≥ 0. Therefore, v ⊤ (A + B)v = 0 holds if and only if v ⊤ Av = 0 and v ⊤ Bv = 0. That is,EQUATION)We note that, due to Proposition 2, the right hand side of this expression is equal to null(A) ∩ null(B), which concludes the proof for (215).Next, note that the column space of a symmetric matrix is the orthogonal complement of its null space. Then, using (215), we haveEQUATION)Since null(A) and null(B) are subspaces of a finite dimensional vector space, we haveEQUATIONThis concludes the proof for (216).In this section, we provide minor modifications to our algorithms that enable constructing rough estimates for threshold upper bounds η * and γ * in practice.Estimating η * : We note that in Algorithm 1, at step t, all matrices investigated must have either rank t or t -1. Therefore, the t -1-th eigenvalues of all of these matrices are non-zero, or equivalently we have ≥ η * . This observation lets us to iteratively build and refine an upper bound on η * during Algorithm 1, which we can use as a surrogate of η * in the subsequent Algorithm 2 and 3. For Algorithm 1 itself, at time t, we can pick k ∈ V t with minimal t-th eigenvalue-this is a threshold-free way of estimating the causal order.Estimating γ * : We first note that definition in (23) depends entirely on the ground truth decoder matrix G and its pseudoinverse. Secondly, we note that Algorithm 3 can be used independently of Algorithm 2 to generate an estimate H of G † up to some estimation error and possible mixing among rows. We can compute the value Guidelines:EQUATION• The answer NA means that the authors have not reviewed the NeurIPS Code of Ethics.• If the authors answer No, they should explain the special circumstances that require a deviation from the Code of Ethics. • The authors should make sure to preserve anonymity (e.g., if there is a special consideration due to laws or regulations in their jurisdiction).The codebase for the experiments can be found at https://github.com/acarturk-e/ finite-sample-linear-crl. 


Decision: The paper provides the first finite sample complexity results for causal representation learning, albeit under some technical conditions.

All reviewers were supportive of accepting the paper.

For a camera-ready version, please take into the account the comments from all reviewers regarding clarifications, for instance: notation/typos, discussion on the assumptions, experiment details and discussion about causal order estimation.

Reviews: review: paper_summary: paper_summary: This paper provides finite sample analysis for causal representation learning with general latent causal models, single-node soft interventions, and linear mixing functions. Sample complexity for identifying graphs up to ancestors and identifying latent variables up to mixing with parent variables have been studied.
strength_weakness: strength_weakness: 1. The considered CRL setting is a little bit restrictive compared to the recent ones. Specifically, the assumption of linear mixing function may restrict its usage in a lot of real-world scenarios. Meanwhile, the indeterminacy is non-trivial, and further process/disentanglement might still be needed. However, it is understandable since sample complexity might be much more difficult to get for more general models.

2. Perhaps it would be better if there was more discussion on the assumptions.
questions: questions: 1. How to make sure that all the imperfections in the experiments are due to score estimates?

2. Are there any intuitive ideas to extend the range of the finite sample analysis to more general settings in CRL?
limitations: limitations: Yes.
paper_summary: paper_summary: This paper provides finite-sample identifiability results for recovering the latent causal graph and the generating latent variables given observations in a high-dimensional space that have been generated by the latent variables by a linear transformation and single-node soft interventional data with interventions on all nodes. Previous work provides identifiability results only in the infinite-sample regime. The algorithms rely on score-based estimators which have been studied in literature. Specifically the correlation matrix of the score differences are shown to be constrained by the graph and the linear transformation matrix. To recover the graph, first a topological order on the variables is estimated. Verifying whether a particular edge exists is checked by a necessary and sufficient condition on the approximate orthogonality of the column space and null space. For the variables recovery, it is sufficient to consider a unit vector in the column space. An existing RKHS-based score estimator is used to obtain the final sample-complexity results.
strength_weakness: strength_weakness: A few weaknesses (with the major one being point 1) follow:

1. A major confusion for me was that the setting of the thresholds was not clear. If the true scores are not known, then how is \eta set from \eta* which depends on the unknown true score estimates. Similarly for \gamma that depends on \gamma^* that depends on the true transformation matrix. What am I missing? 

2. Overall, the paper has multiple instances of notations not defined which make the paper difficult to parse. For example, Definition 2 has P_I, C_{pa}'s diagonal entries undefined? Lemma 2 has the pseudoinverse subscript undefined. 

3. I also think that the experimental section needs basic details. One suggestion would be to cut out the final RKHS sample complexity bounds that are just plugging (29) into Theorem 1 and 2. 

4. The proofs in the paper first derive results for the infinite-sample regime where score estimates are error-free. Finite-sample guarantees are then obtained by bounding the error in the score estimates. It would make the paper stronger if the infinite-sample regime results are compared with existing work. 

5. In Figure 1b, why are there spikes in the recovery rate even for worse score estimate MSE? Does this go against the main message of the paper that worse score estimates should imply poorer graph recovery and latent variable estimation? 

6. What is the reason behind defining the variables recovery as it is defined in Definition 2. Is there a previous instance of such a definition because this doesn't seem standard. Any rationale would help the reader.

Overall, I think the paper is a solid contribution if these questions are clarified. But my current score can't improve given these questions.
questions: questions: Please refer to the weakness section for major issues. A few miscellanous points:
1) I haven't checked the details but it struck me as weird to notice that there is a specific value of \gamma in (20) that verifies whether an edge exists or not. Is there any intuition for this? 
2) Typos - (22) M \subseteq [n], Line 179. 
3) In Figure 1b, is there a reason why even for better score estimates, there is no point for the 10 nodes case?
limitations: limitations: Implicitly addressed in the conclusions.
paper_summary: paper_summary: This paper establishes finite sample guarantees for recovering the latent causal graph under stochastic soft interventions and observations under a unknown linear transformation. Experiments are conducted and code is given.
strength_weakness: strength_weakness: I did not check the proofs in detail but I do not spot any glaring weaknesses.

The experiments are rather small scale (up to 10 nodes).
questions: questions: - What is $P_I$ in equation (12)?
- In practice, do we need to know $\eta^*$ or $\gamma^*$ to invoke your results? How do we know enough samples have been collected?
limitations: limitations: Nil
paper_summary: paper_summary: Interventional causal representation learning aims to recover the causal graph for latent variables and simultaneously recover the latent variables. While identifiability has been established for the infinite sample regime, this is not practical. This paper aims to provide PAC-style bounds on identifiability (of both the graph and latent variables) in the finite-sample regime, assuming stochastic soft interventions and linear transformations from the latent to the observation space. Three prototypical algorithms (causal order estimation, graph estimation, and inverse transformation estimation) are proposed to estimate the graph and the latent variables, which pave the foundation, from which  PAC-style bounds have been established. Numerical assessments are conducted to provide complementary insights into the sample complexity results.
strength_weakness: strength_weakness: The bounds are rather loose, which limits its practical applicability. Apart from that, these bounds do not provide sufficient guidance on designing better algorithms recovering latent variables and graphs. Having said so, I still welcome this type of theoretical papers.
questions: questions: I wonder if causal order estimation is necessary. Since the latent variables have yet to be estimated, their values and semantics are flexible. This implies that any given causal order could be acceptable, as these variables can be learned to fit the specified orders. Some edges would be pruned later, which can be done in the graph estimation stage.

Additionally, I wonder if there are other quantities beyond score differences and the correlation of score differences that could be used to recover the graph and the latent variables.

Furthermore, it would be beneficial if the analysis could assess which variables and which parts of the graph are more prone to error, and to what degree, providing feedback on the estimated/recovered variables and graphs.
limitations: limitations: n.a.

Comments: title: Response to rebuttal
comment: Thank you for the response. I will maintain my positive rating.
The response to the reviewer below addresses my concern about the setting of thresholds. It seems like there still is a need to estimate the true \eta^* and \gamma^* practically. I am satisfied with the other responses and have increased my score accordingly.
Thank you for the response. I will maintain my positive score.
Thanks for the thoughtful question. In practice, we can include routines for estimating $\eta^*$ and $\gamma^*$. These estimates do not have to be highly accurate and even rough estimates suffice to choose reliable thresholds $\eta$ and $\gamma$.  Specifically, based on estimates for $\eta^*$ and $\gamma^*$  we can choose  “safe” thresholds, e.g., one-fourth of the estimates, so that (i) with high probability we satisfy the requirement $\eta < \eta^*$, and (ii) we can avoid collecting excessive samples and compromising the sample complexity bounds. For practical purposes, we can construct such estimates for both $\eta^*$ and $\gamma^*$ as follows.

- **Estimate for $\eta^{*}$:** We note that in Algorithm 1, at step $t$, all matrices investigated must have  either rank $t$ or $t-1$. Therefore, the $t-1$-th eigenvalues of all of these matrices are non-zero, or equivalently we have  $\geq\eta^*$. This observation lets us to iteratively build and refine an upper bound on $\eta^*$ during Algorithm 1, which we can use as a surrogate of $\eta^*$ in the subsequent algorithms 2 and 3. For Algorithm 1 itself, at time $t$, we can pick $k\in{\cal V}_t$ with minimal $t$-th eigenvalue—this is a threshold-free way of estimating the causal order, but might not be amenable to tractable analysis. In order to keep the analysis more amenable to sample complexity analysis, we can use the current statement of Algorithm 1 with thresholded rank statements.

- **Estimate for $\gamma^{*}$:** We first note that $\gamma^*$ definition in eq.(23) depends entirely on the ground truth decoder matrix $\bf G$ and its pseudoinverse. Secondly, we note that Algorithm 3 can be used independently of Algorithm 2 to generate an estimate $\bf H$ of ${\bf G}^\dagger$ up to some estimation error and possible mixing among rows. We can compute the value $\min\_{i\in[n]}\\|{\bf H}\_{i,:}\\|_2\cdot\\|{\bf H}^\dagger\_{:,i}\\|\_2$, which would be equal to $\gamma^*$ had ${\bf H} = {\bf G}^\dagger$, and use it as our estimate of $\gamma^*$.
Thank you for your responses!

I understand from your response that one can select $\eta$ and $\gamma$ to be arbitrarily small for implementation. However, from my understanding of your theoretical claims, this should blow up the required number of samples in order to conclude "something useful", right? Can you give concrete suggestions on how practitioners can use your algorithm while drawing meaningful conclusions?
I'm satisfied with the response and keep the score as 'accept.'