The Theoretical and Experimental Analysis of the Maximal Information Coefficient Approximate Algorithm

2.1 The Definition of MIC

Given two variables X, Y , the scale of the data set D = {(x_i, y_i), i = 1, 2, · · · , n} is n, that is |D| = n. The X, Y axes are divided into x, y grids, respectively. And then an x ×y grid G is obtained. The ratio of data points falling into grids is the approximate probability distribution D|_G. Different grids will result in different probability distribution functions and then different mutual information values are obtained. For the grids with the same numbers of partitions, maximize and normalize the mutual information. The maximal element value of the matrix is the maximal information coefficient (MIC). The specific definition of MIC is as follows.

Definition 1 The data set is D ⊂ R² and the grid is x ×y where x, y are positive integers. The probability distribution is D|_G and then the maximum mutual information is as Formula (1).

In the formula (1), I (D|_G) is the mutual information of the data set D under the gird G and the maximal value in formula (1) is the maximal mutual information on all partitions where the two axes are divided into x, y partitions, respectively.

Definition 2 The entities of the characteristic matrix M (D) of the data set D is defined as formula (2).

The entity of the xth row and the yth column in the matrix M (D) is the maximal mutual information of all x × y grids divided by the minimum value of x, y. That is the entity is the maximal mutual information normalized by the minimum value of x, y. After the normalization step, all entities of the matrix M (D) falls into the domain [0, 1] and the characteristic matrix is symmetric.

Definition 3 The MIC value of the data set D (|D| = n) with two variables is defined as Formula (3).

In Formula (3), the parameter B(n)=n1−ϵ,ϵ∈(0,1).

Taking the computing efficiency into consideration, in order to reduce the calculation workload, only a part, not all, of grids is calculated in the search procedure of calculating MIC. According to the practice experience^[13], the default value B(n) = n⁰.6.

The MIC value of two variables is also falling into the domain [0, 1], because all entities of the characteristic matrix fall into the domain [0, 1]. With the character matrix M (D), besides MIC, other criteria MINE is also summarized in [13].

2.2 The Original Approximate Algorithm Calculating MIC

According to Definitions 1 ∼ 3, the main work of the maximal information coefficient (MIC) is calculating the characteristic matrixM(D) of the data set. That is, based on the fixed number of partitions, maximize the normalized mutual information. In order to calculate the MIC of two variables, reference [13] gives the approximate algorithm calculating MIC of two variables.

In order to improve computing efficiency, besides the definition of MIC, an approximate algorithm calculating MIC of two variables is proposed^[13]. In the approximate algorithm, firstly, one axis is equally partitioned. And then, aiming at maximizing mutual information, the partition of the other axis is obtained via employing the dynamic programming algorithm. The flow chart of the approximate algorithm is shown in Figure 1. The order of the two variable in the data set D is interchanged and then the data set D⊥ is obtained. After that, the maximal mutual information of data sets D,D⊥ of the given numbers of partitions is calculated by the sub algorithm ApproxMaxMI(D, x, y,k̃). Specifically, in Figure 1, given the numbers of partitions of grids, that is the numbers x, y and log min{x, y} are given, the work of calculating the maximal mutual information of the given number of partitions is accomplished by the sub-algorithm ApproxMaxMI(D, x, y,k̃).

Figure 1

The approximate algorithm calculating the MIC of two variables in reference [13]

In the sub-algorithm ApproxMaxMI(D, x, y,k̃), one axis is equally divided into x partitions and the other axis is equally divided into cy partitions. Aiming at maximizing the mutual information, the cy partitions are fused into y partitions by employing the dynamic programming algorithm, and this is the main work of the sub-algorithm OptimizeAxis(D,Q, x). The flow chart of the sub-algorithm OptimizeAxis(D,Q, x) is shown in Figure 2.

Figure 2

The flow chart of the sub algorithm ApproxMaxMI(D, x, y, k)

3.1 The Theoretical Time Complexity Analysis

According to Figures 1 and 2, the approximate algorithm proposed in [13] mainly applies the sub algorithm ApproxMaxM(D, x, y, cx) to calculating the maximal mutual information under the given x, y partitions. Specially, there are four sub-algorithms including Equipartition(D, y), GetClumpsPartition(D,Q), GetSuperClumpsPartition(D,Q,k̃) and OptimizeXAxis(D,Q, x).

The sub-algorithm Equipartition(D, y) equally partitions the data set which has been sorted by the fast sorting algorithm. Then the time complexity of the sub algorithm is equal to that of the fast sorting algorithm, that is O(n log₂n).

The sub-algorithm GetClumpsPartition(D,Q) sorts the data set via employing the classical fast sorting algorithm and then traverses the sorted data set. Then the time complexity of the sub-algorithm GetClumpsPartition(D,Q) is also equal to that of the fast sorting algorithm, that is O(n log₂n).

The sub algorithm GetSuperClumpsPartition(D,Q,k̃) fuses the k̃ = cx partitions into x partitions, where default c = 15. Then the sub algorithm only traverses the k̃ partitions. In other words, the data set is only traversed one time. Then the time complexity of the sub algorithm GetSuperClumpsPartition(D,Q, k̃) is O(n).

The time complexity of the sub algorithm OptimizeXAxis(D,Q, x) is O(k²xy) = O(c²x³y) = O(x²B), where k = cx, xy ≤ B = n^α. Then, if the number of the X partitions is determined, that is X is given, the time complexity of the sub algorithm is O(k²xy) = O(c²x³y) = O(x²B).

According to the original approximate algorithm proposed in [13], the range of x is [2,B/2]. And then the time complexity of whole approximate algorithm proposed in [13] is O(x²B)B = O(x²B²) = O(B⁴) = O(n^4α). If the parameter is the default value α = 0.6, the time complexity of the whole algorithm is n².4, which is relatively high.

3.2 Experimental Analysis

In this paper, the software, developed by Albanese, et al.^[14], for calculating MIC values of two variables is the implementation of the approximate algorithm in [13]. The configuration of the computer used for computation is as follows. Win 7 operation system; CPU: Intel(R) Core(TM) i5-2450, 2.50 GHz; RAM: 4.00 GB. And similar to [13], 9 different type relationships are also adopted in this paper and the specific description is shown in Table 1.

According to the 9 different type bivariate relationships in Table 1, generate the data sets with different scales, and then calculate the MIC values by employing the software of [14] which is the implementation of the original approximate algorithm proposed in [13]. The calculation time is shown in Figure 3.

Figure 3

The calculating time of the bivariate relationships of different types

From Figure 3, the following results can be obtained. For the data sets of the relationships with the same scale, the time of calculating MIC values of random relationships is longest, and with the increasing of the scale, the calculating time of random relationship increases rapidly. In contrast, the computation time of linear relationship is shortest and its growth rate is slowest.

To further study the causes of the varying of the computation time, the computation time of all sub-algorithms and the whole approximate algorithm is compared. And it is found that the majority computation time of the whole algorithm is that of the sub-algorithm OptimizeXAxis(D,Q, x). The comparison of computation time is shown in Figure 4. The parameters are n = 300, α = 0.6, c = 15.

Figure 4

The comparison of the computation time of the sub algorithm OptimizeXAxis(D,Q, x) and the whole algorithm

From Figure 4, the following results can be obtained. The computation time of the random relationship is longer than that of the other type relationships and for every type relationship, the most of the computation time is spent on finding the optimal partitions (maximizing the mutual information under the fixed number of partitions). The results from the experimental aspect are consistent with that from the theoretical aspect. In the theoretical analysis, the time complex of the sub-algorithm OptimizeXAxis(D,Q, x) is O(k²xy) which is relatively higher. In order to study the cause of the varying of different type relationships, the sub-algorithm OptimizeXAxis(D,Q, x) is further analyzed.

In the sub-algorithm OptimizeXAxis(D,Q, x), the number of candidate partitions of different type relationships in the procedure of maximizing mutual information is shown in Table 2. Due to the parameters n = 300, α = 0.6, then the lower integer of n^α is 30. That is the maximal value of the product of the partitions of two axes is 30. Then the product of the first and second column is not larger than 30. The values of the other columns are the number of candidate partitions of different type relationships. Due to the parameter c is equal to 15, the number of candidate partitions is not larger than the product of the value of the second column and the parameter c (c = 15). From Table 2, it can be obtained that the number of partitions of linear relationships is least and the number of the random relationship, which is close to the upper bound of partitions (the product of the parameter c and the number of partitions in the second column), is largest. Due to the difference of partitions, the computation time varies from different type relationships. Among these different type relationships in Table 1, the computation time of linear relationship is shortest and the computation time of random relationship is longest in Figures 3, 4.

Because there are many random bivariate relationships in high dimensional big data and the time of calculating MIC values of random relationships is long, the computation time of detecting bivariate correlations in big data is very long via employing MIC.