変数の独立性を測る指標として有名なものにピアソンの相関係数(Pearson's correlation coefficient)があります。ですがそれはあくまで線形な(直線的な)関係を扱う指標で、たとえば二次関数のようなU字型の関係があるときには「相関がない」と判定してしまいます。このような非線形な関係を取り扱う指標として、MICとHSICが知られているので、そのPythonパッケージを使ってみました。
以下のコードは全て Google Colaboratory で動作確認済みです。
MIC (Maximal Information Coefficient) を計算するPythonパッケージ
次のようにして pip install できます。
!pip install minepy
インストールできたら、次のようにして準備しておきます。
from minepy import MINE
mine = MINE()
HSIC (Hilbert-Schmidt Independence Criterion) を計算するPythonパッケージ
次のように git clone します。
!git clone https://github.com/amber0309/HSIC.git
次のようにして準備しておきます。
from HSIC.HSIC import hsic_gam
一次関数
まずは一次関数について調べてみます。
%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np
X = np.linspace(-1, 1, 51)
Y = X
plt.scatter(X, Y)
ピアソンの相関係数は
np.corrcoef(X, Y)[0, 1]
1.0
MICは
mine.compute_score(X, Y)
mine.mic()
0.9999999999999998
HSSCは、次の値が大きいほど独立性が低く、小さいほど独立性が高いそうです。
testStat, thresh = hsic_gam(X.reshape(len(X), 1), Y.reshape(len(Y), 1))
testStat / thresh
22.610487654066624
一次関数にノイズを加えてみましょう。
import numpy as np
X = np.linspace(-1, 1, 51)
Y = X + np.random.rand(51,)
plt.scatter(X, Y)
ピアソンの相関係数は
np.corrcoef(X, Y)[0, 1]
0.9018483832673769
MICは
mine.compute_score(X, Y)
mine.mic()
0.9997226475394071
HSICは
testStat, thresh = hsic_gam(X.reshape(len(X), 1), Y.reshape(len(Y), 1))
testStat / thresh
13.706487127525389
ノイズの大きさを少しずつ増やしていったときに、ピアソンの相関係数がどう変化するか計算してみましょう。
D = []
C = []
for d in range(100):
X = np.linspace(-1, 1, 51)
Y = X + np.random.rand(51,) * d / 10
D.append(d)
C.append(np.corrcoef(X, Y)[0, 1])
plt.plot(D, C)
同様に、ノイズの大きさを少しずつ増やしていったときに、MICがどう変化するか計算してみましょう。
D = []
C = []
for d in range(100):
X = np.linspace(-1, 1, 51)
Y = X + np.random.rand(51,) * d / 10
D.append(d)
mine.compute_score(X, Y)
C.append(mine.mic())
plt.plot(D, C)
同様に、ノイズの大きさを少しずつ増やしていったときに、HSICがどう変化するか計算してみましょう。
D = []
C = []
for d in range(100):
X = np.linspace(-1, 1, 51)
Y = X + np.random.rand(51,) * d / 10
D.append(d)
testStat, thresh = hsic_gam(X.reshape(len(X), 1), Y.reshape(len(Y), 1))
C.append(testStat / thresh)
plt.plot(D, C)
二次関数
次は二次関数です。
import numpy as np
X = np.linspace(-1, 1, 51)
Y = X**2
plt.scatter(X, Y)
ピアソンの相関係数は
np.corrcoef(X, Y)[0, 1]
-2.3862043230836297e-17
MICは
mine.compute_score(X, Y)
mine.mic()
0.9997226475394071
HSICは
testStat, thresh = hsic_gam(X.reshape(len(X), 1), Y.reshape(len(Y), 1))
testStat / thresh
9.285886928766184
ノイズを加えてみましょう。
X = np.linspace(-1, 1, 51)
Y = X**2 + np.random.rand(51,)
plt.scatter(X, Y)
ピアソンの相関係数は
np.corrcoef(X, Y)[0, 1]
0.04693076570744622
MICは
mine.compute_score(X, Y)
mine.mic()
0.5071787519579662
HSICは
testStat, thresh = hsic_gam(X.reshape(len(X), 1), Y.reshape(len(Y), 1))
testStat / thresh
5.247036645612692
ノイズの大きさを少しずつ増やしていったときに、ピアソンの相関係数がどう変化するか計算してみましょう。
D = []
C = []
for d in range(100):
X = np.linspace(-1, 1, 51)
Y = X**2 + np.random.rand(51,) * d / 10
D.append(d)
C.append(np.corrcoef(X, Y)[0, 1])
plt.plot(D, C)
同様に、ノイズの大きさを少しずつ増やしていったときに、MICがどう変化するか計算してみましょう。
D = []
C = []
for d in range(100):
X = np.linspace(-1, 1, 51)
Y = X**2 + np.random.rand(51,) * d / 10
D.append(d)
mine.compute_score(X, Y)
C.append(mine.mic())
plt.plot(D, C)
同様に、ノイズの大きさを少しずつ増やしていったときに、HSICがどう変化するか計算してみましょう。
D = []
C = []
for d in range(100):
X = np.linspace(-1, 1, 51)
Y = X**2 + np.random.rand(51,) * d / 10
D.append(d)
testStat, thresh = hsic_gam(X.reshape(len(X), 1), Y.reshape(len(Y), 1))
C.append(testStat / thresh)
plt.plot(D, C)
ピアソンの相関係数が二次関数の関係性を検出できていないのに対して、MICもHSICも検出できていることがわかります。また、MICは小さなノイズにはほとんど影響を受けていないだろうことも分かります。
三角関数
三角関数も見てみましょう。
import numpy as np
X = np.linspace(-1, 1, 51)
Y = np.sin((X)*np.pi*2)
plt.scatter(X, Y)
ピアソンの相関係数は
np.corrcoef(X, Y)[0, 1]
-0.37651692742033543
MICは
mine.compute_score(X, Y)
mine.mic()
0.9997226475394071
HSICは
testStat, thresh = hsic_gam(X.reshape(len(X), 1), Y.reshape(len(Y), 1))
testStat / thresh
3.212604441715373
ノイズを加えてみましょう。
import numpy as np
X = np.linspace(-1, 1, 51)
Y = np.sin((X+0.25)*np.pi*2) + np.random.rand(51,)
plt.scatter(X, Y)
ピアソンの相関係数は
np.corrcoef(X, Y)[0, 1]
0.049378342546505014
MICは
mine.compute_score(X, Y)
mine.mic()
0.9997226475394071
HSICは
testStat, thresh = hsic_gam(X.reshape(len(X), 1), Y.reshape(len(Y), 1))
testStat / thresh
1.3135420547013614
ノイズをだんだん大きくしていったときの影響を調べます。ピアソンの相関係数は
D = []
C = []
for d in range(100):
X = np.linspace(-1, 1, 51)
Y = np.sin((X+0.25)*np.pi*2) + np.random.rand(51,) * d / 10
D.append(d)
C.append(np.corrcoef(X, Y)[0, 1])
plt.plot(D, C)
MICは
D = []
C = []
for d in range(100):
X = np.linspace(-1, 1, 51)
Y = np.sin((X+0.25)*np.pi*2) + np.random.rand(51,) * d / 10
D.append(d)
mine.compute_score(X, Y)
C.append(mine.mic())
plt.plot(D, C)
HSICは
D = []
C = []
for d in range(100):
X = np.linspace(-1, 1, 51)
Y = np.sin((X+0.25)*np.pi*2) + np.random.rand(51,) * d / 10
D.append(d)
testStat, thresh = hsic_gam(X.reshape(len(X), 1), Y.reshape(len(Y), 1))
C.append(testStat / thresh)
plt.plot(D, C)
ピアソンの相関係数が、ノイズ付きの三角関数の関係性を検出できなかったのは想定内でしたが、HSICも脱落してしまった感じですね。その中でもMICが関係性を検出できていたのは驚きでした。
シグモイド関数
最後に、シグモイド関数です。
import numpy as np
X = np.linspace(-1, 1, 51)
Y = 1 / (1 + np.exp(-X*10))
plt.scatter(X, Y)
ピアソンの相関係数は
np.corrcoef(X, Y)[0, 1]
0.9354020629807919
MICは
mine.compute_score(X, Y)
mine.mic()
0.9999999999999998
HSICは
testStat, thresh = hsic_gam(X.reshape(len(X), 1), Y.reshape(len(Y), 1))
testStat / thresh
27.433932271008874
ノイズを加えてみましょう。
import numpy as np
X = np.linspace(-1, 1, 51)
Y = 1 / (1 + np.exp(-X*10)) + np.random.rand(51,)
plt.scatter(X, Y)
ピアソンの相関係数は
np.corrcoef(X, Y)[0, 1]
0.7994507640548881
MICは
mine.compute_score(X, Y)
mine.mic()
0.9997226475394071
HSICは
testStat, thresh = hsic_gam(X.reshape(len(X), 1), Y.reshape(len(Y), 1))
testStat / thresh
12.218008194240534
ノイズを増やしていった影響を見てみましょう。ピアソンの相関係数は
D = []
C = []
for d in range(100):
X = np.linspace(-1, 1, 51)
Y = 1 / (1 + np.exp(-X*10)) + np.random.rand(51,) * d / 10
D.append(d)
C.append(np.corrcoef(X, Y)[0, 1])
plt.plot(D, C)
MICは
D = []
C = []
for d in range(100):
X = np.linspace(-1, 1, 51)
Y = 1 / (1 + np.exp(-X*10)) + np.random.rand(51,) * d / 10
D.append(d)
mine.compute_score(X, Y)
C.append(mine.mic())
plt.plot(D, C)
HSICは
D = []
C = []
for d in range(100):
X = np.linspace(-1, 1, 51)
Y = 1 / (1 + np.exp(-X*10)) + np.random.rand(51,) * d / 10
D.append(d)
testStat, thresh = hsic_gam(X.reshape(len(X), 1), Y.reshape(len(Y), 1))
C.append(testStat / thresh)
plt.plot(D, C)
いずれの指標も、関係性を検出できていますが、ノイズに対する強さで言えばMICが最強という感じです。
個人の感想です
MIC MIC にしてやんよ!