The high-throughput technologies, such as microarray (Schena et al., 1995) and RNA-Seq (Wang et al., 2009) in transcriptomic level, generate bunch of data of describing various perspectives of cell state. These data provide unprecedented opportunity to quantify molecular expressions and their relationships. From a systematic perspective, the molecules in a cell orchestrate together to form various integrated and condense network systems of performing comprehensive functions (Liu, 2015). For instance, transcriptional interactions between transcription factor (TF) and target genes are often formulated into gene regulatory network of modeling biological processes (Liu et al., 2014, 2015). Deciphering gene relationships from high-throughput data are crucial to reversely engineer their inner interaction scenarios, as well as profoundly reveal the dysfunctions in certain disorders, such as complex diseases (Liu et al., 2012).

Quantifying the relationship between molecular components becomes fundamental in the new research paradigm from data to knowledge. The data analysis techniques of association support the kind of investigation. Traditionally, when we explore the relationship between two variables, Pearson's correlation coefficient (PCC) is employed to qualify their linear relationship (Zou et al., 2003). From entropy aspects, mutual information (MI) is often used for defining the non-linear relationship between gene variables (Butte and Kohane, 2000). Mathematically, the assumptions underlying these measures are considerable in real applications. Association measures have been developed to meet the requirements of appropriateness and precision in defining relationships from various perspectives.

Detecting gene associations is a fundamental method to reconstruct gene regulatory network from gene expression profiling data (Liu, 2015). Although more integrated methods such as ordinary differential equations are available to model the differential dynamics among genes, the association-based methods are direct, simple, and easy for interpretation as well. With introducing the independence, these measures have been extended to quantify the associations between many genes simultaneously (Stuart et al., 2003). In typical microarray experiments, the gene expression data can often be represented by matrix G,

Where G_ij represents the gene expression value of the i-th gene (1 ≤ i ≤ m) in the j-th experiment (1 ≤ j ≤ n). It is noted that j refers to a sample or a time point with specific phenotype meaning. The association between gene X and gene Y (X, Y ∈ {G₁, G₂, ⋯ , G_m}) is often to indicate their regulatory relationship (Zhang and Horvath, 2005). Let gene expressions be X = (X₁, X₂, …, X_n) and Y = (Y₁, Y₂, …, Y_n). Based on the two vectors, we employ or define an association measure to assess their regulatory strength. Recently, some novel measures besides PCC and MI have been proposed to define the association between two variables (Reshef et al., 2011). It is of great interest to investigate their performances in the reconstruction of gene regulatory network from gene expression data. Figure 1 demonstrates the strategy of inferring gene regulatory network by gene coexpression analysis. Gene regulation, in a particular form of transcriptional regulation, often specifies the regulation from TF to target gene. The quantified gene coexpression evaluates the simultaneous patterns of two gene's redundancy across samples. The expression level of upstream TF's gene is often to approximate its downstream protein product. As shown in Figure 1C, if we set up which ones are TFs by prior knowledge in the gene association network, we can infer a directed gene regulatory network via an undirected association measure.

The coexpression pattern between two genes implies their regulatory aspects. As shown in Figure 1C, it firstly indicates a direct regulatory interaction. In some biological state, gene coexpression exactly responds to the activation or inhibition regulation from a TF to its target gene. The regulation between them is reflected by their highly-related gene expression redundancy. Secondly, gene coexpression is about gene co-regulation. That is to imply the two genes are regulated by the same TF(s) and then they contain highly-related gene expression patterns. Third means that the two genes are functionally-related by participating in the same regulatory circuit or particular signaling pathway. Generally, the dynamic regulations in a cell are inherently embedded with temporal features. Gene regulation is often reflected by time-delayed gene expression patterns from the activation of TF's gene to the downstream target responds (Bar-Joseph et al., 2012). For the simplicity of association measure, the coexpression-based methods are popular in inferring gene regulatory network from gene expression data (Zhang and Horvath, 2005).

In this paper, we provide a comparative study on these available association measures of quantifying gene relationships in regulatory network. Fourteen most-popular association measures or indices will be summarized and compared. Based on some benchmark datasets of gene regulatory network inference challenges, we evaluate their individual performances in the reconstruction of gene regulatory networks. This provides a concise comparison of accuracy and quality in network inference by the association measures. In a case study, we compare the differences of these inferred regulations during the infection of hepatitis C virus on host cells. In data-driven network inference, the characteristics of the association measures in statistics and computations are also analyzed and discussed.

Numerous association measures have been proposed to define the relationship between two random variables. For gene regulations, we collect 14 of them for our assessments of network inference power from data. Table 1 lists the 14 association measures with brief introduction of their statistical assumptions and fundamental properties individually. Some measures are well-known such as PCC, while some become available recently such as maximal information correlation (MIC). For the completeness of introduction and reference, we describe them in details respectively in this section.

For a comparative study of these association measures in inferring gene regulatory relationships, we test these association measures in DREAM3 in silico network challenge datasets (Marbach et al., 2010). In the challenges, gene expression datasets have been generated by some specified network structures. Then, the datasets are open without any information about the network structures. The task is to reconstruct the network structures from the open datasets by developing new inference methods. There are three sizes of networks with 10, 50, and 100 nodes respectively, and multiple datasets for each size (4 for 10-node network, 23 for 50-node network, and 46 for 100-node network). The assessment is to evaluate the consistency between the inferred network and the true network structure (gold standards). Figure 2 illustrate the receiver operating characteristic (ROC) curves of inference performance by these association measures in the 10-node benchmark network. Due to the undirected regulations identified by all these association measures, we omit the regulatory directions when calculating the evaluation metrics of sensitivity (SN), specificity (SP), accuracy (ACC), Matthews correlation coefficient (MCC), F-measure, and area under ROC curve (AUC). Table 2 demonstrates these detailed values of evaluation metrics of these association measures. We find KCCA performs the best in the 14 association measures for inferring 10-node networks and it reaches the AUC of 0.623 ± 0.083 (mean ± standard deviation). Overall, the performances of these methods are comparable with each other in the 10-node network.

For the association measures, it becomes more difficult to achieve high inference performances when the network size becomes bigger from 10, 50 to 100. Although each association measure cannot achieve good inferences for big networks, the performances of them decrease with the same tendency. For 50-node networks, mutual information (MI) achieves the best AUC of 0.569 ± 0.046. Blomqvist's β performs the best for 100-node networks in the inference, while it is not stable for the small-size networks. Figure 3 shows the ranks of their performances according to the mean AUCs in different size of networks individually. From the comparative study, mutual information (MI) performs relatively better with stable ranks for big networks with 50 and 100 nodes. PCC is also stable in the 14 association measures for various sizes of network, as well as KCCA and dCor. This indicates their relative reliability in detecting gene regulatory relationships from expression data. For the other association measures, they accomplish unreliable and unstable regulatory network inferences in the benchmarks.

From the inference performances, we find that most of association-based methods can only achieve limited accuracies in the reconstruction of gene regulatory network from the benchmark datasets, especially for large-size networks. The application scopes of these association measures are mainly determined by the assumptions and characteristics of their definitions listed in Table 1. For instances, PCC is for linear regulatory relationship, MI is for non-linear relationship, KCCA and dCor measure the genuine relationship based on covariance, and the rank-based associations are robust to the noisy and outliers in gene expressions. In practical applications, the selection of suitable association measures could be subjectively determined by research purpose, experimental design, phenotypic condition and data quality. An ensemble and self-adaptive association measures selection strategy is desirable to be proposed for the co-existence of different gene regulatory relationships.

In real microarray data, we perform our comparative study of quantifying gene regulations during hepatitis C virus (HCV) infection on host Huh7 cells. The gene expression data are downloaded from NCBI GEO (accession ID GSE20948) (Edgar et al., 2002). There are 28 samples of 14 HCV infected Huh7 hepatoma cell samples and 14 corresponding mock-infected samples, originally designed three replicates at 6, 12, 18, 24, and 48 h post-infections, respectively. Two samples at 6 h have not been enrolled after quality control. The details can be accessed from Ref. (Blackham et al., 2010). We also download the hepatocellular carcinoma (HCC) gene set from KEGG (Kanehisa and Goto, 2000). The gene set contains 123 genes with 94 genes containing their expression profiles in GSE20948 (Edgar et al., 2002).

For evaluating the inference consistency of these association measures, we calculate the pairwise gene regulatory strengths in the HCC genes by the 14 association measures respectively. In the results of each association measure, the pairs with the top 5% association values are regarded as the identified gene regulations in the context of specific gene expression profiles after HCV infection.

Figure 4 demonstrates the inferred gene coexpression regulatory network in the HCC genes by PCC. There is no information about direction, so we annotate the known human TFs and display them by different color nodes (cyan) with the other genes (green). From Figure 4, we can figure out the regulatory information about positive and negative relationships during HCV infection. As in the former comparisons, we compare the overlapping status of these inferred coexpression relationships by the four association measures with top performances, i.e., Pearson, MI, KCCA and dCor. There exists only one pair of genes (“IFNA1” and “IFNA13”) is identified by the four measures, and the relationship between the two genes can be detected by any of them. Interestingly, Pearson and dCor contain many overlaps (177 regulations). It provides direct evidence that dCor is mainly to extract the linear correlations between genes as that Pearson done in this case study. There are few overlaps (3 regulations) between Pearson and MI, which indicates the linear and non-linear information are inconsistent with each other, and different association measures might identify different gene associations. The selection of suitable association measures is again proved to be very important for inferring gene coexpression regulatory network. The few overlapping regulations also imply the complex and diversity of regulatory relationships underlying gene expressions. More advanced methods beyond association measures are urged for elucidating gene regulatory mechanism from high-throughput data. See Section Discussion for some already available methods.

It is known association is different from causality and correlation does not imply causation (Altman and Krzywinski, 2015). Detecting the causality between genes has been essential in gene regulatory network inference since the availability of high-throughput data (Opgen-Rhein and Strimmer, 2007). Gene association network indicates more general gene-gene relationship than regulation, and gene regulatory network indicates more general gene-gene relationship than causality. The gene causality network, that is to say, the causal regulations between genes are directed in the gene-gene interaction graph with the detailed information of which ones are upstream regulators, and which ones are downstream targets. In the direct regulations, TFs or signal transductors causally affect their target gene expressions. The information flow transits between genes will be revealed if a causal relationship exists. So far, there is no association measure has been defined for describing the causal relationship between genes (Zhang et al., 2014; Zhao et al., 2016), while more advanced methods based on conditional probability, model-based regression and differential equation have been proposed to address the evaluations of causality.

Based on conditional independence, some improved association measures, such as partial correlation coefficient and conditional mutual information, have been proposed to eliminate false positive regulations from gene associations. The original association measures generate the footholds for detecting genuine relationships. Conditioning on another gene or gene set Z, partial correlation measure r_XY·Z between gene X and Y is to access the exact correlation between X and Y and that has no relationship with Z (de la Fuente et al., 2004). It is defined as

Where r refers to PCC. In the similar philosophy of introducing other gene or gene set, the conditional mutual information (CMI; Liang and Wang, 2008) is defined as

Based on CMI and the order of conditioned gene numbers, we proposed a gene regulatory inference method named PCA-CMI (Zhang et al., 2012, 2013), which detect out dedicate associations by removing undirect false positive regulations. For a pair of genes X and Y, Li proposed a conditional coexpression measure named liquid association (LA) between two genes by introducing a third gene Z (Li, 2002). Based on Z, the gene relationship of X and Y is defined as

where n is the sample size. The LA activity determines the functional associations of gene X and Y in the condition of Z.

Currently, the causality between genes is often quantified via Bayesian models (Friedman et al., 2000). According to data, the conditional probability of P(X|Y)=P(Y|X)P(X)P(Y). The probability of gene X conditioned on gene Y, means Y have a causal effect on X because there exists a negative or positive values of the conditional probability. The structured model has been extended and formulated as diagrams using a graphical criterion known as d-separation (Bareinboim and Pearl, 2016). Bayesian network provides a model-based detection of causal regulatory relationships. Gene regulations are then identified from the graphical models (Liu et al., 2013).

Regression and other structured models often extract the effects of regulatory coefficients. The identification of model coefficients determines the global relationship of these individual genes (D'Haeseleer et al., 1999). Specifically, the regression models the response gene as the linear combinations of the other dependent genes, i.e., Y = c₀ + c₁X₁ + c₂X₂ + ⋯ + c_mX_m + ε, m is the number of dependent genes in the regression and ε is the error variable. In generalized linear models, the response gene is changed to θ(Y), and X₁, ⋯ , X_m are replaced by ϕ₁(X₁), ⋯ , ϕ_m(X_m), respectively (Breiman and Friedman, 1985). In the special case of simple linear regression with m = 1, the model is to detect the linear relationship between the response gene and the only one dependent gene. The coefficient of determination denoted by r² is equal to the square of PCC (Altman and Krzywinski, 2016). The coefficient of determination, which represents the proportion of variation due to their linear relationship, generalizes the correlation coefficient for relationships beyond simple linear regression. Often, the regression equations often model the associations between response genes and dependent genes in an inter-coupled system. From a system biology perspective, regression models consider the genes in an integrated manner. Compared to the former pairwise associations, they identify more complicated relationships among genes. After determining the coefficients, the relationships in these genes are quantified correspondingly. How to determine crucial regulators and targets via statistical variable selections techniques, such as lasso (Tibshirani, 1996) and elastic net (Zou and Trevor, 2005), are substantially important.

Similarly, ODE models the derivatives, i.e., dYdt=c0+c1X1+c2X2+⋯+cmXm, and so ODE quantifies the dynamics of the response as a function of the dependents in the system (Wu et al., 2014). The expression change rate of a response gene is modeled by the expressions of dependence genes. The Y might be another dependence gene and thus the system is closed. The system identification is to evaluate the coefficients in the right-hand side of the equation and the coefficient values refer to gene regulatory strengths. When the coefficient is 0, there is no relationship between the responding gene and the depending gene, otherwise the regulatory strength can be represented by positive or negative numeric values.

Compared to association measure, regression model and differential equation model regard gene regulatory network as an integrative system. The gene regulatory network inference is then transformed to a system identification problem of solving the coupled equations. The gene regulation strengths refer to the identified coefficients. From a sequential modeling perspective, the causality between regulators and targets can also be reflected by these system biology techniques.

In machine learning techniques such as clustering (Rui and Wunsch, 2005), there are some metrics have been developed for measuring the association between data points. The distances of Euclidean, cosine, Hamming, Manhattan are often used to measure gene relationships in gene expression clustering (D'Haeseleer, 2005). These distances evaluate the differences including dependences between genes, while these compared association measures focus on quantifying gene relationship such as regulation between genes. In gene expression data analyses of clustering and feature selection, distance metrics provide alternatives to define gene similarities. The distance metrics are not included in the comparative study for their diversity and case-intensity (Santini and Jain, 1999).

In this paper, we summarized and compared the main proximities and metrics for quantifying gene regulatory associations. Written in full, the definitions and descriptions of 14 association measures are summarized and their characteristics with applications in regulatory network inference have been presented. From the benchmark challenge data and real gene expression data, we compared their performances and consistencies in the network inferences. Furthermore, their advantages and limitations are also analyzed and discussed. Currently, developing causality measure is an urgent research topic from driving gene association to regulation causality (Bareinboim and Pearl, 2016). A powerful measure of causality will greatly benefit the discovery of important gene regulations. Moreover, the linear/non-linear regression and differential equation models regard many genes in dynamic systems and the parameters of these models represent the system in details. The model-based gene regulatory network inference methods seem to provide more powerful tools when compared to the association-based methods. However, the association measures contain their flexibility in sense, easy interpretation and large scope of applications.

In conclusion, gene association measures provide fundamental quantifications of detecting gene regulatory relationships from transcriptomic profiling data. The high-throughput technologies advance the measurements of thousands of genes in parallel manners. The association measures effectively accelerate the transformation processes from data to knowledge. Most of the proposed association measures are statistical techniques which focus only on the inter-relationships between genes, and they are very hard to get the causal gene relationships alone. With the improved conditional or joint association measures, such as partial correlation coefficient, conditional mutual information and liquid association, the causality between genes can be partially extracted out from data. The introduction of other genes in evaluating gene regulation provides promising alternatives to grasp the genuine regulations. For an entire system, many genes perform their functions coordinately and cooperatively. So more advanced models are extremely needed to describe the complex system of gene regulations. In such model as ODE, the time-varying regulations are exactly to quantify the gene regulatory interactions with temporal implications. For the model complexity and the data availability, the dynamics underlying the coefficients in regression and ODE will reveal much more complicated regulatory relationships.

ZL conceived and designed the study. ZL wrote the code and analyzed the data. ZL drafted the manuscript.