I have two black and white images and I need to calculate the mutual information.
Image 1 = X Image 2 = Y
I know that the mutual information can be defined as:
MI = entropy(X) + entropy(Y) - JointEntropy(X,Y)
MATLAB already has built-in functions to calculate the entropy but not to calculate the joint entropy. I guess the true question is: How do I calculate the joint entropy of two images?
Here is an example of the images I'd like to find the joint entropy of:
X = 0 0 0 0 0 0 0 0 1 1 0 0 0 0 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 Y = 0 0 0 0 0 0 0 0 0.38 0.82 0.38 0.04 0 0 0.32 0.82 0.68 0.17 0 0 0.04 0.14 0.11 0 0 0 0 0 0 0
Usage: I=mi(A,B), where A and B are equally sized images/signals. Function hist2 (included) is used to determine the joint histogram of the images/signals. Assumptions: 1) 0*log(0)=0, 2) mutual information is obtained on the intersection between the supports of partial histograms.
As such, the joint entropy can be calculated as: jointEntropy = -sum(jointProb1DNoZero. *log2(jointProb1DNoZero));
mi_discrete_cont(x, y, k) computes the mutual information for continuous variable x and discrete (integer values only) y , with k nearest neighbours. mi_cont_cont(x1, x2, k) computes the mutual information for continuous variables x1 and x2 , with k nearest neighbours.
The mutual information can also be calculated as the KL divergence between the joint probability distribution and the product of the marginal probabilities for each variable. — Page 57, Pattern Recognition and Machine Learning, 2006. This can be stated formally as follows: I(X ; Y) = KL(p(X, Y) || p(X) * p(Y))
To calculate the joint entropy, you need to calculate the joint histogram between two images. The joint histogram is essentially the same as a normal 1D histogram but the first dimension logs intensities for the first image and the second dimension logs intensities for the second image. This is very similar to what is commonly referred to as a co-occurrence matrix. At location (i,j)
in the joint histogram, it tells you how many intensity values we have encountered that have intensity i
in the first image and intensity j
in the second image.
What is important is that this logs how many times we have seen this pair of intensities at the same corresponding locations. For example, if we have a joint histogram count of (7,3) = 2
, this means that when we were scanning both images, when we encountered the intensity of 7
, at the same corresponding location in the second image, we encountered the intensity of 3
for a total of 2
times.
Constructing a joint histogram is very simple to do.
256 x 256
matrix (assuming your image is unsigned 8-bit integer) and initialize them to all zeroes. Also, you need to make sure that both of your images are the same size (width and height). 1
. One would be inclined to do this with for
loops, but as it is commonly known, for
loops are notoriously slow and should be avoided if at all possible. However, you can easily do this in MATLAB in the following way without loops. Let's assume that im1
and im2
are the first and second images you want to compare to. What we can do is convert im1
and im2
into vectors. We can then use accumarray
to help us compute the joint histogram. accumarray
is one of the most powerful functions in MATLAB. You can think of it as a miniature MapReduce paradigm. Simply put, each data input has a key and an associated value. The goal of accumarray
is to bin all of the values that belong to the same key and do some operation on all of these values. In our case, the "key" would be the intensity values, and the values themselves are the value of 1
for every intensity value. We would then want to add up all of the values of 1
that map to the same bin, which is exactly how we'd compute a histogram. The default behaviour for accumarray
is to add all of these values. Specifically, the output of accumarray
would be an array where each position computes the sum of all values that mapped to that key. For example, the first position would be the summation of all values that mapped to the key of 1, the second position would be the summation of all values that mapped to the key of 2 and so on.
However, for the joint histogram, you want to figure out which values map to the same intensity pair of (i,j)
, and so the keys here would be a pair of 2D coordinates. As such, any intensities that have an intensity of i
in the first image and j
in the second image in the same spatial location shared between the two images go to the same key. Therefore in the 2D case, the output of accumarray
would be a 2D matrix where each element (i,j)
contains the summation of all values that mapped to key (i,j)
, similar to the 1D case that was mentioned previously which is exactly what we are after.
In other words:
indrow = double(im1(:)) + 1; indcol = double(im2(:)) + 1; %// Should be the same size as indrow jointHistogram = accumarray([indrow indcol], 1); jointProb = jointHistogram / numel(indrow);
With accumarray
, the first input are the keys and the second input are the values. A note with accumarray
is that if each key has the same value, you can simply assign a constant to the second input, which is what I've done and it's 1
. In general, this is an array with the same number of rows as the first input. Also, take special note of the first two lines. There will inevitably be an intensity of 0
in your image, but because MATLAB starts indexing at 1
, we need to offset both arrays by 1
.
Now that we have the joint histogram, it's really simple to calculate the joint entropy. It is similar to the entropy in 1D, except now we are just summing over the entire joint probability matrix. Bear in mind that it will be very likely that your joint histogram will have many 0
entries. We need to make sure that we skip those or the log2
operation will be undefined. Let's get rid of any zero entries now:
indNoZero = jointHistogram ~= 0; jointProb1DNoZero = jointProb(indNoZero);
Take notice that I searched the joint histogram instead of the joint probability matrix. This is because the joint histogram consists of whole numbers while the joint probability matrix will lie between 0
and 1
. Because of the division, I want to avoid comparing any entries in this matrix with 0
due to numerical roundoff and instability. The above will also convert our joint probability matrix into a stacked 1D vector, which is fine.
As such, the joint entropy can be calculated as:
jointEntropy = -sum(jointProb1DNoZero.*log2(jointProb1DNoZero));
If my understanding of calculating entropy for an image in MATLAB is correct, it should calculate the histogram / probability distribution over 256
bins, so you can certainly use that function here with the joint entropy that was just calculated.
So far, we have assumed that the images that you have dealt with have intensities that are integer-valued. What if we have floating point data? accumarray
assumes that you are trying to index into the output array using integers, but we can still certainly accomplish what we want with this small bump in the road. What you would do is simply assign each floating point value in both images to have a unique ID. You would thus use accumarray
with these IDs instead. To facilitate this ID assigning, use unique
- specifically the third output from the function. You would take each of the images, put them into unique
and make these the indices to be input into accumarray
. In other words, do this instead:
[~,~,indrow] = unique(im1(:)); %// Change here [~,~,indcol] = unique(im2(:)); %// Change here %// Same code jointHistogram = accumarray([indrow indcol], 1); jointProb = jointHistogram / numel(indrow); indNoZero = jointHistogram ~= 0; jointProb1DNoZero = jointProb(indNoZero); jointEntropy = -sum(jointProb1DNoZero.*log2(jointProb1DNoZero));
Note that with indrow
and indcol
, we are directly assigning the third output of unique
to these variables and then using the same joint entropy code that we computed earlier. We also don't have to offset the variables by 1 as we did previously because unique
will assign IDs starting at 1.
You can actually calculate the histograms or probability distributions for each image individually using the joint probability matrix. If you wanted to calculate the histograms / probability distributions for the first image, you would simply accumulate all of the columns for each row. To do it for the second image, you would simply accumulate all of the rows for each column. As such, you can do:
histogramImage1 = sum(jointHistogram, 1); histogramImage2 = sum(jointHistogram, 2);
After, you can calculate the entropy of both of these by yourself. To double check, make sure you turn both of these into PDFs, then compute the entropy using the standard equation (like above).
To finally compute Mutual Information, you're going to need the entropy of the two images. You can use MATLAB's built-in entropy
function, but this assumes that there are 256 unique levels. You probably want to apply this for the case of there being N
distinct levels instead of 256, and so you can use what we did above with the joint histogram, then computing the histograms for each image in the aside code above, and then computing the entropy for each image. You would simply repeat the entropy calculation that was used jointly, but apply it to each image individually:
%// Find non-zero elements for first image's histogram indNoZero = histogramImage1 ~= 0; %// Extract them out and get the probabilities prob1NoZero = histogramImage1(indNoZero); prob1NoZero = prob1NoZero / sum(prob1NoZero); %// Compute the entropy entropy1 = -sum(prob1NoZero.*log2(prob1NoZero)); %// Repeat for the second image indNoZero = histogramImage2 ~= 0; prob2NoZero = histogramImage2(indNoZero); prob2NoZero = prob2NoZero / sum(prob2NoZero); entropy2 = -sum(prob2NoZero.*log2(prob2NoZero)); %// Now compute mutual information mutualInformation = entropy1 + entropy2 - jointEntropy;
Hope this helps!
If you love us? You can donate to us via Paypal or buy me a coffee so we can maintain and grow! Thank you!
Donate Us With