Jor*_*rge 23 matlab image image-processing entropy information-theory
我有两个黑白图像,我需要计算互信息.
Image 1 = X
Image 2 = Y
我知道互信息可以定义为:
MI = entropy(X) + entropy(Y) - JointEntropy(X,Y)
MATLAB已经具有内置函数来计算熵但不计算联合熵.我想真正的问题是:我如何计算两幅图像的联合熵？
这是我想要找到的联合熵的图像示例:
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
ray*_*ica 60
要计算关节熵,您需要计算两个图像之间的联合直方图.联合直方图与普通1D直方图基本相同,但第一维记录第一图像的强度,第二维记录第二图像的强度.这与通常称为共生矩阵的非常相似.在(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每个强度值的值.然后,我们将要添加的所有值高达1映射到相同的箱柜,而这就是我们怎么会计算的直方图.默认行为accumarray是添加所有这些值.具体来说,输出accumarray将是一个数组,其中每个位置计算映射到该键的所有值的总和.例如,第一个位置是映射到键1的所有值的总和,第二个位置是映射到键2的所有值的总和,依此类推.
但是,对于联合直方图,您需要确定哪些值映射到相同的强度对(i,j),因此这里的键将是一对2D坐标.这样,在两个图像之间共享的相同空间位置i中具有第一图像和j第二图像中的强度的任何强度都转到相同的键.因此,在2D情况下,输出accumarray将是2D矩阵,其中每个元素(i,j)包含映射到键的所有值的总和(i,j),类似于之前提到的1D情况,这正是我们所追求的.
换一种说法:
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是,如果每个键具有相同的值,您可以简单地为第二个输入分配一个常量,这就是我所做的和它的完成1.通常,这是一个与第一个输入具有相同行数的数组.另外,请特别注意前两行.0图像中不可避免地存在强度,但由于MATLAB开始索引1,我们需要将两个数组都偏移1.
现在我们有联合直方图,计算联合熵非常简单.它类似于1D中的熵,除了现在我们只是对整个联合概率矩阵求和.请记住,您的联合直方图很可能会有很多0条目.我们需要确保跳过这些或log2操作将是未定义的.让我们现在摆脱任何零条目:
indNoZero = jointHistogram ~= 0;
jointProb1DNoZero = jointProb(indNoZero);
请注意,我搜索了联合直方图而不是联合概率矩阵.这是因为联合直方图由整数组成,而联合概率矩阵位于0和之间1.由于划分,我想避免比较此矩阵中的任何条目与0数值舍入和不稳定性.以上也将我们的联合概率矩阵转换为堆叠的1D向量,这很好.
因此,联合熵可以计算为:
jointEntropy = -sum(jointProb1DNoZero.*log2(jointProb1DNoZero));
如果我对在MATLAB中计算图像的熵的理解是正确的,它应该计算在256二进制位上的直方图/概率分布,所以你当然可以在这里使用刚刚计算出的联合熵.
到目前为止,我们假设您处理的图像具有整数值的强度.如果我们有浮点数据怎么办？ 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这些变量的第三个输出分配给我们,然后使用我们之前计算的相同联合熵代码.我们也不必像以前那样unique将变量偏移1,因为它将从1开始分配ID .
实际上,您可以使用联合概率矩阵单独计算每个图像的直方图或概率分布.如果要计算第一个图像的直方图/概率分布,则只需累加每行的所有列.要为第二个图像执行此操作,您只需累积每列的所有行.因此,您可以这样做:
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;
希望这可以帮助!