CausalBench: A Large-scale Benchmark for Network Inference from Single-cell Perturbation Data
Authors: Mathieu Chevalley
Yusuf H Roohani
Arash Mehrjou
Jure Leskovec
Patrick Schwab
Venue: Submitted to ICLR 2023 (2023)
Keywords: causal discovery
gene regulatory networks
gene-gene interaction networks
network inference
single cell RNA sequencing
scRNAseq

Abstract: Mapping biological mechanisms in cellular systems is a fundamental step in early-stage drug discovery that serves to generate hypotheses on what disease-relevant molecular targets may effectively be modulated by pharmacological interventions. With the advent of high-throughput methods for measuring single-cell gene expression under genetic perturbations, we now have effective means for generating evidence for causal gene-gene interactions at scale. However, inferring graphical networks of the size typically encountered in real-world gene-gene interaction networks is difficult in terms of both achieving and evaluating faithfulness to the true underlying causal graph. Moreover, standardised benchmarks for comparing methods for causal discovery in perturbational single-cell data do not yet exist. Here, we introduce CausalBench - a comprehensive benchmark suite for evaluating network inference methods on large-scale perturbational single-cell gene expression data. CausalBench introduces several biologically meaningful performance metrics and operates on two large, curated and openly available benchmark data sets for evaluating methods on the inference of gene regulatory networks from single-cell data generated under perturbations. With real-world datasets consisting of over 200000 training samples under interventions, CausalBench could potentially help facilitate advances in causal network inference by providing what is - to the best of our knowledge - the largest openly available test bed for causal discovery from real-world perturbation data to date.

Under review as a conference paper at ICLR 2023
CAUSALBENCH: A LARGE-SCALE BENCHMARK FOR NETWORK INFERENCE FROM SINGLE-CELL PERTUR- BATION DATA
Anonymous authors Paper under double-blind review
ABSTRACT
Mapping biological mechanisms in cellular systems is fundamental step in early- stage drug discovery that serves to generate hypotheses on what disease-relevant molecular targets may effectively be modulated by pharmacological interventions. With the advent of high-the thggthhruughput methods for measuring single-cell gene ex- pression under genetic perturbations, evidence for causal gene-gene interactions at scale. However, inferring graphi- cal networks of the size typically encountered in real-world gene-gene interaction networks is difficult in terms of both achieving and evaluating faithfulness to the underlying causal graph. Morcover, standardised benchmarks for comparing methods for causal discovery in perturbational single-cell data do not yet exist. Here, we introduce CausalBench a comprehensive benchmark suite for eval- uating network inference methods on large-scale perturbational single-cell gene expression data. CausalBench introduces several biologically meaningful perfor- mance metrics and operates two large, curated and openly available benchmark data sets for evaluating methods on the inference of gene regulatory networks from single-cell data generated under perturbations. With real- world datasets consisting over 200 000 training samples under interventions, CausalBench could poten- tially help facilitate advances in causal network inference by providing what is to the best of our knowledge the largest openly available test bed for causal discovery from real-world perturbation data date.
have effective means for generating
INTRODUCTION
Studying causality in real-world environments is often challenging because uncovering causal rela- tionships generally either requires the ability to intervene and observe outcomes under both inter- ventional and control conditions, or reliance on strong and untestable assumptions that can not be verified from observational data alone (Stone (1993); Pearl (2009); Schwab et al. (2020); Peters cet al. (2017)). In biology, a domain characterised by enormous complexity of the systems studied, establishing causality frequently involves experimentation in controlled in-vitro lab conditions ing appropriate technologies to observe response to intervention, such as for example high-content microscopy (Bray et al. (2016)) and multivariate omics measurements (Bock al. (2016)) High- throughput single-cell methods for observing whole transcriptomics measurements in individual cells under genetic perturbations (Dixit al. (2016); Datlinger al. (2017] 202T)) has recently emerged as promising technology that could theoretically support performing causal inference in cellular systems at the scale of thousands of perturbations per experiment, and therefore holds enor- mous promise in potentially enabling researchers to uncover the intricate wiring diagrams of cellular biology (Yu al. (2004); Chai et al. 2)014; Akers & Murali (202 } Hu al. (2020).
However, while the combination of single-cell perturbational experiments with machine learning holds great promise for causal discovery, making effective use of such datasets is to date still a challenging endeavour due to the general paucity of real-world data under interventions, and the difficulty of establishing causal ground truth datasets to evaluate and compare graphical network inference methods (Neal et al. (2020); Shimoni et al. (2018;; Parikh al. (2022).. In order to progress the causal machine learning field beyond reductionist (semi- semmi-)synthetic experiments towards potential utility in impactful real-world applications, it imperative that the causal machine learning

---

Under review as a conference paper at ICLR 2023
Oheenvtiosal
Intovestioral
Racorabacked Gauel Gragh
Gene
Gene
Gane
Gene
Gasa
Oera
Gume
Network Infernnce Algaritera
Figure I: An overview of causal gene-gene network inference in mixed observational and per- turbational single-cell data. The causal generative process in its unperturbed form is observed in the observational data (left; 10000+ samples in CausalBench) while data under genetic interven- tions (c.g., CRISPR knockouts) are observed in the interventional data (right; 200 000+ samples in CausalBench. Either observational or interventional plus observational data that were sampled from the true causal generative process (bottom distributions) can be used by network inference algorithms (bottom right) infer reconstructed causal graph (top right) that should as closely as possible recapitulate the original underlying functional gene-gene interactions.
research community develop and maintain suitable benchmarks for objectively comparing methods that aim to advance the causal interpretation of real-world interventional datasets.
To facilitate the advancement of machine learning methods in this challenging domain, we in- troduce CausalBench a comprehensive benchmark suite for evaluating network inference meth- ods on perturbational single-cell RNA sequencing data that to the best of our knowledge the largest openly available test bed for causal discovery from real-world perturbation data to date (Figure [T. CausalBench contains meaningful biologically curated set of two large-scale perturbational single-cell RNA sequencing experiments with over 200 000 interventional samples each that are openly available, and integrates numerous baseline implementations of state-of-the-art methods for network inference from single-cell perturbational data. Similar to benchmarks in other domains, e.g. ImageNet in computer vision (Deng et al. 2009 ), we hope CausalBench can help accelerate progress large-scale real-world causal graph inference, and that the methods developed against CausalBench could eventually lead to new therapeutics and deeper understanding of human health through enabling the reconstruction of the functional gene-gene interactome. The source code for CausalBench is openly available at https: /ttpps/////ttbb...// com/ wttp://coo///aanooooo-3333333 43213123 / cauusaaeench bbenh
performance metrics, a
Contributions. Our contributions are as follows:
We introduce CausalBench comprehensive benchmark suite for evaluating network in- ference methods on perturbational single-cell RNA sequencing data consisting of two cu- rated, openly available benchmark datasets with over 200 000 interventional samples each. We introduce a set of meaningful benchmark metrics for evaluating performance including novel statistical metric that leverages single cell perturbational data measure perfor- mance against a larger set of putative gene regulatory relationships than would be possible using observational data alone.
Using CausalBench, we conduct comprehensive experimental evaluation of the perfor- mance of state-of-the-art network inference algorithms in recovering graphical relation- ships for mixed observational and interventional scRNAseq data. We implement relevant state-of-the-ar methods as baselines for network inference from observational and inter- ventional single-cell data.
In addition, we evaluate the performance and scaling characteristics of network inference methods under varying numbers of available training samples and intervention set sizes to establish whether state-of-thh-an network inference algorithms are able to effective use of different scales of intervention and training sample set sizes.
2

---

Under review as a conference paper at ICLR 2023
2
RELATED WORK
Background. Given an observational data distribution, several different causal networks or di- rected acyclic graphs (DAG) could be shown to have generated the data. The causal networks that could equally represent the generative process of an observational data distribution are collectively referred to the Markov Equivalence Class (MEC) of that DAG (Huang et al. (018)). Interven- tional data offers an important tool for limiting the size of the MEC to improve the identifiability of the true underlying causal DAG. In the case of gene expression data, modern gene editing tools such as CRISPR offer a powerful mechanism for performing interventional experiments at scale by altering the expression of specific genes for hundreds of genes at once and observing the resulting interventional distribution across the entire transcriptome (Dixit et al. (2016); 6); Datlinger et al. (2017) 202T)) The ability leverage such interventional experiment data scale could significantly im- prove our ability to uncover underlying causal functional relationships between genes and thereby strengthen our quantitative understanding of biology. Establishing causal links between genes can help implicate genes in biological processes causally involved in disease states and thereby open up new opportunities for therapeutic development (Mehrjou al. (2022): Shifrut al. (2018)).
Network inference in mixed observational and interventional data. Learning network struc- ture from both observational and interventional data presents significant potential in reducing the search space over possible causal graphs. Traditionally this network inference problem has been solved using discrete methods such as permutation-
approaches (Wang et al. (2017); Hauser & Bühhmann] (2012). Recently, several new models have been proposed that can differentiably learn causal structure (Schölkopf et al.] (202T)) However, most of them focus on observational datasets alone. Ke et al. (019)) presented the first differentiable causal learning approach using both obser- vational and interventional data. Lopez et al. (2022) improved the scalability of differentiable causal network discovery for large, high-dimensional datasets by using factor graphs to restrict the search space, and Scherrer al. (2021) introduced an active learning strategy for selecting interventions to optimize differentiable graph inference.
Gene-regulatory network inference. The problem of gene regulatory network inference has been studied extensively in the bioinformatics literature in the case of observational datasets. Early work modeled this problem using a Bayesian network trained on bulk gene expression data (Friedman) et al. (2000)) Subsequent papers approached this as feature ranking problem where machine learning methods such as linear regression (Kamimoto et al.] (2020) or random forests (Huynh-) Thu et al. (2010); Aibar al. (2017)) are used predict the expression of any one gene using the expression of all other genes. However, [Pratapa et al.] (2020) showed that most gene regulatory network inference methods for observational data perform quite poorly when applied to single- cell datasets due to the large size and noisiness of the data. In the case of constructing networks using interventional data, there relatively much lesser work given the recent development of this experimental technology. Dixit et al. (2016) were the first to apply network inference methods to single-cell interventional datasets using linear regression.
In contrast to existing works, we exhaustively test several different network inference approaches for this task using a combination of both observational and interventional data. We present the first standardized benchmark for single-cell causal inference that enables direct comparison of state- of-the-art machine learning methods using a set of quantitative metrics with both biological and statistical motivation. Our contribution is also unique in the significantly larger size of the employed dataset both in terms of number of samples and number of interventions than previous work.
3 METHODOLOGY
In the next section, we briefly introduce the formal framework of SCMs in order to serve as a causal language with which to describe methods, limitations and assumptions in this setting in a formal way, as well as to motivate the quantitative metrics presented subsequently. The pertubational nature of the data necessitate formal statistical language that goes beyond associations and correlations, thus we present and use the causal view as was introduced by Pearl 2009). We will utilise this language to refer back to the setting throughout the rest of the paper.
3

---

Under review as a conference paper at ICLR 2023
Structural Causal Models (SCMs). In this work, we rely the language of SCMs describe causal properties (Peters et al. (2017)) Formally, a SCM M consists of a 4-tuple (U, F.F P(u)), where U is a set of unobserved (latent) variables and X is the set of observed (measured) variables. F is a set of function such that for each X; E X, X; < Ji( U.). U, E U and Pa; E X X, The SCM induces distribution over the observed variables P(x), where the uncertainty stems from the distribution on the unobserved variables P(u). The variable-parent relationships can be represented in a directed graph, where each X; node in the graph, and there is a directed edge between all Pa; to X,, The task of causal discovery can then be described as learning this graph over the variables. The language of causality and in particular the formalism of SCM allows to precisely describe the notion of an intervention in system. In the most general sense, an intervention on a variable X; can be thought as uniformly replacing its structural assignment with a new function X; + f(FFa, U.). For simplicity, we often consider only simple new functional assignment, such as atomic (constant function) and stochastic (the function only depends on Ü.) In this work, we only consider such atomic or stochastic intervention and denote an intervention on X; as (X,). We can then describe the interventional distribution, denoted pxr(x)(x),, as the distribution entailed by the modified SCM, where only the assignment for X, changed and the other assignment are kept similar to the observational distribution. For consistency, we denote the observational distribution as both pD (x) and P(x).
For the rest of this paper, we will utilise the SCM framework to describe methods, assumptions and setting.
Problem Setting. We consider the setting where we are given a dataset of vector samples E Rd, where x; represents the measured expression of gene in a given cell. Modern gene editing tech- nology (e.g., CRISPR) performs direct interventions individual genes to bring their expression towards 0 in case of full on-target efficacy. Given the stochastic nature of cell biology. it is not guaranteed that this results in 0 expression, i.e, the interventions are close to atomic at level O. The goal of a causal model is to learn a causal DAG G, where each node is single gene. The causal DAG G induces a distribution over observed sample P(x), such that:
stochastic even though they
P(x) = II p(a.
The datasets we consider for the benchmark can be thought as consisting of data sampled from (x), as well as (x) for various i. We artificially create different training conditions, namely observational, partially interventional and interventional. The observational setting is re- stricted to samples x ~ (x)... The interventional settings contains the observational settings, as well as interventional data x ~ (x) for all variables i. The partial interventional setting is a middle ground between the above two, where only a fraction of the variables have interventional data associated to them, i.e, where they are the target of the intervention.
Assumptions and challenges. Single-cell data presents idiosyncratic challenges that may break the common assumptions of many existing methods. Apart from the -igensiinnsio in terms of number of variables and large sample size, the distribution of the gene expression present a challenge as it is highly tailed at 0 for some genes (Tracy et al. (2019)). Regarding the underlying true causal graph, biological feedback loops may break the acyclicity assumption underlying the use of DAGs. In addition, different cells sampled from the same batch may not be truly independent as the cells may have interacted and influenced their states. Lastly, cells in scRNAseq experiments are measured at fixed point in time and may therefore have been sampled at various points in their developmental trajectory or cell cycle (Kowalezyk et al. (2013)) making sampling time a potential confounding factor in any analysis of scRNAseq data.
Network inference methods. The benchmarked methods are given data, either only observational or also interventional depending on the setting consisting of the expression of each gene in each cell. For interventional data, the target gene in each cell is also given as input. We do not enforce that the methods need to learn graph on all the variables, and further preprocessing and variables selection is permitted. The only expected output is a list of gene pairs that represent directed edges. No properties of the output network, such as for example acyclicity, are enforced
4

---

Under review as a conference paper at ICLR 2023
either. As such, even though in this first version of the benchmark, we mainly evaluate causal methods, our benchmark can be used to evaluate any graph outputs without restrictions.
Evaluating inferred networks. Contrary to other benchmarks where the true underlying graph is either known, or the data is sampled from a simulated graph, given the nature of inferring causal graphs from real-world experimental data, the true causal graph is unknown in CausalBench as the data is collected from a complex biological process. The lack of a gold standard ground truth network brings about unique challenges in establishing weight of evidence that supports the use of one method over another: In this setting, there is no certitude regarding the existence of a universal causal graph, that would hold across time points, conditions as well cell types. To address this challenge, we aim at developing a collection of synergistic metrics for which higher performance across their set can be assumed to correlate with how close the outputted network is to a true graph. To this end, we implement two types of evaluation: biologically-motie motiivaeed and quantitative statistical evaluation of inferred networks.
Biologically-motivaeed evaluation. The biologically-motivaee evaluation in CausalBench is based on biological datahases of known putatively causal gene-gene interactions. Using these databases of domain knowledge, we can construct putatively true undirected subnetworks to evalu- ate the output networks in the understanding that the discovered edges that are not present in those databases are not necessarily false positives, and undirected edges may not necessarily be direct edge connections in a causal graph, i.e, a functional interaction between two genes may be mediated by one more other genes. The underlying assumption behind the aalogicclly-mooivatttt evaluation approach is that methods that are better will, on average across experimental settings, recall more known gene-gene interactions than those that worse.
Quantitative statistical evaluation. In contrast to the biologically evaluation, the quan- titative statistical evaluation in CausalBench is fully data-driven, cell-specific and prior free. It's also a novel form of network evaluation that is unique to single-cell perturbational data and presents a new approach for approximating ground truth augment the information contatined in biological databases. To conduct the statistical evaluation, we rely on the interventional data at our disposal in perturbational scRNAseq data evaluate edges in the output networks. The main assumption for this evaluation is that if the discovered edge from A to B is true edge denoting functional interaction between the two genes, then perturbing should have a statistically significant effect on the distribution (x((x)) (xg) of values that gene B takes in the transcription profile given enough observed samples, compared to its observational distribution pi (xB) (i.e, compared to control sam- ples where no gene was perturbed). The significance of this change in distribution can be tested with statistical test as well as with distributional distances. One of the challenges in using a statistical ap- proach is that the interventional data may be limited in size in practice, which means that some signal may not be significantly detectable with finite data. Lastly, the quantitative statistical approach can- not differentiate between causal effects from direct edges or from causal paths in the graph. Despite evaluation, by
those limitations, this quantitative evaluation complements the mooloiialll-mmoivatte extending the covered metric set with metrics that are directly derived from the observed data, mak- ing this evaluation data-driven, cell-specific and prior free. Moreover, it presents a novel approach for estimating ground-truth gene regulatory interactions that is uniquely made possible through the size and interventional nature of single-cell perturbational datasets. Full descriptions of both the looogicclll-motivate and quantitative statistical evaluations are presented in Section 4
4 EXPERIMENTAL EVALUATION
Datasets. For the experimental evaluation in CausalBench, we rely on the recently published large- scale seRNAseq dataset by |Replogle et al. (2022). Replogle et al. (2022) performed a genome- scale perturb-seq screen, targeting all expressed genes across millions of human cells with CRISPR perturbations. They generated three datasets, one genome-scale pertubational dataset consisting of 2.5 million K562 cells, and two smaller ones, on the K562 (Andersson et al. (1979) and the RPEI cell lines, targeting only DepMap essential genes (Tsherniak et al. (2017)). We utilize the last two smaller datasets, they focus on putatively important genes, and they allow to fairly compare performance on two distinct cell types. To stabilize evaluation and training, we additionally filter out genes for which there are less than 100 perturbed cells, and hold out 20% of the data for
5

---

Under review as a conference paper at ICLR 2023
dataset
total samples
# observational samples
# gene interventions
Replogle al. 2022)
K562
310 385
10 691
1 158
Reelogle et al. 2022)
RPEI
247 914
11 485
651
Table 1: High-level description of the two large-scale datasets utilised in CausalBench charac- terised by their high numbers of samples and of intervened-upe -unnn variables. Of those samples, after stratification by intervention target (including no target), 20% were kept has held-out data.
evaluation, stratified by intervention target. A summary of the resulting two datasets can be found in T
Metrics. To implement the biologically evaluation, we extract network data from two
widely used open biological databases: CORUM (Giurgiu et al. 2019) and STRING (Von Mering
et al. 2005 2005) Snel et al. 2000; Von Mering et al. 2007 Jensen et al. 2009 Mering et al. [2003 Szklar-
czyk et al. 2010 2015 2016 | |019] Franceschini et al.
CORUM is a repository of experimentally characterized protein complexes from mammalian organisms. The complexes extracted from individual experimental publications, and exclude results from high-throughpe ex- periments. We extract the human protein complexes from the CORUM repository and aggregate them to form a network of genes. STRING is a repository of known and predicted protein-protein interactions. STRING contains both physical (direct) interactions and indirect (functional) inter- actions that we use to create two evaluation networks from STRING: protein-protein interactions (network) and protein-protein interactions (physical). protein-proteii interactions (physical) con- tains only physical interactions, whereas protein- -roetein interactions (network) contains all types of know and predicted interactions. STRING, and in particular string-network can contain less reliable links, as the content of the database is pulled from variety of evidence, such as high-throughput lab experiments, (conserved) co-expression, and text-mining the literature.
To implement the quantitative evaluation, we leverage the hold out observational and interventional data and compute one statistical test and one distribution distance. For each edge A to B, we extract the observational samples for B, and the interventional samples for B where A was targeted. We then perform a two-sided Mann- Whitney rank test (Mann & Whitney 1947 Bucchianico 1999) using the SciPy package (Virtanen et al. 2020) to test the null hyptothesis that the two distributions are equal. We count an edge as true positive if the null hypothesis of the statistical test is rejected. To improve statistical power, given that the distribution are heavy tailed towards 0.0 expression, we trim 90% of the 0.0 values of the observational distribution, and the same fraction of the smallest values of the interventional samples. For the distributional distance, we compute a Wasserstein distance (Ramdas et al. 2017) between the two empirical distributions. We then return the mean Wasserstein distance of all the inferred edges. We do not normalize the computed distance, but we expect this to not have significant effect as the samples always have more or less the same size (exactly the same for the observational samples, and similar for interventional, and least 20 by design),
The advantage compared to a statistical test is that we obtain continuous metric that correlate with the strength of causal effect of parent nodes on children nodes and does not need choosing hyperparameter such the p-value threshold. A higher mean Wasserstein distance would then indicate better performance.
Baselines. We implement representative set of existing state-of-the-art methods for the task of causal discovery from single-cell observation and mixed perturbational data. For the observa- tional setting, we implement PC (named after the inventors, Peter and Clark; constraint-based method) (Spirtes et al. 2000), Greedy Equivalence Search (GES; score-based method) (Chick- ering, 2002), and NOTEARS variants NOTEARS (Linear), NOTEARS (Linear,LII) NOTEARS (MLP) and NOTEARS (MLP.L1) (Zheng et al., 2018 2020). In the interventional setting, we in- cluded Greedy Interventional Equivalence Search (GIES, a score-based method and extension to GES) (Hauser & Bühlmann 2012), and the Differentiable Causal Discovery from Interventional Data (DCDI) variants DCDI-G and DCDI-DSF (continuous optimization based) (Brouillard et al.) 2020). GES and GIES greedily add and remove edges until score the graph is maximized. NOTEARS and DCDI enforce acyclicity via a continuously differentiable constraint, making them suitable for deep learning. We provide broad overview of method classes included below:
PC is one of the most widely used methods in causal inference from observational data that assumes there are no confounders and calculates conditional independence to give asymptotically correct

---

Under review as a conference paper at ICLR 2023
results. It outputs the equivalence class of graphs that conform with the results of the conditional independence tests.
Greedy Equivalence Search (GES) implements a two-phase procedure (Forward and Backward phases that adds and removes edges from the graph) to calculate a score to choose within an equiv- alence class. While GES leverages only observational data, its extension, Greedy Interventional Equivalence Search (GIES), enhances GES by adding a turning phase to allow for the inclusion of interventional data.
NOTEARS formulates the DAG inference problem as a continuous optimization over real-valued matrices that avoids the combinatorial search over acyclie graphs. This is achieved by constructing a smooth function with computable derivatives over the adjacent matrices that vanishes only when the associated graph is acyclic. Various versions of NOTEARS refer to which function approximator is employed (either MLP or Linear) or which regularity term is added to the loss function (e.g., L1 for the sparsity constraint.)
Differentiable Causal Discovery from Interventional Data (DCDI) (Brouillard et al. 2020) lever- ages various types of interventions (perfect, imperfect, unknown), and uses a neural network model to capture conditional densities. DCDI encodes the DAG using a binary adjacency matrix. The intervention matrix is also modeled as binary mask that determines which nodes are the target of intervention. A likelihood- differentiable objective function is formed by using this param- cterization, and subsequently maximized by gradient-base based methods to infer the underlying DAG. DCDI-G assumes Gaussian conditional distributions while DCDI-DSF lifts this assumption by using normalizing flows to capture flexible distributions.
GRNBoost Aiiar et al. 2017 is GRN specific Gradient Boosting tree method, where for every gene, candidate parent gene are ranked based on their predictive power toward the expression profile of the downstream gene. As such, it acts feature selection method toward learning the graph. GRNBoost was identified as one of the best performing GRN method in previous observational data based benchmarks (Pratapa al. 2020
Random (k) is the simplest baseline which outputs a graph from which k nodes are selected at uniformly random without replacement. In the experiments, we tested 1: 100, 1000, and 10000.
Unfortunately, all of the chooe-mentioned methods, with the exception of NOTEARS and GRN- Boost, do not computationally scale to graphs of the size typically encountered in transcriptome- wide gene regulatory network inference. To nonetheless enable meaningful comparison, we pro- pose to partition the variables into smaller subsets and run the subset independently with the final output network being the union of the subnetworks. Acyclicity is guaranteed if each subnetwork is acyelic. The proposed approach breaks the no-Latent-confounder assumption that some methods may make, and it also does not fully leverage the available informa- tion potentially in the data. Methods that do scale to the full graph in a single optimization loop are therefore expected to perform better.
methods on each
For all the above taate-of-the-ar methods, we use publicly available implementations and we imple- ment interfaces between our benchmark suite and the open-source packages. We use causal-learn (Zhag et al.) (022) for PC and GES, the official implementation of the authors for DCDI (Brouillard et al. 22022. , public python implementation for GIES (Gamella 2022). and the official implemen-
lation of the authors for NOTEARS (Zheng et al. 2018 2020).
5 RESULTS AND DISCUSSIONS
Network inference. We here summarize the experimental results of our baselines in both observa- tional and interventional settings, on the quantitative evaluation (statistical test) and the biologically- motivated evaluations (corum, string-network and string physscall............ All results are obtained by train- ing on the full dataset three times with different random seeds. We highlight the trade-off between true positives count and pese rate precision. Indeed, we expect methods to optimize for these two goals, as we want to obtain a high true pesttiveiive--atte precision while maximizing the number of discovered interactions. In the observational setting, we observed this trade-off, some methods can achieve very high true positive rates precision, but with very small output graphs (Figure 2). Also, methods that regularize for sparsity, such as NOTEARS (Linear,L.I) and (MIP.L.)). have
7

---

Under review as a conference paper at ICLR 2023
RPE1 (observational)
KS62 (observational)
True Posiivv Rate [4)
130
Statistical Edge Test
Protein Complexes
100
Statistical Edge Test
Protein Complexes
True Posttvv Rate
50
50
25
True Postive Rate (8)
130
settein-protein
Protein-pracein Interactins (physical)
True Postivv flate (N0
100
Protein-protein
Protein-pratein
75
intercctonnseenn (network)
75
Interactions (network)
interactions (physical)
0
10
10"
10?
10
10
10a
10
10*
10
16a
10F
10
10
True Positives [cound]
True Postives [cound]
Positives [coun]]
True Positives [round]
APE1 (interventional)
KS62 (interventional)
130
100
True Posiive (ate
Statistical Edge Test
Protein Complexes
Thus Pasttive Raze
Statistical Edge Test
Protein Complexes
50
50
25
25
0
True Positive Rate(*)
100
Protein-protein Interactions (network)
Protein-precein meractions (physical)
True Positive Rate]%)
100
Protein-pratein
prctein-protein Interactions (physical)
75
Interactions (neowork)
56
10
10
10
10'
10
10'
10
10
10'
10"
10"
True Positives [cound]
True Postives [cound]
True Poskives [cound]
True Posttives [count]
Figure 2: Performance comparison in terms of Positive Rate (TPR, in %; y-axis) Precision and
True Positives (TP, in counts; x-axis) in correctly identifying edges substantiated by four different biological interaction databases, Statistical Edge Test (ours; top left), protein complexes (CORUM, Giurgiu et al. (2019); top right), protein-protedd interactions from networked sources (STRING Szk-
larczyk et al. 2021 1 bottom left), and ppottin interactions from physical assays (STRING
Szklarczyk et al. (2021); bottom right), between 7 different methods (see legend in top row) using observational (top row) and 6 different methods (see legend in bottom row) using interventional data data (bottom row) in RPE1 (left column) and K562 (right column) cell lines. For each method, we report the results for three independently calculated runs with new random seeds to assess uncer- tainty. Dashed lines in the same color the markers indicate the median performance observed across all randomly seeded runs.
higher true positive rates precision in general, which suggests the importance of sparsity when deal- ing with noisy data. Regarding the quantitative statistical test, antisfactory-resttlle on the RPE1 dataset high precision (around 75%) is only achieved the RPE1 dataset compared to around 50% on the K562 dataset. This could be explained in part by the interventions having greater effects in this cell line, as reported also by Replogle et al. (022)) In the interventional set- ting. results are similar to the observational setting, but with most methods outputting larger graphs with smaller true positive-rateee rates precision. Surprisingly, the difference to the observational setting is not large, especially in the statistical test where access to some of the interventional data should theoretically help. This highlights an opportunity for further method development for causal graph inference in order to be able to fully leverage more interventional data in the future. To offer more principled way to compare the models across all these metrics, we propose a simple and unbiased way to compute ranked scoreboard in Appendix C
Computational performance. All methods were given the same computational resources, which consists of 30 CPU's with 64GB of memory per each. We additionally assign GPU for the DCDI methods. The hyperparameters of each method, such as partition sizes, are chosen such that the running time remains below 30 hours. Partition sizes for each model can be found in Appendix B
8

---

Under review as a conference paper at ICLR 2023
0.30
Moon We
Distence (unitless
0.8
KS62 (observational -varylng dataset size)
0.25
RPE1 (observational varying dataset size)
0.6
020
0.4
0:15
0.2
0.10
0.20
K562 (interventional- varying dataset size)
RPE1 (interventional varying dataset size)
Mian Wasserstein
Distorte (unitess)
0.17
0.15
0.12
0.12
0.10
90 Sedee 1009
0.10
20
Fraction of Dataser [66
Fraction of Dataset [16
Mean We Wescerstin
Distanza [untties]
0.16
K562 (interventional varying intervention set)
0.20
RPE1 (interventional -varying intervention sec)
0.14
0.17
0.15
0.12
0.12
0.10
0.10
20
40
60
20
Fraction of Intervention Set [N]
Fraction of Intervention Set(N) (N)
Figure 3: Performance comparison in terms of Mean Wasserstein Distance (unitless; y-axis) of methods for causal graph inference on observational data (top row; (see legend top right) and 6 methods interventional scRNAseq data (centre row; see legend centre right) when varying the fraction of the full dataset size available for inference (in %; x-axis), and methods on interventional data (bottom row; see legend bottom right) when varying the fraction of the full intervention set used (in %, X-axis). Markers indicate the values observed when running the respective algorithms with one of three random seeds, and colored lines indicate the median value observed across all tested random seeds for a method.
Performance function of sample size. We additionally studied the effect of sample size on the performance of the evaluated state-of-the-art models. We randomly subset the data at different fractions and report the mean Wasserstein distance as explained in the quantitative evaluations part of Section 3 and shown in Figure 3 In the observational setting, the sample size does not seem to have significant effect performance for most methods, indeed having a slightly negative effect for some methods. This suggest that the ability to effectively leverage a large sample size of perturbational data is also open challenge that needs be addressed for high performance. In the interventional setting, a positive impact for larger sample size is observable, especially for the methods that rely deep networks and gradient-based learning such as DCDI, whereas GIES seems to suffer in a large sample setting.
Performance by the fraction of perturbations. Beyond the size of the training set, we also stud- ied the partial interventional setting where only subset of the possible genes to perturb are ex- perimentally targeted. We adapt the fraction of randomly targeted genes from 5% (low ratio of interventions) to 100% (fully interventional). We randomly subset the genes at different fractions, using three different random seeds for each method, and report the mean Wasserstein distance as a measure of quantitative evaluation. We would expect a larger fraction of intervened genes to lead to higher performance, as this should facilitate the identification of the true causal graph. We indeed observe a better performance for DCDI we again observe negative effect of having more data, probably mainly due to the fact that higher fraction of intervened genes implies larger number of training samples in total.
increase the fraction of intervened genes. For GIES,
Limitations. Openly available benchmarks for causal models for large-scale single cell data could potentially accelerate the development of new and effective approaches for uncovering gene reg- ulatory relationships. However, some limitations to this approach remain: firstly, the biological networks used for evaluation do not fully capture ground truth gene regulatory networks, and the reported connection are often biased towards more well studied systems and pathways (Gillis et al. (2014)). True ground truth validation would require prospective and exhaustive interventional wet- lab experiments. However, at present, experiments the scale necessary to exhaustively map gene- gene interactions across the genome are cost prohibitive for all possible edges even for the largest research institutes in the world and even when not considering cell heterogeneity and the confound-

---

Under review as a conference paper at ICLR 2023
ing factor of time. We expect more publicly available experimental validation data aggregate over time, and are looking to incorporate these data into the benchmark in the future. Beyond limitations in data sources used, there are limitations with some of the assumptions in the utilised state-of-the-art models: For instance, feedback loops between genes are a well-known phenomenon in gene regulation (Carthew (2006); [evine & Davidson (2005) that unfortunately presently cannot be represented by existing causal network inference methods.
CONCLUSION
We introduced CausalBench a comprehensive open benchmark for evaluating algorithms for dis- covering gene regulatory networks using large-scale CRISPR-based interventional scRNAseq dataset. CausalBench introduces a set of biologically meaningful performance metrics to compare the graphs proposed by causal inference methods quantitatively using statistical tests, and accord- ing to existing biological knowledge bases on protein-protein interactions and protein complexes. In addition, CausalBench contains implementations of range of state-of-the-art causal discovery algorithms from those that rely on independence tests (e.g., the PC algorithm), score-based methods (e.g.. GIES) and the more recently introduced continuous optimization- based methods (e.g.. differ- ent versions of NOTEARS). CausalBench was designed to lower the barrier of entry for developing causal discovery algorithm for scRNAseq data by fixing mon-model-related components of the eval- uation pipeline and allow users to focus the advancement of causal network discovery methods. With benchmark sample sizes of more than 200 000 interventional samples, CausalBench is built on one of the largest open real-world interventional datasets (Replogle et al. (2022)) ushering in an era which the community of causal machine learning researchers has ready access to real-world large-scale interventional datasets for developing and evaluating causal discovery methods.
REFERENCES
Sara Aibar, Carmen Bravo González-Blas, Thomas Moerman, Vân Anh Huynh-Thu, Hana Imri- chova, Gert Hulselmans, Florian Rambow, Jean-Christophe Marine, Pierre Geurts, Jan Aerts, et al. Scenic: single-cell regulatory network inference and clustering. Nature methods, 14(11): 1083-1086, 2017.
Kyle Akers and TM Murali. Gene regulatory network inference single-cell biology. Current Opinion in Systems Biology, 26:87-97, 2021.
Leif C Andersson, Kenneth Nilsson, and Carl G Gahmberg. K562-a human erythroleukemic cell line. International journal of cancer, 23(2):143- 1979.
Christoph Bock, Matthias Farlik, and Nathan C Sheffield. Multi-omics of single cells: strategies and applications. Trends in biotechnology, 34(8):605-608, 2016.
Mark-Anthony Bray, Shantanu Singh, Han Han, Chadwick T Davis, Blake Borgeson, Cathy Hart- land, Maria Kost-Alimova, Sigrun Gustafsdottir, Christopher C Gibson, and Anne Carpenter. Cell painting, high-content imsge-based based assay for morphological profiling using multiplexed fluorescent dyes. Nature protocols, 1 1777:175711774, 2016.
Philippe Brouillard, Sébastien Lachapelle, Alexandre Lacoste, Simon Lacoste-Julien. and Alexandre Drouin. Differentiable causal discovery from interventional data. Advances in Neural Information Processing Systems, 33:21 865-21877 2020.
Philippe Brouillard, Alexandre Drouin, and Sébastien Lachapelle. Dedi python package, 2022. URL
A Di Bucchianico. Combinatorics, computer algebra, and the Wilcoxon-Mann. Martoo-MMan--hhneeee Whitney test". Jour- nal of Statistical Planning and Inference, 79:349-364, 1999.
Richard W Carthew. Gene regulation by micrornas. Current opinion in genetics development, 16 (2):203-208, 2006.
10

---

Under review as a conference paper at ICLR 2023
Lian En Chai, Swee Kuan Loh, Swee Thing Low, Mohd Saberi Mohamad, Safaai Deris, and Zalmiyah Zakaria. A review on the computational approaches for gene regulatory network con- struction. Computers in biology and medicine, 48:55-65 6555-65,400 2014.
David Maxwell Chickering. Optimal structure identificatin with greedy search. Journal of machine learning research, 3(N0v):507-554 2002.
Paul Datlinger, André Rendeiro, Christian Schmidl. Thomas Krausgruber, Peter Traxler, Johanna Klughammer, Linda C Schuster, Amelie Kuchler, Donat Alpar, and Christoph Bock. Pooled crispr screening with single-cell transcriptome readout. Nature methods, 14(3):297-301, 2017.
Paul Datlinger, André F Rendeiro, Thorina Boenke, Martin Senekowitsch, Thomas Krausgruber, Daniele Barreca, and Christoph Bock. Ultra-bigh-tthroughhuu single-cell ma sequencing and per- turbation screening with combinatorial fluidic indexing. Nature methods, 18(6):635-6422 2021.
Jia Deng, Wei Dong, Richard Socher, Li-Jia Li, Kai Li, and Li Fei-Fei. Imagenet: A large-scale hi- erarchical image database. In 2009 JEEE conference on computer vision and pattern recognition, pp.
Ieee, 2009.
Atray Dixit, Oren Parnas, Biyu Li, Jenny Chen, Charles P Fulco, Livnat Jerby-Amon, Nemanja D) Marjanovic, Danielle Dionne, Tyler Burks, Raktima Raychowdhury, et al. Perturb-seq: dissecting molecular circuits with scalable single-cell ma profiling of pooled genetic screens. cell, 167(7): 1853-1866, 2016.
Andrea Franceschini. Damian Szklarczyk, Sune Frankild, Michael Kuhn, Milan Simonovic, Alexan- der Roth, Jianyi Lin, Pablo Minguez, Peer Bork, Christian Von Mering, et al. String v9. 1: protein- protein interaction networks, with increased coverage and integration. Nucleic acids research, (DI):D808-D8L1 2012.
Andrea Franceschini, Jianyi Lin, Christian von Mering, and Lars Juhl Jensen. Svd-phy: improved prediction of protein functional associations through singular value decomposition of phyloge- netic profiles. Bioinformatics, 32(7):1085-1087, 2016.
Nir Friedman, Michal Linial, Iftach Nachman, and Dana Pe'er. Using bayesian networks to analyze expression data. Journal of computational biology, 703-4):601-600, 2000.
Juan I. Gamella. gies python package, 2022. URL https:::/ttttbb..oo htwp:/thwwwbccmmm/aggaaaaee. laa/
gies
Jesse Gillis, Sara Ballouz, and Paul Pavlidis. Bias tradeoffs in the creation and analysis of protein- protein interaction networks. Journal of proteomics, (00:44-544 2014.
Madalina Giurgiu, Julian Reinhard, Barbara Brauner, Irmtraud Dunger- Dunger-Kaltenbach, Kaltenbach, Gisela Fobo, Goar Frishman, Corinna Montrone, and Andreas Ruepp. Corum: the comprehensive resource of mammalian protein complexes-2019, Nucleic acids research, 47(D1):559--D53, 2019.
Alain Hauser and Peter Bühlmann. Characterization and greedy learning of interventional markov equivalence classes of directed acyclic graphs. The Journal of Machine Learning Research, 13 (1):2409-2464, 2012.
Xinlin Hu, Yaohua Hu, Fanjie Wu, Ricky Wai Tak Leung, and Jing Qin. Integration of single-cell multi-omics for gene regulatory network inference. Computational and Structural Biotechnology Journal, 18: 1e211992-1933 2020.
Biwei Huang, Kun Zhang, Yizhu Lin, Bernhard Schölkopf, and Clark Glymour. Generalized score functions for causal discovery. In Proceedings of the 24th ACM SIGKDD international conference on knowledge discovery data mining, pp. 1551-1560. 2018
Vân Anh Huynh-Thu, Alexandre Irrthum, Louis Wehenkel, and Pierre Geurts. Inferring regulatory networks from expression data using tree-based methods. PloS one, 5(9):e12776, 2010.
Lars Jensen, Michael Kuhn, Manuel Stark, Samuel Chaffron, Chris Creevey, Jean Muller, Tobias Doerks, Philippe Julien, Alexander Roth, Milan Simonovic, et al. String 8-a global view proteins and their functional interactions in 630 organisms. Nucleic acids research, 37(suppl I): D412-D416. 2009.
11

---

Under review as a conference paper at ICLR 2023
Kenji Kamimoto, Christy M Hoffmann, and Samantha A Morris. Celloracle: Dissecting cell identity via network inference and in silico gene perturbation. BioRxiv, 2020.
Nan Rosemary Ke, Olexa Bilaniuk, Anirudh Goyal, Stefan Bauer, Hugo Larochelle, Bernhard Schölkopf. Michael C Mozer, Chris Pal, and Yoshua Bengio. Learning neural causal models from unknown interventions. arXiv preprint arXiv: 1910.0105, 2019.
Monika Kowalezyk, Itay Tirosh, Dirk Heckl, Tata Nageswara Rao, Atray Dixit, Brian Haas, Rebekka K Schneider, Amy J Wagers, Benjamin L Ebert, and Aviv Regev. Single-cell ma-seq reveals changes cell cycle and differentiation programs upon aging of hematopoietic stem cells. Genome research, 25(12):1800-18) 1872, 2015.
Michael Levine and Eric H Davidson. Gene regulatory networks for development. Proceedings of the National Academy of Sciences, (02414)4933-4992, 2005.
Romain Lopez, Jan-Christian Hütter, Jonathan K Pritchard, and Aviv Regev. Large-scale differen- tiable causal discovery of factor graphs. arXiv preprint arXiv:2206.07824 2022.
Henry B Mann and Donald R Whitney. On test of whether one of two random variables is stochas- tically larger than the other. The annals of mathematical statistics, pp. 50 60, 1947.
Arash Mehrjou, Ashkan Soleymani, Andrew Jesson, Pascal Notin, Yarin Gal, Stefan Bauer, and Patrick Schwab. GeneDisco: Benchmark for Experimental Design in Drug Discovery. In International Conference on Learning Representations, 2022.
Christian von Mering, Martijn Huynen, Daniel Jaeggi, Steffen Schmidt, Peer Bork, and Berend Snel. String: a database of predicted functional associations between proteins. Nucleic acids research, 31(1):258-261 2003.
Brady Neal, Chin-Wer Huang, and Sunand Raghupathi. Realcause: Realistic causal inference bench- marking. arXiv preprint arXiv:2011 1.t1:20011.5007,, 2020.
Harsh Parikh, Carlos Varjao. Louise Xu, and Eric Tchetgen. Validating causal inference methods. In International Conference on Machine Learning, pp. 17346-17358, PMLR, 2022.
Judea Pearl. Causality. Cambridge university press, 2009.
Jonas Peters, Dominik Janzing, and Bemhard Schölkopf. Elements of causal inference: foundations and learning algorithms. The MIT Press, 2017.
AdttyaPratapa,
Pratapa, Amogh P Jalihal, Jeffrey N Law, Aditya Bharadwaj, and TM Murali. Benchmarking algorithms for gene regulatory network inference from single-cell transcriptomic data. Nature methods, 17(2):147-15:) 2020.
Aaditya Ramdas, Nicolás García Trillos, and Marco Cuturi. On wasserstein two-sample testing and related families of nonparametric tests. Entropy, 19(2):47, 2017.
Joseph Replogle, Reuben Saunders, Angela Pogson, Jeffrey A Hussmann, Alexander Lenail, Alina Guna, Lauren Mascibroda, Eric J Wagner, Karen Adelman, Gila Lithwick- Yanai, et al. Mapping information-rich genotype-phenotype landscapes with genome-scale perturb-seq. Cell, 2022.
Nino Scherrer, Olexa Bilaniuk, Yashas Annadani, Anirudh Goyal, Patrick Schwab, Bernhard Schölkopf, Michael C Mozer, Yoshua Bengio, Stefan Bauer, and Nan Rosemary Ke. Learning neural causal models with active interventions. arXiv preprint W.Ki::211900244: 2021.
Bernhard Schölkopf, Francesco Locatello, Stefan Bauer, Nan Rosemary Ke, Nal Kalchbrenner. Anirudh Goyal, and Yoshua Bengio. Toward causal representation learning. Proceedings of the IEEE, 10965)3612-634 2021.
Patrick Schwab, Lorenz Linhardt, Stefan Bauer, Joachim M Buhmann, and Walter Karlen. Learn ing Counterfactual Representations for Estimating Individual Dese-Response Curves. In AAAI Conference Artificial Intelligence, 2020.
12

---

Under review as a conference paper at ICLR 2023
Eric Shifrut, Julia Carnevale, Victoria Tobin, Theodore L Roth, Jonathan M Woo, Christina T Bui, Jonathan Li, Morgan E Diolaiti, Alan Ashworth, and Alexander Marson. Genome-wide crispr screens in primary human t cells reveal key regulators of immune function. Cell, 175(7):1958- 1971, 2018. 2018.
Yishai Shimoni, Chen Yanover, Ehud Karavani, and Yaara Goldschmnidt. Benchmarking framework for performance-evaluation of causal inference analysis. arXiv preprint wrrXiv:/822.05000 2018.
Berend Snel, Gerrit Lehmann, Peer Bork, and Martijn A Huynen. String: web-server retrieve and display the repeatedly occurring neighbourhood of gene. Nucleic acids research, 28(18): 3442 3444, 2000.
Peter Spirtes, Clark N Glymour, Richard Scheines, and David Heckerman. Causation, prediction, and search. MIT press, 2000.
Richard Stone. The assumptions which causal inferences rest. Journal of the Royal Statistical Society: Series B (Methodological), 55(2):555-46,, 1993.
Damian Szklarczyk, Andrea Franceschini, Michael Kuhn, Milan Simonovic, Alexander Roth, Pablo Minguez, Tobias Doerks, Manuel Stark, Jean Muller, Peer Bork, al. The string database in 2011: functional interaction networks of proteins, globally integrated and scored. Nucleic acids research, 39(supp..))1551111566 2010.
Damian Szklarczyk, Andrea Franceschini, Stefan Wyder, Kristoffer Forslund, Davide Heller, Jaime Huerta-Cepan Ceeert,CCepa,,.. Milan Simonovic, Alexander Roth, Alberto Santos, Kalliopi P Tsafou, et al. String v10: protein-protein interaction networks, integrated over the tree of life. Nucleic acids research, 43(D1):1447-D445, 2015.
Damian Szklarczyk, John H Morris, Helen Cook, Michael Kuhn, Stefan Wyder, Milan Simonovic, Alberto Santos, Nadezhda Doncheva, Alexander Roth, Peer Bork, et al. The string database in 2017: quality-controlled protein-protein association networks, made broadly accessible. Nucleic acids research, pp. gkw937, 2016.
Damian Szklarczyk, Annika L Gable, David Lyon, Alexander Junge, Stefan Wyder, Jaime Huerta- Cepas, Milan Simonovic, Nadezhda T Doncheva, John H Morris, Peer Bork, al. String vl 1: protein-protein association networks with increased coverage, supporting functional discovery in genome-wide experimental datasets. Nucleic acids research, 47D11)D667-D613. 2019.
Damian Szklarczyk, Annika L Gable, Katerina C Nastou, David Lyon, Rebecca Kirsch, Sampo Pyysalo, Nadezhda T Doncheva, Mare Legeay, Tao Fang, Peer Bork, et al. The string database in 2021: customizable moteinn-poteii networks, and functional characterization of user-uploaded gene/measurement sets. Nucleic acids research, 49(D1):D605-DD61 2021.
Sam Tracy, Guo-Cheng Yuan, and Ruben Dries. Rescue: imputing dropout events in single-cell rna-seueencing data. BMC bioinformatics, 20(1):1-11, 2019.
Aviad Tsherniak, Francisca Vazquez, Phil G Montgomery, Barbara Weir, Gregory Kryukov, Glenn S Cowley, Stanley Gill, William F Harrington, Sasha Pantel, John M Krill-Burger, et al. Defining cancer dependency map. Cell, 17003):364-57,, 2017.
Pauli Virtanen, Ralf Gommers, Travis E. Oliphant, Matt Haberland, Tyler Reddy, David Courna- peau, Evgeni Burovski, Pearu Peterson, Warren Weckesser, Jonathan Bright, Stéfan van der Walt, Matthew Brett, Joshua Wilson, K. Jarrod Millman, Nikolay Mayorov, Andrew R. J. Nel- son, Eric Jones, Robert Kern, Eric Larson, Carey, ilhan Polat, Yu Feng, Eric W. Moore, Jake VanderPlas, Denis Laxalde, Josef Perktold, Robert Cimrman, Ian Henriksen, E. A. Quintero, Charles R. Harris, Anne M. Archibald, Antônio H. Ribeiro, Fabian Pedregosa, Paul van Mul- bregt, and SciPy 1.0 Contributors. SciPy 1.0: Fundamental Algorithms for Scientific Computing in Python. Nature Methods, 17:261-272, 2020. doi:
Christian Von Mering, Lars J Jensen, Berend Snel, Sean Hooper, Markus Krupp, Mathilde Foglierini, Nelly Jouffre, Martijn A Huynen, and Peer Bork. String: known and predicted protein- protein associations, integrated and transferred across organisms. Nucleic acids research, 33 (suppl D):D43-DD437, 2005.
13

---

Under review as a conference paper at ICLR 2023
Christian Von Mering, Lars Jensen, Michael Kuhn, Samuel Chaffron, Tobias Doerks, Beate Krüger, Berend Snel, and Peer Bork. String 7- recent developments the integration and prediction of protein interactions. Nucleic acids research, 35(suppl
2007.
Yuhao Wang, Liam Solus, Karren Yang, and Caroline Uhler. Permutation- based causal inference algorithms with interventions. Advances in Neural Information Processing Systems, 30, 2017.
Jing Yu, Anne Smith, Paul P Wang, Alexander Hartemink, and Erich D Jarvis. Advances to
bayesian network inference for generating causal networks from observational biological data. Bioinformatics, 20( 30(18)33544,303, 2004.
K. Zhang, J. Ramsey, M. Gong, R. Cai, S. Shimizu, P. Spirtes, and C.Glymour. causal-learn python package, 2022. URL httpttttt//gtttttt
com/ cmu phil/ http//aassllaaaaeea
Xun Zheng, Bryon Aragam, Pradeep Ravikumar, and Eric P. Xing. DAGs with NO TEARS: Con- tinuous Optimization for Structure Learning. In Advances in Neural Information Processing Sys- tems, 2018.
Xun Zheng, Chen Dan, Bryon Aragam, Pradeep Ravikumar, and Eric Xing. Learning sparse nonparametrie DAGs. In International Conference Artificial Intelligence and Statislics, 2020.
A BENCHMARK PACKAGE
A.1 INSTALL REQUIREMENTS
The necessary python package requirements to run the benchmark can be found in the requir rementeeentt.omm txt file. They can be installed by running the following command:
pip install
A.2 RUNNING THE BENCHMARK
Here is an example of command model the k562 dataset in the observational regime:
pythen3 causaisceenchappps/aaaa app. py
dataset_name weissmann_k562 \
--data_di -ddaa_directoory /ppth/to//ata//ttt/tt
--training_regime observational \
-model_name model subset_data 1.0
--model_seed 0
Results are written to the folder at /path/to/outputt and processed datasets will be cached htppt//tt//aatee stto/data/storage See the MainApp class for more hyperparameter options, espe- cially in the (partial) interventional setting.
A.3 ADDING A MODEL
New models can easily added to the benchmark suite. The only contract for model to im- plement the returns a fully connected gene-gene graph:
class. Here is simple example where the model
from causalscbench. models. abstract_model import AbstractinferrnceeModel
class FullyConnected [AbstractnferrenceMoolll]]]
def _.init_- sreitt-(self)) (self) -> None:
super init.
def call
14

---

Under review as a conference paper at ICLR 2023
Model name
partition size
PC
30
GES
30
GIES
30
NOTEARS (Linear)
-1
NOTEARS (Linear, LI)
-1
NOTEARS (MLP)
-1
NOTEARS (MLP, L1)
-1
DCDI-DSF
50
DCDI-G
50
Table Partition sizes used for each model. -1 means that the graph was not partitioned.
self inttessiionmmttiix np. array,
interventions: [ist[str]. gene names: List[str]
training_regime: TrainingRegime, seed: int =
-> [upte:
random seed (seed)
edges = set()
for in range (len gene names
a = gene names
for j in range 1, len (gene_nnmes))
b = geme-names[i] edges add b))
edges add ((b, a)) return list (edges)
The new model then only needs to be included in the MainApp model registry.
B PARTITION SIZES
We here recapitulate the partition sizes used to be able to run each methods in Table 2
C MODEL RANKING
We here present a simple unbiased way of ranking the different models based on all the implemented evaluation metrics. We separate the rankings per cell type. First, create preliminary ranking for each evaluation metrics. For metrics based on precision-recall, we compute score by taking the mean of precision and recall (for Statistical, Protein Complexes, PPI network, PPI physical and CHIP-seq network evaluations). The weighting may be modulated. for example depending on the downstream application. Finally, for each model, we take their average rank across the evaluation specific rankings. This ranking thus gives the same weight to each evaluation method. Results are summarized Table 3 for the K562 cell line and in Table 4 for the RPEI cell line. As can be observed, methods able to scale to the full graph (GRNBoost and NOTEARS) perform the best.
D RUN TIME
We here present the average run time for each method in Table 5 for the K562 cell line and in Table 6 for the RPEI cell line
15

---

Under review as a conference paper at ICLR 2023
Model
Average rank
Random (k 100)
10.83
Random (k 1000)
10.5
Random (k= 10000)
8.66
PC
6.66
DCDI-DSF
6.5
NOTEARS (MIP,, L.)
6.0
DCDI-G
5.5
GES
5.16
GIES
4.83
NOTEARS (Linear)
4.16
NOTEARS (Linear, L1)
3.66
NOTEARS (MLP)
3.33
GRNBoost
2.16
Table 3: Ranking of the methods the K562 cell line. Methods able to scale to the full graph (GRNBoost and NOTEARS) perform the best.
Model
Average rank
Random (k 100)
10.5
Random (k 1000)
9,99
Random (k = 10000)
8.66
GIES
8.33
PC
7.66
GES
6.33
NOTEARS (MLP, LI)
5.16
NOTEARS (Linear)
5.0
DCDI-DSF
4.16
NOTEARS (Linear, LI)
3.83
DCDI-G
3.66
NOTEARS (MLP)
2.83
GRNBoost
1.83
Table 4: Ranking of the methods on the RPE1 cell line. Methods able to scale to the full graph (GRNBoost and NOTEARS), well the DCDI methods, perform the best.
Model
Run time (hours)
DCDI-DSF
60.15
DCDI-G
30.49
NOTEARS (MLP)
27.78
NOTEARS (Linear)
14.89
PC
9.58
GES
5.79
GIES
3.97
NOTEARS (Linear, LI)
2.15
NOTEARS (MLP, L1)
0.82
GRNBoost
0.3
Random (k = 10000)
0.0
Random (k 1000)
0.0
Random (k = 100)
0.0
Table 5: Run time in hours for each method on the K562 cell line. observational data only.
GRNBoost was the
16

---

Under review as a conference paper at ICLR 2023
Model
Run time (hours)
DCDI-DSF
37.52
DCDI-G
19.68
PC
13.18
NOTEARS (MLP)
10.9
NOTEARS (Linear)
3.7
GES
3,22
NOTEARS (Linear, LI)
0,96
NOTEARS (MLP, LI)
0.62
GIES
0.58
GRNBoost
0.09
Random (k 10000)
0.0
Random (k = 1000)
0.0
Random (k 100)
0.0
Table 6: Run time hours for each method on the RPEI cell line. GRNBoost was run on the observational data only.
17


Reviews:
[{"review_id": "ojyuTFmbvI", "reviewer": "Reviewer_QjTd", "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": "", "summary_of_the_paper": "The paper aims to establish a dataset / benchmark for gene regulatory network inference from single-cell data.", "summary_of_the_review": "The paper fails in achieving what a good benchmark/dataset paper should demonstrate: clarity, sound assumptions and design choices, and practicality. Reading it did not bring me any closer to being able to quickly, effortlessly evaluate a GRN inference method I may develop.", "strength_and_weaknesses": "The paper is hard to follow for me. From the abstract and introduction, it seems it aims at providing a dataset / benchmark for GRN inference that would be very easy to understand and use by researchers from outside of computational biology or network inference fields - something like DREAM 3/4/5 GRN challenges from a decade ago. However, starting from Section 3, the authors introduce a very specific formal structure (SCM) that also does not necessarily reflect the underlying biology (e.g. the DAG assumption) - it is not clear what purpose this serves. The choice of benchmark evaluation metrics (Section 4) seems somewhat ad-hoc (why STRING and not other network repository, why use just true positives instead of metrics based on precision-recall and ROC used in DREAM). On the other hand, the paper lacks details about how to actually use the dataset in practice, whether there is a clear leaderboard etc. - details that would be essential in establishing CausalBench as the benchmark for the field. It does not help that the github repository mentioned in the introduction is empty.", "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": "1: The contributions are neither significant nor novel.", "empirical_novelty_and_significance": "1: The contributions are neither significant nor novel."}, "extra_text": {}}, {"review_id": "lJ8u_ek0mX1", "reviewer": "Reviewer_2A85", "summary": "", "questions": "", "limitations": "", "rating": 5, "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": "", "summary_of_the_paper": "This paper introduces CausalBench, a framework to benchmark gene regulatory inference methods using perturbational scRNAseq data. Although having access to perturbational scRNA data can help in discovering causal gene networks, generating large-scale perturbational data where every gene undergoes some genetic intervention (knockdown) is very cost prohibitive. There is a recent study that generates genome-wide perturbational scRNA data in only two cell types. The authors of CausalBench have used these data to benchmark some of the well-known network inference methods which use observational and interventional data. For the evaluation, they propose two approaches: biological and quantitative. For the biological evaluation, they use known gene-gene or protein-protein interaction data, and for the quantitative, they use some statistical tests to validate the effects of interventions. ", "summary_of_the_review": "Causal gene network inference (DAGs) are very important in understanding the mechanisms of gene expression. However, learning such DAGs using only observational data is not possible. Recent Perturb-seq data provide helpful interventional data by genetically perturbing the target genes and measuring the expression of other genes. This data empowers the causal DAG learning methods, however, they should be tested and benchmarked via the same metrics to measure their efficacy. \n\nThis paper takes the first step in the right direction by noticing the need for benchmarking causal network inference methods using Perturb-seq data. However, the benchmarking metrics are not comprehensive enough to be certain about the efficacy of such methods. My comments:\n\n- Biological metric does not make so much sense to me. The underlying DAG driving the gene expression programs are very cell-type specific. So, the general gene/protein networks cannot be a good and reasonable metric to evaluate the correctness of the learned DAGs in K562 and RPE1 cells. \n\n- As a qualitative metric, the authors can use held out perturbations and see how well different methods can predict the effect of each perturbation on various genes. This has also been used in Lopez et al. (2022) (cited in the paper). \n\n- Why haven't the authors used the biggest Perturb-seq (2 M cells) in K562 cells as published in (Replogle et al., 2022)? It would be interesting to see how well/fast each method works in a very large dataset.\n\n- Related to the previous comment, I think one another and yet interesting benchmarking metric would be how fast each method could be trained in each dataset. This kind of benchmarking has also been done in Lopez et al. (2022). By growing sequencing technology, it is not so far that we will have access to millions of cells and thousands of perturbations genome-wide, and the methods that could handle this huge computational barrier would be the winners. \n\nMinor commets:\n\nTypos:   - we we observed --> we observed\n             - in not large --> is not large", "strength_and_weaknesses": "Strength:\n\n- This paper tries to take the first step towards benchmarking the causal gene regulatory network inference methods using large-scale perturbational scRNA-seq data, which is interesting and important for causal discovery in biology.\n\nWeakness:\n\n- The idea of having a framework to test the causal network inference methods is great if the benchmarking metrics are powerful and interesting enough. I think the metrics proposed in this paper are not that much interesting. For the biological metrics, as mentioned by the authors, the well-known gene-gene or protein-protein interaction data are not specific for the cell types that they are working. The underlying gene regulation rules are very cell-type specific and the inferred networks could not be reliably evaluated by general networks which have been derived from many cell types. ", "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": "2: The contributions are only marginally significant or novel.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}, {"review_id": "WPuxOLVtXz", "reviewer": "Reviewer_vpJf", "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": "", "summary_of_the_paper": "The authors present CausalBench, a framework to benchmark causal gene regulatory network (GRN) inference methods on perturbational single-cell RNA sequencing data (scRNA-seq). It includes evaluation metrics, baseline implementations of relevant inference methods and access to perturbational scRNA-seq data. With CausalBench the authors evaluate the ability of recovering “silver-standard” ground truth networks and how they are affected by sample size.", "summary_of_the_review": "The author's contributions to the GRN evaluation field are somewhat new but some aspects already exist in previous work. Moreover, the authors are missing some previous contributions, limiting the added cumulative value of it. On the other hand, the description of the problem, its assumptions and the description of the methods is overall clear and of good quality.\n", "strength_and_weaknesses": "Strengths:\n- Authors are aware of and nicely describe the assumptions and limitations of this evaluation task and data.\n- Including both biological and statistical metrics for evaluation of method performance.\n- Evaluation frameworks allow to easily test if novel methods improve over baselines, motivating the development of improved methods.\n\nWeaknesses:\n- Authors propose novel evaluation metrics but did not implement classical ones such as the comparison against ChIP-seq derived networks, see (Pratapa et al. 2020).\n- Authors miss to evaluate classic methods, such as SCENIC (Aibar et al. 2017), that have been shown to be the top performers in other benchmarks (Pratapa et al. 2020). Although these methods are not causal, it would be beneficial to assess if the modeling of causality actually improves GRN inference.\n- While the authors provide different metrics, it is still not clear which methods perform better than others consistently. A “consensus” score based on rankings would improve interpretation. Moreover, there is no discussion about which methods overperform the others and why.\n- Although the authors are aware of the limitation of perturbation experiments, no quality control assessment is performed to test whether the perturbation actually worked. Samples where the perturbed gene still shows high levels of gene expression should be removed or at least accounted for.\n- The partition of the data due to scalability issues raises concerns since, as the authors already mention, it breaks the no-latent-confounder assumption. Authors chose partition sizes such that the running time remains below 30 hours but do not mention anywhere the actual numbers for each method. In order to comply with this rule, some methods might use fractions of the feature space that are too small to generate valuable graphs.\n- In figure 3, authors mention “significant” increase but do not perform any statistical test to corroborate it.\n\n", "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": "2: The contributions are only marginally significant or novel.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}, {"review_id": "TXYbM4jnIy", "reviewer": "Reviewer_WpyP", "summary": "", "questions": "", "limitations": "", "rating": 5, "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": "", "summary_of_the_paper": "The authors of this paper propose a suite of tools for benchmarking network inference on both observational and interventional data single-cell gene expression data.  They run several network inference / causal discovery methods for recovering graphical relationships from scRNAseq data and evaluate the performance of the methods using existing putatative knowledge about specific gene-gene interactions as well as by comparing performance of observational and interventional data.  ", "summary_of_the_review": "Overall, I think the authors take on a challenging and important task of collating multiple datasets and implementing many causal discovery methods on these datasets.  Overall, I think the paper lacks clarity and justification around the statistical benchmarks employed.  As such, it's unclear how much value is actually added in the suite of benchmarks proposed.", "strength_and_weaknesses": "- I think the paper tackles an important problem and provides a great service by colating these datasets and metrics.  \n- My main complaint is that the paper advertises it's contribution as providing a benchmark for network inference on large scale scRNA seq data although they admit that there are no gold standard datasets available from which to actually benchmark these methods.  This is a noted limitation, but my concern is that this is a really big limitation.  It's really hard to know whether the results are biologically meaningful, or really just noisy and reflect the difficulty of identifying useful benchmarks.\n- Along these lines, I think the \"Metrics\" section needs a lot more clarification.  The authors refer to equality tests on \"the two distributions\" as well as Wasserstein distance on \"the two empirical distributions\".  What are these distributions of? There could be a lot more detail here regarding the data / statistics used to create these distributions and how the tests were calculated.\n- While I understand some of the difficulties, it's hard to know how much to read into the TPR vs TP counts graphs.  The authors need to explicitly discuss the potential for false positives / false positive rates.  Are there any ways to identify putative false positives that can be used as a pseudo-benchmark?\n\n", "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": "2: The contributions are only marginally significant or novel.", "empirical_novelty_and_significance": "2: The contributions are only marginally significant or novel."}, "extra_text": {}}]

Author rebuttal:

Decision: Reject