MI.py 1011 B

123456789101112131415161718192021222324252627282930313233343536
  1. # -*- coding: utf-8 -*-
  2. """
  3. Created on Wed Nov 10 23:33:33 2021
  4. @author: arefk
  5. """
  6. from __future__ import print_function # print('me') instead of print 'me'
  7. from __future__ import division # 1/2 == 0.5, not 0
  8. # - import common modules
  9. import numpy as np # the Python array package
  10. def mutualInfo(Im1,Im2):
  11. t1_slice = Im1
  12. t2_slice = Im2
  13. hist_2d, x_edges, y_edges = np.histogram2d(t1_slice.ravel(),t2_slice.ravel(),bins=20)
  14. hist_2d_log = np.zeros(hist_2d.shape)
  15. non_zeros = hist_2d != 0
  16. hist_2d_log[non_zeros] = np.log(hist_2d[non_zeros])
  17. pxy = hist_2d / float(np.sum(hist_2d))
  18. px = np.sum(pxy, axis=1) # marginal for x over y
  19. py = np.sum(pxy, axis=0) # marginal for y over x
  20. px_py = px[:, None] * py[None, :] # Broadcast to multiply marginals
  21. # Now we can do the calculation using the pxy, px_py 2D arrays
  22. nzs = pxy > 0 # Only non-zero pxy values contribute to the sum
  23. MI = np.sum(pxy[nzs] * np.log(pxy[nzs] / px_py[nzs]))
  24. return MI