MIC - 最大信息系数

MIC

前言

皮尔逊相关系数即我们通常说的(线性)相关系数，是用来反映两个变量线性相关程度的统计量，变化范围为-1到1。

• 系数的值为1意味着X和Y可以很好的由直线方程来进行描述，所有的数据点都很好的落在一条直线上，且随着X的增加而增加。
• 系数的值为−1意味着所有的数据点都落在直线上，且随着X的增加而减少。系数的值为0意味着两个变量之间没有线性关系。

斯皮尔曼相关系数又称秩相关系数，是衡量两个变量的依赖性的非参数指标。它利用两变量的秩次大小作线性相关分析，对原始变量的分布不作要求，属于非参数统计方法，适用范围要广些。对于服从皮尔逊相关系数的数据也可以计算斯皮尔曼相关系数，但统计效能要低一些。

互信息(Mutual Information)是信息论里一种有用的信息度量，它可以看成是一个随机变量中包含的关于另一个随机变量的信息量，或者说是一个随机变量由于已知另一个随机变量而减少的不肯定性。

​ 设两个随机变量$（X,Y）$的联合分布为 $p(x,y)$，边缘分布分别为 $p(x)$、$p(y)$，互信息 $I(X;Y)$ 是联合分布与边缘分布的相对熵。计算公式为： $I(X;Y)={\sum }_{x\in X}{\sum }_{y\in Y}p(x,y)lo{g}_2\frac{p(x,y)}{p(x)p(y)}$

MIC介绍

MIC（Maximal Information Coefficient）最大信息系数，是检测变量之间非线性相关性的新方法。它的值域在 0 和 1 之间，值越高表示相关性越强。

• 广泛性（普适性），不管是什么函数关系，都可以识别。线性和非线性函数关系，对这个系数而言是一样的。
• 公平性，是指在样本量足够大时能为不同类型单噪声程度相似的相关关系给出相近的系数。

注意，MIC的计算理念是针对大数据集提出来的，是针对大数据之间的计算的。只有在大数据下，MIC的计算结果才近似等于真实值，才相对可靠。对于小数据集的MIC计算没有意义，因为误差会非常大。

常用系数对比

	适用范围	是否标准化	计算复杂度	鲁棒性
皮尔逊系数	线性数据	是	低	低
斯皮尔曼系数	线性数据、简单单调非线性数据	是	低	中等
K近邻	线性数据、非线性数据	否	高	高
MIC	线性数据、非线性数据	是	低	高

从中可以看出，MIC是常用系数中的较优的一种方案，适用范围广、复杂度低、鲁棒性高。

---

在机器学习中，想把互信息直接用于特征选择其实不是很方便：

1. 它不属于度量方式，也没有办法归一化，在不同数据及上的结果无法做比较
2. 对于连续变量的计算不是很方便（X和Y都是集合，x，y都是离散的取值），通常变量需要先离散化，而互信息的结果对离散化的方式很敏感

而最大信息系数 MIC 克服了这两个问题。它首先寻找一种最优的离散化方式，然后把互信息取值转换成一种度量方式，取值区间在[0，1]。

MIC简单计算方式：

1. 给定 i、j，对XY构成的散点图进行 i 列 j 行网格化，并求出最大互信息值。$I(x;y)\approx I[X;Y]$
2. 对最大互信息值进行归一化。 $mic(x;y)=ma{x}_{a\ast b<B}\frac{I(x;y)}{lo{g}_{2}min(a,b)}$
3. 选择不同 i、j 尺度下互信息的最大值作为MIC值

（用 a,b 表示在 x,y 方向上的划分格子的个数，设置a*b<B，其中 B 的大小为数据量的0.6次方左右，这是一个经验值）

---

MIC库

MIC算法在C、C++、Python中都有对应的实现，分别在头文件<mine.h>、<cppmine.h>和 库 <minepy>中

下面我们重点介绍 Python API

class minepy.MINE(alpha=0.6, c=15, est=“mic_approx”)

基于信息的最大非参数探索

参数：

• alpha：float类型，取值范围为（0 ，1.0 ] 或 > = 4，如果alpha的取值范围在（0,1]之内，那么B的取值范围为（N^αlpha，4）其中n是样本的数目。如果alpha的取值范围是> = 4， alpha直接定义B参数。如果alpha高于样本数 n，则它将被限制为n，因此 B = min（alpha，n）
• c：float类型，取值必须大于0，确定每个分区中的列数比列数多了多少。默认值为15，这意味着当尝试在x轴上绘制x网格线时，算法将以最多15 * x块开始。
• est：估算器，参数为“mic_approx”或者“mic_e”。使用est =” mic_approx”，将计算原始MINE统计信息；使用est =” mic_e”，将计算等特征矩阵，并且 mic() 和 tic() 方法将分别返回 MIC_e 和 TIC_e 值。

方法：

• compute_score(x, y)：计算(equi)特征矩阵(即最大标准化互信息得分)
• mic()：返回最大的信息系数(MIC或MIC_e)
• mas()：返回最大不对称分数（MAS）
• mev()：返回最大边缘值（MEV）
• mcn(eps=0)：返回eps> = 0的最小单元数（MCN）
• gmic(p=-1)：返回广义最大信息系数（GMIC）
• tic(norm=false)：返回总信息系数（TIC或TIC_e）。如果norm == True TIC将在[0，1]中标准化

好用的功能函数

minepy.pstats（X，alpha = 0.6，c = 15，est =“ mic_approx”）

计算变量之间的成对统计量（MIC和归一化的TIC）

对于每个统计量，矩阵的上三角按行存储（压缩矩阵）。如果m是变量的数量，则对于i <j <m，（行）i和j之间的统计信息存储在

k = m * i-i *（i + 1）/ 2-i-1 + j中。向量的长度为n = m *（m-1）/ 2。

参数：

• X：由n个变量和m个样本组成的n×m数组

• alpha、c、est 与 minepy.MINE 一致

返回值：

• mic：长度为n *（n-1）/ 2的压缩MIC统计矩阵
• tic：长度为n *（n-1）/ 2的压缩归一化TIC统计矩阵

minepy.cstats（X，Y，alpha = 0.6，c = 15，est =“ mic_approx” ）

计算两个变量集合中的每一对之间的统计信息（MIC和归一化的TIC）

如果n和m分别是X和Y中的变量数，则（行）i（对于X）和j（对于Y）之间的统计信息存储在mic [i，j]和tic [i，j]中。

参数：

• X：由n个变量和m个样本组成的n×m数组

• Y：由p个变量和m个样本组成的p x m数组

• alpha、c、est 与 minepy.MINE 一致

返回值：

• mic：MIC统计矩阵
• tic：标准化的TIC统计矩阵

Python实例

安装 Python 的 Minepy 库请参考这篇文章 教你在windows10环境下如何安装minepy并成功运行！

API给出了三个 MIC 使用的例子（位于examples/python_example.py中）

例1（简化版）

import numpy as np
from minepy import MINE

# np.linspace(start, end, step)主要用来创建等差数列，返回在start到end之间间隔均匀的数据
x = np.linspace(0, 1, 1000)
# np.sin(x) 对矩阵x中每个元素取正弦，np.pi代表π
y = np.sin(10 * np.pi * x) + x

mine = MINE(alpha=0.6, c=15)
mine.compute_score(x, y)

print("Without noise:")
# 得到MIC的值
print("MIC", mine.mic())

np.random.seed(0)
# numpy.random.uniform(low,high,size)从一个均匀分布[low,high)中随机采样，size为输出样本数目
y += np.random.uniform(-1, 1, x.shape[0])  # 添加一些噪声，x.shape[0]是数组的行数
# x、y之间的最大归一化互信息得分
mine.compute_score(x, y)

print("With noise:")
print("MIC", mine.mic())

例2

from __future__ import division
import numpy as np
import matplotlib.pyplot as plt
from minepy import MINE

# np.around(a, decimals=0) 对输入数组 a 的元素返回四舍五入后的值
# decimals 保留小数点后 n 位。 默认值为0。
# np.corrcoef()返回皮尔逊相关系数
def 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=7)
    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 ax

# np.dot()处理一维数组则得到两数组的内积，处理矩阵则得到矩阵乘积
def 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]
    # enumerate()用于将一个可遍历的数据对象(如列表、字符串)组合为一个索引序列，同时列出数据和数据下标，一般用在for 循环当中。
    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]]
    # np.random.multivariate_normal方法用于根据实际情况生成一个多元正态分布矩阵
    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) / 3
    mysubplot(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))
    # np.concatenate((a1, a2, ...), axis=0)函数。能够一次完成多个数组的拼接
    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.savefig('MIC.png')
plt.show()

从图中我们可以看出

• 在数据线性相关时，皮尔逊系数和 MIC 都有较好的表现
• 在数据非线性相关时，MIC 表现仍然出色

例3

import numpy as np
from minepy import pstats, cstats

np.random.seed(0)

# 随机建立X矩阵，8个变量，320个样本，样本值在[0,1)之间，可以近似理解为8行320列
X = np.random.rand(8, 320)

# 建立Y矩阵，4个变量，320个样本
Y = np.random.rand(4, 320)

# 计算 X 中样本间的两两统计量 MIC_e 和 归一化 TIC_e
# B=9, c=5
mic_p, tic_p = pstats(X, alpha=9, c=5, est="mic_e")

# 计算 X 和 Y 中每对样本之间的统计量
mic_c, tic_c = cstats(X, Y, alpha=9, c=5, est="mic_e")

print("normalized TIC_e (X):")
print(tic_p)
print("MIC_e (X vs. Y):")
print(mic_c)

MIC缺陷

然而，“天下没有免费的午餐”（no free lunch），机器学习领域的同学应当都知道这句话。

虽然 MIC 可以有效地计算出非线性数据之间的相关性，但是这并不意味着 MIC 是万能的。MIC的统计能力也遭到了一些学者的质疑，通过实验发现，当零假设不成立的时候，MIC的计算结果会受到非常大的影响。

零假设（null hypothesis），统计学术语，又称原假设，指进行统计检验时预先建立的假设。 零假设成立时，有关统计量应服从已知的某种概率分布。

当统计量的计算值落入否定域时，可知发生了小概率事件，应否定原假设。

Reshef等人精心打造了一些数据集，将MIC的统计能力与标准皮尔逊相关和距离相关进行了对比，发现MIC在有些数据集上的相关性不如皮尔逊强，线性相关的数据集尤为严重。（这并不意味着MIC的统计能力就不好，而是在有些数据集上存在这样的问题）

他们得出一个结论，如果将MIC用于大规模的探索性分析，可能会产生不少的错误，MIC的 “公平性” 并不是真正的公平。这也在提示我们，具体问题需要具体分析。

以下是文章正文（部分）和实验结果的截图，有兴趣的小伙伴可以再深入阅读和思考。

文章链接 http://statweb.stanford.edu/~tibs/reshef/comment.pdf

参考文章

最大信息系数详解

MIC：最大信息系数

minepy - Maximal Information-based Nonparametric Exploration

特征选择 (feature_selection)