AnAlysis specific to biological context (Fig. 1c) . The notion of a causal link is fundamentally distinct from a correlational link (Fig. 1d) . Causal network inference is profoundly challenging 15, 16 , and many methods for inferring regulatory networks connect correlated, or mutually dependent, nodes that might not have any causal relationship. Some approaches (e.g., directed acyclic graphs [17] [18] [19] ) can in principle be used to infer causal relationships, but their success can be guaranteed only under strong assumptions 15, 20 that are almost certainly violated in biological settings. This is due to many limitations-some possibly fundamentalin our ability to observe and perturb biological systems. These observations imply that it is essential to undertake careful empirical assessment in order to learn whether computational methods can provide causal insights in specific biological settings. Network inference methods are often assessed using data simulated from a known causal network structure (a so-called gold-standard network 5, 17 ). Such studies (and their synthetic biology counterparts 21 ) are convenient and useful, but at the same time they are limited because it is difficult to truly mimic specific biological systems of interest. Inferred networks are often compared to the literature, but for the purpose of learning novel, potentially context-specific, regulatory relationships, this is an inherently limited approach, and experimental validation of network inference methods has remained limited 9, 10, 19, 22 . With the support of the Heritage Provider Network (HPN), we developed the HPN-DREAM challenge to assess the ability to learn causal networks and predict molecular time-course data. The Dialogue for Reverse Engineering Assessment and Methods (DREAM) project 23 (http://dreamchallenges.org) has run several challenges focused on network inference 22, [24] [25] [26] [27] . Here we focused on causal signaling networks in human cancer cell lines. Protein assays were carried out using reverse-phase protein lysate arrays 28, 29 (RPPAs) that included functional phosphorylated proteins. The HPN-DREAM challenge comprised three sub-challenges. Sub-challenge 1 was to infer causal signaling networks using protein time-course data. To focus on networks specific to genetic and epigenetic background, the task spanned 32 different contexts, each defined by a combination of cell line and stimulus, and each with its own training and test data. The test data were used to assess the causal validity of inferred networks, as described below. A companion in silico data task also focused on causal networks but by design did not allow the use of known biology. Sub-challenge 2 was to predict phosphoprotein time-course data under perturbation. This sub-challenge comprised both an experimental data task and an in silico data task, and the same training data sets were used as in sub-challenge 1. Sub-challenge 3 was to develop methods to visualize these complex, multidimensional data sets. Across all sub-challenges, the scientific community contributed 178 submissions. In the network inference sub-challenge we found that several submissions achieved statistically significant results, providing substantive evidence that causal network inference may be feasible in a complex, mammalian setting (we discuss a number of relevant caveats below). The use of pre-existing biological knowledge (e.g., from online databases) seemed to be broadly beneficial. However, FunChisq, a method that did not incorporate any known biology whatsoever, was not only the top performer in the in silico data task but also highly ranked in the experimental data task. Challenge data, submissions and code have been made available as a community resource through the Synapse platform 30 , which was used to run the challenge (https://www.synapse.org/HPN_ DREAM_Network_Challenge; methods applied in the challenge are described in Supplementary Notes 1-3). For the experimental data network inference task, participants were provided with RPPA phosphoprotein data from four breast cancer cell lines under eight ligand stimulus conditions 31 . The 32 (cell line, stimulus) combinations each defined a biological context. Data for each context comprised time courses for ~45 phosphoproteins (Supplementary Table 1 ). The training data included time courses obtained under three kinase inhibitors and a control (dimethyl sulfoxide (DMSO)) ( Fig. 2a ; details of the experimental design, protocol, quality control and pre-processing can be found in the Online Methods). The data set is also available in an interactive online platform (http://dream8.dibsbiotech.com) that uses the Biowheel design developed by the winning team of the visualization sub-challenge. Participants were tasked with using the training data to learn causal networks specific to each of the 32 contexts. Networks had to comprise nodes corresponding to each phosphoprotein with directed edges between the nodes. The edges were required to have weights indicating the strength of evidence in favor of each If node A causally influences node B via measured node C, the causal network should contain edges from A to C and from C to B, but not from A to B (top). However, if node C is not measured (and is not part of the network), the causal network should contain an edge from A to B (bottom). Note that in both cases inhibition of node A will lead to a change in node B. (c) Causal edges may depend on biological context; for example, a causal edge from A to B appears in context 1, but not in context 2 (lines in graphs are as defined in a). (d) Correlation and causation. Nodes A and B are correlated owing to regulation by the same node (C), but in this example no sequence of mechanistic events links A to B, and thus inhibition of A does not change the abundance of B (lines in bottom right graph are as defined in a). Therefore, despite the correlation, there is no causal edge from A to B.  possible edge, but they did not need to indicate sign (i.e., whether activating or inhibitory). For the companion in silico data task, participants were provided with data generated from a nonlinear differential equation model of signaling 12 . The task was designed to mirror some of the key features of the experimental setup, and participants were asked to infer a single directed, weighted network (Online Methods and Supplementary Fig. 1 ). Whereas the experimental data task tested both data-driven learning and use of known biology, the in silico data task focused exclusively on the former, and for that reason node labels (i.e., protein names in the underlying model) were anonymized. An incorrect causal network can score very well on standard statistical assessments of goodness of fit or predictive ability; for example, two nodes that are highly correlated but not causally linked (Fig. 1d ) may predict each other well. For the experimental data task, we therefore developed a procedure that leveraged interventional data to assess inferred networks in a causal sense. The key idea was to assess the extent to which causal relationships encoded in inferred networks agreed with test data obtained under an unseen intervention (Fig. 2a) . Specifically, for a given context c, we identified the set of nodes that showed salient changes under a test inhibitor (here an mTOR inhibitor) relative to the DMSO-treated control ( Fig. 2b and Online Methods). These nodes can be regarded as descendants of the inhibitor target (mTOR) in the underlying causal network for context c. We denote this gold-standard descendant set by D c GS (Supplementary Fig. 2 ). Note that D c GS may include both downstream nodes and those influenced via feedback loops within the experimental time frame. We emphasize that these 'gold-standard' sets are derived from (held-out) experimental data and should not be regarded as representing a fully definitive ground truth. For each submitted context-specific network, we computed a predicted set of mTOR descendants (D c pred ) and compared it with D c GS to obtain an area under the receiver operating characteristic curve (AUROC) score (Fig. 2c) . Teams were ranked in each of the 32 contexts by AUROC score, and the mean rank across contexts was used to provide an overall score and final ranking (Online Methods, Fig. 2d and Supplementary Fig. 3a) . We tested the robustness of the rankings using a subsampling strategy (Online Methods). In addition to mean ranks, we used mean AUROC scores (across the contexts) in the analyses described below; these scores complement the mean ranks by giving information on the absolute level of performance, and the two metrics are highly correlated (Supplementary Fig. 3c ). For the in silico data task, the true causal network was known (Online Methods and Supplementary Fig. 4 ) and was used to obtain an AUROC score for each participant that determined the final rankings (Supplementary Fig. 3b ). An alternative scoring metric to AUROC is the area under the precision-recall curve (AUPR), which is often used when there is an imbalance between the number of positives and negatives in the gold standard 32 . Some of our settings were imbalanced, and we therefore compared rankings based on the AUROC and AUPR, which showed reasonable agreement  Train (three inhibitors and DMSO)  Across the 32 contexts included in the experimental data network inference task ( Fig. 3a) , a mean of 11.8 teams (s.d. = 7.3; Supplementary Fig. 7 includes a full set of counts by context) achieved statistically significant AUROC scores (FDR < 0.05; multiple testing correction performed within each context with respect to the number of teams; Online Methods). For the in silico data task, the top 14 teams achieved significant AUROC scores ( Supplementary  Fig. 3b ). The fact that several teams achieved significant scores with respect to causal performance metrics suggests that causal network inference may be feasible in this setting. Supplementary Table 2 presents a summary of submissions. Scores on the experimental data and in silico data network inference tasks were modestly correlated (r = 0.35, P = 0.011) but were better correlated when only teams that did not use prior information were compared (r = 0.68, P = 0.002; Fig. 3b and Supplementary Note 4). To identify teams that performed well across both tasks, we averaged ranks for experimental and in silico data tasks ( Fig. 3b and Supplementary Fig. 3d ). To test the notion of 'crowdsourcing' 22, 27, 33, 34 for causal network inference, we combined inferred networks across all teams (Online Methods, Fig. 3c and Supplementary Fig. 8a ). For the experimental data task, this ensemble or aggregate submission slightly outperformed the highest-ranked individual submission (mean AUROCs of 0.80 and 0.78, respectively), and for the in silico data task it ranked within the top five (AUROC of 0.67). Combinations of as few as 25% of randomly chosen submissions performed well on average (mean AUROCs of 0.72 and 0.64 for experimental and in silico data tasks, respectively; Fig. 3d and Supplementary Fig. 8b ). Methodological details were provided by 41 of the 80 participating teams (Supplementary Note 1), allowing us to classify submissions (Fig. 3e,  The best-scoring method for the experimental data task, "PropheticGranger with heat diffusion prior, " by Team1, used a prior network created by averaging similarity matrices. The matrices were obtained via simulated heat diffusion applied to links derived from the Pathway Commons database 35 . The prior network was then coupled with an L1-penalized regression approach that considered not only past but also future time points (a detailed description is presented in Supplementary Note 1). The best scoring approach for the in silico data network inference task, and the most consistent performer across both data types, was the FunChisq method by Team7 (Supplementary Note 1) . This approach used a novel functional χ 2 test to examine functional dependencies among the variables and did not use any biological prior information. Before FunChisq was applied, the abundance of each protein was discretized via the Ckmeans.1d.dp method 36 , with the number of discretization levels selected using the Bayesian information criterion on a Gaussian mixture model. incorporating pre-existing biological knowledge On average, teams that used prior biological information outperformed those that did not ( Fig. 4a ; one-sided rank-sum test, P = 0.032). The submission ranked second used only a prior network and did not use the protein data. However, use of a prior network did not guarantee good performance, with mean AUROC scores ranging from 0.49 to 0.78 for teams using a prior network. Interestingly, the same prior network that was itself ranked second was used in both the top-performing submission and the submission ranked 43rd, the difference being the approach used to analyze the experimental data. Conversely, not using a prior network did not necessarily result in poor performance; mean AUROC scores ranged from 0.49 to 0.71 for teams not using a prior network. The top-performing teams using prior networks in the experimental data task did not perform as well in the in silico data task (Fig. 3b) . To further investigate the influence of known biology, we combined submitted prior networks to form an aggregate prior network (Online Methods). This outperformed the individual prior networks and had a score similar to that of the aggregate submission described above (mean AUROC of 0.79). We combined the aggregate prior network with each of the two top methods (PropheticGranger and FunChisq) in varying proportions (Fig. 4b) . Combining FunChisq with the aggregate prior improved upon the aggregate prior alone (this was not the case for PropheticGranger). Finally, we considered three-way combinations of PropheticGranger, FunChisq and the aggregate prior; the highest-scoring combination consisted of 20% PropheticGranger, 50% FunChisq and 30% aggregate prior (mean AUROC of 0.82; Supplementary Fig. 9 ). We set the combination weights by optimizing performance on the test data; we note that because additional test data were not available, we could not rigorously assess the combination analyses.  The overall score in the experimental data task was an average over all contexts; to gain additional insight, we further investigated performance by context. In line with their good overall performance, aggregate submission and prior performed well relative to individual submissions in most contexts (Fig. 4c) . The aggregate prior network performed particularly well for cell line MCF7 but less well for BT549, supporting the notion that biological contexts differ in the extent to which they agree with known biology. The aggregate submission offered the greatest improvements over the aggregate prior in settings where the aggregate prior performed less well, suggesting that combining data-driven learning with known biology might offer the most utility in noncanonical settings. The context-specific aggregate submission networks (see Fig. 5a for an example) provided crowdsourced signaling hypotheses. Comparing the aggregate submission with the aggregate prior network helped to highlight potentially novel edges; we have provided a list of context-specific edges with their associated scores as a resource (Supplementary Table 3) . Dimensionality reduction suggested that differences between cell lines were more prominent than those between stimuli for a given cell line ( Fig. 5b and Online Methods), in line with the notion that (epi)genetic background has a key role in determining network architecture. In the time-course prediction sub-challenge, participants predicted phosphoprotein time courses obtained under interventions not seen in the training data (Online Methods). We assessed predictions by direct comparison with the test data using root-mean-square (r.m.s.) error (Online Methods and Supplementary Note 6), focusing on predictive ability rather than causal validity. Supplementary Table 4 and Supplementary Note 2 present team scores and descriptions of submissions. Testing the robustness of team ranks gave two top performers for the experimental data task and a single top performer for the in silico data task (Online Methods). The two top performers for the experimental data task took different approaches. Team42 (ranked second) simply calculated averages of values in the training data. Team10 (ranked third) used a truncated singular value decomposition to estimate parameters in a regression model. This method also ranked highly for the in silico data task and was the most consistent performer across both data types. Team44, the top-ranked team, was not eligible to be named as a top performer because of an incomplete submission (Supplementary Note 7) , but their approach also consisted of calculating averages. The good performance of averaging may be explained to some degree by a shortcoming with the r.m.s. error metric used here (Supplementary Fig. 10) . Team34, the top performer for the in silico data task, used a model informed by networks learned in the network inference sub-challenge. This suggests that networks can also have a useful role in purely predictive analyses. A total of 14 teams submitted visualizations that were made available to the HPN-DREAM Consortium members, who then voted for their favorite (Online Methods). The winning entry, Biowheel, is designed to enhance the visualization of time-course protein data and aid in their interpretation (Supplementary Note 3) . The data associated with a cell line are plotted to depict protein-abundance levels by color, as in a heat map, but are displayed as a ring, or wheel. Time is plotted along the radial axis and increases from the center outward. The interactive tool provides a way to mine data by displaying data subsets in various ways.  Inferring molecular networks remains a key open problem in computational biology. This study was motivated by the view that empirical assessment will be essential in catalyzing the development of effective methods for causal network inference. Such methods will be needed to systemically link molecular networks to the phenotypes they influence. Although causal network inference may fail for many theoretical and practical reasons, our results, obtained via a large-scale, community effort with blinded assessment, suggest that the task may be feasible in complex mammalian settings. By "feasible, " we mean capable of reaching a performance level significantly better than that achieved by chance, and this was accomplished by a number of submissions, including approaches that did not use any prior information. Our assessment approach focuses on causal validity and is general enough to be applicable in a variety of settings, such as gene regulatory or metabolic networks. However, it is important to take note of several caveats. First, the procedure relies on the specificity of test inhibitors. However, if the inhibitor were highly nonspecific, it would probably not be possible to achieve good results or for a prior network to perform well, because the predictions themselves are based on assumed specificity. In addition, data suggest that the mTOR inhibitor used here is reasonably specific 37 . Second, the procedure used only one of the inhibitors for testing, whereas rankings could be changed by the inclusion of additional inhibitors. Hill et al. 31 used a cross-validation-type scheme that iterated over inhibitors. Such an approach, although more comprehensive, is not possible in a 'live' challenge setting, as training and test data must be fixed at the outset. Third, the procedure does not distinguish between direct and indirect causal effects. Finally, all downstream nodes were weighted equally, regardless of whether they were context specific. Metrics that better emphasize context-specific effects will be an important avenue for future research and would probably shed further light on the utility of priors (which are not usually context specific). We also emphasize that further work is needed to clarify the theoretical properties of the score used here with respect to capturing agreement with the (unknown) ground truth. Several submissions used novel methods or incorporated novel adaptations of existing methods (Supplementary Tables 2 and 4) . Notably, the best-performing team for the network inference in silico data task developed a novel procedure (FunChisq) that also performed well on the experimental data task without use of prior information, increasing confidence in its robustness. Indeed, the ability to make such comparisons is a key benefit of running experimental and in silico challenges in parallel. Although some approaches performed well on one data type only (Fig. 3b) , the overall positive correlation between experimental and in silico scores is striking given that they were based on different data and assessment metrics. Teams that did not use prior information were relatively well correlated (Fig. 3b) , suggesting that good performers among these teams on the in silico data task could perform competitively on experimental data if their methods were extended to incorporate known biology. The observation that prior information alone performs well reflects the fact that much is already known about signaling in cancer cells and suggests that causal networks are not entirely 'rewired' in those cells. However, our analysis revealed contexts that deviate from known biology; such deviations are likely to be particularly important for understanding disease-specific dysregulation and therapeutic heterogeneity. Furthermore, it is possible that the literature is biased toward cancer, and for that reason priors based on the literature may be less effective in other disease settings. We anticipate that in the future a combination of known biology and data-driven learning will be important in elucidating networks in specific disease states. A previous DREAM challenge also focused on signaling networks in cancer 26 . However, the scoring metric was predictive rather than causal (r.m.s. error between predicted and test data points) with a penalty related to sparseness of the inferred network. Our assessment approach shares similarities with other approaches in the literature, including those used by Maathuis et al. 38 , who focused on inferring networks from static observational data, and Olsen et al. 39 , who used a different scoring metric, considering predicted downstream targets in close network proximity to the inhibited node. It remains unclear to what extent the ranking of specific submitted methods could be generalized to different data types and biological processes. In our view, it is still too early to say whether there could emerge broadly effective 'out-of-the-box' methods for causal network inference analogous to methods used for some tasks in statistics and machine learning. Given the complexity of causal learning and the wide range of application-specific factors, we recommend that at the present time network inference efforts should whenever possible include some interventional data and that suitable scores, such as those described in this paper, be used for empirical assessment in the setting of interest. AnAlysis of the challenge. D.E.C. and Y.Z. performed analyses to compare top-performing approaches submitted for the network inference sub-challenge. D.E.C., A.S., E.O.P., C.K.W., K.G., A.B. and J.M.S. designed the top-performing approach in the experimental data network inference task. Y.Z., H.W. and M.S. designed the approach that performed best in the in silico data network inference task and was the highest ranked across both experimental and in silico data network inference tasks. F.Z. and Y.G. developed an algorithm that was a top performer in the experimental data time-course prediction task and was also the highest ranked across both experimental and in silico data time-course prediction tasks.  The authors declare no competing financial interests. Brief. Bioinform. 15, 195-211 (2014
Of particular note is the DREAM network inference challenge where priorbased methods took the top two positions in an independent evaluation REF .
Inferring causal molecular networks: empirical assessment through a community-based effort
Nature Methods
Nature Methods
Medicine Biology