1. 概念

1.1 MIC

1.2 互信息

2. MIC的优点

3. 算法原理

3.1 MIC公式原理

3.2 MIC计算步骤

（1）计算最大互信息值

（2）对最大的互信息值进行归一化

（3）选择不同尺度下互信息的最大值作为MIC值

4. 代码实现

6. 其他方法

1. 概念

1.1 MIC

MIC，即（Maximal Information Coefficient）最大信息系数，属于Maximal Information-based Nonparametric Exploration (MINE) 最大的基于信息的非参数性探索，用于衡量两个变量X和Y之间的关联程度，线性或非线性的强度，常用于机器学习的特征选择。

MIC相较于Mutual Information（MI）互信息而言有更高的准确度，是一种优秀的数据关联性的计算方式。

1.2 互信息

互信息（Mutual Information）是信息论里一种有用的信息度量，它可以看成是一个随机变量中包含的关于另一个随机变量的信息量，或者说是一个随机变量由于已知另一个随机变量而减少的不肯定性。这个已经是机器学习中老生常谈的内容了，如果想不起来，请参考百度百科-互信息。

2. MIC的优点

算法对比表：

根据 MIC 的性质，MIC 具有普适性、公平性和对称性。

所谓普适性，是指在样本量足够大(包含了样本的大部分信息)时，能够捕获各种各样的有趣的关联，而不限定于特定的函数类型(如线性函数、指数函数或周期函数)，或者说能均衡覆盖所有的函数关系。一般变量之间的复杂关系不仅仅是通过单独一个函数就能够建模的，而是需要叠加函数来表现。对于普适性较好的函数，不同类型的关联关系其起点应当是接近的。而且是接近于一的。

所谓公平性，是指在样本量足够大时能为不同类型单噪声程度相似的相关关系给出相近的系数。例如，对于一个充满相同噪声的线性关系和一个正弦关系，一个好的评价算法应该给出相同或相近的相关系数。而对于公平性较好的比较方法，随着噪音的增加，不同类型关联关系函数变化应当是相近的。

3. 算法原理

3.1 MIC公式原理

MIC基本原理会利用到互信息概念，互信息的概念使用以下方程来说明：

$p(x,y)$ 为变量x和y之间的联合概率，一般情况下联合概率计算相对来说比较麻烦。

MIC的想法是针对两个变量之间的关系，将其离散在二维空间中，并且使用散点图来表示，将当前二维空间在 x,y 方向分别划分为一定的区间数，然后查看当前的散点在各个方格中落入的情况，这就是联合概率的计算，这样就解决了在互信息中的联合概率难求的问题。下面的公式给出 mic 的计算公式：

上式中 a,b 是在 x,y 方向上的划分格子的个数，本质上就是网格分布，B 是变量，在原作者的论文当中提到 B 的大小设置是数据量的 0.6 次方左右。

3.2 MIC计算步骤

MIC计算分为三个步骤：

给定i、j，对X、Y构成的散点图进行i列j行网格化，并求出最大的互信息值
对最大的互信息值进行归一化
选择不同尺度下互信息的最大值作为MIC值

（1）计算最大互信息值

给定i、j，对X、Y构成的散点图进行i列j行网格化，并求出最大的互信息值，如下图所示：

值得注意的是，给定i和j后，可以得出多种不同的网格化方案。我们需要从这些不同的网格化方案中找到使互信息最大的网格化方案。举个例子，假设i=2，j=2，则可能有以下红、黄、绿三种网格化方案（其实更多，这里只是随便挑三种方案作说明），分别计算每个网格化方案对应的互信息值，找出使互信息值最大的网格化方案。

那么，给定了某个网格化方案后，如何计算其对应的互信息值呢？这里以上图中红色的网格化方案为例进行说明。红色网格化方案将所有数据点分为四个区域：左上，右上，左下，右下。每个区域对应的数据点数量为1，4，4，1。将数据点数归一化得到四个区域的数据点频率，分别为0.1，0.4，0.4，0.1。也就是说，此时，X有两种取值：左和右，Y有两种取值：上和下。P(X=左,Y=上)=0.1，P(X=右,Y=上)=0.4，P(X=左,Y=下)=0.4，P(X=右,Y=下)=0.1。并且，P(X=左)=0.5，P(X=右)=0.5，P(Y=上)=0.5，P(Y=下)=0.5。根据互信息计算公式，得到X和Y在这种分区下的互信息为：

计算完红色网格方案的互信息后，以此类推，再分别计算出黄、绿等网络方案的互信息，最后比较哪种方案得到的互信息值最大，即为该（i，j）取值下的最大的互信息值。

（2）对最大的互信息值进行归一化

将得到的最大互信息除以 $log(min(X,Y))$ ，即为归一化，这个与互信息原公式有关。

（3）选择不同尺度下互信息的最大值作为MIC值

上面是在给定i和j的情况下计算M(X,Y,D,i,j)。第三步就是给定很多(i，j)值，计算每一种情况下M(X,Y,D,i,j)的值，将所有M(X,Y,D,i,j)中的最大值作为MIC值。注意的是，这里的(i，j)是有条件的。当然，B(n)这个值可以自己定，这里是别人做实验认为效果最好的值。

4. 代码实现

在Python中的minepy类库中实现了MIC算法，具体使用如下。除此之外值得一提的是，minepy含有很多其他系数，有兴趣的话也可以研究一下。

参数解释：

alpha：float数据类型，取值范围为（0 ，1.0 ] 或 > = 4，如果alpha的取值范围在（0,1]之内，那么B的取值范围为（N ^α，4）其中n是样本的数目。如果alpha的取值范围是> = 4， alpha直接定义B参数。如果alpha高于样本数（n），则它将被限制为n，因此B的取值实际上是个分段函数，具体公式为：B = min（alpha，n）。
c（float 取值范围为大于）） - 确定比每个分区中的列多多个块。默认值为15，这意味着当尝试在x轴上绘制x网格线时，算法将以最多15 * x个团块开始。

代码一：展示直接使用MIC

import numpy as np
from minepy import MINEx = np.linspace(0, 1, 1000)
y = np.sin(10 * np.pi * x) + x
mine = MINE(alpha=0.6, c=15)
mine.compute_score(x, y)print("Without noise:")
print("MIC", mine.mic())np.random.seed(0)
y += np.random.uniform(-1, 1, x.shape[0])  # add some noise
mine.compute_score(x, y)print("With noise:")
print("MIC", mine.mic())

结果为：

Without noise: MIC 1.0000000000000002

With noise: MIC 0.5057166934173714

代码二：展示如何在sklearn的单变量选择方法中使用该函数。

from minepy import MINE
from numpy import array
from sklearn.datasets import load_iris
from sklearn.feature_selection import SelectKBestirisdata = load_iris()def mic(x, y):m = MINE()m.compute_score(x, y)return (m.mic(), 0.5)#选择 K 个最好的特征，返回特征选择后的数据
irisdata_new =  SelectKBest(lambda X, Y: tuple(map(tuple,array(list(map(lambda x:mic(x, Y), X.T))).T)),k=3).fit_transform(irisdata.data, irisdata.target)print(irisdata.data.shape, irisdata_new.shape)

(150, 4) (150, 3)

代码三：

from __future__ import division
import numpy as np
import matplotlib.pyplot as plt
from minepy import MINEdef mysubplot(x, y, numRows, numCols, plotNum,xlim=(-4, 4), ylim=(-4, 4)):r = np.around(np.corrcoef(x, y)[0, 1], 1)mine = MINE(alpha=0.6, c=15)mine.compute_score(x, y)mic = np.around(mine.mic(), 1)ax = plt.subplot(numRows, numCols, plotNum,xlim=xlim, ylim=ylim)ax.set_title('Pearson r=%.1f\nMIC=%.1f' % (r, mic),fontsize=10)ax.set_frame_on(False)ax.axes.get_xaxis().set_visible(False)ax.axes.get_yaxis().set_visible(False)ax.plot(x, y, ',')ax.set_xticks([])ax.set_yticks([])return axdef rotation(xy, t):return np.dot(xy, [[np.cos(t), -np.sin(t)],[np.sin(t), np.cos(t)]])def mvnormal(n=1000):cors = [1.0, 0.8, 0.4, 0.0, -0.4, -0.8, -1.0]for i, cor in enumerate(cors):cov = [[1, cor],[cor, 1]]xy = np.random.multivariate_normal([0, 0], cov, n)mysubplot(xy[:, 0], xy[:, 1], 3, 7, i+1)def rotnormal(n=1000):ts = [0, np.pi/12, np.pi/6, np.pi/4, np.pi/2-np.pi/6,np.pi/2-np.pi/12, np.pi/2]cov = [[1, 1],[1, 1]]xy = np.random.multivariate_normal([0, 0], cov, n)for i, t in enumerate(ts):xy_r = rotation(xy, t)mysubplot(xy_r[:, 0], xy_r[:, 1], 3, 7, i+8)def others(n=1000):x = np.random.uniform(-1, 1, n)y = 4*(x**2-0.5)**2 + np.random.uniform(-1, 1, n)/3mysubplot(x, y, 3, 7, 15, (-1, 1), (-1/3, 1+1/3))y = np.random.uniform(-1, 1, n)xy = np.concatenate((x.reshape(-1, 1), y.reshape(-1, 1)), axis=1)xy = rotation(xy, -np.pi/8)lim = np.sqrt(2+np.sqrt(2)) / np.sqrt(2)mysubplot(xy[:, 0], xy[:, 1], 3, 7, 16, (-lim, lim), (-lim, lim))xy = rotation(xy, -np.pi/8)lim = np.sqrt(2)mysubplot(xy[:, 0], xy[:, 1], 3, 7, 17, (-lim, lim), (-lim, lim))y = 2*x**2 + np.random.uniform(-1, 1, n)mysubplot(x, y, 3, 7, 18, (-1, 1), (-1, 3))y = (x**2 + np.random.uniform(0, 0.5, n)) * \np.array([-1, 1])[np.random.random_integers(0, 1, size=n)]mysubplot(x, y, 3, 7, 19, (-1.5, 1.5), (-1.5, 1.5))y = np.cos(x * np.pi) + np.random.uniform(0, 1/8, n)x = np.sin(x * np.pi) + np.random.uniform(0, 1/8, n)mysubplot(x, y, 3, 7, 20, (-1.5, 1.5), (-1.5, 1.5))xy1 = np.random.multivariate_normal([3, 3], [[1, 0], [0, 1]], int(n/4))xy2 = np.random.multivariate_normal([-3, 3], [[1, 0], [0, 1]], int(n/4))xy3 = np.random.multivariate_normal([-3, -3], [[1, 0], [0, 1]], int(n/4))xy4 = np.random.multivariate_normal([3, -3], [[1, 0], [0, 1]], int(n/4))xy = np.concatenate((xy1, xy2, xy3, xy4), axis=0)mysubplot(xy[:, 0], xy[:, 1], 3, 7, 21, (-7, 7), (-7, 7))plt.figure(facecolor='white')
mvnormal(n=800)
rotnormal(n=200)
others(n=800)
plt.tight_layout()
plt.show()

结果：