Learning Mixtures of Unknown Causal Interventions

Authors: abhinav  kumar
kirankumar  shiragur
caroline  uhler
Venue: NeurIPS_cc_2024_Conference
Published: 2024-05-13 06:51:24

Abstract: The ability to conduct interventions plays a pivotal role in learning causal relationships among variables, thus facilitating applications across diverse scientific disciplines such as genomics, economics, and machine learning. However, in many instances within these applications, the process of generating interventional data is subject to noise: rather than data being sampled directly from the intended interventional distribution, interventions often yield data sampled from a blend of both intended and unintended interventional distributions.

We consider the fundamental challenge of disentangling mixed interventional and observational data within linear Structural Equation Models (SEMs) with Gaussian additive noise without the knowledge of the true causal graph. We demonstrate that conducting interventions, whether do or soft, yields distributions with sufficient diversity and properties conducive to efficiently recovering each component within the mixture. Furthermore, we establish that the sample complexity required to disentangle mixed data inversely correlates with the extent of change induced by an intervention in the equations governing the affected variable values. As a result, the causal graph can be identified up to its interventional Markov Equivalence Class, similar to scenarios where no noise influences the generation of interventional data. We further support our theoretical findings by conducting simulations wherein we perform causal discovery from such mixed data.

Main paper: 
 introduction: 
 Interventions are experiments that can help us understand the mechanisms governing complex systems and also modify these systems to achieve desired outcomes [1, 5, 21] . For example, in causal discovery, interventions are used to infer causal relationships between variables of interest, which has applications in various fields such as biology [5, 12] , economics [8] , and psychology [13, 21] .Given their extensive applications, there has been significant research in developing methods and experimental design strategies to conduct interventions under different scenarios [32, 24] . Despite significant efforts to develop sophisticated experimental techniques to perform interventions, they often encounter noise [4, 6, 30] . For example, the CRISPR technology, extensively used to perform gene perturbations (or interventions), is known to have off-target effects, meaning the interventions do not always occur on the intended genes [6, 30] . Consequently, in many applications, performing interventions generates data from a mixture of intended and unintended interventional distributions. Analyzing such mixed data directly can lead to incorrect conclusions, adversely affecting downstream applications. Hence, it is essential to disentangle the mixture and recover the components corresponding to each individual intervention for further use in downstream tasks like causal discovery.In our work, we formally address the challenge of disentangling mixtures of unknown interventional and observational distributions within the framework of linear structural equation models (Linear-SEM) with additive Gaussian noise. Given iid samples from a mixture with a fixed number of components as input, we present an efficient algorithm that can learn each individual component. Our results are applicable to both do and more general soft interventions.We chose to study our problem in the Linear-SEM with additive Gaussian noise framework for its fundamental importance in the causal discovery literature. Shimizu et al. [25] showed that observational data is sufficient for learning the underlying causal graph when the data-generating process is a Linear-SEM with additive non-Gaussian noise with no latent confounders. However, in the same setting with Gaussian noise, the causal graph is only identifiable up to its Markov Equivalence Class (MEC) [18, 25] . Thus, performing interventions (possibly noisy) is necessary to identify the causal graph, making it an interesting framework for our problem.Our contributions First, we show that given samples from a mixture of unknown interventions within the framework of Linear-SEM with additive Gaussian noise, there exists an efficient algorithm to uniquely recover the individual components of the mixture. The sample complexity of our procedure scales polynomially with the dimensionality of the problem and inversely polynomially with the accuracy parameter and the magnitude of changes induced by each intervention. Our findings indicate that the recovery error for each individual interventional distribution approaches zero as the number of samples increases. Therefore, in the infinite sample regime, we can recover the true interventional distributions, even when the targets of the interventions are unknown. Second, if the input distributions satisfy a strong interventional faithfulness assumption (as defined in Squires et al. [26] ), we can utilize the results from [26] to identify the targets of the interventions, thereby enabling causal discovery using these accurately recovered interventional distributions. Finally, we conduct a simulation study to validate our theoretical findings. We show that as sample size increases, one can recover the mixture parameters, identify the unknown intervention targets, and learn the underlying causal graph with high accuracy. 
prior work: 
 Mixture of DAGs and Interventions. There has been a lot of interest in understanding the mixture distribution arising from a collection of directed acyclic graphs (DAGs) [27, 23, 29, 14, 9, 28] . Saeed et al. [23] studied the distribution arising from a mixture of multiple DAGs with common topological order, thus a generalization of our problem. Their method identifies the variables whose conditional distribution across DAGs varies. However, there is no theoretical guarantee for the identifiability of the mixture's components. Like us, Kumar and Sinha [14] also studied the mixture arising from interventions on a causal graph and gives an algorithm to identify the mixture component. However, they assume knowledge of the correct topological order and sample access to the observational distribution. Thiesson et al. [28] also studied the problem of learning a mixture of DAGs using the Expectation-Maximization (EM) framework. However, there is no theoretical guarantee about the identifiability of individual components.Learning Mixture of Gaussians. Learning a mixture of Gaussians is a heavily studied problem [2, 7, 15, 11] . There exist efficient algorithms both in terms of runtime and sample complexity for a fixed number of components in the mixture [11, 15, 2] . Ge et al. [7] gave an efficient algorithm for the case when the number of components is almost √ n where n is the dimension of the variables. However, this method only guarantees identifiability in the perturbative setting, where the true parameters are randomly perturbed before the samples are generated. 
notation and preliminaries: 
 Notation. We use the upper-case letter X to denote a random variable and lower-case "x" to denote the value taken by the random variable X. Let, the uppercase bold-face letter X denote a set of random variables and the lowercase bold-face letter x denote the corresponding value taken by X. Let the conditional probability function P(X = x|Y = y) be denoted by P(x|y). We use the calligraphic letter S to denote a set and |S| to denote the cardinality of the set S. Let [n] denote the set of natural numbers {1, . . . , n}. For any vector v, we use the notation [v] j to denote it's j th entry and for any matrix M we use [M ] i,j to denote the entry in the i th row and j th column of M . We use R + to denote positive scalars.Structural Equation Model (SEM) Following Definition 7.1.1 in Pearl [18] , let a causal model (or SEM) be defined by a 3-tuple M = ⟨V , U , F⟩, where V = {V 1 , . . . , V n } denotes the set of observed (endogenous) variables and U = {U 1 , . . . , U n } denotes the set of unobserved (exogenous) variables that represent noise, anomalies or assumptions. Next, F denotes a set of n functions {f 1 , . . . , f n }, each describing the causal relationships between the random variables having the form:v i = f i (pa(V i ), u i ),where U i ⊆ U and pa(V i ) ⊆ V are such that the associated causal graph (defined next) is acyclic.A causal graph G M 1 is a directed acyclic graph (DAG), where the nodes are the variables V and the edges U with edges pointing from pa(V i ) to V i for all i ∈ [n].Linear-SEM (with causal sufficiency) In this work, we study a special class of such causal models (Gaussian Linear-SEMs) where the function class of each f i is restricted to be linear and of the formv i = f i (pa(V i ), u i ) = vj ∈pa(Vi) α ij v j + u i ,where α ij ̸ = 0, ∀V j ∈ pa(V i ). The causal Sufficiency assumption states that U i = {U i }, i.e., U i is the only exogenous variable that causally affects the endogenous variable V i . This is equivalent to the absence of any latent confounder (Chapter 9 in [20] ). In our work, we consider causally-sufficient Linear-SEMs; with a slight abuse of nomeclature, we will call them Linear-SEMs. The functional relationship between the exogenous and endogenous variables is deterministic, and the system's stochasticity comes from a probability distribution over the exogenous noise variables U . Thus, the probability distribution over the exogenous variable P(U ) defines a probability distribution over the endogenous variable P(V ). Without loss of generality, let the nodes {V 1 , . . . , V n } of the underlying causal graph be topologically ordered. Then, we can equivalently write the above set of equations as:EQUATIONwhere A, with A ij = α ij ,contains the causal effects between the endogenous variables. Thus, the matrix A, hereafter described as the adjacency matrix, characterizes the causal relationships between the endogenous variables (V ), where A ij ̸ = 0 denotes an edge between the variable V i and V j in G.We further specialize the exogenous variable u i (henceforth referred to as noise variable) to be Gaussian with mean µ i and variance σ i , i.e. u i ∼ N (µ i , σ i ). Thus, the joint distribution of the exogenous variables is given by a multivariate Gaussian distribution u ∼ N (µ, D) where [µ] i = µ i and the covariance is given by a diagonal matrix D with [D] ii = σ i . Thus, the endogenous variables also follow a multivariate Gaussian distribution with P(v) = N (m, S), where m ≜ Bµ i S ≜ BDB T and B ≜ (I -A) -1 . Causal discovery aims to identify the unknown adjacency matrix A given observational or other auxiliary data.Interventions. Following Definition 7.1.2 from Pearl [18] , the new causal model describing the interventional distribution, where the variables in a set I are set to a particular value, is given byM I = ⟨U , V, F I ⟩, where F I = {f j : V j / ∈ I} ∪ {f ′ i : V i ∈ I}and the functional relationship of every node V i ∈ I with their parents and corresponding exogenous variable U i is changed from f i to f ′ i . In particular, the functional relationship of node V i is changed toEQUATION1 We will drop the subscript M when it is clear from context where u i ′ ∼ N (µ i ′ , σ i ′ ). Such interventions are broadly referred to as "soft". Several other kinds of interventions are also defined in the literature, e.g., do, uncertain, soft etc. [3] . We consider three different types of widely studied specializations of soft interventions in our work:(1) shift: the mean of the noise distribution is shifted by a particular value, i.e., µ i ′ = µ i + κ for some κ ∈ R, and everything else remains the same, i.e., σ ′ i = σ i anda ij ′ = a ij , ∀j ∈ [n].(2) stochastic do (henceforth referred as stochastic): where all the incoming edges from parents are broken, i.e., α ij ′ = 0, andu i ′ ∼ N (µ i ′ , σ i ′ ).(3) do: in addition to breaking all incoming edges, i.e., α ij ′ = 0, we also set the variance of the noise distribution to 0 and the mean to any value of choice, i.e., u i ∼ N (µ i ′ , 0).In this work, we consider soft interventions where only one node is intervened at a time, i.e., |I| = 1. Thus after a soft intervention on node V i , the adjacency matrix is modified such thatA i ≜ A -e i (a i -a ′ i ) T = A -e i c T i 2 , where c T i ≜ (a i -a ′ i ) T , a T i is the i th row of matrix A and a ′ Ti is the new row after intervention such that [a i ] k = 0, ∀k ≥ i, and e i is the unit vector with entry 1 at the i th position and 0 otherwise. Thus, the linear SEM from Eq. 1 is:EQUATIONwhereu i ∼ N (µ i , D i ), µ i = µ + γ i e i for some γ i ∈ R, D i = D -δ i e i e T i where δ i ≜ (σ i -σ ′ i), is a diagonal matrix with the i th diagonal entry as σ ′ i and rest is same as D. Thus, the interventional distribution of the endogenous variables is also a multivariate Gaussian distribution, i.e.,P i (V ) = N (m i , S i ) where m i ≜ B i µ i , S i ≜ B i D i B T i and B i ≜ (I -A + e i c T i ) -1 .Question: For each theoretical result, does the paper provide the full set of assumptions and a complete (and correct) proof? Answer: [Yes]• It should be clear whether the error bar is the standard deviation or the standard error of the mean. • It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have a 96% CI, if the hypothesis of Normality of errors is not verified. • For asymmetric distributions, the authors should be careful not to show in tables or figures symmetric error bars that would yield results that are out of range (e.g. negative error rates). • If error bars are reported in tables or plots, The authors should explain in the text how they were calculated and reference the corresponding figures or tables in the text. 
problem formulation and main results: 
 In §4.1, we begin by formulating the problem of learning the mixture of interventions and then state our main result on the identifiability of parameters of the mixture. As a consequence of our identifiability result, under an interventional faithfulness assumption (Squires et al. [26] ), we show in §4.2 that the underlying true causal graph can be identified up to its I-MEC using a mixture of unknown interventions, thereby obtaining the same identifiability results as in the unmixed setting.We begin by formally defining the mixture of interventions over Linear-SEM with additive Gaussian noise and then state the main result of our paper -a mixture of interventions can be uniquely identified under a mild assumption discussed below.Definition 4.1 (Mixture of Soft Atomic Interventions). Let M = ⟨V , U , F⟩ be an unknown Gaussian Linear SEM where the distribution of the endogenous variables is given by P(V ) (see §3). Let I = {i 1 , . . . , i K }, hereafter referred to as intervention target set, be a set of unknown soft atomic interventions where each i k generates a new interventional distribution P i (V ). Then, the mixture of soft atomic intervention is defined as:EQUATIONwhere π i k ∈ R + , henceforth referred to as mixing weight, is a positive scalar such thati k ∈I π i k = 1, i k ∈ [n] ∪ {0}where n = |V | is the number of endogenous variables. We also allow i k = 0, which denotes the setting when none of the nodes is intervened, i.e., P 0 (V ) ≜ P(V ). Using Eq. 1 and 3, a mixture defined over a linear SEM with additive Gaussian noise is a mixture of Gaussians with parametersθ = {(m i k , S i k , π i k )} i k ∈I where P i k (V ) = N (m i k , S i k ).Having defined a mixture of interventions, we then aim to answer the following questions: (1) Does there exist an algorithm that can uniquely identify the parameters (θ) of the mixture of interventions under an infinite sample limit? (2) What is the run time and sample complexity of such an algorithm? It is immediate that if the intervention doesn't change the causal mechanism in any way, then the interventional distribution is equal to the observational distribution, and we would not be able to distinguish between them. This discussion suggests that it is necessary to put an additional constraint on the interventions performed. Below, we formally state the assumption that will ensure this and that it is sufficient for the identifiability of mixture distribution.Assumption 4.1 (Effective Intervention). Let P i (V ) = N (m i , S i ) be an interventional distribution after intervening on node v i , wherem i = B i µ i = B i (µ + γ k e i ), B i = (I -A + e i c T i ) and c i = (a i -a ′ i ), S k = B i D i B T i and D i = D -δ i e i e Ti (see atomic intervention paragraph in §3). Then, at least one of the following holds:γ i ̸ = 0 or ∥c i ∥ ̸ = 0 or δ i ̸ = 0.Now, we are ready to state the main result of our work that will help us answer the above questions.For an exact expression of sample complexity and runtime, see Lemma 5.1. i k ∈I ∥m i k -mρ(i k ) ∥ 2 + ∥S i k -Ŝρ(i k ) ∥ 2 + |π i k -πρ(i k ) | 2 ≤ ϵ 2for some permutation ρ : {1, 2, . . . , |I|} → {1, 2, . . . , |I|} and arbitrarily small ϵ > 0.Theorem 4.1 helps us separate the mixture of interventions P mix (V ) and provides us with the parameters {(m i , S i )} i∈I of the distribution of all the components in the mixture. However, it does not reveal which nodes were intervened, corresponding to the different components recovered from the mixture. There has been recent progress in performing causal discovery with a disentangled set of unknown interventional distributions [26, 16] . Specifically, Squires et al. [26] proposes an algorithm (UT-IGSP) that greedily searches over the space of permutations to determine the I-MEC and the unknown intervention target of each component. UT-IGSP is consistent, i.e., it will output the correct I-MEC as the sample size goes to infinity. Thus, combining Theorem 4.1 with the UT-IGSP algorithm implies that, as sample size goes to infinity, we can recover the underlying causal graph up to its I-MEC given a mixture of interventions over a Linear-SEM with additive Gaussian noise:Corollary 4.1.1 (Mixture-MEC). Given samples from a mixture of interventions P mix (V ) over a linear-SEM with additive Gaussian noise and samples from the observational distribution P(V ), there exists a consistent algorithm that will identify the I-MEC of the underlying causal graph under the I-faithfulness assumption (defined in Squires et al. [26] ) (and restated in §A.1).Proof. The proof follows from the identifiability of the parameters of the mixture distribution (Theorem 4.1) and the consistency of UT-IGSP given by Squires et al. [26] .Remark. The I-faithfulness assumption imposes certain restrictions on both observational and interventional distributions. However, as noted by Squires et al. [26] , in the case of Linear Gaussian distributions, the set of distributions excluded by this assumption is of measure zero. This is because the Linear Gaussian distributions, defined by a matrix A, that do not meet this assumption are subject to multiple polynomial constraints of the form p(A) = 0. It is a well-known result that for a random matrix A, the set of matrices that satisfy such polynomial equalities has measure zero [17] .5 Proof Sketch of Theorem 4.1Here we provide an overview of the proof of Theorem 4.1. Definition 4.1 tells us that the mixture of interventions defined over a Linear-SEM with additive Gaussian noise is a mixture of Gaussians. Learning mixtures of Gaussian is well-studied in the literature [2, 15, 7] . Since most of these approaches require some form of separability between the distribution of components or the parameters of the distributions, in the following lemma we will first show that the covariance matrix and the mean of any interventional or observational distribution taken pairwise is well-separated when Assumption 4.1 holds. In particular, this seperation of parameters will ensure that the Gaussian mixture can be uniquely identified using the results from Belkin and Sinha [2] .Lemma 5.1. [Parameter Separation] Let P 0 (V ) denote the observational distribution of a linear SEM with additive Gaussian noise "M" (see §3) with "n" endogenous variables. For some i, j ∈[n] ∪ {0}, let P i (V ) = N (m i , S i ) and P j (V ) = N (m j , S j ) be two interventional distributions (observational if one of "i" or "j" =0). Then the separation between covariance S i and S j and mean m i and m j is lower bounded by:∥S i -S j ∥ 2 F + ∥m i -m j ∥ 2 F ≥ f (B, D) ∥c i ∥ 2 + ∥c j ∥ 2 + h(B, D, µ) γ 2 i + γ 2 j + g(B) |δ i | min |δ i |, λ min (D) + |δ j | min |δ j |, λ min (D) ,where after intervention on node k ∈ {i, j}, ∥c k ∥ is the norm of the perturbation (or change) in the k th row of the adjacency matrix, γ k is the perturbation in the mean and |δ k | is the perturbation in the variance of the noise distribution of node k (as defined in the Atomic Intervention paragraph of §3). Also, B = (I -A) -1 , and D and µ are the covariance matrix and mean of the noise distribution in M, respectively. Furthermore, "f ", "g", and "h" are positive valued polynomial functions of B, µ and smallest eigenvalue of D (see the proof in §A.2 for the exact expressions).Remark. If one of the distributions is observational, i.e., say i = 0, then the above bound holds with ∥c i ∥ 2 = 0, γ i = 0, and|δ i | = 0.Remark. Different types of interventions will allow us to change certain parameters in the above bound. Let the exogenous noise variables u i ∼ N (µ i , σ i ) when not intervened. Now, if we intervene on node v i , then the new noise variable has distributionu ′ i ∼ N (µ i + γ i , σ ′ i )given as follows for different intervention types:(1) do intervention:∥c i ∥ ̸ = 0 if v i is not a root node, |δ i | ̸ = 0 if σ i ̸ = 0, and γ i ̸ = 0 if the value of node v i is set to any value other than µ i .(2) stochastic intervention: ∥c i ∥ ̸ = 0 if v i is not a root node. (3) shift intervention: γ i ̸ = 0. (4) soft intervention is the most general case, and nothing is guaranteed to be non-zero.Next, we restate a definition from Belkin and Sinha [2] that defines the radius of identifiability (R(θ)) of a probability distribution. If R(θ) > 0, this implies that we can uniquely identify the distribution. Definition 5.1. Let P θ , θ ∈ Θ, be a family of probability distributions. For each θ we define the radius of identifiabiity R(θ) as the supremum of the following set:{r > 0 | ∀θ 1 ̸ = θ 2 , (∥θ 1 -θ∥ < r, ∥θ 2 -θ∥ < r) =⇒ (P θ1 ̸ = P θ2 )}. In other words, R(θ) is the largest number, such that the open ball of radius R(θ) around θ intersected with Θ is an identifiable (sub) family of the probability distribution. If no such ball exists, R(θ) = 0.Next, we restate a result from [2] adapted to our setting, which shows that there exists an efficient algorithm for disentangling a mixture of Gaussians as long as the parameters are separated, which will ensure that the radius of identifiability R(θ) > 0.Theorem 5.2 (Theorem 3.1 in Belkin and Sinha [2] ). Let P mix (V ) be a mixture of Gaussians with parameters θ = {(m 1 , S 1 , π 1 ), . . . , (m |I| , S |I| , π |I| )} ∈ Θ where Θ is the set of parameters within a ball of radius Q. Then, there exists an algorithm which given ϵ > 0 and 0 < δ < 1 and poly n, max( 1 ϵ , 1 R(Θ) ), 1 δ , Q samples from P mix (V ), with probability greater than (1 -δ), outputs a parameter vector θ = θ = ( m1 , Ŝ1 , π1 ), . . . , ( m|I| , Ŝ|I| , π|I| ) ∈ Θ such that there exists a permutation ρ : {1, . . . , |I|} → {1, . . . , |I|} satisfying:i k ∈I ∥m i k -mρ(i k ) ∥ 2 + ∥S i k -Ŝρ(i k ) ∥ 2 + |π i k -πρ(i k ) | 2 ≤ ϵ 2 ,where the radius of identifiability R(θ) is lower bounded by:R(θ) 2 ≥ min 1 4 min i̸ =j (∥m i -m j ∥ 2 + ∥S i -S j ∥ 2 ), min i π i .Algorithm 1: Mixture-UTIGSP input :mixed dataset (D mix ), observational data (D obs ), number of nodes (n), cutoff ratio (τ ) output :Mixture Distribution Parameters (θ), Intervention Targets (I), causal graph ( Ĝ)1. Esstimate θ k ≜ (μ 1 , Ŝ1 ), . . . (μ k , Ŝk ) = GaussianMixtureModel D mix ,k for each possible number of component in the mixture i.e k ∈ [n + 1]. Define Θ ≜ {θ 1 , . . . , θ n+1 }be the set of estimated parameters and L = {l 1 , . . . , l n+1 } be the log-likelihood of the mixture data corresponding to the models with a different number of components. 2. To estimate the number of components in the mixture (k * ) iterate over k = (n + 1) to 2:(a) stop where the relative change in the likelihood increases above a cutoff ratio i.e|l k -l k-1 | l k > τ (b) k * = k if the stopping criteria is met otherwise k * = 1. 3. I, Ĝ = UT-IGSP D obs , θ k * [26]return Θ, I, Ĝ Our Theorem 4.1 along with Assumption 4.1 states that for every pair i, j∈ ([n] ∪ {0}) ⊗2 we have ∥S i -S j ∥ 2 F + ∥m i -m j ∥ 2 F > 0.Also, by construction of the mixture of interventions (in Definition 4.1), we have π i > 0, ∀i ∈ [n] ∪ {0}. This implies that the radius of convergence R(θ) > 0 and thus the parameters of the mixture of interventions P(V ) can be identified uniquely given samples from the mixture distribution with sample size inversely proportional to R(θ).6 Empirical Results 
experiment on simulated datasets: 
 Proposition 4.1.1 establishes that given samples from the mixture distribution, one can identify the underlying causal graph up to its I-MEC. To learn the causal graph, we first disentangle the mixture. Theorem 4.1 and 5.2 show that the sample complexity of our mixture disentangling algorithm is inversely proportional to various parameters of the underlying system and the intervention parameters (γ, |δ| and ∥c∥). Our simulation study further validates our theoretical results and characterizes the end-to-end performance of identifying the causal graph with such mixture data and its dependence on the above-mentioned parameters.Simulation Setup We consider data generated from a Linear-SEM with additive Gaussian noise, x = (I-A) -1 η (see §3), with n endogenous variables and corresponding exogenous (noise) variables. Here η ∼ N (0, D), where the noise covariance matrix is diagonal with entries D = diag(σ 2 , . . . , σ 2 ) and σ = 1 unless otherwise specified. A is the (lower-triangular) weighted adjacency matrix whose weights are sampled in the range [-1, -0.5] ∪ [0.5, 1] bounded away from 0. Let G * denote the causal graph corresponding to this linear SEM with edge i → j ⇔ A ji > 0. By sampling from the resulting multivariate Gaussian distribution, we obtain observational data. Next, for each causal graph, we generate separate interventional data by intervening on a given set of nodes in the graph one at a time (atomic interventions), which is again a Gaussian distribution but with different parameters (see Atomic Intervention paragraph in §3). We experiment with two settings:(1) all: where we perform an atomic intervention on all nodes in the graph.(2) half : where we perform an atomic intervention on a randomly selected half of the nodes. Then, the mixed data is generated by pooling all the individual atomic interventions and observational data with equal proportions into a single dataset. The decision to use equal proportions of samples from all components is solely intended to simplify the design choices for the experiment setup. In our experiment, we vary the total number of samples in the mixed dataset as N ∈ {2 10 , 2 11 , . . . , 2 17 }. In particular, we perform two kinds of atomic interventions: "do" and "stochastic" (see Interventions paragraph in §3 for a formal definition). The initial noise distribution for all the nodes is univariate Gaussian distribution N (0, 1). In our experiments for do interventions, instead of setting the final variance of noise distribution to 0, we set it to a very small value of 10 -9 for numerical stability. Unless otherwise specified, we perform 10 runs for each experimental setting and plot the 0.05 and 0.95 quantiles. See B for additional details on the experimental setup. 2 0 4 8  4 0 9 6  8 1 9 2  1 6 3 8 4 3 2 7 6 8 6 5 5 3 6 1 3 1 0 7 2 2 6 2 1 4 4 5 2 4 2 8 8 1 0 4 8 5 7 1 0 2 4  2 0 4 8  4 0 9 6  8 1 9 2  1 6 3 8 4 3 2 7 6 8 6 5 5 3 6 1 3 1 0 7 2 2 6 2 1 4 4 5 2 4 2 8 8 1 0 4 8 5 7 num intervention = all (↑) 2 4  2 0 4 8  4 0 9 6  8 1 9 2  1 6 3 8 4 3 2 7 6 8 6 5 5 3 6 1 3 1 0 7 2 2 6 2 1 4 4 5 2 4 2 8 8 1 0 4 8 5 7 1 0(f) num intervention = half (↓)Figure 1 : Performance of Alg. 1 as we vary sample size and number of nodes: The first row (a-c) shows the performance when the mixed data contains atomic intervention on all the nodes and observational data. The second row (d-f) shows the performance when the number of atomic interventions (chosen randomly) in mixed data is taken to be half of the number of nodes along with observational data. The column shows different evaluation metrics, i.e., Parameter Estimation Error, Average Jaccard Similarity, and SHD. The symbols (↑) represent higher is better, and (↑) represents the opposite (see Evaluation metric paragraph in §6). In summary, performance improves for both cases as the number of samples increases. However, the graph with more nodes requires a larger sample to perform similarly. For a detailed discussion, see §6.1.Given the mixed data generated from the underlying true causal graph G * , the goal is to estimate the underlying causal graph Ĝ. We break down the task into two steps. First, we disentangle the mixture data and identify the parameters of the individual interventional and/or observational distributions. Our theoretical result (Theorem 4.1) uses [2] for identifiability of the parameters of a mixture of Gaussians. Since they only show the existence of such an algorithm, we use the standard sklearn python package [19] that implements an EM algorithm to estimate the parameters of the mixture. Importantly, our experiment doesn't require prior knowledge about the number of components (k) in the mixture. We train separate mixture models varying the number of components. Then we select the optimal number of components using the log-likelihood curve of the fitted Gaussian mixture model using a simple thresholding heuristic (see step 2 in Alg. 1). We leave the exploration of better model selection criteria for future work. For all our simulation experiments, unless otherwise specified, we use a cutoff threshold of 0.07, chosen arbitrarily. In §B, we experiment with different values of this threshold and show that Mixture-UTIGSP is robust to this choice. The intervention targets present in the mixture are still unknown at this step. Next, we provide the estimated mixture parameters to an existing causal discovery algorithm with unknown intervention targets (UT-IGSP [26] ), which internally estimates the unknown intervention targets and outputs one of the possible graphs from the estimated I-MEC. We assume that observational data is given as an input to the UT-IGSP algorithm. The proposed algorithm is provided in Alg. 1. See B for our hyperparameter choice and other experimental details.Evaluation Metrics We evaluate the performance of Mixture-UTIGSP (Alg. 1) on three metrics:Parameter Estimation Error: This metric measures the least absolute error between the estimated parameters (mean and covariance matrix defining each individual distribution) after the first step of Mixture-UTIGSP matched with the ground truth parameters considering all possible matchings between the components averaged over all runs. See §B.5 for details.Average Jaccard Similarity (JS): We measure the average Jaccard similarity between the estimated intervention target and the corresponding ground truth (atomic) intervention target. We use the matching between the estimated and ground truth components found while calculating the parameter estimation error averaged over all runs. See §B.5 for details.Average Structural Hamming Distance (SHD): Given the estimated and ground truth graphs, we compute the SHD between the two graphs averaged over all runs.Results Fig. 1 shows the performance of Alg. 1 in the two settings, "all" in the first row and "half" in the second row (see Simulation setup above). The first column shows the performance of the first step of Alg. 1 where the mixture parameters are identified. We observe that parameter estimation error decreases as the number of samples increases in both settings. As expected, larger graphs require a larger sample size to perform similarly to smaller-sized graphs within each setting.Step 1 of our Alg. 1 only recovers the parameters ({(m i , S i )} k i=1 ) of the components present in the mixture distribution. In step 2 of our Alg. 1, we call UT-IGSP [26] that identifies the individual intervention targets from the estimated distribution parameters and also returns a causal graph from the estimated I-MEC. Fig. 1b and 1e show the average Jaccard Similarity between the ground truth and the estimated intervention targets. The colored lines denote experiments on graphs with different numbers of nodes. The corresponding dotted lines show the oracle performance of UT-IGSP when the separated ground truth mixture distributions were given as input. As expected, in both the settings (Fig. 1b and 1e ), the oracle version performs much better compared to its non-oracle counterpart for small sample sizes (2 10 to 2 14 ) but performs similarly as sample size increases.Finally, in the third column (Fig. 1c and 1f ), we calculate the SHD between the estimated causal graph Ĝ and the ground truth causal graph G * . The SHD of the graph estimated by Mixture-UTIGSP and the oracle version are similar for different node settings and all sample sizes. This suggests that small errors in the estimation of the parameters of the mixture distribution don't affect the estimation of the underlying causal graph.Additional Experiments: In §B, we provide details on the experimental setup and additional results. In Fig. 2 , we plot two additional metrics for the simulation experiments. The first metric is the number of estimated components in the mixture and the second metric is the the error in estimation of the mixing coefficient. In Fig. 3 , we study the sensitivity of the cutoff ratio used by Mixture-UTIGSP to select the number of components in the mixture. Next, in Fig. 4 , we evaluate the performance of our Alg. 1 as we vary the density, i.e., the expected number of edges in the graph. In Fig. 6 , we show how the sparsity of the graph and other intervention parameters like the value of the new mean and variance of the noise distribution after intervention affects the performance of Alg. 1.We evaluate our method on the Protein Signaling dataset [22] to demonstrate real-world applicability. The dataset is collected from flow cytometry measurement of 11 phosphorylated proteins and phospholipids and is widely used in causal discovery literature [31, 26] . The dataset consists of 5846 measurements with different experimental conditions and perturbations. Following Wang et al. [31] , we define the subset of the dataset as observational, where only the receptor enzymes were perturbed in the experiment. Next, we select other 5 subsets of the dataset where a signaling protein is also perturbed in addition to the receptor enzyme. The observational dataset consists of 1755 samples, and the 5 interventional datasets have 911, 723, 810, 799, and 848 samples, respectively. The mixed dataset is created by merging all the observational and interventional datasets.The total number of nodes in the underlying causal graph is 11. Thus, the maximum number of possible components in the mixture is 12 (11 single-node interventional distribution and one observational). In the mixture dataset described above, we have 6 components (1 observational and 5 interventional). The second column in Table 1 shows that Mixture-UTIGSP recovers 4 components close to the ground truth 6 when the cutoff ratio is 0.01 (step 2 of Alg. 1). Next, we give the disentangled dataset from the first step of our algorithm to identify the unknown target. Though the Jaccard similarity of the recovered target is not very high (average of 0.03 shown in the last row of the third column, where the maximum value is 1.0), it is similar to that of the oracle performance of UT-IGSP when the disentangled ground truth mixture distributions were given as input. This shows that it is difficult to identify the correct intervention targets even with correctly disentangled data. Also, the SHD between the recovered graph and the widely accepted ground truth graph for Mixture-UTIGSP (ours) and UT-IGSP (oracle) is very close. Overall, at a lower cutoff ratio, the performance of Mixture-UTIGSP is close to the Oracle UT-IGSP algorithm. Unlike the simulation case (see Fig. 3 ), Mixture-UTIGSP's performance is sensitive to the choice of the cutoff ratio on this dataset. In Fig. 5 , we plot the ground truth graph curated by the domain expert alongside the estimated causal graph for visualization. 
conclusion: 
 We studied the problem of learning the mixture distribution generated from observational and/or multiple unknown interventional distributions generated from the underlying causal graph. We show that the parameters of the mixture distribution can be uniquely identified under the mild assumption that ensures that any intervention changes the distribution of the observed variables. As a consequence of our identifiability result, under an interventional faithfulness assumption (Squires et al. [26] ), we show that the underlying true causal graph can be identified up to its I-MEC based on a mixture of unknown interventions, thereby obtaining the same identifiability results as in the unmixed setting. Finally, we conduct a simulation study to validate our findings empirically. We demonstrate that as the sample size increases, we obtain parameter estimates of the mixture distribution that are closer to the ground truth and, as a result, we eventually recover the correct underlying causal graph. 
limitations and future work: 
 Since our work is the first to study the problem of a mixture of causal interventions without assuming knowledge of causal graphs, we have restricted our attention to one particular family of causal models-Linear-SEM with additive Gaussian noise. In the future, it would be interesting to study this problem for a more general family of causal models. Further, our work uses Belkin and Sinha [2] for identifying the parameters of a mixture of Gaussians, which assumes that the number of components in the mixture is fixed. Recent progress in [7] gives an efficient algorithm for recovering the parameters when the number of components is almost √ n, where n is the number of variables. However, they only work in perturbative settings, and proving the result for non-perturbative settings is out of the scope of the current paper. Finally, to identify the parameters of the mixture distribution in our empirical study, we use heuristics to estimate the number of components. It would be interesting to explore other methods to automatically select the number of components.Question: For each experiment, does the paper provide sufficient information on the computer resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?Answer: [Yes]Justification: See §B.Guidelines:• The answer NA means that the paper does not include experiments.• The paper should indicate the type of compute workers CPU or GPU, internal cluster, or cloud provider, including relevant memory and storage. • The paper should provide the amount of compute required for each of the individual experimental runs as well as estimate the total compute. • The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary or failed experiments that didn't make it into the paper). 
code of ethics: 
 Question: Does the research conducted in the paper conform, in every respect, with the NeurIPS Code of Ethics https://neurips.cc/public/EthicsGuidelines? 
broader impacts: 
 Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?Answer: [NA] Justification: We theoretically study the problem in a simple setting that might not be directly applicable in real-world scenarios.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: We don't release any data or models.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: To the best of our knowledge, we have cited all the prior work and any package or library used in our empirical study.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. 
 
 Appendix: A.1 I-faithfulness in UT-IGSP algorithm.Here, we restate the assumption from Squires et al. [26] that is needed for consistency of the UT-IGSP algorithm:Assumption A.1 (I-faithfulness assumption). Let I be a list of intervention targets. The set of distributions {f obs } ∪ {f I } I I is I-faithful with respect to a DAG G if f obs is faithful with respect to G and for any I k ∈ I and disjoint A, C ⊆ [p], we have that(A ⊥ ⊥ ζ k |C ∪ ζ I\{k} ) G I if and only if f k (x A |x C ) = f obs (x A |x C ).A.2 Proof of Lemma 5.1 Lemma 5.1. [Parameter Separation] Let P 0 (V ) denote the observational distribution of a linear SEM with additive Gaussian noise "M" (see §3) with "n" endogenous variables. For some i, j ∈ [n] ∪ {0}, let P i (V ) = N (m i , S i ) and P j (V ) = N (m j , S j ) be two interventional distributions (observational if one of "i" or "j" =0). Then the separation between covariance S i and S j and mean m i and m j is lower bounded by:∥S i -S j ∥ 2 F + ∥m i -m j ∥ 2 F ≥ f (B, D) ∥c i ∥ 2 + ∥c j ∥ 2 + h(B, D, µ) γ 2 i + γ 2 j + g(B) |δ i | min |δ i |, λ min (D) + |δ j | min |δ j |, λ min (D) ,where after intervention on node k ∈ {i, j}, ∥c k ∥ is the norm of the perturbation (or change) in the k th row of the adjacency matrix, γ k is the perturbation in the mean and |δ k | is the perturbation in the variance of the noise distribution of node k (as defined in the Atomic Intervention paragraph of §3). Also, B = (I -A) -1 , and D and µ are the covariance matrix and mean of the noise distribution in M, respectively. Furthermore, "f ", "g", and "h" are positive valued polynomial functions of B, µ and smallest eigenvalue of D (see the proof in §A.2 for the exact expressions).Proof. First, we state a lemma that will give us a lower bound on the separation between the covariance matrix of two intervention distributions (or one could be observational). The proof is given in §A.3.Lemma A.1 (Minimum Covariance Separation). Let S i and S j be the covariance matrix of the Gaussian distribution corresponding to intervention on node V i and V j in the causal graph. Let, without loss of generality, V i be topologically greater than V j in the causal graph. Let B = (I -A) -1 , D = diag(σ 1 , . . . , σ n ) be the covariance matrix of the noise distribution, let S = BDB T be the covariance matrix of the observed distribution P (V ), and let c i = a i -a ′ i be the soft intervention performed on node V i (see §3 for definitions). Then we have:∥S i -S j ∥ 2 F ≥ f (B, D) ∥c i ∥ 2 + ∥c j ∥ 2 + g(B) |δ i | min |δ i |, λ min (D) + |δ j | min |δ j |, λ min (D) .If one of the covariance matrices is from the observational distribution, i.e., say S j = S then the above lower bounds still holds with ∥c j ∥ = 0 and δ j = 0.Using the above Lemma A.1 and substituting the explicit form of f (B, D) from Eq. 19 we obtain:EQUATIONNext, we state another lemma that will give us a lower bound on the separation of the mean of two interventional distributions (or one of them could be observational). The proof of the lemma below is given in §A.5.Let m i and m j denote the mean of the Gaussian distribution corresponding to the intervention on node V i and V j in the causal graph. Then we have:∥m i -m j ∥ 2 F ≥                γ 2 i 4∥B -1 ∥ 2 F + γ 2 j 4∥B -1 ∥ 2 F , ψ + i , ψ + j are active γ 2 i 4∥B -1 ∥ 2 F , ψ + i , ψ - j are active and γ 2 j 4∥Bµ∥ 2 ≤ ∥c j ∥ 2 γ 2 j 4∥B -1 ∥ 2 F , ψ - i , ψ + j are active andγ 2 i 4∥Bµ∥ 2 ≤ ∥c i ∥ 2 0,ψ - i , ψ - j are active andγ 2 i 4∥Bµ∥ 2 ≤ ∥c i ∥ 2 , γ 2 j 4∥Bµ∥ 2 ≤ ∥c j ∥ 2 where ψ + i ≜ |c T i Bµ| < |γi| 2 or |c T i Bµ| > 3|γi| 2 , ψ - i ≜ |γi| 2 ≤ |c T i Bµ| ≤ 3|γi|2 and similarly for ψ + j and ψ - j . If one of the means say m j is from the observational distribution, i.e., m j = m = Bµ, then setting γ j = 0 above will give the appropriate bounds and only case 2 and 4 are applicable. Now, From Case 1 of the above Lemma A.2 (Eq. 42, ψ + i , ψ + j is active) and above equation we have:EQUATIONFrom Case 2 of the above Lemma A.2 (Eq. 42, ψ + i , ψ - j is active) and using γ 2 j 4∥Bµ∥ 2 ≤ ∥c j ∥ 2 we have:EQUATIONFrom Case 3 of the above Lemma A.2 (Eq. 42, ψ - i , ψ + j is active) and usingEQUATIONSimilarly, for Case 4 and usingγ 2 i 4∥Bµ∥ 2 ≤ ∥c i ∥ 2 and γ 2 j 4∥Bµ∥ 2 ≤ ∥c j ∥ 2we have:EQUATIONNext, combining Eq. 6, 7, 8, 9 we have:EQUATIONFor the case when one of the distributions is observational, say w.l.o.g. S j = S and m j = m, then the same analysis holds since c j = 0 and δ j = 0 (from Lemma A.1) and only the analyses of case 2 and case 4 are applicable (from Lemma A.2), thereby completing the proof.Proof. Without loss of generality, throughout our analysis, we will assume that the endogenous variables V 1 , . . . , V n are topologically ordered based on the underlying causal graph. Thus, the corresponding adjacency matrix A is lower triangular. The value of ∥S i -S j ∥ 2 F will be unaffected by any permutation of the node order in the matrix A as shown next. Let Ã = P AP T be the adjacency matrix when the nodes are permuted by the permutation matrix P where P P T = I. Also, let the corresponding permuted covariance matrix be S = B D BT where B = (I -Ã) -1 = P (I -A) -1 P T = P BP T =⇒ S = P BP T P DP T P B T P T = P BDB T P T = P SP T . Similarly, we have Si = P S i P T and Sj = P S j P T . Now we have:EQUATIONsince the permutation matrix only permutes the row and column in the above equation, the Frobenius norm remains the same.First, we state a lemma that characterizes the covariance matrix S i of the interventional distribution.The proof of this lemma can be found in §A.4.Lemma A.3 (Covariance Matrix Update). Let P i (V ) = N (m i , S i ) be an interventional distribution and let the endogenous nodes V 1 , . . . V n be topologically ordered based on the underlying causal graph, then we have:S i = S -δ i r i r T i , for root node B i DB T i -δ i r i r T i , otherwise where r i = Be i , δ i = σ i -σ ′ i , B i = B -Be i c T i B, B = (I -A) -1, S is the covariance matrix of the observational distribution and D is the covariance matrix of the observational noise distribution (see Atomic Intervention paragraph in §3).S i = B i DB T i -δ i v i v T i .Substituting the value of B i from the above lemma again we get:EQUATIONwhereu i = B T c i , v i = Be i , η i = u T i Du i = c T i Sc i .Next, multiplying the LHS and RHS of the above equation with B -1 and B -T we get:EQUATIONThe absolute value of (i, k) th (k ∈ [n]) entry of this matrix is given by:EQUATIONwhere in the second case c T i Be i = 0 and D kk = σ k is the k th diagonal entry of the matrix D. Similarly, the (j, k) th entry of this matrix is given by:EQUATIONSimilarly, the covariance matrix of node V j and j ̸ = i is given by: S j = S -BDB T c j e T j B T -Be j c T j BDB T + η j Be j e T j B T -δ j Be j e T j B T B -1 (S j -S)B -T = η j e j e T j -DB T c j e T j -e j c T j BD -δ j e j e T j .(16) Thus the absolute value of the (i, k) th entry of the above matrix is given by:EQUATIONsince we had assumed without loss of generality that V j ≺ V i =⇒ e i B T c j = 0. Similarly, the (j, k) th entry of this matrix is given by:EQUATIONwhere in the second case e T j B T c j = 0. Thus, taking the difference of the i th and j th row of both the matrix (from Eq. 14, 15, 17 and 18 and using the fact that c T i Be i = 0, c T j Be j = 0 and c T j Be i = 0 since without loss of generality we have V j ≺ V i ) we obtain the following lower bound:EQUATION) after substituting η i = c T i Sc i ≥ λ min (S)∥c i ∥ 2 , η j = c T j Sc j ≥ λ min (S)∥c j ∥ 2 and using Lemma A.4 that gives us a lower bound for the second and third term. Thus, we obtain:EQUATION)where λ min (S) ≥ λ min (D)λ 2 min (B) = λ min (D) (∵ λ min (ST ) ≥ λ min (S)λ min (T ) and λ min (B) = λ min (B T ) = 1 since it is a lower triangular matrix where all diagonal entries are 1) and ∥ST ∥ F ≤ ∥S∥ F ∥T ∥ F . For the case when one of the covariance matrices is S corresponding to no intervention and the other is S i corresponding to intervention on node V i , then from Eq. 14 and the following similar analysis we have:EQUATIONthereby completing our proof of this lemma.Lemma A.4. Given λ ≥ 0 and η ≥ 0, the following inequality holds true:EQUATIONProof. Case 1:δ < 0 Let T ≜ λη + (η -δ) 2 , then we have T ≥ λη + δ 2 ≥ λη + δ 2 /4.Case 2: δ > 0 and [0 ≤ η < δ/2 or η > 3δ/2] Again in this case T ≥ λη + δ 2 /4.Case 3: δ > 0 and δ/2 ≤ η ≤ 3δ/2 In this case, we have T ≥ λη. Also, we have η > δ/2 =⇒ T ≥ λδ/2 which together implies:EQUATIONThus, combining the above three cases we have the following lower bound on the value of T :EQUATIONEQUATIONwhich completes the proof.Proof. From Eq. 3, we haveS i = B i D i B T i where B i = (I -A + e i c T i ) -1 and D i = D -(σ i - σ ′ i )e i e T i ≜ D -δ i e i e T i. Since e i c T i is a rank-1 update to the (I-A) matrix, using the Sherman-Morrison identity we obtain:EQUATIONwhere r i ≜ Bei di , scalar d i ≜ 1 + c T i Be i = 1 and q i ≜ B T c i . The scalar d i = 1 since q T i e i = c T i Be i = 0 given by the following lemma:Lemma A.5. Let A and B be a lower triangular matrix and c i be a vector such that c i = a i -a′ i and [c i ] t = 0, ∀t ≥ i (see Atomic Interventions paragraph of §3 for definition), then q T i e i = c T i Be i = 0. Also, if c i ̸ = 0, then q i = B T c i ̸ = 0.The proof of the above lemma can be found below after the proof of the current lemma. Next, S i is given by:EQUATION)If the intervened node V i is one of root nodes of the underlying unknown causal graph, then 25) , and thus:c i = 0 =⇒ q i = B T c i = 0 =⇒ B i = B (from Eq.EQUATION)where S is the covariance matrix of the observational distribution. Otherwise, for non-root nodes, term 1 in the above equation is:EQUATION)where w i = BDq i and scalar η i = q T i Dq i . Next simplifying term 2 in Eq. 26 we get:EQUATION)since q T i e i = c T i Be i = 0 from Lemma A.5 and r i ≜ Be i (from Eq. 25). Thus the covariance matrix of the endogenous variables in the intervened distribution is given by:EQUATIONthereby completing our proof.Proof of Lemma A.5. Showing q i ̸ = 0: We are given that B is a lower triangular matrix with nonzero diagonal entries. Also,c i = a i -a ′ i ̸ = 0 and [c i ] t = 0, ∀t ≥ i =⇒ ∃t s.t.[c i ] t ̸ = 0 and let t * be that last index where c i is non-zero. Now lets look at the[B T c i ] t * = c t * • [B T ] t * ,t * ̸ = 0 since B T is a upper triangular matrix and [c i ] t * ̸ = 0.i e i = 0: Be i is the i th column of the lower triangular matrixB =⇒ [Be i ] t = 0, ∀t < i =⇒ c T i Be i = 0 since [c i ] t = 0, ∀t ≥ i. A.5 Proof of Lemma A.2Proof. Similar to the proof of Lemma A.1, without loss of generality, we will assume that the endogenous variables are topologically ordered based on the underlying causal graph such that the adjacency matrix A is lower triangular. The permutation matrix will only permute the rows of the mean vectors m i and m j ; hence, there will be no effect on the value of ∥m i -m j ∥ 2 F . Now, let the mean of the Gaussian distribution corresponding to intervention on node V i be given by (see Atomic Intervention paragraph in §3 for definition):EQUATIONwhereB i = (I -A i ) -1 = B -Be i c T i B (from Lemma A.3), µ i is the new mean vector for the noise distribution, µ is the mean vector of the observational noise distribution and γ i e i is the update to the mean of the noise distribution when intervening on node V i . Similarly, the mean corresponding to the intervened distribution on node V j is given by:EQUATIONLooking at the i th entry of the vector B -1 m i , we get:EQUATIONand the i th entry of the vector B -1 m j is given by:EQUATIONSimilarly, the j th entry of the vector B -1 m i is given by:EQUATIONand the i th entry of the vector B -1 m j is given by:EQUATIONThus, the difference in the mean of both distributions can be lower bounded by:EQUATIONNow based on the value of |γ| i and |c T i Bµ|, we have:EQUATIONwhere ψ + i and ψ - i are used to denote two different mutually exclusive different events for ease of exposition later. In the case when ψ - i is true, we have:EQUATIONSimilarly, we have:EQUATIONand in the event when ψ - j is true, we have:EQUATIONCombining Eq. 38 and 40 and using ∥B -1 (m i -m j )∥ F ≤ ∥B -1 ∥ F ∥m i -m j ∥ F , we have:∥m i -m j ∥ 2 F ≥              γ 2 i 4∥B -1 ∥ 2 F + γ 2 j 4∥B -1 ∥ 2 F , ψ + i , ψ + j are active γ 2 i 4∥B -1 ∥ 2 F , ψ + i , ψ - j are active, Eq. 41 holds γ 2 j 4∥B -1 ∥ 2 F, ψ - i , ψ + j are active, Eq. 39 holds 0, ψ - i , ψ - j are active, Eq. 39 and 41 holds.(42)For the case when one of the mean say m j is observational then:EQUATIONThus the i th and j th entry of m j is e T i µ and e T j µ, respectively. Thus using Eq. 33 and 35 we get:EQUATION)Then again, following a similar analysis as in Eq. 38 and 39 we have the following cases:EQUATIONwhich is equivalent to say that only case 2 and 4 are applicable in Eq. 42. This completes the proof.B.1 Experimental Setup Discussion -Simulation Random Graph Generation: As mentioned in the Simulation setup in §6, we randomly generate the adjacency matrix A of the causal graphs used to simulate the mixture distribution. All the weights in the adjacency matrix are sampled in the range [-1, -0.5] ∪ [0.5, 1] bounded away from zero. This gives us a complete graph of all nodes. Thus, to introduce sparsity in the graph, we only keep an edge with probability ζ by setting the corresponding value in the adjacency matrix to zero if the edge is removed. Unless otherwise specified, we set ζ = 0.8 for all our experiments.Step 1 of Alg. 1: We use the GaussianMixture class from the scikit-learn Python package to disentangle the components of the mixture. For all our experiments, we use the default tol = 10 -3 used by GaussianMixture to decide on convergence of the underlying EM algorithm.Step 2 of Alg. 1: We run the UT-IGSP algorithm in Step 2 of our Alg. 1 with the standard parameter mentioned in the documentation. Specifically, we use α = 10 -3 for both MemoizedCITester and MemoizedInvarianceTester functions used by UT-IGSP.Specific to results in Appendix Unless otherwise specified, for all the results in the appendix we run all the experiments for 5 random settings. We ignore the error bars in the appendix for clarity in exposition. Also for all the experiments, the half-intervention setting i.e. where the mixture contains intervention on half of the randomly selected nodes is considered.Dataset: In the interventional data the signaling protein is also perturbed along with the receptor enzymes. The different perturbed signaling proteins (along with the number of samples corresponding) are: Akt (911), PKC (723), PIP2 (810), Mek (799), PIP3 (848). The observational data contained 1755 samples so overall 5846 samples were used for our experiments. For details see Wang et al. [31] and Sachs et al. [22] .The source code to all the experiments can be found in the following GitHub repository: https://github.com/BigBang0072/mixture_mecWe use an internal cluster of CPUs to run all our experiments. We run 10 random runs for each of the experimental configurations and report the mean (as points) and 5 th and 95 th quantiles as error bars for all our experiments in the main paper, and we report only mean the mean for the experiments in the appendix to declutter the figures. For the all setting, the actual number of components corresponding to the system with nodes 4,6 and 8 are 5,7,9 respectively (one intervention on each node + one observational distribution). We observe that Mixture-UTIGSP is able to correctly estimate the number of components even with a small number of samples. Similarly, for half setting, the actual number of components corresponding to the system with nodes 4,6 and 8 are 3,4,5 respectively (intervention on half of node and one observational distribution). Even for this case Mixture-UTIGSP is able to correctly estimate the number of components. The second column shows the error in the estimation of the mixing coefficient (π i 's, see Definition 4.1). For both cases, we observe that the error in the estimation of the mixing coefficient goes to zero as the sample size increases.Parameter Estimation Error Given the mixture distribution, the first step of our Alg. 1 estimates the parameters of every interventional and observational distribution present in the mixture (mixing weight π i , mean m i and covariance matrix S i ). Alg. 1 returns the maximum number of possible components, i.e. k = n + 1, by default. For all our experiments, we used this default value and left searching over the number of components for future work. To calculate the parameter estimation error, we first find the best match of the estimated parameters with the ground truth parameters based on the minimum error between the mean and covariance matrix. More specifically, let k * be the ground truth number of components in the mixture. Then, we iterate over all possible k * sized subsets of the estimated parameters and choose the one with the smallest absolute error sum between the mean and covariance matrix. Formally: Figure 3 : Performance of Alg. 1 as we change the cutoff ratio used for automatic component selection: We consider graphs with 6 nodes in this experiment with half intervention setting. In step 2 of Mixture-UTIGSP, we select the number of components using the log-likelihood curve. We scan the curve starting from the mixture model with the largest number of components to the smallest and stop where the relative change in the likelihood increases above a cutoff ratio (to select the elbow point of the curve). The cutoff ratio in the algorithm is chosen to be an arbitrary number close to zero.Parameter Estimation Error ≜ min k * i=1 |m i -mρ(i) |+|S i -Ŝρ(i) | : ρ ∈ perm([n+1]) (46)Here we compare the performance of Mixture-UTIGSP on all three metrics for the half setting of Fig. 1 as we vary the cutoff ratio. We observe that for the cutoff ratio close to zero i.e. 0.01, 0.15,0.3 the performance remains similar showing that the model selection criteria are robust to the selected cutoff ratio. The number of nodes where perm([n + 1]) represents all possible permutations of the indices [0, 1, . . . , n] and ρ * is the corresponding permutation with the minimum error.Average Jaccard Similarity (JS) In step 2 of our Alg. 1, UT-IGSP estimates the unknown interventional targets for each of the individual components disentangled in Step 1. We use the same matching (ρ * ) found in the parameter estimation error step (as mentioned above) to calculate the Jaccard similarity between the estimated and ground truth intervention target for that component. Formally:EQUATION)where t i is the ground truth intervention target set and tj is the estimated target set. For the case when both t i = ti = ϕ, then JS(t i , ti ) ≜ 1. The mixture data contains atomic interventions on all nodes as well as observational data (half setting as described in the results in §6). The column shows different evaluation metrics, i.e., Parameter Estimation Error, Average Jaccard Similarity, and SHD (see Evaluation metric paragraph in §6). In this experiment, we vary the density of the underlying causal graph by keeping the edges in a fully connected graph with a fixed probability, labeled as density in the legend of the above plots (see random graph generation paragraph in §B.1 for details). The maximum possible density is 1, i.e., the probability of keeping an edge is 1, corresponding to a fully connected graph, and the lowest possible density is 0. We observe that as the density of the graph increases, we require more samples to achieve similar performance to less dense graphs on all three metrics. Our Theorem 4.1 shows that the sample complexity required for estimating the parameters of the mixture is proportional to the norm of the adjacency matrix ∥A∥ and as the density of the graph increases ∥A∥ increases. Thus, as the density increases, we require more samples to achieve a similar performance in estimating the parameters of the mixture, as seen in the parameter estimation error plot above. [31] . 5b shows the graph estimated by our Mixture-UTIGSP and 5c is the graph estimated by oracle UT-IGSP when they are given the ground truth disentangled mixture. The blue colored arrow in 1b and 1c shows the correctly recovered edges in the domain expert graph. Green shows the edges with the same skeleton in the domain expert graph but in a reversed direction. The red shows the edges that are incorrectly added in the estimated graph. We observe that Mixture-UTIGSP correctly identifies two more edges (PKA->ERK and PKA-> Akt) as compared to an oracle which could be due to randomness in the UTIGSP algorithm. For this estimation, the best-performing cutoff of 0.01 was selected (see Table 1 ). (a) Varying the mean of the noise distribution (|γi|): The initial mean of the noise distribution is 0.0 for all the nodes. Upon intervening on a node to generate the interventional distribution, we change the mean of the noise to a different value. The variance of the noise distribution of the intervened node is kept the same as the initial distribution i.e. 1.0. From Theorem 4.1, we expect that as the new mean increases further from the initial mean = 0.0, the parameter estimation error should be lower for a given sample size and lead to better performance in intervention target estimation and causal discovery. As expected, the setting with the smallest change in the mean of the noise distribution (blue curve) has the worst performance. The case when the new noise mean is 10.0 (orange curve) is unusual where we see an unexpected increase in parameter estimation error. ) Varying the variance of the noise distribution (|δi|): The initial variance of the noise distribution is 1.0 for all nodes, and we change the variance of the new noise distribution upon intervention in this experiment. The mean of the noise distribution of the intervened distribution is kept the same as the initial distribution i.e. 0.0. From Theorem 4.1, we see that sample complexity to recover the parameters of the mixture distribution is inversely proportional to the change in the noise variance |δi|. Thus, we expect that the performance of Alg. 1 should improve as the new noise variance moves away from the initial noise variance 1.0. We can see in the above plot that the performance of the green curve (δi = 0) is worst in terms of the Jaccard similarity and SHD of the recovered graph, validating our expectation. Parameter estimation cannot be directly compared since as the variance increases, the norm of the covariance matrix increases, and thus, the overall error in the estimation error increases. We observe that compared to changing the mean (Fig. 6a ) increasing the variance gives slightly lower performance gains. Guidelines:• The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those are not discussed in the paper. • The authors are encouraged to create a separate "Limitations" section in their paper.• The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification, asymptotic approximations only holding locally). The authors should reflect on how these assumptions might be violated in practice and what the implications would be. • The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on a few datasets or with a few runs. In general, empirical results often depend on implicit assumptions, which should be articulated. • The authors should reflect on the factors that influence the performance of the approach.For example, a facial recognition algorithm may perform poorly when image resolution is low or images are taken in low lighting. Or a speech-to-text system might not be used reliably to provide closed captions for online lectures because it fails to handle technical jargon. • The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size. • If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and fairness. • While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, a worse outcome might be that reviewers discover limitations that aren't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play an important role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations.Justification: Yes, we do.Guidelines:• 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).38th Conference on Neural Information Processing Systems (NeurIPS 2024).We have used the subscript "i" notation in Ai to denote the adjacency matrix of the intervened distribution. We will use a similar subscript notation to represent other variables related to intervened distributions 


Decision: The reviewers are recommending acceptance after the rebuttal, although some had reservations about the contribution being on the minor end since the method heavily draws from disentangling a mixture of Gaussians. The additive noise assumption is also seen as a weakness of the analysis. However, I believe the sample-complexity results are much needed in the causal discovery literature and might justify the needed ANM assumption. Reviewer SSWb gave an 8 with a confidence of 5, but I do not think their review content justifies this high score/confidence pair. Looking at the other scores, this paper is a bit above borderline in its current form in my opinion.

Reviews: review: paper_summary: paper_summary: This paper considers a setup where an intervention results in obtaining iid data from a mixture of multiple interventional data and observational data and the goal is to learn the mixing weights and the resulting interventional distributions. Particular focus is on the linear SEM setting with Gaussian noise where interventions are allowed to be soft with some restrictions. The goal is to then learn the mixing weights and the interventional distribution (all multivariate Gaussians) parameters without knowledge of the causal graph. Given the large amount of literature on learning mixtures of Gaussians, it only remains to be shown that the individual components are well-separated in terms of changes in the interventional distributions.
strength_weakness: strength_weakness: While I understand that this is the first step in learning mixtures of interventions, I believe that it has limited novelty both conceptually and technically since it largely follows from existing work on learning Gaussian mixtures. The considered setup also assumes unconfoundedness on top of linear SEMs and Gaussian noise. Few more detailed questions follow in the next section.

The experimental evaluation is also limited to simulation data which don't completely endorse the validation of the algorithm. See next section on suggestions/questions. 

The writing can be more careful. There are typos in the main assumption (4.1) and there's some confusing notation about the high probability delta and the difference in variance delta_i.
questions: questions: 1. Theorem 4.1 is an existence theorem and in the experimental section, EM is used to learn the parameters. Given that the polynomial time guarantee includes information about A, how would we claim to not requiring knowledge of the causal graph?
2. While misspecifying k for the 'half' interventions case is indeed a feature, is there an improvement if k is correctly specified. Currently, it's not completely clear to me that the error goes down in (d). 
3. Is there any intuition for why when SHD and Jaccard similarity metrics improve, the parameter estimation error still does not? See nodes = 8 in half interventions case.
limitations: limitations: The authors have addressed limitations in a separate section.
paper_summary: paper_summary: The paper proposes linear structural equation models with additive Gaussian noise to address the challenge of disentangling mixed interventional and observational data. The problem is highly relevant to real-world applications with mixed data.
strength_weakness: strength_weakness: * The use of linear structural equation models is well-established in recent literature. Implementing SEMs with unknown interventions is not a novel approach.

* The identifiability guarantees rely on the assumption of soft interventions. This assumption may restrict the broader applicability of the model.

* The proposed method lacks sufficient experimental support. And it seems that the proposed method is primarily validated through simulations in the experiments, is it possible to provide real-world data examples?
questions: questions: See Weaknesses.
limitations: limitations: * The method focuses on linear structural equation models with additive Gaussian noise only, and the theoretical guarantees rely on the linear- SEM assumptions, which may limit the method's applicability in practice.
paper_summary: paper_summary: The paper proposes a method to disentangle the observational and interventional data under linear SEM with Gaussian noise.
strength_weakness: strength_weakness: The experiments are not comprehensive enough to thoroughly assess the effectiveness of the proposed method. More experiments on real-world datasets are required. This is especially necessary since the problem setup for the mixture of interventions is a bit restrictive. Also, it would be interesting to analyze how the method empirically performs on datasets that are not generated from Linear SEM or do not have Gaussian noises.
questions: questions: 1) Could you please elaborate on how realistic the problem setting is? Especially concerning the noisy intervention setting [6, 28], how compatible is the paper's problem setup with these real-world scenarios?

2) Could you please elaborate on the connection between your proposed method and existing literature on learning mixture of Gaussians?
limitations: limitations: Yes
paper_summary: paper_summary: This paper explores the challenge of disentangling and identifying causal relationships in situations where interventional data is noisy and mixed with both intended and unintended effects. It focuses on applying interventions within linear SEMs with Gaussian noise without prior knowledge of the true causal graph.  It presents an efficient algorithm that can learn and separate individual components of a mixture of unknown interventions, allowing for the recovery of true interventional distributions even when intervention targets are unknown.
strength_weakness: strength_weakness: The methods developed in the paper are specific to linear SEMs with Gaussian noise, which might limit their applicability in scenarios where these assumptions do not hold. Non-linear relationships or non-Gaussian noise structures, which are common in many real-world datasets, may not be adequately addressed by the proposed approach.
questions: questions: 1. I suggest improving readability by properly describing insights related to theorems, even though some of the proofs focus on matrix computations.
2. In the experimental setup, what would be the effect of different variances in Gaussian noise?
3.There might be some symbols used incorrectly:
In line 104, "the edge U" should be checked.
In line 128, the definition of M_I needs verification.
In line 134, the letter \(\kappa\) is used when mentioning shift interventions, but \(\gamma\) is used below. I know that do-interventions also change \(\gamma\), but I suggest keeping the notation consistent.
In line 279, "A_{ij} > 0" should possibly be "|A_{ij}| > 0".
limitations: limitations: The author has explained the limitations of the study. This paper has no possible negative social impact.

Comments: title: Final response
Response to rebuttal
Response to the Author
comment: We thank the reviewer for the question regarding the generalizability of automatic component selection due to manual thresholding. Manual thresholding was mainly chosen due to the short time constraint for the rebuttal. However, in the final revision, we will try to include other well-known methods for model selection, like the BIC criterion, that would generalize to other datasets.
Thanks for the example and the clarification on the automatic component selection. Regarding the latter, I am not sure if this is then generalizable beyond the examples considered? I have increased my score accordingly.
We thank the reviewer for their response and for acknowledging that the experimental evaluation is thorough now. Below, we answer the concerns and questions in detail:

>**Limited Novelty**

Here is an example where there is no separation between the parameters in the mixture.
Let the system consist of three nodes $X, Y$ and $Z$, the corresponding SEM is defined as:

$X= N(0,1)$
 
$Y=X$ and 

$Z=X-Y+N(0,1)$


Let the mixture distribution consist of two components: a) *observational* and b) *stochastic do* intervention on node $Z$. Here, we keep the noise distribution the same as before post-intervention, i.e., we set $Z = N(0,1)$. We can show that the mean and covariance of both components are the same. **Thus, there is zero parameter separation in spite of the adjacency matrix being different for both components**.

We agree that the above constructed example violates the faithfulness assumption. Thus, we needed a careful analysis to capture such nuances and other complexities in the parameter separation calculation. Our lower bound of parameter separation in Lemma 5.1 captures this i.e., parameter separation is zero for the above example since $f(B, D)=\lambda_{min}^{2}(D)/4||I-A||_{F}^{4} = 0$  since the noise covariance matrix $D$ has one of the entries as zero and the other parameters $\delta_i=\delta_j=\gamma_i=\gamma_j=0$ (see Appendix A2 for exact expression of $f(B, D)$). 

Also, it is somewhat intuitive that for non-degenerate cases, the parameters will be separated, but the exact dependence on the parameters of intervention was not clear. In Lemma 5.1, we show that the lower bound of the parameter separation is a polynomial function of the parameters of intervention ($||c_i||$, $\delta_i$ and $\gamma_i$’s). This also enables us to show the **polynomial sample complexity of the existence algorithm** in Theorem 5.2. 





\
\
>**Threshold for Automatic Component Selection:** 

Currently, we manually set the percentage threshold to a fixed value (7% in our simulation experiments, selected arbitrarily). As we increase the number of components in the mixture model the log-likelihood of the model is expected to increase due to increasing model capacity. However, we expect that as we reach the correct number of components any further increase in the number of components will most likely lead to minimal change in the log-likelihood. Thus, we chose to keep the threshold a small number close to zero.
I thank the authors for the detailed responses, which address some of my concerns, I will maintain my score.
As the deadline approaches, we would greatly appreciate hearing from you. Please let us know if you need any further clarifications.
Thanks for the response and apologies for my delayed response. 
Regarding limited novelty, I am not convinced yet mainly because I don't see non-trivial cases where separation wouldn't hold. Perhaps it would be more constructive to consider a non-trivial example of intervention in your model where separation does not hold? I understand this request comes late. 
Regarding Question 1 - Yes, I meant to write sample complexity. This is convincing since you're not really claiming polynomial sample complexity but just identifiability so the mismatch in the theorem statement and experimental algorithm is ok. 

Experimental Evaluation - 
The authors have responded with a thorough evaluation including experiments on a real-world dataset. A question regarding the automatic component selection algorithm - how are the cutoff percentages decided? 

I still find the novelty issue a minus but the experimental evaluation is now thorough for me to increase my score.
Thank you for your response. The authors have partially addressed my concerns, thus I will keep my rating.
Thanks for your response, it addresses my concerns. I will keep my rating.