要计算联合熵，需要计算两幅图像之间的联合直方图。联合直方图与正常的一维直方图基本相同，但第一维记录第一张图像的强度，第二维记录第二张图像的强度。这与通常所说的co-occurrence matrix 非常相似。在联合直方图中的位置(i,j)，它告诉您我们遇到了多少强度值，在第一张图像中强度为i，在第二张图像中强度为j。
重要的是，它记录了我们在在相同的相应位置看到这对强度的次数。例如，如果我们的联合直方图计数为(7,3) = 2，这意味着当我们扫描两张图像时，当我们遇到7 的强度时，在第二张图像的相同对应位置，我们遇到的强度为3 共2 次。
构建联合直方图非常简单。
- 首先，创建一个
256 x 256 矩阵（假设您的图像是无符号8 位整数）并将它们初始化为全零。此外，您需要确保两个图像的大小（宽度和高度）相同。
- 完成此操作后，查看每个图像的第一个像素，我们将其表示为左上角。具体来说，查看该位置的第一张和第二张图像的强度。第一张图像的强度将作为行，而第二张图像的强度将作为列。
- 在矩阵中找到这个位置并将矩阵中的这个位置增加
1。
- 对图像中的其余位置重复此操作。
- 完成后，将所有条目除以任一图像中的元素总数（记住它们的大小应该相同）。这将为我们提供两幅图像之间的联合概率分布。
人们倾向于使用for 循环来执行此操作，但众所周知，for 循环非常缓慢，应尽可能避免使用。但是，您可以在 MATLAB 中通过以下方式轻松做到这一点没有循环。假设im1 和im2 是您要比较的第一张和第二张图片。我们可以做的是将im1 和im2 转换为向量。然后我们可以使用accumarray 来帮助我们计算联合直方图。 accumarray 是 MATLAB 中最强大的函数之一。您可以将其视为一个微型 MapReduce 范例。简单地说，每个数据输入都有一个键和一个关联的值。 accumarray 的目标是对属于同一键的所有值进行分箱，并对所有这些值进行一些操作。在我们的例子中，“关键”是强度值，而值本身就是每个强度值的1 值。然后，我们希望将映射到同一 bin 的 1 的所有值添加，这正是我们计算直方图的方式。 accumarray 的默认行为是添加所有这些值。具体来说，accumarray 的输出将是一个数组，其中每个位置计算映射到该键的所有值的总和。例如，第一个位置是映射到键 1 的所有值的总和，第二个位置是映射到键 2 的所有值的总和，依此类推。
但是，对于联合直方图，您想弄清楚哪些值映射到同一强度对 (i,j)，因此这里的键将是一对 2D 坐标。因此，在第一张图像中具有i 强度和在第二张图像中具有j 强度的任何强度在两个图像之间共享的相同空间位置 转到相同的键。因此，在二维情况下，accumarray 的输出将是一个二维矩阵，其中每个元素 (i,j) 包含映射到键 (i,j) 的所有值的总和，类似于前面提到的一维情况，这正是我们追求。
换句话说：
indrow = double(im1(:)) + 1;
indcol = double(im2(:)) + 1; %// Should be the same size as indrow
jointHistogram = accumarray([indrow indcol], 1);
jointProb = jointHistogram / numel(indrow);
对于accumarray，第一个输入是键，第二个输入是值。 accumarray 的一个注释是，如果每个键都有 same 值，您可以简单地为第二个输入分配一个常量，这就是我所做的，它是 1。通常，这是一个与第一个输入具有相同行数的数组。另外，请特别注意前两行。您的图像中不可避免地会有0 的强度，但由于 MATLAB 从1 开始索引，我们需要将两个数组偏移1。
现在我们有了联合直方图，计算联合熵真的很简单。它类似于 1D 中的熵，只是现在我们只是对整个联合概率矩阵求和。请记住，您的联合直方图很可能有许多 0 条目。我们需要确保跳过这些，否则log2 操作将未定义。现在让我们摆脱任何零条目：
indNoZero = jointHistogram ~= 0;
jointProb1DNoZero = jointProb(indNoZero);
请注意，我搜索的是联合直方图而不是联合概率矩阵。这是因为联合直方图由整数组成，而联合概率矩阵将位于0 和1 之间。由于除法，由于数值舍入和不稳定性，我想避免将此矩阵中的任何条目与0 进行比较。以上还将我们的联合概率矩阵转换为堆叠的一维向量，这很好。
因此，联合熵可以计算为：
jointEntropy = -sum(jointProb1DNoZero.*log2(jointProb1DNoZero));
如果我对在 MATLAB 中计算图像熵的理解是正确的，它应该计算256 bins 上的直方图/概率分布，因此您当然可以在此处使用该函数和刚刚计算的联合熵。
如果我们有浮点数据怎么办？
到目前为止，我们假设您处理的图像具有整数值的强度。如果我们有浮点数据怎么办？ accumarray 假设您正在尝试使用整数对输出数组进行索引，但是我们仍然可以通过这个小小的障碍来完成我们想要的。您要做的就是简单地为两个图像中的每个浮点值分配一个唯一 ID。因此，您可以将accumarray 与这些ID 一起使用。为了便于分配此 ID，请使用 unique - 特别是函数的第三个输出。您将获取每张图像，将它们放入unique 并将这些索引输入到accumarray。换句话说，改为这样做：
[~,~,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));
请注意，对于indrow 和indcol，我们直接将unique 的第三个输出分配给这些变量，然后使用我们之前计算的相同联合熵代码。我们也不必像以前那样将变量偏移 1，因为unique 将分配 ID从 1 开始。
一边
您实际上可以使用联合概率矩阵单独计算每张图像的直方图或概率分布。如果您想计算第一张图像的直方图/概率分布，您只需累积每行的所有列。要为第二个图像执行此操作，您只需累积每列的所有行。因此，您可以这样做：
histogramImage1 = sum(jointHistogram, 1);
histogramImage2 = sum(jointHistogram, 2);
之后，您可以自己计算这两者的熵。要仔细检查，请确保将这两个文件都转换为 PDF，然后使用标准方程计算熵（如上）。
我如何最终计算互信息？
要最终计算互信息，您将需要两张图像的熵。您可以使用 MATLAB 的内置 entropy 函数，但这假设有 256 个唯一级别。您可能希望将其应用于存在 N 不同级别而不是 256 的情况，因此您可以使用我们在上面对联合直方图所做的操作，然后在上面的代码中计算每个图像的直方图，然后计算每个图像的熵。您只需重复联合使用的熵计算，但将其分别应用于每个图像：
%// 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;
希望这会有所帮助！