ITK vs PyTorch calculation for Mutual Information

Beginner Questions
1

Hi All,

I’m calculating the mutual information metric value between two 3D images using ITK and PyTorch and getting two different values for the same bin size. The sample values are (ITK and Pytorch):

  1. 0.5314091942961909, 0.5432233810424805
  2. 0.578950776982905, 0.46988388895988464
  3. 0.576229949241947, 0.4693763554096222

Please refer below for the functions of ITK (5.3.0) and PyTorch (2.0.1) versions:

# Computes similarity between two images using ITK
def get_metric(fixed, moving):
    fixed = convert_image_type(fixed, itk.UC, itk.F)
    moving = convert_image_type(moving, itk.UC, itk.F)
    if similarity_metric == "NCC":
        metric = itk.CorrelationImageToImageMetricv4[
            fixed, moving
        ].New()
    elif similarity_metric == "MI":
        metric = itk.MattesMutualInformationImageToImageMetricv4[
            fixed, moving
        ].New()
        metric.SetNumberOfHistogramBins(no_of_bins_MI)
    elif similarity_metric == "MSE":
        metric = itk.MeanSquaresImageToImageMetricv4[
            fixed, moving
        ].New()
    transform = itk.Euler3DTransform[itk.D].New()
    interpolator = itk.LinearInterpolateImageFunction[type(
        fixed), itk.D].New()
    metric.SetFixedImage(fixed)
    metric.SetMovingImage(moving)
    metric.SetTransform(transform)
    metric.SetMovingInterpolator(interpolator)
    metric.Initialize()
    return metric

# Computes MI similarity between two images using PyTorch
def get_mi_metric(image_1, image_2):   
    bins = no_of_bins_MI
    joint_hist = torch.histc(image_1 * bins + image_2, bins=bins*bins, min=0, max=bins*bins-1)
    joint_hist = joint_hist.view(bins, bins)
    
    joint_prob = joint_hist / torch.sum(joint_hist)
    
    x_hist = torch.sum(joint_hist, dim=1)
    y_hist = torch.sum(joint_hist, dim=0)
    
    x_prob = x_hist / torch.sum(x_hist)
    y_prob = y_hist / torch.sum(y_hist)
    
    mutual_info = 0.0
    for i in range(bins):
        for j in range(bins):
            if joint_prob[i, j] > 0:
                mutual_info += joint_prob[i, j] * torch.log(joint_prob[i, j] / (x_prob[i] * y_prob[j] + 1e-10) + 1e-10)
    
    return -1 * mutual_info

Could you please point out why there is a difference between these two values?

Thank you!