FCause: Flow-based Causal Discovery
Authors: Tomas Geffner
Emre Kiciman
Angus Lamb
Martin Kukla
Miltiadis Allamanis
Cheng Zhang
Venue: ICLR 2022 Submitted (2022)

Abstract: Current causal discovery methods either fail to scale, model only limited forms of functional relationships, or cannot handle missing values. This limits their reliability and applicability. We propose FCause, a new flow-based causal discovery method that addresses these drawbacks.  Our method is scalable to both high dimensional as well as large volume of data, is able to model complex nonlinear relationships between variables, and can perform causal discovery under partially observed data. Furthermore, our formulation generalizes existing continuous optimization based causal discovery methods, providing a unified view of such models. We perform an extensive empirical evaluation, and show that FCause achieves state of the art results in several causal discovery benchmarks under different conditions reflecting real-world application needs.

Under review as a conference paper at ICLR 2022
FCAUSE: FLOW-BASED CAUSAL DISCOVERY
Anonymous authors
Paper under double-blind review
ABSTRACT
Current causal discovery methods either fail to scale, model only limited forms of functional relationships, cannot handle missing values. This limits their reliability and applicability. We propose FCause a new flow-based causal dis- covery method that addresses these drawbacks. Our method is scalable both high-dimensional as well as large volume of data, is able to model complex non- linear relationships between variables, and can perform causal discovery under partially observed data. Furthermore, our formulation generalizes existing con- tinuous optimization-base
causal discovery methods, providing a unified view of such models. We perform an extensive empirical evaluation, and show that FCause achieves state of the art results in several causal discovery benchmarks under different conditions reflecting real-world application needs.
INTRODUCTION
Understanding causal relationships between variables crucial in many applications, including bi-
ology (Koller& & Friedman 20099 Sachs et al. 2005),, economics (Varian 2016 Cunningham 2020). and healthcare (Tu et al. 2019). In addition, such information can also be used to advance other ma- chine learning domains, such as fairness (Kusner al., [017 Chiappa 2019). privacy (Tople al.]
2020 Chandrasckaran et al. (021) and robustness (Arjovsky et al. 2019 Zhang et al. 2020 Zheng
al. 2021). In practice, however, we often do not have a priori knowledge of causal relationships. We can gain this knowledge experimentally, e.g., through randomized control trials (Hariton & Lo- cascio 2018), but experiments are sometimes cost-prohile prohibiiieee or unethical. In such cases, causal discovery, the task of learning the causal relationships from existing data is essential. While causal discovery methods have been studied for decades, classical approaches tend have scaling issues and rely on domain-specifie assumptions.
Recent advances causal discovery use continuous optimization methods to overcome the scal- ability issues in classical approaches (Zheng et al.] 2018 2020). They frame the combinatorial optimization problem of learning a directed acyclic graph (DAG) as a smooth optimization problem under constraints, allowing the of efficient algorithms for continuous optimization. In addition to scalability advantages, this represents a first step in bridging the gap between deep learning and causal discovery. However, continuous pptimizattionbbsse causal discovery remains experimental. These methods cannot handle datasets with missing values and suffer from simple distributional assumptions, making performance sensitive changes in the dataset's scale (Reisach et al. 2021 Kaiser & Sipos 2021). These limitations directly impact these methods' reliability and applicability.
In contrast, when the causal graph-the structure but not functional forms of causal relationships-is already known, causal-aware (Bayesian) deep learning methods have many application advantages. One recent framework, Carefl (Khemakhem al., 202I), fits an autoregressive flow model using the variables' causal ordering for the autoregressive transformations. This yields model that can be used to answer complex causal queries, such as counterfactuals. Its performance for end-to- end causal discovery, however, is limited by its reliance classical causal discovery methods for inferring causal structures, the cost of fitting several flow-based models, and an inability to handle datasets with missing values.
In this work, we show how continuous optimization methods can be combined with flow-based models to address these deficiencies. We propose FCause, efficient flow-based causal discovery method able to capture complex nonlinear relations between variables, robust to data scaling, and able to handle datasets with missing values (missing at random). FCause uses powerful flow-based

---

Under review as a conference paper at ICLR 2022
models to capture complex relationships between variables (Khemakhem et al. 2021). performs efficient causal discovery using ideas from continuous optimization bassod causal discovery methods (Zheng et al. 2018). and handles partially partially-observed datasets using amortized approximate inference to estimate distributions over missing values (Kingma & Welling 2013).
In |ection we perform extensive empirical evaluations of FCause, and observe that it achieves state-of-the-art results on multiple causal discovery benchmarks, and that its performance robust with respect to changes in the scaling of the data (e.g., standardized data or not). In addition, FCauseremains competitive when used on datasets with 30% of missing values, even comparing against baselines used on the corresponding fully-observed dataset.
Finally, we present a unified formulation of several continuous optimization-bare causal discovery algorithms based on flow-based models (Khemakhem al.) 2021). This unified perspective facil- itates the development of general new techniques to improve these methods, and allows a simple comparison between methods, shedding light on the benefits and drawbacks of each one.
RELATED WORK
Approaches for causal discovery from observational data be classified into three broad groups: constraint-based, score-based and functional causal models (Glymour et al. (019). Constraint-based methods employ conditional independence tests discover the underlying causal structure (Spirtes) & Glymour 1991; Spirtes et al. 2000). Score-based methods find a causal structure by optimizing a score function (Chickering & Meek] 2015 Chickering, 2020p. Functional causal models represent each variable as a function of its direct causes and some noise term (Shimizu et al., 2006 Hoyer al., 2008). Classical causal discovery methods in these groups often struggle to scale to high dimensions, as the combinatorial space of possible causal orderings among variables grows super- exponentially (Peters et al., 2017). Thus, scalable causal discovery algorithms tend to reduce the search space through reliance on domain-specific assumptions.
Recently, Zheng et al. 2018) introduced Notears, a new family of score-based method based on continuous optimization. Notears uses novel algebraic characterization of directed acyclic graphs (DAG), which allows an gquality-constraingg rameters and adjacency relationships between variables. Notears inspired the development of other methods, Notears- MLP and Notears-Sob (Zheng et al. 2020), Grandag (Lachapelle et al.) 2019), and DAG-GNN (Yu et al. 2019), which extend the original formulation to model nonlinear rel relation- ships between variables. The methods' main benefits are its scalability and simplicity, a consequence of the fact that standard numerical solvers can be used to solve the resulting optimization problem. However, they have often been observed to be sensitive to different data scalings (Kaiser & Sipos 2021]. and cannot handle missing values.
optimization problem to jointly learn the model pa-
Additionally, normalizing flows (Rezende & Mohamed, 2015) have been used to build causal aware models (Cai et al., 2018 Khemakhem et al.] (021). These are based on the fact that the variables' causal ordering can be used for the transformations used by autoregressive flows (Kingma et al.) 2016 Huang et al. 2018) These methods are able to model complex nonlinear relationships be- tween variables. However, they rely on prior traditional methods or local search algorithms for causal discovery. Specifically, Self (Cai et al., 2018) uses a hill climbing procedure, and Carefl (Khemakhem et al. 202T) proposes use constraint-basee method (e.g. PC) to find the graphs skeleton with as many oriented edges as possible and to fit several flow models to determine the ori- entation of the remaining edges. Apart from aforementioned methods, causal aware deep learning models have, in general, shown better properties regarding generalization and robustness (Arjovsky
al. 2019 Zhang al. 2020 Kyono et al. 2020
Tople et al.) 2020)
PRELIMINARIES
In this section, we describe the representations of causal relationships (SEM and DAGs) and ex- plain how causal discovery can be formulated as an optimization task. We bring these components together when introduce FCause in section 4
While the formal connection to restricted type of autoregressive flows was proposed by (Khemakhem 202T,, some previous methods (Mooij al. 2011 Cai al. 2018) use closely related ideas.
2

---

Under review as a conference paper at ICLR 2022
Causal graphical models (Peters al. 2017) are widely used to represent causal relationships between variables. Each variable is assigned a node in a graph (typically a DAG), and causal re- lationships represented directed edges between nodes. These models are encoded using the graph's 's adjacency matrix A E (0, 1 1]dxd where d is the number of variables, or by explicitly stat- ing the set of parents for each variables. If each of the edges is additionally assigned weight, weighted adjacency matrix is used W E R instead. Here, zero entries indicate the absence of an edge, and non-zero entries indicate the presence of an edge with the corresponding weight.
Structural Equation Models (SEM) are commonly used to describe a causal system. They char- acterize the value of variable as function of its causal parents and some external noise. Let x = [z1..... Fd] represent d causally related variables, A E [0,1] the binary adjacency matrix between variables, and = [2,.... 2d] pairwise independent random noise. SEMs model each variable as p()) 2,), where fi a scalar function and T(p()) {az Az.3 = 1} represents the set of parents of 2; according to A.
Additive Noise Models (ANM) represent specific type of SEM for which the noise is additive:
= f((paa()) +
i == d....,d.
Equivalently, given a fixed adjacency matrix A, these d equations can be expressed jointly as
x = fa(x) +
where fa( outputs d-dimensional vector whose i-th component given by (()))) This is one of the most common functional forms in causal discovery. A crucial question regarding ANMs involves the identifiability of the underlying causal graph under observational data. Can the true adjacency matrix be recovered in the limit of infinite data? Hoyee et al. (2008) showed that it is possible except when the underlying system a combination of linear functions with Gaussian noise. In this work, assume that the underlying functions nonlinear without any specific assumption on their functional form. In this setting, the causal graph can be identified.
A new algebraic characterization of DAGs was recently proposed by Zheng et al. 2018). They showed that, given a weighted adjacency matrix W E Rddd the quantity
h(W)- = tr (e (eWow) -
(3)
is non-negative, and zero if and only if W represents a DAG. Since its introduction, this algebraic characterization of DAGs has been widely used to frame causal discovery problems with ANMs as continuous optimization tasks (Zheng et al.) 2018: 2020) [ Lachapelle et al., 2019]. All these methods propose to train a model's parameters & by maximizing score subject
the constraint
h(W (0)) = O, where W represents the weighted adjacency matrix as function of the parameters 0. Due to the non-convexity of the set {0: : h(W (0) 0}. these methods use h(W (0) penalty term added to the loss, whose weight is increased as optimization proceeds.
4 FCAUSE: FLOW-BASED CAUSAL DISCOVERY
We present FCause, a novel flow-based causal discovery method for nonlinear additive noise mod- els. Our method learns binary adjacency matrix jointly with the model parameters, able to cap- ture complex nonlinear relationships between variables, is stable to data re- re-scalings, and can handle datasets with missing values. We present the method for fully-observed datasets in section 4.]] and its extension to datasets with missing values in section 4.2 For the latter we assume that values missing completely at random or missing at random, meaning that the missingness pattern either has cause has fully observed causes, which is common setting (Rubin 1976 Stekhoven & Bühlmann 2012 Ma et al. 2018 Mattei & Frellsen 2019 Strobl et al. 2018).
4.1 FCAUSE WITH FULLY-OBSERVED DATA
FCause takes a Bayesian approach to causal discovery (Heckerman et al. 1999). We model the adjacency matrix jointly with the observations as
A) = ((A) II Pa (A))
(4)
3

---

Under review as a conference paper at ICLR 2022
We propose to fit the model parameters 8, and to use the posterior distribution over binary adjacency matrices pa (A|z² .."l",,,,..,)) to characterize the causal structure learned by the model.
There are two challenges to applying this approach: (i) The posterior distribution over A is in- tractable, and (ii) maximum likelihood cannot be used to fit the model parameters, due the the pres- ence of the latent variable A. We overcome these using variational inference (Jordan et al. 1999] Blei et al. 2017: Zhang et al. 2018). We define a variational distribution (s(A) to approximate the intractable posterior (A(A ...... x ). and use to build the evidence lower bound (ELBO), a lower bound on the marginal likelihood that can be used as surrogate objective. This yields
ELBO Φ) (g)
logp(A)
II [Po(=") 11 (ga) S log ,N),
(5)
where (9s)) represents the entropy of the distribution qa- (Details for the derivation in appendix
Modeling and optimization details Following Carefl, we set pa (x* to an autoregressive flow with base distribution pz and transformation p(A) to some prior over adjacency matrices that places zero mass on matrices that do not represent DAGs (needed to guarantee invertibility of the transformation, ef. appendix [E). We get FCause's final objective by replacing po(x" () in eq. (5) using the change of variable formula from random variables (Kingma et al.. 2016). This yields
9A(z)
(z) (ef. equations] | and 2 and
ELBO(B, 6) ( log log
II p= (x" fa(x":0)) det (x(x";0) 0)1 + H(o)),
(6)
where (a;;; is the Jacobian of the transformation used, (a 0) x fx( We set the base noise distribution p= to be factorized Gaussian with mean and learnable variances, the variational distribution (A) to be the product of independent Bernoulli distributions, one for each entry in the matrix A and use the Gumbel-softmas to get stochastic estimates of the ELBO's
gradient (Maddison al. 2016 Jang al. 2016)
Interestingly, the objective from eq. (6) may be simplified by noting that the bacobian-determinant term is always one for matrices A that represent DAG (see e.g. (Mooij et al., (011). We formalize this in lemma [ and include a proof in appendix D). Thus, since the prior over A only allows DAGs, this term does need be computed.
Lemma 1. Let A represent binary adjacency matrix, fA fA: Rdd + Rdd a function whose i-th output only depends on the parents of x;, and J (2) the Jacobian of ga(x) = x 1a--JA JA(I). the adjacency matrix A represents a DAG, then det JA (z)
The choice for fA : Rd - Rdd must satisfy the adjacency relations specified by A. That is, if AJ. 0, then the function f(( (i.e., the i-th component of the output of fa(z)) must satisfy (f((z)/0x,, O. Inspired by Graph Neural Networks (Hamilton) (020), we propose a flexible parameterization that satisfies this by setting
f((x) = ha Σ Ass (5
where zi and h; (i 1,.... d) are multi-layer perceptrons. A naive implementation of this requires training 2d neural networks. To avoid this, we propose to parameterize these functions as h₂(-) = h(,,,)) -) and 9,(+) = g(u,,)), where u, is an d-dimensional trainable embedding. This simple idea reduces the number of MLPs needed from 2d to just 2.
Finally, mentioned above, the prior p(A) must assign zero probability matrices that do not represent a DAG. While such a prior can be obtained in practice, it leads to a poorly conditioned
We the parameterization proposed by Lippe et al. (021)) which uses two parameters per each edge, one to model the existence and the other for the orientation.
'hhe use of the Gumbel-softmal to learn hinary adjacency matrix was also previously used by Ng et al.
2019a) and Brouillard al. 2020).
4

---

Under review as a conference paper at ICLR 2022
optimization problem. Therefore, following Zheng et al. (2018), "soft" alternative instead and gradually anneal it towards the "hard" DAG-forcing prior during training. This soft prior is
p(A) & exp (-eellvee(A)]] - - ((A)),
where h(-) is the non-DAG penalty from eq. (3) (h(A) 0 for non DAGs and exactly zero for DAGs). The first term favours sparsity in the adjacency matrix A, and the latter two favour DAGs. We anneal this soft prior by increasing the weights p and a optimization proceeds. We include details on the full optimization procedure used in appendix A
Note that the use of this soft prior means that only the latter stages of training does the soft prior approach the hard DAG-enforcing constraint. While this means that in the imitial training stages the assumptions of lemma do not hold, we find that in latter stages all sampled matrices A are DAGs, satisfying lemma I assumptions and providing strong empirical performance.
4.2 FCAUSE WITH PARTIALLY-OBSERVER DATA
We introduce FCause for parteelly-observed datasets. Missing values naturally present in many domains (e.g. online education (Wang et al. (02I)). Thus, developing methods to handle them effectively is crucial. Traditional deletion-based methods fail in many applications since, after dele- tion, very few data points are left. In contrast, FCause handles missing values efficiently, by using approximate inference estimate distribution over missing values.
We use r, to denote the observed components of sample x", x to denote the unobserved compo- nents, and po (x%, x"]) to denote the probability of the sample obtained by combining the values in x² and z/ As in the formulation without missing values, exact maximum likelihood is intractable. Thus, we propose to introduce a distribution approximate the true posterior over the latent vari- ables and to optimize the ELBO instead. However, in this case the latent variables include not only the adjacency matrix also all unobserved values in the dataset. We define the variational approximation
go. go. (A, 2 "., I (A) II 9u. (x
(9)
which leads to the objective
ELBO o, v/) Eggo. log A) II po (x aB + H((s)) + E((g.))
(10)
This is the objective maximized by FCause when the data has missing values. As explained before, compute p(=" [A) via the change of variable formula for random variables, and (A) to be the product of independent Bernoulli distributions. Finally, we set qv. (x:") [=") to be a factorized Gaussian, and reparameterization get unbiased gradients with respect to s//x-
For efficiency, we use an amortization network (Kingma & Welling 2013) to avoid training N independent Gaussian approximations. To do so, train single neural network with parameters vb. which receives as input the zero imputed sample x' concatenated with its "missingness" mask- which indicates indices of components missing for that sample-and outputs the mean and scale of q This allows us to train a single set of parameters v which are shared across samples in the dataset, instead of N set of independent parameters, {vi,.... 3/ }. one for each sample.
5 UNIFIED FLOW-BASED FORMULATION AND ANALYSIS
This section introduces a simple analysis showing that, similarly to FCause, most causal discovery methods based on continuous optimization can be framed from a probabilistic perspective as fitting flow. The benefits of this unified perspective are twofold. First, it allows simple comparison between methods, and sheds light on the different assumptions used by each one, their benefits and drawbacks. Second, it simplifies the development of new tools to improve these methods, since any improvements to one them can be casily mapped to the others by framing them this unified framework (e.g. our extension handle missing values can be easily integrated with Notears).
5

---

Under review as a conference paper at ICLR 2022
The connection between causal discovery methods based on continuous optimization and flow-based models uses the concept of a weighted adjacency matrix W (0) E IR linked to function f(z; 0) R + R.. Loosely speaking, these matrices can be seen as characterizing how likely is each output of f(x; 0) to depend on each component of the input x. For instance, W W = 0 indicates that (r is completely independent of z.j. Such adjacency matrices be constructed efficiently for a wide range of parameterizations for f. such as multi layer perceptrons and weighted combinations of nonlinear functions (Zheng et al. 2018 2020). We give precise details on their definition and construction in appendix
Lemma 2. Let f(x;; 6) Rd - Rdd be a 0-parameterized function with weighted adjacency matrix (0)) E Rdxd, Given a dataset {xt,..... ....,.,,,,,, fittin a flow with the transformation = x=f((;;) base distribution pz and a hard acyclicity constraint on W (0) is equivalent to solving
max
((x;; Loggpa - ((x;0)
s.L.
(() = -0,
(11)
where i(.) is the algebraic characterization of DAGs from eq.
Proof. The acyclicity constraint is enforeed by constraining the optimization domain to θ {0 h0) = 0}. Then, the maximum likelihood objective can be written as
Σ (og))("") Σ log pz px f(x"; + log det d(x" -xxx:
(12)
dx
T2
Σ log p. (x" f()
(lemma
(13)
where the first equality we the change of variable formula, valid because the transformation 2 = g(x; 0) = x - f(x; 0) is invertible for any 0 E e.
Lemma the main building block the formulation of continuous optimization- causal discovery methods from probabilistic perspective as fitting flow models. This is simply because the objective used by each of the methods can be exactly recovered from eq. (TT) with specific
choices for fa; ()) and pz-
Notears (Zheng et al. 2018) uses standard Gaussian for Ps and a linear transformation for [ar 0). (See appendix for more details.)
Notears-MLP Zheng et al. 2020) uses standard Gaussian for Ps and d independent multi-layer perceptrons, one for each component of fxx f).
Notears-Sob (Zheng al. 2020) uses standard Gaussian for pa and weighted linear combina- tion of nonlinear basis functions.
GAE (gg et al., 2019b) uses a standard Gaussian for ps and a GNN for (xx
Grandag (Lahaapelle et al., 2019) uses factorized Gaussian with mean zero and learnable scales for Ps and d independent multi layer perceptrons, one for each component of f(x, ().
Golem Ng et al. 2020). This is a linear method whose original formulation was already in a prob- abilistic perspective, using the linear transformation for [(;; 0).
In summary, recently proposed causal discovery methods based on continuous optimization can be formulated from probabilistic perspective as fitting flow with different constraints, transforma- tions, and base distributions. From these three components, the base distribution used by the flow sometimes overlooked but plays an important role. It has been observed, both theoretically and empirically, that methods that use standard Gaussian (e.g. Notears, Notears-MLPP Noteurs-Sob) may fail to recover the true underlying causal graph simple cases, and that their performance may change significantly when re-scaling variables in the dataset (c.g. standardizing data not) (Loh) & Bühlmann] 2014; Kaiser & Sipos 2021 Reisach et al.) (021). These issues can be overcome by setting the flow's 's base distribution to the actual noise distribution used to generate the data. While this distribution is typically unknown, methods such as FCause, Grandag and Golem take a step in this direction by optimizing some of its parameters, instead of fixing all of them arbitrary values.
Finally, as mentioned above, this unified formulation simplifies the development of new tools to improve these methods. For instance, in this work we introduce method to handle datasets with

---

Under review as a conference paper at ICLR 2022
missing values for FCause. Turns out that this idea is readily applicable any of the causal discov- ery methods described in this section. This can be seen by formulating them as probabilistic models po(x), and replacing the likelihood objective used for fully observed datasets with the ELBO from eq. (TO) (and possibly removing the random variable A from the objective).
EXPERIMENTS
This section presents results that empirically validate FCause for causal discovery. All simulations were performed using NVIDIA V100 GPUs.
Datasets. We consider synthetic, pseudo-real, and real data. For the synthetic data we follow the approach from Lachapelle al. al.19) and [Zheng et al. (2020). We sample a DAG following two different random graph models, Erdos-Renyi (ER) and scale-free (SF), For each graph we simulate 3; (i(pp)()) + 2i, where fi sampled from a Gaussian process and ,, (N of), with Ø; (/(02,1)) 1). We consider three different dimensionalities d E {16, 31, 64). and for each d we consider two possible number of edges E 4d]. Each resulting dataset is identified as ER e) or SF(d,e), e). All datasets have n = 1000 training samples. (In appendix F we show results in the exact same setting but with higher levels of noise, achieved by sampling
(().,,2).
For the pseudo-real dataset we consider data generated with the SynTReN generator (Van den Bul- cke et al. (006), which creates synthetic transcriptional regulatory networks and produces simu- lated gene expression data that approximates experimental data. We use the datasets generated by Lachapelle al. (2019) (dimension 20), and take n 400 for training. Finally, for the real dataset we use the widely used dataset that measures the level of several proteins in human cells from (aachs al. (005). We use training set with = 800 observational samples of dimension d 11.
Baselines. We compare against PC (Kliisch & Bühlman 2007p, (linear) Notears (Zheng et al. 2018), the nonlinear variants Notears-MLP and Notears-Sob (Zheng et al. 2020). Grandag [Lachapelle et al. 2019) (without the preliminary neighborhood search), GES (Chickering 2002). and ICALiNGAM (Shimizu al.] (0066). GES obtained from the GES package, and all other baselines were obtained from the gcastle package.
Causality metries used. Following common practice (Glymour al., 2019 Tu et 2019) we report adjacency and orientation metrics (recall, precision, and Fl score), and causal accuracy (Claassen Heskes 2012). For FCause, which returns distribution over adjacency matrices, we report the expected values of these metrics estimated using samples A ~ (A)).
Avg Rank (lower is better) Based on
Adjacency F1
Orientation F1
Causal Accuracy
FCause
1.4 + 0.7
1.3 + 1.0
1.4 + 1.3
Matears-MLPP
3.5 + 1.4
3.6
4.3 0.8
Notears-Sob
4.7
0.9
4.0
GES
2.6
3.4
2.1
PC
3.0 1.2
1.
3.0
Grandag
7.4 0.8
6.5
6.8
ICALiNGAM
6.4
6.4
6.8
Notears
6.9 0.8
7.2
0.6
Table I: FCause ranks highest among all baselines. Best results in bold. The table shows the methods average rank across datasets for different metrics (lower is better, is the best possible, and 8 is the worst).
In all cases standardize data and repeat all simulations for four different random seeds. (For the synthetic and pseudo-real datasets this includes generating four datasets, one for each seed.) Table presents high level summary of the results, by ranking the methods according to different metries. can be observed that FCause ranks highest, followed by Notears-ME MLLP. Notears- Sob, GES and PC (all four perform similarly), and finally followed by Grandag, ICALiNGAM and Notears.

---

Under review as a conference paper at ICLR 2022
FCause
Notears- Noteare-Sob
Grandag Notears
ICALINGAM
FCause (PO)
PC
0.8
Orientation F1 Adjacency
0.6
0.4
0.2
0.8
0.6
0.4
0.2
Causal Accuracy
0.8
0.4
0.2
0.0
ER(16,16)
ER(16.64) 64)
(R(32,32)
ER 128)
ER(64.64)
ER(64.256) 256)
SF(16, 16)
SF(16.64) 64)
5F(32,32)
SF 128)
SF(64.64) 64)
SF 256)
Protein
SynTReN
Figure 1: FCause achieves better results than the baselines in all three metries shown. The legend "FCause (PO)" corresponds to running FCause with 30% of the training data missing com- pletely at random. For readability, we highlight the FCause and FCause (PO) results by connecting them with soft lines. The figure shows mean results across four different random seeds. Results including standard deviations are included in appendix G
Fig. presents more detailed results. It shows the adjacency Fl score, orientation F1 score, and causal accuracy achieved by each method in each of the datasets considered. (We report full results, including standard deviations and all remaining metrics precision and recall metrics for adjacency and orientation in appendix G Again, it can be observed that FCause tends perform best across most datasets, followed by Notears-MLP, Notears-Sol GES and PC.
The figure also shows the results achieved by FCause of the values missing. This is identified as FCause (PO). As expected, the method performance gets slightly worse in the presence of missing values. However, it can be observed that it is still very competitive, and that on many datasets its performance is better than that of the baselines, even when these are ran using the dataset with no missing values. This shows the effectiveness of FCause's way of dealing with missing values.
the exact same datasets with 30%
Results using datasets generated with higher level of noise in are shown in fig. 2 (appendix F. As expected, as a consequence of the higher level of noise, all methods perform slightly worse However, the main conclusions remain unchanged.
Stability of FCause with respect to data scaling Since FCause learns the scales of the base distribution, we expect its performance to be robust to different re-scaling of the variables in the dataset (cf. section 5). This is in contrast to other methods, such Notears-MLP and Notears-Sob, who fix the base distribution to a standard Gaussian, and whose performance has been observed to change when re-scaling variables (Kaiser& & Sipos 2021 Reisach al. 2021). We consider simple experiment to verify this empirically. We ran FCause, Notears-MLP and Notears- on all datasets with and without data standardization, and compare the results obtained in both cases.
Table 2 shows the performance achieved by each method when ran with and without data standard- ization. It can be observed that FCause's mean performance is essentially the same in both cases.
8

---

Under review as a conference paper at ICLR 2022
Standard data
Adjacency
Orientation
Causal Accuracy
Recall
Prec
F1
Recall
Prec
F1
True
0.57
0.77
0.61
0.49
0.65
0.53
0.50
FCause
False
0.58
0.76
0.62
0.50
0.68
0.55
0.51
True False
0.37
0.85
0.49
0.27
0.62
0.35
0.27
Notears-MLP
0.41
0.80
0.45
0.31
0.69
0.37
0.31
True 0.34
0.85
0.44
0.27
0.71
0.37
0.27
Notears-Sob
False
0.49
0.51
0.42
0.38
0.45
0.35
0.38
Table 2: FCause performance is robust with respect to data standardization. The table shows the mean performance across datasets achieved by each method with and without data standardization.
Metric used to compute relative difference
Adjacency
Orientation
Causal Accuracy
Recall
Prec
F1
Recall
Prec
F1
FCause 0.04
0.10
0.08
0.10
0.13
0.13
0.10
Notears-MLP 0.25
0.12
0.11
0.29
0.23
0.18
0.30
Notears-Sob 0.32
0.51
0.14
0.33
0.51
0.19
0.33
Table 3: The relative difference between the performance of FCause with and without stan- dardized data lowest across all metrics. Best results in bold. The table shows, for each metric, the relative difference between the performances achieved by each method with and without data standardization (lower is better). The table reports mean values across datasets.
The cannol be said for Notears-MLP and Noteary-Sob, whose performances vary (sometimes significantly) depending on whether the data is standardized or not.
To quantify this variation, table shows the frelative difference between the performances achieved by each method with and without standardizing the data. For each method, metric and dataset, the relative difference between the performance obtained with and without data standardization is computed as RD b) 2/a (a + b). The table reports mean values across datasets. It can be observed that FCause achieves the lowest relative difference across all metrics by a large margin, followed by Notears-MLP. and finally by Notears-Sob.
Finally, we note that FCause can used for missing value imputation, using the mean of go (Ira (xa) to impute missing values. We show imputation results in appendix H where it can be observed that FCause performs similarly to mice (Van Buuren & Shoothuis-Oudshoo (Stckhoven & Bühlmann 2012), and outperforms PVAE (Ma et 2018).
2011 and missforest
DISCUSSION AND FUTURE WORK
We proposed FCause, a scalable causal discovery method able to model complex nonlinear rela- tionships between variables that can handle datasets with missing values. Our results indicate that FCause could be impactful in real-world applications, where the functional form for the underlying model is not known, data normalization is out of the practitioner's control, and there may be missing values in the datasets. FCause has shown to yield good results in all these scenarios.
As future work, we identify two possible paths for further improvement. One involves the use of more powerful density approximators to model the noise variables. While we use Gaussian with learnable scales, other options, such as one dimensional flows, could be tried. The other improve- ment involves extending the method to models more general than ANMs, such as the post nonlinear (Zhang & Hyvarinen 2012) or affine [Khemakhem et al. 2021). These extensions should be pos- sible, since all transformations involves in these kinds of models are invertible. Interestingly, once developed, these extensions could be readily applied to other causal discovery methods based on continuous optimization, thanks to the unifying framework presented in section 5

---

Under review as a conference paper at ICLR 2022
REFERENCES
Martin Arjovsky, Léon Bottou, Ishaan Guirajani, and David Lopez-Paz. Invariant risk minimization. arXiv preprint arXiv:/907 000999 2019.
David Blei, Alp Kucukelbir, and Jon D McAuliffe. Variational inference: review for statisti- cians. Journal of the American statistical Association, 1126518)859-8877 2017.
Philippe Brouillard, Séhastien Lachapelle, Alexandre Lacoste, Simon Lacoste-Julien, and Alexan- dre Drouin. Differentiable causal discovery from interventional data. arXiv preprint arXiv:2007.0175: 015544 2020.
Ruichu Cai, Jie Qiao, Zhenjie Zhang, and Zhifeng Hao. Self: structural equational likelihood frame- work for causal discovery. In Thirty-second AAAI conference on artificial intelligence, 2018.
Varun Chandrasekaran, Darren Edge, Somesh Jha, Amit Sharma, Cheng Zhang, and Shruti Tople. Causally constrained data synthesis for private data release. arXiv preprint arXiv:2105.
2021.
Silvia Chiappa. Path-specific counterfactual faireess. In Proceedings of the AAAI Conference on Artificial Intelligence, volume 33, pp. 7801-7808, 2019.
David Maxwell Chickering. Optimal structure identification with greedy search. Journal of machine learning research, 5077545 2002.
David Maxwell Chickering and Christopher Meek. Selective greedy equivalence search: Find- ing optimal bayesian networks using polynomial number of score evaluations. arXiv preprint arXiv:1506.022 113, 2015.
Max Chickering.
Statistically efficient greedy equivalence search. In Conference on Uncertainty in Artificial Intelligence, pp. 241 249. PMLR, 2020.
Tom Claassen and Tom Heskes. A bayesian approach to constraint based causal inference. arXiv preprint arXiv:/2/0.4800 2012.
Scott Cunningham. Causal inference. The Mixtape, 1, 2020.
Clark Glymour, Kun Zhang, and Peter Spirtes. Review of causal discovery methods based on graph- ical models. Frontiers in genetics, 10:524, 2019.
William Hamilton. Graph representation learning. Synthesis Lectures on Artifical Intelligence and Machine Learning, 14(3):1-159, 2020. 2020.
Eduardo Hariton and Joseph Locascio. Randomised controlled the gold standard for effec-
tiveness research. BJOG: an international journal of obstetrics and gynaecology, 125(13):1716,
2018.
David Heckerman, Christopher Meek, and Gregory Cooper. A bayesian approach to causal discov- ery. Compulation, causation, and discovery, 19: 19:141-166, 1999.
Patrik O Hoyer, Dominik Janzing, Joris Mooij, Jonas Peters, Bernhard Schölkopf, et al. Nonlinear causal discovery with additive noise models. In NIPS, volume 21, pp. 689-696. Citeseer, 2008.
Chin-Wei Huang, David Krueger, Alexandre Lacoste, and Aaron Courville. Neural autoregressive flows. In International Conference on Machine Learning, pp. 2078-2087 PMLR, 2018.
Eric Jang, Shixiang Gu, and Ben Poole. Categorical reparametterization with gumbel-softmax arXiv preprint xrXiv:1611.011 2016.
Michael Jordan, Zoubin Ghahramani, Tommi S Jaakkola, and Lawrence K Saul. An introduction to variational methods for graphical models. Machine learning, 37(2):183-233 1999.
Marcus Kaiser and Maksim Sipos. Unsuitability of notears for causal graph discovery. arXiv preprint arXiv:2/04 05441, 2021.
10

---

Under review as a conference paper at ICLR 2022
Markus Kalisch and Peter Bühlman. Estimating nigh-dimensional directed acyclic graphs with the pe-algorithm. Journal of Machine Learning Research, 8(3), 2007.
Ilyes Khemakhem, Ricardo Monti, Robert Leech, and Aapo Hyvarinen. Causal autoregressive flows. In International Conference Artificial Intelligence and Statistics, pp. 3520-3528. PMLR, 2021.
Diederik Kingma and Jimmy Ba. Adam: A method for stochastic optimization. arXiv preprint
arXiv:1412.6988 2014.
Diederik P Kingma and Max Welling 1312.66 2013.
Auto-encoding variational bayes. arXiv preprint
Durk P Kingma, Tim Salimans, Rafal Jozefowicz, Xi Chen, Ilya Sutskever, and Max Welling. Im- proved variational inference with inverse autoregressive flow. Advances in neural information processing systems, 29:4743-4755 4751, 2016.
Daphne Koller and Nir Friedman. Probabilistic graphical models: principles and techniques. MIT press, 2009.
Matt Kusner, Joshua Loftus, Chris Russell, and Ricardo Silva. Counterfactual faimess. In Pro- ceedings of the 31st International Conference on Neural Information Processing Systems, pp. 4069-4079, 2017.
Trent Kyono. Yao Zhang, and Mihacla van der Schaar. Castle: regularization via auxiliary causal graph discovery. arXiv preprint 1xXiv:2009.13/800 2020.
Sébastien Lachapelle, Philippe Brouillard, Tristan Deleu, and Simon Lacoste-Julien. Gradient-based neural dag learning. arXiv preprint (:Xiv:9966002220 2019.
Phillip Lippe. Taco Cohen, and Efstratios Gavves. Efficient neural causal discovery without acyclic- ity constraints. arXiv preprint arXiv:2107.10488, 2021.
Po-Ling Loh and Peter Bühlmann. High-dimensional learning of linear causal networks via inverse covariance estimation. The Journal of Machine Learning Research, 15(1):3065-33 -30000 2014.
Chao Ma. Sebastian Tschiatschek, Konstantina Palla, José Miguel Hernández-Looato Sebastian Nowozin, and Cheng Zhang. Eddi: Efficient dynamic discovery of high-value information with partial vae. arXiv preprint arXiv: 11142, 2018.
Chris J Maddison, Andriy Mnih, and Yee Whye Teh. The concrete distribution: A continuous relaxation of discrete random variables. arXiv preprint arXiv:1611 1.007/2. 2016.
Alexad-Ae Mattei and Jes Frellsen. Miwae: Deep generative modelling and imputation of incomplete data sets. In International Conference on Machine Learning, pp. 4413-4423, PMLR, 2019.
Joris M Mooij. Dominik Janzing, Tom Heskes, and Bernhard Schölkopf. On causal discovery with cyclic additive noise model. 2011.
AS Nemirovsky. Optimization ii numerical methods for nonlinear continuous optimization. 1999.
Ignavier Ng, Zhuangyan Fang, Shengyu Zhu, Zhitang Chen, and Jun Wang. Masked gradient-based causal structure learning. arXiv preprint arXiv:1910.08527, 2019a.
Ignavier Ng, Shengyu Zhu, Zhitang Chen, and Zhuangyan Fang. A graph autoencoder approach to causal structure learning. arXiv preprint 19Xiv:1911.07422) 2019b.
Ignavier Ng, AmirEmad Ghassami, and Kun Zhang. On the role of sparsity and dag constraints for learning linear dags. arXiv preprint 20Xiv:2006.100 2020.
Jonas Peters, Dominik Janzing, and Bernhard Schölkopf. Elements of causal inference: foundations and learning algorithms. The MIT Press, 2017.
Alexander G Reisach. Christof Seiler, and Sebastian Weichwald. Beware of the simulated dag! varsortability in additive noise models. arXiv preprint arXiv:2102.1 13647, 2021.
11

---

Under review as a conference paper at ICLR 2022
Danilo Rezende and Shakir Mohamed. Variational inference with normalizing flows. In Interna- tional conference on machine learning, pp. 1530-1538. PMLR, 2015.
Donald B Rubin. Inference and missing data. Biometrika, 63(3):581-592, 1976.
Karen Sachs, Omar Perez, Dana Pe'er, Douglas A Lauffenburger, and Garry P Nolan. Causal protein-signaling networks derived from multiparameter single-cell data. Science, 308(5721): 523-529, 2005.
Shohei Shimizu, Patrik O Hoyer, Aapo Hyvärinen, Antti Kerminen, and Michael Jordan. A linear non-gaussian acyclic model for causal discovery. Journal of Machine Learning Research, 7(10), 2006.
Peter Spirtes and Clark Glymour. An algorithm for fast recovery of sparse causal graphs. Social science computer review, 9(1):62-72, 1991.
Peter Spirtes, Clark N Glymour, Richard Scheines, and David Heckerman. Causation, prediction, and search. MIT press, 2000.
Daniel Stekhoven and Peter Bühlmann. Missforest Miisffrestt-no-teee -----ooo--aaaa-ee non-parametric missing value imputation for mixed-type data. Bioinformatics, 28(1): 12-118, 2012.
Eric V Strobl, Shyam Visweswaran, and Peter Spirtes. Fast causal inference with non-random missingness by test-wise deletion. International journal of data science and analytics, 6(1):47- 62, 2018.
Shruti Tople, Amit Sharma, and Aditya Nori. Alleviating privacy attacks via causal learning. In International Conference Machine Learning, pp. 9537-9547. PMLR, 2020.
Ruibo Tu, Cheng Zhang, Paul Ackermann, Karthika Mohan, Hedvig Kjellström, and Kun Zhang. Causal discovery in the presence of missing data. In The 22nd International Conference on Artificial Intelligence and Statistics, pp. 1762-1770. PMLR, 2019.
Stef Van Buuren and Karin Groathuis-Oudssoomm mice: Multivariate imputation by chained equa- tions in Journal of statistical software, 45(1):1-67, 2011.
Tim Van den Bulcke, Koenraad Van Leemput, Bart Naudts, Piet Remortel, Hongwu Ma, Alain Verschoren, Bart De Moor, and Kathleen Marchal. Syntren: a generator of synthetic gene ex- pression data for design and analysis of structure learning algorithms. BMC bioinformatics, 7(1): 1-12, 2006.
Hal R Varian. Causal inference in economics and marketing. Proceedings of the National Academy of Sciences, 113(27):731--7155, 2016.
Zichao Wang, Sebastian Tschiatschek, Simon Woodhead, José Miguel Hernández-Lotd Lobato, Si- mon Peyton Jones, Richard G Baraniuk, and Cheng Zhang. Educational question mining at scale: Prediction, analysis and personalization. In Proceedings of the AAAI Conference on Artificial Intelligence, volume 35, pp. 15669-15677, 2021.
Yue Yu, Jie Chen, Tian Gao, and Mo Yu. Dag-gnn: Dag structure learning with graph neural networks. In International Conference on Machine Learning, pp. 7154-7163. PMLR, 2019.
Cheng Zhang, Judith Bütepage, Hedvig Kjellstrom, and Stephan Mandt. Advances in variational inference. IEEE transactions on pattern analysis and machine intelligence, (1( 2026, 2018.
Cheng Zhang, Kun Zhang, and Yingzhen Li. A causal view on robustness of neural networks. Advances in Neural Information Processing Systems, 33, 2020.
Kun Zhang and Aapo Hyvarinen. On the identifiability of the post-nonlinear causal model. arXiv preprint 1rXiv:1005.2500 2012.
Xiangyu Zheng, Xinwei Sun, Wei Chen, and Tie- Yan Liu. Causally invariant predictor with shift- robustness. arXiv preprint arXiv:2107.01892 2021.
12

---

Under review as a conference paper at ICLR 2022
Xun Zheng, Bryon Aragam, Pradeep Ravikumar, and Eric P Xing. Dags with no tears: Continuous optimization for structure learning. arXiv preprint arXiv:1803.0143: 2018
Xun Zheng, Chen Dan, Bryon Aragam, Pradeep Ravikumar, and Eric Xing, Learning sparse non- parametric dags. In International Conference Artificial Intelligence and Statistics, pp. 3414- 3425. PMLR, 2020.
13

---

Under review as a conference paper at ICLR 2022
A OPTIMIZATION DETAILS
As mentioned in the main text, we gradually increase the values of p and a as optimization pro- ceeds, that non DAGs are heavily penalized. Inspired by Notears, do this with method that resembles the updates used by the augmented Lagrangian procedure for optimization (Nemirovsky) 1999). The optimization process interleaves of p and a for a certain number of steps; and (ii) Update the values of the penalty parameters p and a. The whole optimization process involves running the sequence (i)-(ii) until convergence, until the maximum allowed number of optimization steps is reached.
steps: (i) Optimize the objective for fixed values
Step (i). Optimizing the objective for some fixed values of p and a using Adam (Kingma & Ba 2014). We optimize the objective for a maximum of 500 epochs or until convergence, whichever happens first (we assume convergence if the loss does not improve during 100 epochs. If so, we move to step (ii)). We use Adam, initialized with step-size of 0.01. During training, we reduce the step-size by a factor 10 the training loss does not improve for 50 steps. We do this a maximum of two times. If we reach the condition third time, we do not decrease the step-size and assume optimization converged, and move to step (ii).
0. At the beginning of step (i) measure the dag-penalty Pi = (a)h((A.. Then, we run step (i) as explained above. At the beginning of step (ii) we measure the dag-penalty again, P3 = Eg and update
Iterating (i)-(ii). We initialize p 1 and
If P2 < 0.65 Pt. we leave p unchanged a Pz. Otherwise, if P2 2 0.65 P1, leave unchanged and update We repeat the sequence (i)-(ii) for maximum 25 steps or until convergence (measured as a or p reaching some max value), whichever happens first. The approach is summarized in algorithms and 2
(A)). (A))).())
10 p.
Other details. We A (the regularization factor favouring sparsity).
Algorithm Outer training loop (iterating (i) and (ii)
Require: Sparsity factor Xa temperature for Gumbel-softmax т, dataset D = {x"J,
Jn=1, maximum
values Pmaxx: mmax
Initialize (p, a.) (1,0)
Initialize each component of qu (A) to a Bernouli with p 0.1
Initialize Niters
while Pmax & a < Omma & Nitors < 25 do
Estimate DAG-penale P1
1 Σ m=1 h(Am), with A; ~ (A) M
Run inner optimization loop for and p
D See algorithm
Estimate DAG-penalty P2 = M
1 Em h(Am), with A, ~ (A)
if P2 < 0.65 x P then
a = a + pP2
else. is odd = 10 x
end AfNiters Niters +
end while
B ELBO DERIVATION
The goal of maximum likelihood involves maximizing the likelihood of the observed variables. For FCause (with fully observed datasets) this corresponds to the -og-migginal likelihood
(x (.l"...,,')) - log P(A) II po po(x" (A).
(14)
Computing the marginalization from the equation above is intractable, for moderately low dimensions, since the number of terms in the sum grows exponentially with the size of (which grows quadratically with the problem's dimensionality).
14

---

Under review as a conference paper at ICLR 2022
Algorithm 2 Inner training loop (step (i))
Require:
Sparsity factor Aę temperature for Gumbel-softmax T, dataset D = {x" J Initialize Nepochs
Initialize step-size 7) = 0.01
while Not converged Nopochs < 500 do
for X E D) do
Batch of size 100
Sample A (A) with gumbel-softmax trick with temperature т Compute ELBO from eq. (6) using sample A and (=(:;) for each XB
Update & and using Vo.eEBBO with Adam update and step-size 7) if Condition to reduce step-size then
17 = 1 10 end if end for N Nepochs N opochs end while
Variational inference proposes to use distribution (x to build the ELBO, lower bound of the objective from eq. (4). as follows:
log po ..... x N) log Σ ((A) II (a(x"IA)
(15)
log Σ (A) ((A) II,, (o(x" (A)
(16)
(A)
((A) II. (a (A)
log (A)(A)
(17)
Tol (A)
log p(A) II,, (o(x" (A)
2 (g)(A)
(Jensen's inequality)
(18)
(A)
(g)(A)
log pogp(A)
II 12
Po (x*1A) + (1a))
(19)
=EBOO 0).
(20)
where we used that H (4) log (A) is the entropy of the distribution 9a. Interestingly, the distribution that maximizes the ELBO is exactly the one that minimizes the KL-divergence
between the approximation and the true posterior, (T(go (A)) [] po (p(((pp (Aa= (ar.....x%)) (see, e.g.,
(Blei
et al. 2017)). This is why 46 be used posterior approximation.
C DETAILS ON THE UNIFIED FLOW-BASED FORMULATION
Details weighted adjacency matrices These are used by Notears, Notears-MLP. Notears-Sob and Grandag, among other methods. We will present the formulation from (Zheng et al. 2020) since
is the more general one. Simply put, given function f(r, Rd R" Zheng al. (020))
propose to compute each entry of the weighted adjacency ((0) as
W(0,,, = 10,
(21)
where f. denotes the partial derivative of fe(x; 9) with respect to zj. Clearly, the output of fa 0) does not depend the value of zy, then we get (0)),, 0.
While this is a useful general notion to build these weighted adjacency matrices, it is not easy to apply in practice. Zheng al. (2020) propose the case where each J. is given by a multi layer perceptron (MLP). or as a weighted combination of nonlinear basis functions. We include the details for MLPs here. In this case, each fa given by an independent MLP. Let Q represent the weights of the first layer of ft. If the j-th column of Q(1) is all zeros, then the output of S1 cannot depend on Tj. Using this insight Zheng al. (2020)
specific alternatives for few cases, including
15

---

Under review as a conference paper at ICLR 2022
propose to build the weighted adjacency matrix
5 I| colmm("
(22)
Details for the Notears case We show full details for the flow formulation of (linear) Notears. The details for the other methods follow this closely, Notears (linear) finds the optimal weighted adjacency W by solving
N 1
min
W Wx"ll
s.t.
h(W) = 0,
(23)
W WEREX
where h(W) = O exactly characterizes DAGS (eq. (3) As proposition 1 states, the exact same optimization problem is obtained when fitting a flow under acyclicity constraints.
Proposition 1. Fitting a flow using a standard Gaussian as base distribution pz, the transforma- tion Wx, and a hard constraint enforcing acyelicity on W results in exactly the same optimization problem as the one in eq. (23) used by Notears (Zheng et al. 2018).
Proof. We enforce the acyclicity constraint on W by restricting the optimization domain to {W: h(W) 0). Then, the maximum likelihood objective be written as
Σ log p(x")
Σ 1ogps("" Wz") + log det(I W)l
(change of variable)
(24)
= Σ log P2(x" -Wx")
(lemma
(25)
Σ
||"" W
(p2 is std Gaussian),
(26)
where eq. (24 is valid because the transformation : == =g(x) = =x Wa is invertible for any W E 0.
We that we recover exactly the Notears objective from eq. 223
D PROOF OF LEMMA 1
We split the proof in several simple steps.
1. ga(x) = x - fa(r) has Jacobian-determinant det - Ja())) where Ja(z) is the Jacobian of fz(
(x))) has non-zero entries exactly in the positions where A is non-zero. Therefore, it retains the DAG structure.
3. Matrices with DAG structure are nilpotent (i.e., all eigenvalues are zero). Thus, Ja(x) can be factorized as (x) QUQ", where Q unitary and U is strictly upper triangular (Schur factorization).
Finally, det (T - fa(x)) = det - QUQ') = det(T - -U)- =1
E DAGS AND INVERTIBILITY
As mentioned section if A is not DAG the transformation (x = r. fa()) may not be invertible. A simple example suffices to show this. Consider the two-variable scenario with
A 0 1
0
and fa(x) (1(22)) where (1(2)) = x2 and f2(x1) F1- Then, computing fA(Z) we get
σ T2
2/ = Z2 T1 Clearly, the function 2 = (x) = I - fx(x) is not invertible, since for any pair (x), 22) satisfying Ji 32 we get 21 22 0.
16

---

Under review as a conference paper at ICLR 2022
F RESULTS WITH HIGHER LEVEL OF NOISE
This section shows results in the exact same setting as the one described in section 6 with the differ- ence being that synthetic datasets generated with higher level of noise (the real and pseudo-real datasets are exactly the same as before). This is achieved by sampling the noise scale d; ~ (4(022, 2) (instead of o; (0.2, 1)). Results are shown in fig. 2 As expected, as consequence of the higher level of noise, all methods perform slightly worse. However, the main conclusion remains un- changed, FCause outperforms other methods and is still competitive even when ran on the datasets with 30% of missing values.
Notears-MLLF P
Grandag Notears
ICALINGAM
FCause (PO)
Noteera-Sob
Adjacency F1
0.8
0.6
0.4
Orientation F1
0.8
0.8
0.4
0.2
Causal Accuracy
0.8
0.6
0.4
0.2
0.0
ER(16.16)
ER(16,84)
32(32.32)
ER(32, 128)
ER(64,84)
E5(64.250)
SF(16,16) 16)
SF(16,84)
9F(32,32)
SF(82,188)
SF(04 64)
SF(64,256)
Protein
SynTReN
Figure 2: FCause achieves better results than the baselines in all three metrics shown. Results obtained with the dataset with higher levels noise. The legend "Flow ViCause (PO)" corresponds to rumning FCause with 30% of the training data missing completely at random. For readability, highlight the FCause and FCause (PO) results by connecting them with soft lines. For clarity the figure shows mean results across four different random seeds.
G FULL RESULTS WITH STANDARD DEVIATIONS
Results for all metrics, including standard deviations, are shown in fig. 3 can be observed that some methods achieve higher precision values than FCause for some datasets. Looking at the sults detail, this is not surprising. The methods that achieve better precision values always have extremely low recall. Meaning that they detect very few edges, which tend be correct. This leads to quite high precision values, but very low recall, Fl and causal accuracy. That is, their overall performance is quite poor, since they fail detect most edges. (e.g. Suppose that problem has 100 edges. A method that detects just one edge, and that edge is correct, gets precision of one, but a recall and Fl of almost zero). On the other hand, while FCause tends achieve slightly lower precision values, its overall performance is significantly better.
H MISSING VALUES IMPUTATION
Datasets. We use the same datasets described section 6, all standardized. We remove 30% of the values in the training randomly. The testing sets composed of n. 500 samples for the
17

---

Under review as a conference paper at ICLR 2022
Grandag Notears
CCLINGAM
FCause
Nolears-Sob
PC
Adjacency F1 vdjacency procision Adjacency recall
1.00
0.75
025
0.00
1.00
0.75
0.50
025
0.75
0.50
0.25
0.00
Orientation recall
0.75
0.25
0.00
Orientation precision
1.0
0.5
Orientation F1
0.75
0.50
0.25
0.00
Causal Accuracy
0.50
0.25
ER(16,16)
84(16,84)
3R(32,32)
ER(32, 128)
E4(64.84)
ER(G4 256)
SF(16,16) 16)
54(16.54)
SF(82.32)
9F(22118)
SF(64 64)
SF(64.250 256)
Protein
SynTReN
Figure 3: Results for all metrics including standard deviations computed over runs with four different random seeds. For readability, we highlight the FCause and FCause (PO) results by connecting them with soft lines.
synthetic datasets, = 100 samples for the pseudo-real data obtained with the SynTReN generator, and 53 samples for the protein cells dataset.
18

---

Under review as a conference paper at ICLR 2022
Baselines. We compare FCause against mean imputing, mice (Van Buuren & Groothuis- Oudshoom 2011, missforest (Stekhoven & Bühlmann 2012) and PVAE (Ma et al., 2018). We also include as baseline FCause with an adjacency matrix fixed to strictly upper triangular. We call this method
UTFlow.
Metrics. We report the normalized Root Mean Squared Error (NRMSE) computed on the test set. All results are computed as the mean across four different random seeds.
Faause UTFFiw
Mice
Meen impuling
0.17
NRMSE
0.16
0.1
0.15
18(18.18)
ER(16.64)
ER 32)
ER(32 128)
ERY84 84)
ER 2080
SF(18, 16)
9F(16.64)
SF(32 32)
SF(32, 128)
SF/64 64)
SF 258)
0.08
Protein
SynTRaN
Figure 4: Imputation results. All results are obtained as the mean over four runs with different random seeds. For readability, highlight the FCause results by connecting them with soft lines. For SynTRen, mean imputing achieves a NRMSE of 0.29, not shown in the plot.
Results are shown in fig. 4 For clarity, all values reported the mean across four runs with different random seeds. There is no single method that dominates across all datasets. Overall, it can be observed that Missforest achieves the best performance, followed by FCause and UTFlow, and finally Mice. PVAE, while better than mean imputing, performs worse than other methods in most datasets. Perhaps surprisingly, UTFlow and FCause achieve very similar performances.
19


Reviews:
[{"review_id": "xBqQXyYepr-", "reviewer": "Reviewer_2kV5", "summary": "", "questions": "", "limitations": "", "rating": 3, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "The paper builds up on Carefl (Khemakhem et al 2021), which fits an autoregressive flow model using the variables' causal ordering, extending it to use a flow-based model for inferring the causal structure. \n\nOne of the issues I have is that it is supposed to provide a more scalable alternative to better understood methods, e.g. constraint based methods like PC or FCI. Both of these methods can easily scale up to hundreds of variables (also since they do not exhaustively search over the space of all DAGs but work on equivalence classes) and provide theoretical guarantees. Moreover it does not provide any identifiability result, showing that the learnt graph is indeed causal. Given the approximated posterior, this cannot be directly derived from standard score-based methods.", "summary_of_the_paper": "The paper proposes a flow-based causal discovery method that is supposed to be scalable (in the experiments d=16,32,64) and able to impute missing values. ", "summary_of_the_review": "The paper seems to be a relatively incremental contribution, that is not embedded in the related work on causal discovery and does not provide any theoretical guarantee that one would learn a causal graph. ", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "2: Several of the paper’s claims are incorrect or not well-supported.", "technical_novelty_and_significance": "2: The contributions are only marginally significant or novel.", "empirical_novelty_and_significance": "2: The contributions are only marginally significant or novel."}, "extra_text": {}}, {"review_id": "RM5ZSpq6Omc", "reviewer": "Reviewer_q8GH", "summary": "", "questions": "", "limitations": "", "rating": 3, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "Strengths:\n\n•\tFirst paper to combine VB with NOTEARS\n\n•\tEmpirically, it does well in simulations\n\n•\tIt can handle missing values and allows nonlinear functional relationships\n\nWeaknesses:\n\n•\tThe motivation of using DAG is unclear. DAG is a special case of directed cyclic graphs (DCG). One important advantage of DAG over DCG is its simple factorization, which leads to more efficient computation than DCG. With the NOTEARS framework, this advantage seems largely lost because during the training, at least in the early phase, the graphs may be not be acyclic. Then why one even makes an often unrealistic assumption of acyclicity to begin with? In a typical run of the proposed algorithm, is the DAG factorization ever exploited? If it is exploited in iterations when the graph is acyclic, then typically in how many iterations the graph is acyclic? In summary, acyclicity is, in my opinion, a modeling restriction, not a feature. Is the proposed algorithm faster than the same algorithm without the h(A) penalty (which is both simpler computationally and more general in terms of modeling)? \n\n•\tThe proposed approach seems doing very well in various simulations but doesn’t do well in the real data, which questions the real practical advantage of the proposed method. Perhaps including more real data analyses/comparison can help support it. \n\n•\tIt would be helpful for readers if computation/algorithm/implementation details are given. \n\n•\tNot clear why Carefl is relevant to this paper.  Maybe I missed something, but the only place Carefl is mentioned in the method description is around Equation (6). But it looks to me just a standard definition of a nonlinear SEM? See e.g., Equation (4) in “On Causal Discovery with Cyclic Additive Noise Models”. \n\n•\tScalability is claimed but not tested: the largest dimension investigated is 64. \n\n•\tAlthough the proposed method can handle missing data, missing at random is a rather strong assumption. Under this assumption, practically any Bayesian methods can easily impute the missing values. \n\n•\tIs Lemma 1 a well-known result for any DAG? The Jacobian can be transformed to an upper or lower triangular matrix with unit diagonals by permuting the nodes and hence the determinant is trivially 1. I might have missed something here. \n\n•\tWhat is the rationale of imposing prior on A but not on theta. This choice loses some uncertainty in addition to that lost by using Variational Bayes. ", "summary_of_the_paper": "The paper proposes a flow-based causal discovery method to learn nonlinear causal DAGs with potential missing values. The proposed approach is based on variational Bayes and continuous optimization approach. ", "summary_of_the_review": "In summary, while the paper has some interesting idea of combining variational Bayes algorithm with NOTEARS (continuous optimization) for nonlinear DAG learning, the significance of the contribution is not yet strongly supported theoretically or empirically. ", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "2: Several of the paper’s claims are incorrect or not well-supported.", "technical_novelty_and_significance": "2: The contributions are only marginally significant or novel.", "empirical_novelty_and_significance": "2: The contributions are only marginally significant or novel."}, "extra_text": {}}, {"review_id": "NTclBs180ZG", "reviewer": "Reviewer_Du7p", "summary": "", "questions": "", "limitations": "", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "The authors consider a fundamental problem of learning directed causal graphs with nonlinear continuous data.  They propose a scale causal discovery method based on the flow, variational inference, GNN, and Notears.  Their method can handle high dimensional and missing data.   \n \nHowever, there are some concerns / suggestions:\n\n1. It's not clear to me whether the proposed optimization formulation (Equation 6 or 10) can always render the highest score? I notice that it does not restrict the distribution of noise z. Is there any theorem to support it? Please clarify it.\n2. I strongly recommend that the authors should give a summary, such as an algorithm, to show the optimization process.  \nRegarding the unified flow-based formulation (lemma 2). \n3. As far as I know, in the additive noise setting, the log-likelihood of observed variables can be converted into the log-likelihood over the noise estimations, such as in the following paper. This is not surprising to me though the author provides another strategy to prove it.\n \nCai R, Qiao J, Zhang Z, et al. Self: structural equational likelihood framework for causal discovery[C]//Thirty-second AAAI conference on artificial intelligence. 2018.\n \n4.  The authors should explain why the metrics of the other methods cannot perform well. For example,  Notear perform incredibly bad, why?  \n5. More experiments are needed, report based on only four data sets is unconvincing.   \n\nminor typos:\n\nEquation 3: no symbol description for tr,e.\n \n \nIf the authors sufficiently address the mentioned concerns, I am happy to change my assessment.", "summary_of_the_paper": "This  paper  focuses on causal discovery,  especially  on learning DAG using continuous optimization methods. \\\nIn the nonlinear additive noise setting, the authors propose a flow-based method, called FCause, to optimize the log-likelihood of the observed variables.  ", "summary_of_the_review": "The authors address an important and hard problem in learning DAG. The method is a novel and interesting contribution to the continuous optimization literature \n", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "3: Some of the paper’s claims have minor issues. A few statements are not well-supported, or require small changes to be made correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "Not applicable"}, "extra_text": {}}, {"review_id": "3aXOIhsLJBy", "reviewer": "Reviewer_m3am", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "### Pros\n- The paper is well written and easy to follow.\n- The paper develops a unified fiew of existing continous optimization methods for learning DAGs, which provide insights into how those methods relate to one another.\n\n### Cons\n- In terms of methodology, the idea of normalizing flow has been explored by [1]. I would encourage the authors to include explanation about the difference between their work and [1].\n- Also, the methods of using gumbel-softmax and binary matrix have been explored by [1, 2], which may be worth mentioning and comparing to.\n- Regarding the experiments, it would be better to include more nonlinear/nonparametric baselines, such as GES with generalized score [3], PC with kernel-based test [4], DAG-GNN [5], CAM [6], and [2].\n- The results of GraN-DAG seem to be much worse than those reported in the original paper. Did the authors use the preliminary neighborhood selection and pruning steps (Appendix A.3 in [7])?\n- It would be better to include the structural intervention distance [8] that measures the accuracy of interventional queries and could be more informative.\n- The proposed $f_i(x)$ with graph neural network seems to be restrictive, since it assumes a causal additive model (i.e. elementwise nonlinear relationship $g_j$) with a \"post-nonlinear\" operation $h_i$. However, it is surprising that this function works well with the Gaussian process dataset considered in the experiments, although the functional assumption is not met. Could the authors provide explanation on this?\n\n### Other comments\n- $h$ is used for different purposes, i.e. neural network and algebraic characterization of DAGs.\n- I am not sure if it is appropriate to claim about causality without further conditions. As the authors know, previous works [5, 7, 9] focus on learning DAGs instead of claiming about causal discovery. Perhaps could the authors include some explanation on this?\n- The graph neural networks used in Section 4.1 seem to be similar to [5, 10] which fall under the proposed unfied view; the authors could include comparison with them in Section 5.\n\n### References\n1. Differentiable Causal Discovery from Interventional Data, 2020.\n2. Masked Gradient-Based Causal Structure Learning, 2020.\n3. Generalized Score Functions for Causal Discovery, 2018.\n4. Kernel-based Conditional Independence Test and Application in Causal Discovery, 2011.\n5. DAG-GNN: DAG Structure Learning with Graph Neural Networks, 2019.\n6. CAM: Causal additive models, high-dimensional order search and penalized regression, 2014.\n7. Gradient-Based Neural DAG Learning, 2020.\n8. Structural Intervention Distance (SID) for Evaluating Causal Graphs, 2015.\n9. DAGs with NO TEARS: Continuous Optimization for Structure Learning, 2018.\n10. A Graph Autoencoder Approach to Causal Structure Learning, 2019.", "summary_of_the_paper": "This paper proposes a general flow-based approach to learn DAGs from data which provides a unified view of existing continuous optimization methods for structure learning. As a side benefit, the authors demonstrate that the proposed method could naturally be modified to handle missing data. The authors provide empirical studies to show that their proposed method outperforms the other baselines.", "summary_of_the_review": "The paper develops a unified fiew of existing continous optimization methods for learning DAGs, which provide insights into how those methods relate to one another. However, I have some concerns about the experiments, and am willing to increase my score if the authors could address my comments.", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "3: Some of the paper’s claims have minor issues. A few statements are not well-supported, or require small changes to be made correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}, {"review_id": "0DDuCj91LvK", "reviewer": "Reviewer_DG4z", "summary": "", "questions": "", "limitations": "", "rating": 8, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "Strengths:\n\n- The method presented in the paper is scalable and empirically gives good results in various settings.\n\n- All aspects of the method are founded on theoretical developments in many earlier papers.\n\n- The developments are also used to state a framework that makes it much easier to state the relations between different algorithms in prior work.\n\nConcern:\n\n- notears turned out to have a particular weakness to different noise variances. Could you discuss possible weaknesses that your proposed method might have?\n\nMinor points:\n\n- Final sentence in section 2, \"better properties\" - better than what?\n\n- page 3: double \"to\"\n\n- page 5: double \"to\"\n\n- page 9: \"practitioners\" -> \"practitioner's\"", "summary_of_the_paper": "This paper proposes to combine the continuous optimization-based causal discovery approach from notears with flow-based function learning. An extension is given to data missing (completely) at random. Also, a unifying framework is presented facilitating the comparison and exchange of ideas between different continuous optimization-based causal discovery methods, viewing all in terms of flows.", "summary_of_the_review": "I think this is a strong paper, both for the main result as well as for the unifying framework in section 5. My only concern is that the proposed method has many \"moving parts\" making it hard to see how it might or might not succeed in particular situations. But the empirical results are very strong, so that I recommend acceptance.", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "4: All of the claims and statements are well-supported and correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "Not applicable"}, "extra_text": {}}]

Author rebuttal:

Meta-review: This paper proposed a flow-based approach FCause to Bayesian causal discovery that is scalable, flexible, and adaptive to missing data. Reviewers were split on this paper and could not reach a consensus during the discussion, and no reviewer pushed for acceptance. After taking a closer look myself, I agree with several of the reviewers that while the core ideas here are interesting and novel, there remain too many unresolved issues that require another round of revision.

I encourage the authors to carefully take in account the reviewers' comments and re-submit this promising work to another ML venue.
Decision: Reject