Number of bins when computing mutual information

Question

I want to quantify the relationship between two variables, A and B, using mutual information. The way to compute it is by binning the observations (see example Python code below). However, what factors determines what number of bins is reasonable? I need the computation to be fast so I cannot simply use a lot of bins to be on the safe side.

from sklearn.metrics import mutual_info_score

def calc_MI(x, y, bins):
    c_xy = np.histogram2d(x, y, bins)[0]
    mi = mutual_info_score(None, None, contingency=c_xy)
    return mi

Answer 1

There is no best number of bins to estimate mutual information (MI) with histograms. The best way is to choose it via cross-validation if you can, or to rely on a rule of thumb. This the reason why many other estimators of MI which are not based on histograms have been proposed.

The number of bins will depend to the total number of data points $n$. You should try to avoid too many bins to avoid estimation errors for the joint distribution between the two variables. You should also avoid too few bins to be able to capture the relationship between the two variables. Given that np.histogram2d(x, y, D) generates a 2D histogram with D equal width bins for both x and y I would personally choose: $$ D = \lfloor \sqrt{n/5} \rfloor$$ In this case on average for two uniformly distributed random variables you will have at least $5$ points for each cell of the histogram: $$ \frac{n}{D_X D_Y} \geq 5 \Rightarrow \frac{n}{D^2} \geq 5 \Rightarrow D^2 \leq n/5 \Rightarrow D = \lfloor \sqrt{n/5} \rfloor$$ This is one possible choice that simulates the adaptive partitioning approach proposed in (Cellucci, 2005). The latter approach is often used to estimate MI to infer genetic networks: e.g. in MIDER.

If you have lots of data points $n$ and no missing values you should not worry too much about finding the best number of bins; e.g. if $n = 100,000$. If this is not the case, you might consider to correct MI for finite samples. (Steuer et al., 2002) discusses some correction for MI for the task of genetic network inference.

---

Estimating the number of bins for a histogram is an old problem. You might be interested in this talk by Lauritz Dieckman about estimating the number of bins for MI. This talk is based on a chapter in Mike X Cohen's book about neural time-series.

You might choose $D_X$ and $D_Y$ independently and use the rule of thumb used for the estimating the number of bins in 1D histograms.

Freedman-Diaconis' rule (no assumption on the distribution): $$D_X = \lceil \frac{\max{X} - \min{X}}{2 \cdot \mbox{IQR} \cdot n^{-1/3}} \rceil$$ where $\mbox{IQR}$ is the difference between the 75-quantile and the 25-quantile. Look at this related question in SE.

Scott's rule (normality assumption): $$D_X = \lceil \frac{\max{X} - \min{X}}{3.5 \cdot s_X \cdot n^{-1/3}} \rceil$$ where $s_X$ is the standard deviation for $X$.

Sturges' rule (might underestimate the number of bins but good for large $n$): $$D_X = \lceil 1 + \log_2{n} \rceil$$

---

It is difficult to correctly estimate MI with histograms. You might then choose a different estimator:

• Kraskov's $k$NN estimator, which is a bit less sensitive to parameter choice: $k = 4$ or $k = 6$ nearest neighbours is often used as default. Paper: (Kraskov, 2003)
• Estimation of MI with Kernels (Moon, 1995).

There are lots of packages for estimating MI:

• Non-Parametric Entropy Estimation Toolbox for Python. site.
• Information-dynamics toolkit in Java but available also for Python. site.
• ITE toolbox in Matlab. site.

Answer 2

The following binning rules should be added to Simone's list, which have become even more commonplace:

Given that mutual information is the sum of marginal entropies adjusted by their joint entropy, $$I(X,Y) = H(X) + H(Y) - H(X,Y) $$

The optimal binning rule for marginal entropy $H(X)$, as well as $H(Y)$, found by Hacine-Gharbi et al. (2012) is

$$B_X = round\bigg(\frac{\xi}{6} + \frac{2}{3\xi} + \frac{1}{3} \bigg) $$ where $\xi = \big( 8 + 324N + 12 \sqrt{36N + 729N^2}\big)^{\frac{1}{3}} $

while the optimal binning rule for joint entropy $H(X,Y)$ according to Hacine-Gharbi and Ravier (2018) is

$$B_X = B_Y = round\Bigg[ \frac{1}{\sqrt{2}} \Bigg(1 + \sqrt{1+\frac{24N}{1-\rho^2}} \Bigg)^{\frac{1}{2}} \Bigg] $$

Applying these binning rules when measuring the individual terms of $I(X,Y)=H(X)+H(Y)−H(X,Y)$, you should have an optimally binned low-bias estimator of mutual information.

Hacine-Gharbi, A., and P. Ravier (2018): “A Binning Formula of Bi-histogram
for Joint Entropy Estimation Using Mean Square Error Minimization.”
Pattern Recognition Letters, Vol. 101, pp. 21–28.

Hacine-Gharbi, A., P. Ravier, R. Harba, and T. Mohamadi (2012): “Low Bias
Histogram-Based Estimation of Mutual Information for Feature Selection.”
Pattern Recognition Letters, Vol. 33, pp. 1302–8.

Answer 3

I prefer minepy to get and estimate of mutual information in python.

You can see the implementation details of the package here, and an example code here. For the sake of easier reference, I copy paste the example and it's output here:

import numpy as np
from minepy import MINE

def print_stats(mine):
    print "MIC", mine.mic()
    print "MAS", mine.mas()
    print "MEV", mine.mev()
    print "MCN (eps=0)", mine.mcn(0)
    print "MCN (eps=1-MIC)", mine.mcn_general()

x = 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_stats(mine)
print

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_stats(mine)

Which gives this as output:

Without noise:
MIC 1.0
MAS 0.726071574374
MEV 1.0
MCN (eps=0) 4.58496250072
MCN (eps=1-MIC) 4.58496250072

With noise:
MIC 0.505716693417
MAS 0.365399904262
MEV 0.505716693417
MCN (eps=0) 5.95419631039
MCN (eps=1-MIC) 3.80735492206

My experience is that the results are sensitive to alpha, and the default value .6 is a reasonable one. However, on my real data alpha=.3 is much faster and the estimated mutual informations have a really high correlation with the case that alpha=.6. So in case you're using MI to select the ones with a high MI, you can simply use a smaller alpha and use the highest values as a replacement with a good accuracy.