BACKGROUND
the discovery of molecular pathways that drive diseases and vital cellular functions is a fundamental activity of biomedical research.
unraveling disease pathways is a major component in the efforts to develop new therapies that will effectively fight deadly diseases.
furthermore, knowing pathways significantly facilitates the design of personalized medicine modalities for diagnosis, prognosis, and management of diseases.
the discovery of pathways is a challenging problem and its solution to a large extent relies on the identification of causal molecular interactions in genomics data.
by causal molecular interactions or relations we mean interactions of molecular variables that match the notion of randomized controlled experiment, which is the de facto standard for assessing causation in the general sciences and biomedicine  <cit> .
assume that a hypothetical experimenter can change the distribution of a variable x .
we say that x is a cause of y  and denote this by x→y if the probability distribution of y changes for some experimental manipulation of x.
causal molecular interactions can be discovered using randomized experiments such as interference with rna ; however such experiments are often costly, infeasible, or unethical.
fortunately, over the last  <dig> years many algorithms that infer causal interactions from observational data have been developed  <cit>  and some of them have been adopted to the high dimensionalities of modern genomics data  <cit> .
outside of biomedicine, two nobel prizes have recently been awarded in  <dig> and  <dig> for methods which seek to discover causal relations from non-experimental data  <cit> .
in our prior work we evaluated the ability of state-of-the-art causal discovery algorithms to de-novo identify unoriented edges in genome-scale regulatory networks  <cit> , which represent causal interactions between transcription factors and their target genes without distinguishing the mechanistic role of the involved molecular variables .
we deliberately avoided performing causal orientation of the discovered unoriented edges  because this problem has previously been deemed worst-case unsolvable in observational data using existing algorithms  <cit>  due to the statistical indistinguishability of causal models in the same markov equivalence class.
for example, causal models x→y and x←y have generally been assumed to be statistically indistinguishable given only observational data for x and y.
over the last  <dig> years researchers have developed a new class of methods that, given observational data for just two causally related variables x and y, aim to determine which variable is the cause and which one is the effect   <cit> .
these causal orientation methods aim to solve the problem of statistical indistinguishability of graphs in the markov equivalence class by exploiting asymmetries in the shapes of the conditional probability densities and without requiring randomized experiments.
these methods could have significant implications for the field of causal discovery because they can orient unoriented edges that are discoverable by other established techniques, e.g.
generalized local learning  or local-to-global learning   <cit> .
therefore, at face value, these causal orientation methods have the potential to reduce the number of, and in some cases even eliminate, randomized experiments needed for causal orientation of edges in the markov equivalence class of graphs and make the causal model fully identifiable from observational data alone.
as promising as these new causal orientation methods are, they have not been previously applied in genomics, where the data is usually noisy and the sample sizes are relatively small compared to prior test applications of these methods  <cit> .
in this paper we report results of an extensive study of recent causal orientation techniques in the genomics domain by  testing  <dig> methods/variants to distinguish cause  from effect  in  <dig>  causal interactions and  conducting sensitivity analyses with respect to noise and sample size for the best-performing methods.
in addition,  we introduce a new ensemble technique for causal orientation that is shown to be more accurate than any best individual causal orientation method in the tested data.
the results of this study can serve as a foundation for further development of causal orientation techniques for genomics data and establishing a context for wide applications in molecular or biomedical research.
methods
causal orientation methods
as mentioned above, the purpose of the tested causal orientation methods is to separate cause from effect given data for just two variables x and y that have a causal relation .
typically these methods are not designed to be used to causally orient pairs of variables that only have univariate association/correlation due to the possibility of confounding.
for example, in the majority of distributions, a causal structure x ← t → y implies that x and y are associated even though they are not causally affecting each other.
therefore, the presence of association in a pair of variables x and y is, in general, a necessary but not sufficient condition to be eligible for causal orientation.
a rigorous approach would involve first using correct causal discovery methods  to identify unoriented edges that denote the existence of unconfounded causal relations and then applying causal orientation techniques to orient unoriented edges.
while each causal orientation method has its own principles and sufficient assumptions that are outlined in table  <dig> , most of these techniques are based on the idea that the factorization of the joint probability distribution p into pp yields a simpler representation than the factorization into pp.
one can furthermore show that if the marginal probability distribution of the cause p is independent of the causal mechanism p, then the factorization pp has lower complexity than the factorization pp.
given two causally related variables x and y, estimating the complexity of the two different factorizations of p or determining independence between marginal and conditional distributions can thus provide the basis for causal orientation techniques.
in practice, however, it is difficult to directly test independence between p and p or estimate  their complexity; hence the methods typically use simplifying assumptions or rely on approximate formulations.
• x and e <dig> are independent;
• f is non-linear, or one of x and e is non-gaussian;
• probability densities are strictly positive;
• x and e <dig> are independent;
• either f <dig> or e <dig> is gaussian;
• f is continuous and invertible;
• x and e are independent;
• e is gaussian;
• x and e are independent;
• e is gaussian.
• x and e <dig> are independent;
the last column indicates whether a method is sound, i.e.
it can provably orient a causal structure under its sufficient assumptions.
because causal orientation methodologies are fairly new and not completely characterized, it is possible that proofs of correctness will become available for gpi-mml, anm-mml, gpi, and anm-gauss.
all methods implicitly assume that there are no feedback loops.
the noise term in the models is denoted by small "e".
the majority of tested causal orientation methods  output two scores indicating likelihood of the forward causal model  and the backward one .
other tested methods  output two p-values indicating significance of the forward and backward causal models.
in order to make all methods directly comparable to each other, we decided to force them to make causal orientation decisions for all tested causal interactions.
this was achieved by comparing scores or p-values of the forward and backward causal models and selecting the most likely orientation.
this approach follows the practices of previously published applications of causal orientation methods by their inventors  <cit> .
we also note that an alternative approach for the anm, pnl, and gpi methods is to select a model  that is statistically significant at a given alpha level.
the latter approach can sometimes improve accuracy of the causal orientation method while reducing the fraction of causally oriented edges  <cit> .
while results for this approach are not central to this manuscript, we report them in the additional file  <dig>  we also note that the main findings for the primary approach are for the most part consistent with the alternative approach.
finally, prior to application of the causal orientation methods, we standardized the data to mean zero and standard deviation one.
gold standard construction and observational data
the primary challenges in evaluating causal orientation methods for genomics applications are  limited availability of known gold standards of causal molecular interactions and  limited sample sizes of the available observational data.
to overcome these challenges we focused on transcription factor-gene regulatory interactions that have been discovered on the genome wide level and experimentally verified in model organisms  <cit>  and more recently in human cell lines  <cit> .
therefore, the gold standards in this work contain tuples of genes  with orientation x→y, where x is a transcription factor and y is its target gene.
we used the following four gold standards:  interactions of the notch <dig> transcription factor and its target genes in human t-cell acute lymphoblastic leukaemia ;  interactions of the rela transcription factor and its target genes in human t-cell acute lymphoblastic leukaemia ;  interactions of  <dig> transcription factors and their target genes in escherichia coli ; and  interactions of  <dig> transcription factors and their target genes in saccharomyces cerevisiae .
we used microarray gene expression data from the public domain for orientation of causal relations in each gold standard.
the summary statistics of gold standards and corresponding microarray gene expression datasets are given in table  <dig> and table  <dig>  and details of gold standard creation are provided below.
task
ecoli
<dig>
yeast
<dig>
notch1
553
rela
931
"tf" stands for "transcription factor".
statistically significant associations were determined using fisher's z-test at 5% fdr in microarray gene expression data .
ecoli
yeast
notch1
rela
only t-all samples were selected for notch <dig> and rela in order to match cell population used for creation of the respective gold standard.
once each of the gold standards was constructed, we removed interactions without statistically significant associations  according to fisher's z-test  <cit>  at 5% fdr  <cit> .
we performed this filtering because presence of association is a necessary condition for detecting causal relations in most practical settings.
all gold standards and microarray gene expression datasets are available for download from http://www.nyuinformatics.org/downloads/supplements/causalorientation.
creation of notch <dig> and rela gold standards: these gold standards contain genes that are directly downstream of a particular transcription factor  and are functionally regulated by it.
the gold standards were obtained using the method described in  <cit> .
functional gene expression data was first used to identify genes that are downstream  of a particular transcription factor.
the samples in such data are randomized to either 'experiment'  or 'control' treatment.
all genes that are differentially expressed between 'experiment' and 'control' treatments are expected to be downstream of the transcription factor.
we have used a t-test with α =  <dig>  to identify such genes.
genome-wide binding data  was then employed to identify direct binding targets of each transcription factor.
specifically, for each studied transcription factor we obtained the set of genes with detected promoter region-transcription factor binding according to the primary study that generated binding data.
we note that using genome-wide binding data by itself is insufficient to find downstream functional targets of a transcription factor, because many binding sites can be non-functional  <cit> .
therefore, the final step in gold standard creation required overlapping the list of direct binding targets  with the list of downstream functional targets .
knowledge gained by integration of data from these two sources is believed to provide high confidence that a given transcription factor directly regulates a particular gene  <cit> .
also, integration of data from two different sources contributes to the reduction of false positives in the resulting gold standards.
creation of gold standard for yeast and ecoli: these gold standards contain genes that bind to and are likely to be regulated by the known transcription factors in saccharomyces cerevisiae and escherichia coli.
the saccharomyces cerevisiae  gold standard was built by identifying the promoter sequences that are both bound by transcription factors according to chip-on-chip data at  <dig>  alpha level and conserved within  <dig> related species in the saccharomyces genus  <cit> .
binding information is essential because transcription factors must first bind to a gene to induce or suppress expression, while conservation information is important because true-positive transcription factor-dna interactions are often conserved within a genus.
the escherichia coli  gold standard was constructed from regulondb , a manually curated database of regulatory interactions obtained mainly through a literature search  <cit> .
chip-qpcr data has shown regulondb to be approximately 85% complete  <cit> .
evidence for each regulatory interaction in regulondb is classified as "strong" or "weak", depending on the type of experiment used to predict the interaction.
for example, binding of a transcription factor to a promoter is considered strong evidence, whereas gene-expression based computational predictions are considered weak evidence.
for the purposes of our study, we created a conservative gold standard of only strong interactions.
the gold standards yeast and ecoli contain relations of the type "transcription factor → gene" and "transcription factor → transcription factor" .
we decided to simplify the setting of our evaluation when we assess whether the inferred causal orientation x→y or x←y is correct, and restricted attention to only interactions of the type "transcription factor → gene".
this results in minimizing the number of cases with feedback that can be represented by causal edges in both directions.
note that it is not currently possible to comprehensively apply this filtering step to notch <dig> and rela gold standards because the transcription factors are not well known in human cells.
performance metrics and statistical significance testing
two metrics were used to assess performance of causal orientation algorithms.
the first metric is accuracy which is the percentage of causal interactions that have been oriented correctly.
a method that orients all causal interactions in the gold standard as "transcription factor → gene" would achieve an accuracy of 1; a method that orients all interactions as "gene → transcription factor" would achieve an accuracy of 0; and a method that flips a fair coin to make every orientation decision would on average achieve an accuracy of  <dig> .
the second metric is area under roc curve , which is known to be more discriminative than accuracy because it takes into account the confidence of orientation decisions  <cit> .
the roc curve is the plot of sensitivity versus 1-specificity for a range of threshold values on the difference between the scores/p-values of the forward and backward causal models  <cit> .
auc ranges from  <dig> to  <dig>  where auc =  <dig> corresponds to perfectly correct prediction of causal orientation, auc =  <dig>  corresponds to prediction by chance, and auc =  <dig> corresponds to completely incorrect prediction of causal orientation.
computation of sensitivity/specificity and auc requires a response variable with both positive and negative labels.
we created such a response variable by representing gold standard edges  in the following two equivalent ways: 50% of the edges were represented as "transcription factor → gene" and the other 50% were represented as "gene ← transcription factor".
the edges "→" were labeled as positives and "←" were labeled as negatives.
this process is illustrated in table 4; in particular note that the direction of causality always points from transcription factor to gene.
auc was then computed according to the formula given in  <cit> , with the difference in scores/p-values serving as a predictor.
note that the auc can also be interpreted as the probability that the difference between scores/p-values of the forward and backward causal models for a randomly chosen positive instance is higher than the difference between scores/p-values for a randomly chosen negative instance.
since each of the four gold standards has a large number of edges , the variance in auc due to different choices of edges for negative and positive labels is minimal and typically smaller than  <dig>  auc, as estimated by computation of auc for  <dig>  random choices of positive and negative labels.
a fragment of the gold standard is shown in a).
the edges always point from a transcription factor  to its target gene.
50% of the edges are represented as "transcription factor → gene" and the other 50% as "gene ← transcription factor" in b).
this constructs a response variable with positives corresponding to "→" edges  and negatives corresponding to "←" edges .
given values of the above two performance metrics , we need to assess their statistical significance, i.e.
find out if the causal orientation is better than by chance.
notice that our gold standards are such that many causal edges are not independent because they share the same transcription factor.
that is why we chose to apply an exact statistical testing procedure that can accurately estimate the variance of orientation by chance in our setting  <cit> .
a schematic illustration of the statistical testing procedure is given in figure  <dig>  first we compute auc using real gene expression data .
then we replace the real data with random data from the normal distribution with mean  <dig> and standard deviation  <dig>  and compute auc for the same gold standard as used with the real data.
this step is repeated with  <dig>  different randomly generated datasets .
finally, we compare aucs from the null distribution to the auc obtained in the real data and compute a p-value that corresponds to the proportion of random datasets that yield higher aucs than the real data .
a downside of the above statistical testing procedure is that it is computationally expensive and requires running each causal orientation method  <dig> , <dig> times  in order to assess its significance in all  <dig> gold standards used in our study.
to make this analysis feasible, we assessed the statistical significance of only the two best performing methods  and utilized the asclepius compute cluster at the center for health informatics and bioinformatics  at new york university langone medical center .
methodology for sensitivity analyses
in order to study sensitivity to sample size , we sample without replacement from the original gene expression data nested subsets of size  <dig>   <dig>   <dig>  ..., n, where n is the number of samples in the dataset.
specifically, the subset of size  <dig> is included in the subset of size  <dig>  which is in turn included in the subset of size  <dig>  and so on.
we then run the causal orientation algorithms on each subset and compute performance.
this process is repeated with different sampled nested subsets, and mean performance and variance are estimated over all runs.
for the notch <dig> and rela gold standards, we used  <dig> subsets of each size, while for the more computationally expensive yeast and ecoli gold standards we used  <dig> subsets of each size.
for the sensitivity analysis to noise, we add a certain proportion  of random gaussian noise to the gene expression data for both transcription factors and their target genes, run causal orientation methods in the noisy data, report their performance, and repeat the entire process to assess variance .
denoting by x the transcription factor and by y its target gene, the noisy transcription factor x′ and gene y′ are defined as follows: x′ =  · x + p · n and y′ =  · y + p · n, where n is a normally distributed random variable with mean m and standard deviation s, and mx, my and sx, sy are means and standard deviations of x and y in the original data .
we use the following proportions of noise :  <dig> ,  <dig> ,  <dig> , ...,  <dig> ,  <dig> ,  <dig> .
a new ensemble method for causal orientation
as an enhancement to using individual causal orientation techniques, we introduce ensemble causal orientation models that combine decisions of all available individual causal orientation methods in order to produce a more powerful predictor of causal orientation.
these methodologies are popular in the field of supervised learning, where non-random weak learners are often combined to produce a more accurate predictor  <cit> .
the use of ensemble modelling for causal orientation is motivated by our empirical observations that there is no single causal orientation methodology that performs perfectly , many causal orientation methods appear to perform different than chance, and causal orientation methods often make errors in orienting different edges.
in this study we experimented with a straightforward approach to ensemble modelling, where we train a logistic regression model  <cit>  on predictions of all  <dig> tested causal orientation methods .
for the response variable, we follow the same approach as for the computation of auc which was described above in the subsection on performance metrics.
namely, 50% of edges are represented as "transcription factor → gene" and the other 50% as "gene ← transcription factor".
then the response variable is constructed by labelling "→" edges as positives and "←" edges as negatives.
as with every supervised learning procedure that is trained and tested using the same dataset, it is essential to split the available data into non-overlapping training and testing sets, whereby the training set is used to fit a learning model and the testing set is used to estimate its performance  <cit> .
for each gold standard we used 30% of the causal interactions  for training and the remaining 70% for testing.
we decided to use a training set that was smaller than the testing set so that our study resembles a possible practical application where only a small portion of the gold standard is known and the rest is to be discovered.
the predictions of the ensemble model in each gold standard are compared with the predictions of the best individual causal orientation technique in the same testing set with 70% of the data .
finally, in addition to exploring holdout validation performance of the ensemble models, we trained and tested the ensemble models on different gold standards.
in practice this approach can be justified if the data distributions in the gold standards used for training and testing of ensemble models are similar.
it also resembles a practical situation when a gold standard is known in a previously studied dataset but is not known in a new but distributionally similar one.
RESULTS
evaluating causal orientation methods with the accuracy metric
the causal orientation accuracy values are given in table  <dig> for  <dig> causal orientation methods  and  <dig> gold standards.
the performance ranks of methods with accuracies higher than  <dig>  are given in table  <dig>
for each gold standard  dark orange cells correspond to methods that have high values of accuracy, while white cells correspond to methods that have low values of accuracy.
accuracies higher than  <dig>  are shown in bold.
ranks were computed only for the methods with accuracies greater than  <dig> .
the lower the rank, the better the accuracy of the causal orientation method for the given gold standard.
the computation of rank took into account statistical variability, e.g.
accuracies  <dig>  and  <dig>  obtained by the two igci methods in the ecoli gold standard are statistically indistinguishable; that is why they have the same rank value.
as can be seen, igci gaussian/entropy and igci gaussian/integral methods achieve the highest accuracies in each of the four gold standards.
in general, the other causal orientation methods perform worse, and some methods  consistently prefer wrong decisions and have accuracies lower than  <dig> .
interestingly, if we consider the best performing method  with the average rank  <dig> , its results are statistically significant at alpha =  <dig>  according to the exact test  only for the ecoli gold standard .
the second best performing method  with the average rank  <dig>  achieves significance in two out of four gold standards  at alpha =  <dig> .
the reason why the igci gaussian/integral method achieves significance in more gold standards than the best performing technique igci gaussian/entropy is the small variance of the former method.
the detailed statistical significance results including null distributions are given in the additional file  <dig> in figure s <dig>  and figure s <dig> .
it is worth noting that statistical significance was achieved by neither of the two best performing methods in notch <dig> and rela gold standards that have only  <dig> transcription factor.
this was primarily due to a large variance of causal orientation accuracies in the null distribution .
on the other hand, if we join these two gold standards into one with  <dig> transcription factors, both methods achieve statistical significance at alpha =  <dig>  .
the superior and often statistically significant performance of the two igci methods compared to other techniques was a surprising finding that we did not expect theoretically.
igci assumes a noise free model  that is unrealistic in genomics data.
on the other hand, other methods that have a priori more realistic assumptions perform worse.
we hypothesize that sufficient assumptions of the igci methods are too strict in practice and can be mitigated in many ways that are currently not well understood.
evaluating causal orientation methods with the auc metric
the causal orientation auc values are given in table  <dig> for  <dig> causal orientation methods  and  <dig> gold standards.
the performance ranks of methods with aucs higher than  <dig>  are given in table  <dig>
for each gold standard  dark orange cells correspond to methods that have high values of auc, while white cells correspond to methods that have low values of auc.
aucs higher than  <dig>  are shown in bold.
ranks were computed only for the methods with aucs greater than  <dig> .
the lower the rank, the better the auc of the causal orientation method for the given gold standard.
the computation of rank took into account statistical variability, e.g.
the aucs of  <dig>  and  <dig>  obtained by the two igci methods in the ecoli gold standard are statistically indistinguishable; that is why they have the same rank value.
similarly to the accuracy results, igci gaussian/entropy and igci gaussian/integral methods achieve the highest aucs in each of the four gold standards.
other causal orientation methods perform worse, and some methods  consistently prefer wrong decisions and have aucs lower than  <dig> .
the statistical significance analysis of igci gaussian/entropy and igci gaussian/integral is described in detail in the additional file  <dig> in figure s <dig>  and figure s <dig> .
in summary, both methods achieve statistical significance of causal orientations  in ecoli and yeast, only igci gaussian/integral achieves significance in rela, and none of the two methods achieves significance in notch <dig>  on the other hand, similarly to results for the accuracy metric, both methods achieve statistically significant results in the joined notch <dig> and rela gold standard with  <dig> transcription factors.
sensitivity analysis to noise
the results of sensitivity analysis to noise for the two best performing methods  are given in figures  <dig>   <dig>   <dig>   <dig> for notch <dig>  rela, ecoli, and yeast gold standards, respectively.
in all gold standards except for yeast, the accuracy of the methods decreases with increasing noise proportion.
on the other hand, in yeast gold standard the performance of causal orientation methods significantly increases when a small amount of noise is added, and then gradually decreases for higher proportions of noise.
the additional file  <dig> provides a detailed analysis of this phenomenon.
whereas in notch <dig> and rela gold standards it takes only 5-10% of noise to make the results statistically indistinguishable from orientation by chance, in yeast and ecoli gold standards the methods can tolerate much higher proportions of noise and still produce statistically significant results.
this can be attributed to a larger number of transcription factors in yeast and ecoli gold standards, as well as larger sample sizes in the corresponding datasets which both decrease the variability of the results.
a decrease in performance upon the addition of noise is theoretically expected since igci assumes a noise-free model, and the addition of gaussian noise violates its sufficient assumptions.
also, as can be seen in the figures, the igci gaussian/integral method has lower variance than the igci gaussian/entropy method.
the above results are consistent with our prior findings and statistical significance testing by the exact test .
sensitivity analysis to sample size
the results of sensitivity analysis to sample size for the two best performing methods  are given in figures  <dig>   <dig>   <dig>   <dig> for notch <dig>  rela, ecoli, and yeast gold standards, respectively.
in all gold standards except for yeast, the accuracy of the methods decreases as the sample size gets smaller.
on the other hand, in yeast gold standard the performance of causal orientation methods slightly increases by reducing the sample size and then gradually decreases for smaller sample sizes.
the additional file  <dig> provides a detailed analysis of this phenomenon.
whereas in notch <dig> and rela gold standards results become statistically indistinguishable from orientation by chance when the sample size is <80- <dig>  in yeast and ecoli gold standards the methods yield statistically significant results for smaller sample sizes.
this can be attributed to a larger number of transcription factors in yeast and ecoli gold standards which decreases variability of the results.
also, as can be seen in the figures, the igci gaussian/integral method has lower variance than the igci gaussian/entropy method.
the above results are consistent with our prior findings in statistical significance testing by the exact test .
ensemble causal orientation
for each gold standard, table  <dig> compares the auc achieved by the best individual causal orientation method to the auc achieved by the ensemble method, which combines the predictions of all  <dig> methods using a logistic regression model.
a detailed description of the ensemble modeling methodology is given in the methods section.
as can be seen, ensemble causal orientation achieves higher values of auc than any individual causal orientation method in all four gold standards.
it is worthwhile to highlight the magnitude of the improvement in the yeast gold standard: the ensemble approach improves performance over the best individual causal orientation method  by  <dig>  auc.
table  <dig> provides coefficients for the ensemble logistic regression model in the yeast gold standard.
bold values correspond to coefficients that are statistically significant at  <dig>  alpha level.
the model preserves its performance  when it is trained/tested using only the  <dig> causal orientation methods that have statistically significant coefficients.
therefore the improvement in auc in the yeast gold standard can be attributed to effectively combining the igci and anm-mml causal orientation predictions by the logistic regression ensemble model.
bold p-values indicate a statistically significant performance improvement by using an ensemble causal orientation.
the p-values were obtained from delong's test for auc comparison  <cit> .
bold values correspond to coefficients that are statistically significant at  <dig>  alpha level.
we note that due to multicollinearity among the igci uniform methods and among the igci gaussian methods, care must be taken when interpreting the logistic regression coefficients  <cit> .
the above results were obtained by holdout validation where we used different portions of the same gold standard for training and testing ensemble models.
we also experimented with training and testing ensemble models on different gold standards.
first we experimented with the rela and notch <dig> gold standards that were derived from the same organism and phenotype, and thus are likely to be distributionally similar and support cross-gold standard application of the ensemble model.
we find that an ensemble logistic regression model trained on rela obtains auc =  <dig>  when tested on notch <dig>  and likewise an ensemble model trained on notch <dig> obtains auc =  <dig>  when tested on rela.
both these results significantly improve performance over the best individual causal orientation method  in both notch <dig> and rela gold standards .
in addition, we experimented with the yeast and ecoli gold standards which originate from different organisms and thus are unlikely to be distributionally similar; for this reason they a priori question cross-gold standard application of the ensemble model.
indeed, our results confirm this expectation: an ensemble logistic regression model trained on yeast performs with auc =  <dig>  in ecoli, and an ensemble model trained in ecoli performs with auc =  <dig>  in yeast.
neither of these results improves the best individual causal orientation method in the respective gold standard.
this suggests that the success of cross-gold standard application of ensemble models is grounded on similarity of the underlying distributions.
discussion
this work represents the first comprehensive effort to evaluate performance and furthermore enhance the recently introduced causal orientation methods  <cit>  in genomics data.
one of the main challenges is the limited availability of gold standards of causal molecular interactions.
that is why we have focused on regulatory interactions between transcriptions factors and their target genes that have been recently identified on a genome-wide level in model organisms and human cell lines.
these interactions have a well-defined causal directionality  and can be readily used for an evaluation study such as ours.
however, it is possible that some edges in the gold standards  have causal relationships in both directions due to feedback mechanisms.
since the signal in the direction from a transcription factor to its target gene is expected to be stronger than in the opposite direction , we are implicitly assuming that in such cases a causal orientation method would prefer the direction from transcription factor to gene.
if this assumption is not true, this does not invalidate results of the methods  but provides additional explanations as to why some methods prefer the opposite causal direction.
even though the choice of the gold standard with transcription factor-gene regulatory interactions enables this study, its practical relevance may be limited in the organisms/settings where all transcription factors have already been identified.
that is why we plan to work on extending this evaluation to other types of causal molecular interactions, for example in cellular protein signaling networks  <cit> .
in this study we have implicitly assumed that unoriented edges  are given by an oracle and we have evaluated performance of only causal orientation methods.
however, in practical tasks one typically has to both discover and orient edges.
although we have previously evaluated methods for discovery of unoriented edges  <cit> , it will be interesting to assess the performance of the two classes of methods  when they are applied together.
finally, we think that a fruitful area of research will be to extend this study by comparison with classical causal orientation techniques that output markov equivalence classes of graphs  and thus, in general, can orient only a subset of edges in the graph  <cit> .
CONCLUSIONS
in this paper we have taken a first step toward practical use of recent causal orientation techniques in the genomics domain.
first of all, we report results of an extensive study of causal orientation methods in genomics data that utilized  <dig> methods/variants to distinguish cause  from effect  in  <dig>  causal interactions.
we have found that igci gaussian methods  <cit>  often accurately infer directionality of the causal interaction, and they outperform other causal orientation techniques.
in addition, we have performed sensitivity analyses that allow us to empirically establish the minimal requirements for the sample size and maximal level of noise that can be tolerated by the best performing causal orientation techniques.
second, we described a novel ensemble technique for causal orientation that combines decisions of individual causal orientation methods to provide a more powerful predictor of causal directionality.
the proposed ensemble method was found to be more accurate than any best individual causal orientation method in the tested data.
in summary, our results suggest that causal orientation methods have significant potential to facilitate reconstruction of molecular pathways by minimizing the number of required randomized experiments to find causal directionality and by avoiding experiments that are infeasible and/or unethical.
competing interests
the authors declare that they have no competing interests.
authors' contributions
as conceived the study.
all authors participated in the design of the study.
mh and nil performed all computational experiments.
as and mh drafted the manuscript.
all authors edited, read, and approved the final manuscript.
supplementary material
additional file 1
this file contains  brief description of causal orientation algorithms;  results of causal orientation methods anm, pnl, and gpi obtained by assessing statistical significance of the forward and backward causal models;  detailed results of significance testing of igci gaussian/entropy and gaussian/integral methods;  explanation of performance increase due to adding small amount of noise or reducing the sample size in yeast gold standard.
click here for file
acknowledgements
the authors would like to acknowledge dominik janzing and joris m. mooij, who contributed to the development of the majority of causal orientation methods used in this study, and thank them for providing  software implementations of causal orientation algorithms,  help with stating assumptions of the tested methods,  ideas about statistical significance testing approach, and  feedback on other aspects of the manuscript and, in particular, interpretation of the results.
the authors are also grateful to efstratios efstathiadis and eric peskin for the help with providing access and running experiments on the high performance computing facility.
finally, the authors would like to thank ioannis aifantis for providing experimental data for notch <dig> that was used for the development of the corresponding gold standard.
the empirical evaluation was supported in part by the grants 1r01lm011179-01a <dig> from the national library of medicine and 1ul1rr <dig> from the national center for research resources, national institutes of health.
this article has been published as part of bmc genomics volume  <dig> supplement  <dig>  2012: proceedings of the international conference on intelligent biology and medicine : genomics.
the full contents of the supplement are available online at http://www.biomedcentral.com/bmcgenomics/supplements/13/s <dig>

 the discovery of molecular pathways that drive diseases and vital cellular functions is a fundamental activity of biomedical research. unraveling disease pathways is a major component in the efforts to develop new therapies that will effectively fight deadly diseases. furthermore, knowing pathways significantly facilitates the design of personalized medicine modalities for diagnosis, prognosis, and management of diseases. the discovery of pathways is a challenging problem and its solution to a large extent relies on the identification of causal molecular interactions in genomics data. by causal molecular interactions or relations we mean interactions of molecular variables that match the notion of randomized controlled experiment, which is the de facto standard for assessing causation in the general sciences and biomedicine. assume that a hypothetical experimenter can change the distribution of a variable x. we say that x is a cause of y and denote this by x→y if the probability distribution of y changes for some experimental manipulation of x. causal molecular interactions can be discovered using randomized experiments such as interference with rna; however such experiments are often costly, infeasible, or unethical. fortunately, over the last years many algorithms that infer causal interactions from observational data have been developed and some of them have been adopted to the high dimensionalities of modern genomics data. outside of biomedicine, two nobel prizes have recently been awarded in and for methods which seek to discover causal relations from non-experimental data. in our prior work we evaluated the ability of state-of-the-art causal discovery algorithms to de-novo identify unoriented edges in genome-scale regulatory networks, which represent causal interactions between transcription factors and their target genes without distinguishing the mechanistic role of the involved molecular variables. we deliberately avoided performing causal orientation of the discovered unoriented edges because this problem has previously been deemed worst-case unsolvable in observational data using existing algorithms due to the statistical indistinguishability of causal models in the same markov equivalence class. for example, causal models x→y and x←y have generally been assumed to be statistically indistinguishable given only observational data for x and y. over the last years researchers have developed a new class of methods that, given observational data for just two causally related variables x and y, aim to determine which variable is the cause and which one is the effect. these causal orientation methods aim to solve the problem of statistical indistinguishability of graphs in the markov equivalence class by exploiting asymmetries in the shapes of the conditional probability densities and without requiring randomized experiments. these methods could have significant implications for the field of causal discovery because they can orient unoriented edges that are discoverable by other established techniques, e.g. generalized local learning or local-to-global learning. therefore, at face value, these causal orientation methods have the potential to reduce the number of, and in some cases even eliminate, randomized experiments needed for causal orientation of edges in the markov equivalence class of graphs and make the causal model fully identifiable from observational data alone. as promising as these new causal orientation methods are, they have not been previously applied in genomics, where the data is usually noisy and the sample sizes are relatively small compared to prior test applications of these methods. in this paper we report results of an extensive study of recent causal orientation techniques in the genomics domain by testing methods/variants to distinguish cause from effect in causal interactions and conducting sensitivity analyses with respect to noise and sample size for the best-performing methods. in addition, we introduce a new ensemble technique for causal orientation that is shown to be more accurate than any best individual causal orientation method in the tested data. the results of this study can serve as a foundation for further development of causal orientation techniques for genomics data and establishing a context for wide applications in molecular or biomedical research. methods causal orientation methods as mentioned above, the purpose of the tested causal orientation methods is to separate cause from effect given data for just two variables x and y that have a causal relation. typically these methods are not designed to be used to causally orient pairs of variables that only have univariate association/correlation due to the possibility of confounding. for example, in the majority of distributions, a causal structure x ← t → y implies that x and y are associated even though they are not causally affecting each other. therefore, the presence of association in a pair of variables x and y is, in general, a necessary but not sufficient condition to be eligible for causal orientation. a rigorous approach would involve first using correct causal discovery methods to identify unoriented edges that denote the existence of unconfounded causal relations and then applying causal orientation techniques to orient unoriented edges. while each causal orientation method has its own principles and sufficient assumptions that are outlined in table, most of these techniques are based on the idea that the factorization of the joint probability distribution p into pp yields a simpler representation than the factorization into pp. one can furthermore show that if the marginal probability distribution of the cause p is independent of the causal mechanism p, then the factorization pp has lower complexity than the factorization pp. given two causally related variables x and y, estimating the complexity of the two different factorizations of p or determining independence between marginal and conditional distributions can thus provide the basis for causal orientation techniques. in practice, however, it is difficult to directly test independence between p and p or estimate their complexity; hence the methods typically use simplifying assumptions or rely on approximate formulations. • x and e are independent; • f is non-linear, or one of x and e is non-gaussian; • probability densities are strictly positive; • x and e are independent; • either f or e is gaussian; • f is continuous and invertible; • x and e are independent; • e is gaussian; • x and e are independent; • e is gaussian. • x and e are independent; the last column indicates whether a method is sound, i.e. it can provably orient a causal structure under its sufficient assumptions. because causal orientation methodologies are fairly new and not completely characterized, it is possible that proofs of correctness will become available for gpi-mml, anm-mml, gpi, and anm-gauss. all methods implicitly assume that there are no feedback loops. the noise term in the models is denoted by small "e". the majority of tested causal orientation methods output two scores indicating likelihood of the forward causal model and the backward one. other tested methods output two p-values indicating significance of the forward and backward causal models. in order to make all methods directly comparable to each other, we decided to force them to make causal orientation decisions for all tested causal interactions. this was achieved by comparing scores or p-values of the forward and backward causal models and selecting the most likely orientation. this approach follows the practices of previously published applications of causal orientation methods by their inventors. we also note that an alternative approach for the anm, pnl, and gpi methods is to select a model that is statistically significant at a given alpha level. the latter approach can sometimes improve accuracy of the causal orientation method while reducing the fraction of causally oriented edges. while results for this approach are not central to this manuscript, we report them in the additional file we also note that the main findings for the primary approach are for the most part consistent with the alternative approach. finally, prior to application of the causal orientation methods, we standardized the data to mean zero and standard deviation one. gold standard construction and observational data the primary challenges in evaluating causal orientation methods for genomics applications are limited availability of known gold standards of causal molecular interactions and limited sample sizes of the available observational data. to overcome these challenges we focused on transcription factor-gene regulatory interactions that have been discovered on the genome wide level and experimentally verified in model organisms and more recently in human cell lines. therefore, the gold standards in this work contain tuples of genes with orientation x→y, where x is a transcription factor and y is its target gene. we used the following four gold standards: interactions of the notch transcription factor and its target genes in human t-cell acute lymphoblastic leukaemia; interactions of the rela transcription factor and its target genes in human t-cell acute lymphoblastic leukaemia; interactions of transcription factors and their target genes in escherichia coli; and interactions of transcription factors and their target genes in saccharomyces cerevisiae. we used microarray gene expression data from the public domain for orientation of causal relations in each gold standard. the summary statistics of gold standards and corresponding microarray gene expression datasets are given in table and table and details of gold standard creation are provided below. task ecoli yeast notch1 553 rela 931 "tf" stands for "transcription factor". statistically significant associations were determined using fisher's z-test at 5% fdr in microarray gene expression data. ecoli yeast notch1 rela only t-all samples were selected for notch and rela in order to match cell population used for creation of the respective gold standard. once each of the gold standards was constructed, we removed interactions without statistically significant associations according to fisher's z-test at 5% fdr. we performed this filtering because presence of association is a necessary condition for detecting causal relations in most practical settings. all gold standards and microarray gene expression datasets are available for download from http://www.nyuinformatics.org/downloads/supplements/causalorientation. creation of notch and rela gold standards: these gold standards contain genes that are directly downstream of a particular transcription factor and are functionally regulated by it. the gold standards were obtained using the method described in. functional gene expression data was first used to identify genes that are downstream of a particular transcription factor. the samples in such data are randomized to either'experiment' or'control' treatment. all genes that are differentially expressed between'experiment' and'control' treatments are expected to be downstream of the transcription factor. we have used a t-test with α = to identify such genes. genome-wide binding data was then employed to identify direct binding targets of each transcription factor. specifically, for each studied transcription factor we obtained the set of genes with detected promoter region-transcription factor binding according to the primary study that generated binding data. we note that using genome-wide binding data by itself is insufficient to find downstream functional targets of a transcription factor, because many binding sites can be non-functional. therefore, the final step in gold standard creation required overlapping the list of direct binding targets with the list of downstream functional targets. knowledge gained by integration of data from these two sources is believed to provide high confidence that a given transcription factor directly regulates a particular gene. also, integration of data from two different sources contributes to the reduction of false positives in the resulting gold standards. creation of gold standard for yeast and ecoli: these gold standards contain genes that bind to and are likely to be regulated by the known transcription factors in saccharomyces cerevisiae and escherichia coli. the saccharomyces cerevisiae gold standard was built by identifying the promoter sequences that are both bound by transcription factors according to chip-on-chip data at alpha level and conserved within related species in the saccharomyces genus. binding information is essential because transcription factors must first bind to a gene to induce or suppress expression, while conservation information is important because true-positive transcription factor-dna interactions are often conserved within a genus. the escherichia coli gold standard was constructed from regulondb, a manually curated database of regulatory interactions obtained mainly through a literature search. chip-qpcr data has shown regulondb to be approximately 85% complete. evidence for each regulatory interaction in regulondb is classified as "strong" or "weak", depending on the type of experiment used to predict the interaction. for example, binding of a transcription factor to a promoter is considered strong evidence, whereas gene-expression based computational predictions are considered weak evidence. for the purposes of our study, we created a conservative gold standard of only strong interactions. the gold standards yeast and ecoli contain relations of the type "transcription factor → gene" and "transcription factor → transcription factor". we decided to simplify the setting of our evaluation when we assess whether the inferred causal orientation x→y or x←y is correct, and restricted attention to only interactions of the type "transcription factor → gene". this results in minimizing the number of cases with feedback that can be represented by causal edges in both directions. note that it is not currently possible to comprehensively apply this filtering step to notch and rela gold standards because the transcription factors are not well known in human cells. performance metrics and statistical significance testing two metrics were used to assess performance of causal orientation algorithms. the first metric is accuracy which is the percentage of causal interactions that have been oriented correctly. a method that orients all causal interactions in the gold standard as "transcription factor → gene" would achieve an accuracy of 1; a method that orients all interactions as "gene → transcription factor" would achieve an accuracy of 0; and a method that flips a fair coin to make every orientation decision would on average achieve an accuracy of. the second metric is area under roc curve, which is known to be more discriminative than accuracy because it takes into account the confidence of orientation decisions. the roc curve is the plot of sensitivity versus 1-specificity for a range of threshold values on the difference between the scores/p-values of the forward and backward causal models. auc ranges from to where auc = corresponds to perfectly correct prediction of causal orientation, auc = corresponds to prediction by chance, and auc = corresponds to completely incorrect prediction of causal orientation. computation of sensitivity/specificity and auc requires a response variable with both positive and negative labels. we created such a response variable by representing gold standard edges in the following two equivalent ways: 50% of the edges were represented as "transcription factor → gene" and the other 50% were represented as "gene ← transcription factor". the edges "→" were labeled as positives and "←" were labeled as negatives. this process is illustrated in table 4; in particular note that the direction of causality always points from transcription factor to gene. auc was then computed according to the formula given in, with the difference in scores/p-values serving as a predictor. note that the auc can also be interpreted as the probability that the difference between scores/p-values of the forward and backward causal models for a randomly chosen positive instance is higher than the difference between scores/p-values for a randomly chosen negative instance. since each of the four gold standards has a large number of edges, the variance in auc due to different choices of edges for negative and positive labels is minimal and typically smaller than auc, as estimated by computation of auc for random choices of positive and negative labels. a fragment of the gold standard is shown in a). the edges always point from a transcription factor to its target gene. 50% of the edges are represented as "transcription factor → gene" and the other 50% as "gene ← transcription factor" in b). this constructs a response variable with positives corresponding to "→" edges and negatives corresponding to "←" edges. given values of the above two performance metrics, we need to assess their statistical significance, i.e. find out if the causal orientation is better than by chance. notice that our gold standards are such that many causal edges are not independent because they share the same transcription factor. that is why we chose to apply an exact statistical testing procedure that can accurately estimate the variance of orientation by chance in our setting. a schematic illustration of the statistical testing procedure is given in figure first we compute auc using real gene expression data. then we replace the real data with random data from the normal distribution with mean and standard deviation and compute auc for the same gold standard as used with the real data. this step is repeated with different randomly generated datasets. finally, we compare aucs from the null distribution to the auc obtained in the real data and compute a p-value that corresponds to the proportion of random datasets that yield higher aucs than the real data. a downside of the above statistical testing procedure is that it is computationally expensive and requires running each causal orientation method, times in order to assess its significance in all gold standards used in our study. to make this analysis feasible, we assessed the statistical significance of only the two best performing methods and utilized the asclepius compute cluster at the center for health informatics and bioinformatics at new york university langone medical center. methodology for sensitivity analyses in order to study sensitivity to sample size, we sample without replacement from the original gene expression data nested subsets of size..., n, where n is the number of samples in the dataset. specifically, the subset of size is included in the subset of size which is in turn included in the subset of size and so on. we then run the causal orientation algorithms on each subset and compute performance. this process is repeated with different sampled nested subsets, and mean performance and variance are estimated over all runs. for the notch and rela gold standards, we used subsets of each size, while for the more computationally expensive yeast and ecoli gold standards we used subsets of each size. for the sensitivity analysis to noise, we add a certain proportion of random gaussian noise to the gene expression data for both transcription factors and their target genes, run causal orientation methods in the noisy data, report their performance, and repeat the entire process to assess variance. denoting by x the transcription factor and by y its target gene, the noisy transcription factor x′ and gene y′ are defined as follows: x′ = · x + p · n and y′ = · y + p · n, where n is a normally distributed random variable with mean m and standard deviation s, and mx, my and sx, sy are means and standard deviations of x and y in the original data. we use the following proportions of noise:,,,...,,,. a new ensemble method for causal orientation as an enhancement to using individual causal orientation techniques, we introduce ensemble causal orientation models that combine decisions of all available individual causal orientation methods in order to produce a more powerful predictor of causal orientation. these methodologies are popular in the field of supervised learning, where non-random weak learners are often combined to produce a more accurate predictor. the use of ensemble modelling for causal orientation is motivated by our empirical observations that there is no single causal orientation methodology that performs perfectly, many causal orientation methods appear to perform different than chance, and causal orientation methods often make errors in orienting different edges. in this study we experimented with a straightforward approach to ensemble modelling, where we train a logistic regression model on predictions of all tested causal orientation methods. for the response variable, we follow the same approach as for the computation of auc which was described above in the subsection on performance metrics. namely, 50% of edges are represented as "transcription factor → gene" and the other 50% as "gene ← transcription factor". then the response variable is constructed by labelling "→" edges as positives and "←" edges as negatives. as with every supervised learning procedure that is trained and tested using the same dataset, it is essential to split the available data into non-overlapping training and testing sets, whereby the training set is used to fit a learning model and the testing set is used to estimate its performance. for each gold standard we used 30% of the causal interactions for training and the remaining 70% for testing. we decided to use a training set that was smaller than the testing set so that our study resembles a possible practical application where only a small portion of the gold standard is known and the rest is to be discovered. the predictions of the ensemble model in each gold standard are compared with the predictions of the best individual causal orientation technique in the same testing set with 70% of the data. finally, in addition to exploring holdout validation performance of the ensemble models, we trained and tested the ensemble models on different gold standards. in practice this approach can be justified if the data distributions in the gold standards used for training and testing of ensemble models are similar. it also resembles a practical situation when a gold standard is known in a previously studied dataset but is not known in a new but distributionally similar one.  evaluating causal orientation methods with the accuracy metric the causal orientation accuracy values are given in table for causal orientation methods and gold standards. the performance ranks of methods with accuracies higher than are given in table for each gold standard dark orange cells correspond to methods that have high values of accuracy, while white cells correspond to methods that have low values of accuracy. accuracies higher than are shown in bold. ranks were computed only for the methods with accuracies greater than. the lower the rank, the better the accuracy of the causal orientation method for the given gold standard. the computation of rank took into account statistical variability, e.g. accuracies and obtained by the two igci methods in the ecoli gold standard are statistically indistinguishable; that is why they have the same rank value. as can be seen, igci gaussian/entropy and igci gaussian/integral methods achieve the highest accuracies in each of the four gold standards. in general, the other causal orientation methods perform worse, and some methods consistently prefer wrong decisions and have accuracies lower than. interestingly, if we consider the best performing method with the average rank, its results are statistically significant at alpha = according to the exact test only for the ecoli gold standard. the second best performing method with the average rank achieves significance in two out of four gold standards at alpha =. the reason why the igci gaussian/integral method achieves significance in more gold standards than the best performing technique igci gaussian/entropy is the small variance of the former method. the detailed statistical significance results including null distributions are given in the additional file in figure s and figure s. it is worth noting that statistical significance was achieved by neither of the two best performing methods in notch and rela gold standards that have only transcription factor. this was primarily due to a large variance of causal orientation accuracies in the null distribution. on the other hand, if we join these two gold standards into one with transcription factors, both methods achieve statistical significance at alpha =. the superior and often statistically significant performance of the two igci methods compared to other techniques was a surprising finding that we did not expect theoretically. igci assumes a noise free model that is unrealistic in genomics data. on the other hand, other methods that have a priori more realistic assumptions perform worse. we hypothesize that sufficient assumptions of the igci methods are too strict in practice and can be mitigated in many ways that are currently not well understood. evaluating causal orientation methods with the auc metric the causal orientation auc values are given in table for causal orientation methods and gold standards. the performance ranks of methods with aucs higher than are given in table for each gold standard dark orange cells correspond to methods that have high values of auc, while white cells correspond to methods that have low values of auc. aucs higher than are shown in bold. ranks were computed only for the methods with aucs greater than. the lower the rank, the better the auc of the causal orientation method for the given gold standard. the computation of rank took into account statistical variability, e.g. the aucs of and obtained by the two igci methods in the ecoli gold standard are statistically indistinguishable; that is why they have the same rank value. similarly to the accuracy results, igci gaussian/entropy and igci gaussian/integral methods achieve the highest aucs in each of the four gold standards. other causal orientation methods perform worse, and some methods consistently prefer wrong decisions and have aucs lower than. the statistical significance analysis of igci gaussian/entropy and igci gaussian/integral is described in detail in the additional file in figure s and figure s. in summary, both methods achieve statistical significance of causal orientations in ecoli and yeast, only igci gaussian/integral achieves significance in rela, and none of the two methods achieves significance in notch on the other hand, similarly to results for the accuracy metric, both methods achieve statistically significant results in the joined notch and rela gold standard with transcription factors. sensitivity analysis to noise the results of sensitivity analysis to noise for the two best performing methods are given in figures for notch rela, ecoli, and yeast gold standards, respectively. in all gold standards except for yeast, the accuracy of the methods decreases with increasing noise proportion. on the other hand, in yeast gold standard the performance of causal orientation methods significantly increases when a small amount of noise is added, and then gradually decreases for higher proportions of noise. the additional file provides a detailed analysis of this phenomenon. whereas in notch and rela gold standards it takes only 5-10% of noise to make the results statistically indistinguishable from orientation by chance, in yeast and ecoli gold standards the methods can tolerate much higher proportions of noise and still produce statistically significant results. this can be attributed to a larger number of transcription factors in yeast and ecoli gold standards, as well as larger sample sizes in the corresponding datasets which both decrease the variability of the results. a decrease in performance upon the addition of noise is theoretically expected since igci assumes a noise-free model, and the addition of gaussian noise violates its sufficient assumptions. also, as can be seen in the figures, the igci gaussian/integral method has lower variance than the igci gaussian/entropy method. the above results are consistent with our prior findings and statistical significance testing by the exact test. sensitivity analysis to sample size the results of sensitivity analysis to sample size for the two best performing methods are given in figures for notch rela, ecoli, and yeast gold standards, respectively. in all gold standards except for yeast, the accuracy of the methods decreases as the sample size gets smaller. on the other hand, in yeast gold standard the performance of causal orientation methods slightly increases by reducing the sample size and then gradually decreases for smaller sample sizes. the additional file provides a detailed analysis of this phenomenon. whereas in notch and rela gold standards results become statistically indistinguishable from orientation by chance when the sample size is <80- in yeast and ecoli gold standards the methods yield statistically significant results for smaller sample sizes. this can be attributed to a larger number of transcription factors in yeast and ecoli gold standards which decreases variability of the results. also, as can be seen in the figures, the igci gaussian/integral method has lower variance than the igci gaussian/entropy method. the above results are consistent with our prior findings in statistical significance testing by the exact test. ensemble causal orientation for each gold standard, table compares the auc achieved by the best individual causal orientation method to the auc achieved by the ensemble method, which combines the predictions of all methods using a logistic regression model. a detailed description of the ensemble modeling methodology is given in the methods section. as can be seen, ensemble causal orientation achieves higher values of auc than any individual causal orientation method in all four gold standards. it is worthwhile to highlight the magnitude of the improvement in the yeast gold standard: the ensemble approach improves performance over the best individual causal orientation method by auc. table provides coefficients for the ensemble logistic regression model in the yeast gold standard. bold values correspond to coefficients that are statistically significant at alpha level. the model preserves its performance when it is trained/tested using only the causal orientation methods that have statistically significant coefficients. therefore the improvement in auc in the yeast gold standard can be attributed to effectively combining the igci and anm-mml causal orientation predictions by the logistic regression ensemble model. bold p-values indicate a statistically significant performance improvement by using an ensemble causal orientation. the p-values were obtained from delong's test for auc comparison. bold values correspond to coefficients that are statistically significant at alpha level. we note that due to multicollinearity among the igci uniform methods and among the igci gaussian methods, care must be taken when interpreting the logistic regression coefficients. the above results were obtained by holdout validation where we used different portions of the same gold standard for training and testing ensemble models. we also experimented with training and testing ensemble models on different gold standards. first we experimented with the rela and notch gold standards that were derived from the same organism and phenotype, and thus are likely to be distributionally similar and support cross-gold standard application of the ensemble model. we find that an ensemble logistic regression model trained on rela obtains auc = when tested on notch and likewise an ensemble model trained on notch obtains auc = when tested on rela. both these results significantly improve performance over the best individual causal orientation method in both notch and rela gold standards. in addition, we experimented with the yeast and ecoli gold standards which originate from different organisms and thus are unlikely to be distributionally similar; for this reason they a priori question cross-gold standard application of the ensemble model. indeed, our results confirm this expectation: an ensemble logistic regression model trained on yeast performs with auc = in ecoli, and an ensemble model trained in ecoli performs with auc = in yeast. neither of these results improves the best individual causal orientation method in the respective gold standard. this suggests that the success of cross-gold standard application of ensemble models is grounded on similarity of the underlying distributions. discussion this work represents the first comprehensive effort to evaluate performance and furthermore enhance the recently introduced causal orientation methods in genomics data. one of the main challenges is the limited availability of gold standards of causal molecular interactions. that is why we have focused on regulatory interactions between transcriptions factors and their target genes that have been recently identified on a genome-wide level in model organisms and human cell lines. these interactions have a well-defined causal directionality and can be readily used for an evaluation study such as ours. however, it is possible that some edges in the gold standards have causal relationships in both directions due to feedback mechanisms. since the signal in the direction from a transcription factor to its target gene is expected to be stronger than in the opposite direction, we are implicitly assuming that in such cases a causal orientation method would prefer the direction from transcription factor to gene. if this assumption is not true, this does not invalidate results of the methods but provides additional explanations as to why some methods prefer the opposite causal direction. even though the choice of the gold standard with transcription factor-gene regulatory interactions enables this study, its practical relevance may be limited in the organisms/settings where all transcription factors have already been identified. that is why we plan to work on extending this evaluation to other types of causal molecular interactions, for example in cellular protein signaling networks. in this study we have implicitly assumed that unoriented edges are given by an oracle and we have evaluated performance of only causal orientation methods. however, in practical tasks one typically has to both discover and orient edges. although we have previously evaluated methods for discovery of unoriented edges, it will be interesting to assess the performance of the two classes of methods when they are applied together. finally, we think that a fruitful area of research will be to extend this study by comparison with classical causal orientation techniques that output markov equivalence classes of graphs and thus, in general, can orient only a subset of edges in the graph.  in this paper we have taken a first step toward practical use of recent causal orientation techniques in the genomics domain. first of all, we report results of an extensive study of causal orientation methods in genomics data that utilized methods/variants to distinguish cause from effect in causal interactions. we have found that igci gaussian methods often accurately infer directionality of the causal interaction, and they outperform other causal orientation techniques. in addition, we have performed sensitivity analyses that allow us to empirically establish the minimal requirements for the sample size and maximal level of noise that can be tolerated by the best performing causal orientation techniques. second, we described a novel ensemble technique for causal orientation that combines decisions of individual causal orientation methods to provide a more powerful predictor of causal directionality. the proposed ensemble method was found to be more accurate than any best individual causal orientation method in the tested data. in summary, our results suggest that causal orientation methods have significant potential to facilitate reconstruction of molecular pathways by minimizing the number of required randomized experiments to find causal directionality and by avoiding experiments that are infeasible and/or unethical. competing interests the authors declare that they have no competing interests. authors' contributions as conceived the study. all authors participated in the design of the study. mh and nil performed all computational experiments. as and mh drafted the manuscript. all authors edited, read, and approved the final manuscript. supplementary material additional file 1 this file contains brief description of causal orientation algorithms; results of causal orientation methods anm, pnl, and gpi obtained by assessing statistical significance of the forward and backward causal models; detailed results of significance testing of igci gaussian/entropy and gaussian/integral methods; explanation of performance increase due to adding small amount of noise or reducing the sample size in yeast gold standard. click here for file acknowledgements the authors would like to acknowledge dominik janzing and joris m. mooij, who contributed to the development of the majority of causal orientation methods used in this study, and thank them for providing software implementations of causal orientation algorithms, help with stating assumptions of the tested methods, ideas about statistical significance testing approach, and feedback on other aspects of the manuscript and, in particular, interpretation of the results. the authors are also grateful to efstratios efstathiadis and eric peskin for the help with providing access and running experiments on the high performance computing facility. finally, the authors would like to thank ioannis aifantis for providing experimental data for notch that was used for the development of the corresponding gold standard. the empirical evaluation was supported in part by the grants 1r01lm011179-01a from the national library of medicine and 1ul1rr from the national center for research resources, national institutes of health. this article has been published as part of bmc genomics volume supplement 2012: proceedings of the international conference on intelligent biology and medicine: genomics. the full contents of the supplement are available online at http://www.biomedcentral.com/bmcgenomics/supplements/13/s 
the discovery of molecular pathways is a challenging problem and its solution relies on the identification of causal molecular interactions in genomics data. our results suggest that these methods have the potential to facilitate reconstruction of molecular pathways by minimizing the number of required randomized experiments to find causal directionality and by avoiding experiments that are infeasible and/or unethical. we have found that a particular family of causal orientation methods is often able to accurately infer directionality of causal interactions, and that these methods usually outperform other causal orientation techniques. fortunately, algorithms that infer causal interactions from observational data have been in development for decades, predominantly in the quantitative sciences, and many of them have recently been applied to genomics data. while these algorithms can infer unoriented causal interactions between involved molecular variables, causally orienting all inferred molecular interactions was assumed to be an unsolvable problem until recently. we also introduced a novel ensemble technique for causal orientation that combines decisions of individual causal orientation methods. the ensemble method was found to be more accurate than any best individual causal orientation method in the tested data. this work represents a first step towards establishing context for practical use of causal orientation methods in the genomics domain. we have found that some causal orientation methodologies yield accurate predictions of causal orientation in genomics data, and we have improved on this capability with a novel ensemble method. causal molecular interactions can be discovered using randomized experiments; however such experiments are often costly, infeasible, or unethical. in this work, we use transcription factor-target gene regulatory interactions in three different organisms to evaluate a new family of methods that, given observational data for just two causally related variables, can determine which one is the cause and which one is the effect. 22- april the international conference on intelligent biology and medicine nashville, tn, usa
 the discovery of molecular pathways is a challenging problem and its solution relies on the identification of causal molecular interactions in genomics data. causal molecular interactions can be discovered using randomized experiments; however such experiments are often costly, infeasible, or unethical. fortunately, algorithms that infer causal interactions from observational data have been in development for decades, predominantly in the quantitative sciences, and many of them have recently been applied to genomics data. while these algorithms can infer unoriented causal interactions between involved molecular variables, causally orienting all inferred molecular interactions was assumed to be an unsolvable problem until recently. in this work, we use transcription factor-target gene regulatory interactions in three different organisms to evaluate a new family of methods that, given observational data for just two causally related variables, can determine which one is the cause and which one is the effect.  we have found that a particular family of causal orientation methods is often able to accurately infer directionality of causal interactions, and that these methods usually outperform other causal orientation techniques. we also introduced a novel ensemble technique for causal orientation that combines decisions of individual causal orientation methods. the ensemble method was found to be more accurate than any best individual causal orientation method in the tested data.  this work represents a first step towards establishing context for practical use of causal orientation methods in the genomics domain. we have found that some causal orientation methodologies yield accurate predictions of causal orientation in genomics data, and we have improved on this capability with a novel ensemble method. our results suggest that these methods have the potential to facilitate reconstruction of molecular pathways by minimizing the number of required randomized experiments to find causal directionality and by avoiding experiments that are infeasible and/or unethical. 22- april the international conference on intelligent biology and medicine nashville, tn, usa 