{'fraction_non_alphanumeric': 0.03361142307075245, 'fraction_numerical': 0.011992619926199263, 'mean_word_length': 4.576604786401252, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Causal inference using deep neural networks

Causal inference from observation data is a core problem in many scientific fields. Here we present a general supervised deep learning framework that infers causal interactions by transforming the input vectors to an image-like representation for every pair of inputs. Given a training dataset we first construct a normalized empirical probability density distribution (NEPDF) matrix. We then train a convolutional neural network (CNN) on NEPDFs for causality predictions. We tested the method on several different simulated and real world data and compared it to prior methods for causal inference. As we show, the method is general, can efficiently handle very large datasets and improves upon prior methods.


Introduction
An important problem in many scientific areas is the inference of causal interactions between two variables from observational data. Examples include causal inference between proteins in biological studies [37,22], studies related to climate and global warming [30] and studies of economic implications of various potential interventions [8]. Given its central role, several supervised and unsupervised methods have been developed for predicting such causal interactions. Unsupervised methods utilize various correlation metrics [6], regression analysis [8] or attempt to learn joint probability distributions using directed graphical models [18,35]. Supervised methods often start with a few labeled samples and then use various classifiers or their modifications to predict causal interactions in test data [9,20,23]. While some recent work attempted to use deep learning techniques for causal inference [7,13], in many cases the way the observational data is obtained and constructed does not preserve any local information. Consider, for example, weather related data from several locations. For each location we may have several different values (avg temp, wind, precipitation etc.) but it is not clear in what order they should be provided for the network to fully utilize dependency between values that are in close proximity in the input space. A similar case occurs when genes in one cell are profiled. In such cases, the ordering of the input vectors is often arbitrary.
Here we extend regression methods for causality inference by enabling the use of deep neural networks (NN) for causal inference. For this, we extend a method first applied to inferring interactions in gene expression data [37]. In contrast to gene expression data, where the main question is undirected interaction prediction, causality inference is non symmetric. To infer causal interactions we reformat the input and obtain an image-like representation for the data. For this, we first convert the joint observations for each pair of samples to a normalized discrete probability matrix, where each entry represents the joint density of a set of observations. This preserves the local relationships between the values observed for each of the samples. Next, a CNN is used to classify based on these distributions. An overview of the method is illustrated in Figure 2.
We tested our method on a series of simulated network structures including "V structure" (X → Z ← Y ), "chain structure" (X → Z → Y ), and "reverse structure" (X ← Y → Z), and also applied it to four real world datasets in various fields. In all cases our method improved on prior methods suggested for causal inference.

Related work
The gold standards for causal inference are randomized controlled experiments in which X is perturbed and its impact on Y is studied [36,12]. However, for many real world datasets such experiments are either impossible to perform (for example, manipulating genes within an individual) or very challenging [16].
Another popular strategy for causal inference, which is often applied to time-series data, uses regression analysis. One of the most popular methods for such analysis is Granger causality. Granger causality inference learns two models Y t (X t ) and X t (Y t ) where t < t and compares between them to determine which better explains the other variable [8,17,5]. Some methods extend Granger causality using supervised training. For example, SIGC attempts to classify Granger causality results based on difference between conditional distributions observed in training data [3]. Ganger causality, and other regression based causal inference methods can utilize either linear or non-linear regression [17,31,28,10]. However, the specific type of regression model (as opposed to the parameters) is a user defined issue and is not learned from the data.
A more global strategy attempts to learn joint efficient probability distributions based on directed edges. These methods mainly rely on variants of probabilistic graphical model (PGM) [35]. Unlike pairwise analysis methods, PGM methods can leverage information from other variables when determining the specific pairwise relationship between two variables. However, despite being directed, PGMs are not meant to provide causal relationship inference and in many cases edge directions in such networks do not imply causality [25].
Other methods focus primarily on causality inference for paired variables. One idea for inferring casual relationships between such variables is based on the independent additive noise assumption. Under such model we assume that if X → Y , then the observed values will follow the equation: Y = f (X) + noise [31,28,39]. However, how to find a suitable f and noise function form is a challenge for such methods [10]. Despite several of the methods that we will discuss as comparisons, an extension of such notion is the TiMINo model, which uses both additive noise and time information to infer causality from time-course data [27].
When labeled data is available, any one of several supervised methods can be used to perform such inference. One recent example is RCC which uses kernel mean embedding to construct features from sample pairs to train random forest classifiers [20]. Another is based on neural networks for embedding samples in lower dimension before performing causal inference [21]. See [14] for a recent comprehensive review of causal inference methods.

Motivation
For pair-wise causality inference tasks, most prior methods rely on specific assumptions about the type of noise (additive, multiplicative etc.) or the type of regression relationship expected (linear f (X), nonlinear f (X), Gaussian and other noise distributions [14]). While such approaches work well for certain types of data, they may not be general enough to accommodate datasets that display mixtures of these noise or regression characteristics. For example, consider the causal relationships displayed in Figure 1 which are all part of the Kaggle causal-effect pair dataset [9]. While some of these relationships are indeed linear, others are not and no specific regression model is likely to be able to accomodate all the different types of relationships observed in this data.  [9]. Given paired samples for one variable pair ( X, Y ) a normalized empirical probability distribution is constructed.
To enable the use of general strategy for causal inference we developed a new method that replaces the regression component of current methods with a general CNN that can be used to learn both linear and non linear relationships depending on the input data. Specifically, consider the Granger's causality formulation [8], which learns a model for multiple time-series variables X 1 , X 2 , ...X n with time-lag p by setting : X j (t) = n i=1 p a=1 β ji,a X i (t − a) + j (t), where β is the matrix of coefficients modeling the effect of time series X i on the time series target X j . In this approach we need to set both the time lag and the types of coefficients we learn. Imn our model we replace this with a general CNN that can be trained to determine both an appropriate lag and the specific type of coefficients needed to accurately infer casuality as we discuss below.

Distribution classification
Our first goal is to obtain a representation of the input data that can enable the use of convolutional neural networks for this task. Many causal problems use as input a collection samples. For example, in biology, when attempting to learn gene-gene interactions we may have samples with the joint expression of these genes in thousands of cells. Such problems attempt to infer causality from an observation set and a small number of labeled pairs: where n i is the number of observations for the ith pair (x i , y i ), N is the number of pairs and L i ∈ {1, −1, 0} is the label for the causal interactions such that L i = 1 if x i → y i and L i = −1 if y i → x i . L i could also be 0 meaning the two variables are independent. One way to infer causal interactions from this data is by learning a mapping function: f : R 2×ni → L where L is label set, and n i is the number of observations for each pair.
We aim to develop an efficient method that can learn casual relationships without relying on a pre-determined statistical model. If we assume (x ij , y ij ) then the problem can be converted to distribution classification problem [19,20,29]: f : P → L where P is a joint PDF family and P i ∈ P. Note that in many cases the pairwise values are indeed independent between the samples and so this assumption holds for many real world datasets.

NEPDF generation
Given a set of observations for two entities in our data {(x ij , y ij )} ni j=1 ), L i , the first step is to estimate the pair's joint PDF P i . P i is often estimated using a 2D kernel density estimation [33]: where I is the indicator function.
Given samples of the pair, the joint PDF can thus be be estimated as To enable the use of joint density P (x i , y i ) by a CNN we need to convert it to a 2D histogram. We thus obtain a 2D joint empirical PDF of size of K × K by setting: Where x o (x o+1 ) is the minimal (maximal) x value for the oth bin of x, y p (y p+1 ) is the minimal (maximal) y value for the pth bin of y, and 1 ≤ o, p ≤ K. In this paper we uniformly divided the observation ranges either in the original or in log space, though other ways that result in non uniform bins are also possible. We next normalize the EP DF i and obtain a Normalized EP DF i between [0, 1].
Our representation of each pair by a N EP DF converts the inputs to an image-like matrix of size of K × K. This provides a number of advantages. First, such representation preserves locality (unlike the original vector space [37]). Second, since the NEPDF matrix size is fixed, for each pair the training and test run time is fixed regardless of the number of observations. Finally, a very common problem in dealing with scientific dataset is the missing of observations, but the construction of N EP DF can easily handle such cases without the need to add additional modifications to the observations.
When generating a train/test dataset for N EP DF , for each positive pair (x, y), the model will also learn a N EP DF for (y, x) which will be assigned a label of -1. Note that the N EP DF , of (x, y), is the transpose matrix of (y, x). CNNs are known to learn localized features and thus will assign different values for such rotations. While this is a regarded as a problem for other applications [32], we can take advantage of such property when learning the networks by using both the original and rotated versions of the N EP DF with opposite labels.

CNN classifier
We next use the constructed N EP DF to train a CNN based on VGGNet [34]. The CNN consists of q stacked convolutional layers of size 3 × 3, where q is the power of 2, and interleaved maxpooling layer of size of 2 × 2. The final layer is a logistic regression for binary label task and softmax classifier for multilabel task, so that the prediction is with a crossentropy loss function where, L i,l = 1 if sample i belongs to lth class, y i,l is the prediction probability that sample i belongs to lth class by CNN. Stochastic gradient descent (SGD) is used to train the model. Note that the size of N EP DF is an important hyper-parameter for the CNN. As a result size can be optimized by cross validation. We provide an implementation of N EP DF and CNN classifier in the github. 2 .

Extension to multivariate case
For multivariate cases, for example, drawing the pathway of RNA/DNA flows or determining the causal graphs between multiple cause and effect data, the CNN method can be used in two ways. The CNN method can be used in a general PC-algorithm [35] framework, and provide causal direction for pairwise data beyond V-structure cases. Another way of mining three-variable relationships is by construction of 3D N EP DF s.

Experiments
To test the usefulness of our CNN method, we conduct series of experiments on both simulated and real data, and compare our algorithm to other state of the art methods including the baseline regression method (BivariateFit), Randomized Causation Coefficient RCC [20], Neural Causation Coefficient NCC [21] ,Additive Noise Model ANM [1], Information Geometric Causal Inference method IGCI [11],Regression Error based Causal Inference RECI [2], and Jarfo [6].

Numerical Simulations
A Structure Equation Model [26] for an acyclic graph G(V, E) specifies for all i ∈ V a function for the random variable where PA i denotes the parents of X i and N i , i ∈ V are jointly independent. In the bivariate case of X and Y , the causal relationship can be expressed according to the function Y = f (x) + n [31,39,10]. Such additive function can be further extended to the time series model [27]: For each i from a finite set of variables V , let (X i t ) t∈N denotes the time series. Thus, the goal is to acquire a summary time graph, with the time series components as vertices and an edge between X i and X j if there is an edge from X i t−k to X j t for some k. As a simplification to the time series model, we performed simulated experiments using three random variables (X, Y, Z) with a skeleton network X − Y − Z. We selected the three most common causal structures for such variables in directed acyclic graphs (DAGS): , and X ← Y → Z (Reverse V structure). For each of these causal structures we selected specific functions (see below) and repeated the simulation 10, 000 times with different random seeds. Then, we divided the generated causal pairs into training, validation and testing, and obtained the score using a five-fold cross validation analysis.
Since many of the methods we compared to formulate the problem as a binary classification task, we followed the approach presented in the ChaLearn's Challenge [20] to evaluate performance. Specifically, we trained two classifiers for causal/reverse-causal and dependent/independent predictions. Since these classifiers can lead to contradictory results (for example y causal and y ind ), the final prediction is computed as: y pred = y ind · (2 · y causal − 1).

V structure
We used the following non-Gaussian additive equations to construct the scenario where two variables share a common confounder with one time-lag. The equations were repeated 1, 000 times starting with a random seed to obtain 1000 values for our three variables (X, Y and Z).
where X noise,t , Z noise,t and Y noise,t are generated from N (µ, σ). At each step, the µ and σ are sampled from the uniform distributions U[0, 10] and U[0, 50], respectively. We experimented with different settings of α and β as we show in Table 1. Larger α and β put more weight on the previous values whereas smaller ones put more weight on the influence of noise.
We used the simulation data to construct six 2D N EP DF s with ground truth labels: (X, Y, 1), (Y, X, −1), (X, Z, 0), (Z, X, 0), (Z, Y, 1) and (Y, Z, −1). Here the labels of 1 and −1 stand for causal/ reverse-causal relationship, while 0 means that the two variables are independent. Table 1 shows average AUROC for one label versus the other two for all methods we tested. As can be seen, all methods perform less satisfactorily when the data is more noisy (lower α and β). However, when we increase the weight of past observation our CNN approach greatly outperforms all other methods in terms of AUC score.

Chain Structure
For the chain structure, 1, 000 time-step data pairs were simulated using the following non-linear equations with non-Gaussian noise:   . We again tested different settings for α, β and γ. Table 2 shows the average AUROC for all methods. In all cases, CNN achieved an AUROC score above 0.9, and it can successfully detect the independence between X and Y variable when Z was the intermediate node. For this setting CNN achieved much better performance than all methods we compared to for all parameters we tested.

Tübingen dataset
The Tübingen dataset consists of 108 heterogeneous, hand-crafted X −Y pairs from real world causal discoveries, including physical phenomenons, social and economical effects, biological observations, and etc [24]. Since this dataset is fairly small, it is not possible to generate sufficient N EP DF s for CNN to train only from this data. Instead, we followed [21] and constructed a larger dataset by adding heteroscedastic additive noise to obtain additional pairs. We used this to construct 15, 000 synthetic samples, where each sample S i included m i ∼ U[100, 1000] pairs of {x ij , y ij } mi j=1 . A similar X − Y pair generation was used in our experiment [21], where each of the cause x ij was drawm from a k i -mixture Gaussian distribution, with each Gaussain sampled with From this we constructed N EP DF with size of 16 × 16 for both (x ij , y ij ) (label −1) and (y ij , x ij ) (label 1). For CNN training, we applied five-fold cross validation and achieved an average test accuracy of 0.913 on the synthetic causal pairs. In the directional classification task for Tübingen dataset, we calculated a weighted accuracy score to avoid excess similarity in the testing data, as suggested by [24]. The classification accuracy of our CNN approach was 0.784, which is comparable to performance reported for state-of-the-art methods NCC (0.79) [21] and RCC(0.75) [20].

Gene relationship inference
While causal inference is challenging in all scientific areas, following decades of research biologist have been able to map several causal interactions between genes using knockout experiments. Current knowledge about such interactions is often encoded in the form of genetic pathways in which directed edges connect causal genes [4]. We thus tested our method using data from the Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway database [15] which compiles a set of experimentally validated gene pathways graphs where each node represents one gene and each edge represents directed or undirected interaction. For our analysis we focused on the directed gene pairs for the positive labels. As for the observation data, we used the Encyclopedia of DNA Elements (ENCODE) [38] gene expression dataset that contains 249 samples, each of which has more than forty thousand genes.
We compared our method with the same set of prior methods we used for the simulated data. Given a known gene pair gene1 → gene2, a N EP DF of (gene1, gene2) with label 1 and a N EP DF of (gene2, gene1) with label -1 were generated. Table 4 shows the comparison result with three-fold cross validation. For RCC, several combinations of hyperparameters were used, and we reported the series of hyperparameters with the best performance. As can be seen, our CNN method outperformed all prior methods for gene causality prediction. Interestingly, unlike for the simulation data we observed that only the supervised learning methods were able to perform well on this task whereas all the unsupervised methods, including ANM, RECI, IGCI, BivariateFIT and CDS do not perform well. This may indicate that for real world data (including gene regulation network inference) the strong assumptions used by unsupervised methods about the statistical model of dependency may lead to inaccurate results.

ChaLearn's Challenge Data(Kaggle Competition)
The Kaggle causal-effect pairs competition collected data from different scientific areas [9]. It provided more than 15, 000 variable pairs, each of which has thousands of observation values and a target label vector. As with our experiments, there were labels of 1 and −1 meaning the causal and anti-causal relationships, 0 indicating no causal relationship. The competition also provided a test dataset of 4, 050 pairs with observation values but without target labels. On the test data, our single neural network model achieved a mean test test bidirectional AUC of 0.74, which is comparable with the score obtained by RCC (0.74) [20]. Note that to run our method we only need a single CNN network trained on NEPDF inputs. In contrast, RCC uses a large engineered feature vector and an ensemble learning strategy based on tens of classifiers. While 0.74 is less than the winning score Figure 3: Results for CNN, PC, MI, and GC for stock sector prediction. While all CNN, PC and MI have shown increases in AUROC score than simple random guess and GC, CNN is also able to obtain a larger AUROC.
(0.82) of the winning method, Jarfo. However, Jarfo was heavily tailored to the data and so may not generalize to other causality prediction tasks. Table 5 shows the comparison result we compare with other supervised and unsupervised baselines.

Stock sector classification
While our main focus is on causal inference, our method can be trained on any type of dependency label a user supplies. To illustrate a potential use for such training we used our method to predict relationships between pairs of stocks. Specifically, we trained a model to predict if two stocks are from the same sector based on their price information. For this we downloaded data from Yahoo Finance covering the years between 2014 and 2019, and containing both price and sector information for all stocks in the S&P 500 index. This dataset contained 1259 observed prices for each stock. Data was processed to compute the daily changing percentage of stock prices which was used as an input to construct N EP DF s.We divided the stock into five equal size folds for cross-validation analysis. For training, we assigned a label of 1 if stocks x i and x j belong to the same sector and 0 otherwise.
Most of the methods we compared to above focus on predicting the direction of an interaction and cannot be used to predict new interactions. Thus, unlike our method they cannot be used to infer interactions between stocks. We thus compared our results with the most popular statistical measure, Pearson correlation(PC) and mutual information(MI), as well as Granger's causality(GC), which is a popular method to determine the relationships in time-series data. For Granger's causality, we use the p-values as the final prediction. Figure 3 presents the comparison between the methods we tested. As can be seen, our method achieved a higher AUROC when compared to the other two methods. While PC and MI performed well, CNN still improved on MI for this task as well (AUROC of 0.945 vs. MI(0.908) and PC(0.918)).

Conclusion and Future Work
In this paper, we presented an efficient method to learn casual relationships that does not rely on an explicit statistical model. Our method converts general pair-wise observation data to normalized empirical probability density distribution (N EP DF ) and uses this representation to train a convolutional neural network to predict causal or other relationships between input samples. The use of NEPDFs enables our method to utilize methods that work very well for image classification problems for the task of causal inference. Experiments on simulated and real data demonstrates that the method is general and works well on a diverse set of input data types.
There are several ways in which we can further improve the methods and its results. First, improving the N EP DF construction using other density estimation methods can lead to better representation and results. Second, we would like to extend pairwise inference that is currently supported to learn interactions and causality between a larger set of variables. Finally, we would like to use the method and representation to infer patterns that are indicative of causal relationships and to interpret these to better understand the underlying mechanisms that lead to the causal relationship observed.

Broader Impact
Here we present a general pair-wise variable relationship inference method based on probability estimation and deep learning technologies without relying on an explicit statistical model. Such method would benefit a range of scientific fields that are interested in pair-wise variable relationship exploration while lack known models, for example, gene regulation in biology, stock relationship in economy, climate changes impact among different regions and so on.