信息熵介绍及其在图像处理中的应用与实现

Nov 3,2018   7042 words   26 min


熵这个概念源自于热力学,意义是体系混乱程度的度量,熵越大越混乱。1948年,香农将统计物理中熵的概念引申到信道通信的过程中,从而开创了“信息论”这门学科。香农定义的“熵”又被称为“香农熵”或“信息熵”。话说我第一次听到香农这个名字的时候有种说不上来的感觉(接地气?和农业有关?),它其实是英文Shannon的音译。 熵的英文是Entropy,如果你之前玩过崩坏3或看过漫画,应该知道里面有一个和天命对着干的组织就叫逆熵,英文是Anti-Entropy。 傲娇双马尾的特斯拉和冷冷帅气的爱茵。从字面翻译就是逆混乱或阻止混乱,因为我们的奥拓主教总是利用崩坏企图制造混乱达到他的目的,所以逆熵的人要阻止混乱。有点扯远了。。。

本篇博客主要介绍简单直观的介绍一下信息熵的概念,然后介绍其在图像处理中的应用。

1.熵的直观理解

如何理解信息熵呢?

(1)信息

首先得要从信息说起,什么是信息,信息我的个人理解就是对某件不确定事情的描述(可能不是很准确)。 香农给信息的定义是——信息是用来消除随机不确定性的东西。 例如,我告诉你明天太阳一定会升起。这句话对你来说就不是“信息”,是废话。因为只要地球不爆炸、宇宙不毁灭,太阳每天都会升起。 而如果我告诉你明天会下雨,这句话就是很有用的信息了。因为并不是每天都会下雨,存在着下或不下两种不确定性。

(2)信息熵

通过这个例子你可以发现,一句话中所包含的信息量是不同的。太阳会升起这是常识,对你而言信息量几乎为0。但下不下雨则不是的。 但我们能直观感受到这两句话对你的影响不同,但是这只是感性上的,我们无法从数学上精确比较。 一个事件的信息量可以认为是这个事件发生的概率的负对数信息熵与所有可能性有关,每个可能事件的发生都有个概率。 信息熵就是平均而言发生一个事件我们得到的信息量大小。所以数学上,信息熵其实是信息量的期望。 而信息熵就是这样一个指标,用于评价衡量一段信息中信息量的多少。 信息熵的定义如下。设有一个离散的随机变量X,它有m种可能情况,每种情况发生的概率为Pi,i=1,2,3…,m,则它的信息熵H(X)为

\[H(X)=-\sum_{i=1}^{m}p_{i}(x)logp_{i}(x)\]

这样,通过这个简单的公式,我们便能精确计算出信息量的多少了。信息熵有三条性质:

  • 单调性,即发生概率越高的事件,其所携带的信息熵越低
  • 非负性,即信息熵不能为负
  • 累加性,即多随机事件同时发生存在的总不确定性的量度是可以表示为各事件不确定性的量度的和

还是以上面的两句话为例。对于太阳会升起这件事情,其发生的概率可以认为就是1。从信息论角度,认为这句话没有消除任何不确定性。那么可以算出它的信息熵H为

\[H = -1 \times log(1)=-1\times 0 = 0\]

也就是说对于一件确定的事情,它包含的信息量为0!同样一句话,你们再重复一遍,等于你们也有责。。。额,错了,应该说人类的本质是复读机。 那我们再看看“明天会下雨”这句话。这件事无非就是两种情况,就是下雨或者不下雨。这里log的底取2,所以信息熵H如下

\[H = \left (- \frac{1}{2}\times log_{2}(\frac{1}{2}) \right )+\left ( - \frac{1}{2}\times log_{2}(\frac{1}{2}) \right )=-log_{2}(\frac{1}{2})=log_{2}2=1\]

这句话的信息熵为1,也即刚好可以用1bit的数据表达这个事件。因为在计算机中一般为2进制,只有0和1两种情况,所以一般在信息论中信息熵计算一般以2为底,具体计算采用换底公式即可。

通过以上例子可以发现一个事件其不确定性越高,熵越大,所包含的信息量越大。但需要注意的是,信息量越多(信息熵越大),不见得对你越有用,有些时候你希望信息熵小一点。

举个简单的生活中的例子,还是与天气有关,一个人说明天可能会下雨,这句话的信息熵为1(以2为底,共两种情况,概率都为1/2),而如果另一个人说明天有90%可能性会下雨,这句话的信息熵为0.47(以2为底,共两种情况,概率为0.1和0.9)。虽然后一个人的话信息量(熵)小了,但显然对于你而言是更有参考价值的。因为下不下雨本来就是1/2的概率,第一个人并没有降低系统任何的混乱程度。在这种场景下,你当然是希望系统的混乱程度越小越好,最好混乱程度降为0,也就是说可以肯定地告诉你明天下或者不下雨,以便于你更好地决策。因此说的更直白一些就是一句话说的越肯定,其包含的信息量越小。因为它所描述的内容的不确定性降低了。一句话说的越是云里雾里,打太极,它的信息量越大,因为可能性很多。仔细琢磨一下其实是个挺有意思的结论。

再举一个例子。例如有100个抽奖球,编号为1到100,其中只有一个球有奖,你可以指定抽几号球。首先我们说抽球的不确定性是100(可以理解为抽一次有100种可能的结果)。

如果你对抽球的人说:“我也不知道几号是中奖球”,那么抽球的人只能凭感觉随机选择一个号码,一共有100种可能性(不确定性)。这句话对抽球人而言显然是没有任何价值的信息,他抽一次中奖的概率还是1/100。但你说的这句话信息量(熵)是非常大的,约为6.64(以2为底,共100种可能情况,每种情况概率都为1/100)。这句话对于抽奖人而言就像是放了个巨大的烟雾弹。

如果你对抽球的人说:“中奖球在60-70号之间”,那么抽球的人肯定会在60到70号之间选(假设抽奖人完全相信你说的话),一共就只有10种可能性(不确定性)。这句话的信息量约为3.32(以2为底,共10种情况,每种情况概率都为1/10)。你这句话描述内容的可能性(熵)减少了,从而信息量降低了。但这对抽奖者而言就是个很好的消息,中奖概率由1/100变成了1/10,一下缩小了范围。

而如果你直接对抽球的人说:“偷偷告诉你中奖球是66号”,那么抽球的人直接就会选择66号(假设抽奖人完全相信你说的话),一共只有1种可能性(就是抽奖者选择66号球)。那么你说的这句话信息量为0(以2为底,共1种情况,概率都为1)。因为不存在不确定性了,系统的混乱程度(熵)降为0了。但此时抽奖人中奖概率就是1。

作为抽奖者而言,自然是希望系统的混乱度(熵)越低越好,因为这样中奖的概率才会更高。而作为举办方,则自然是希望系统的混乱度(熵)越高越好,这样才能真真假假,虚虚实实,可能大家都抽不到奖,这样就不用兑现奖品或者内部消化了。

这个例子同样可以从抽奖者的角度运用信息熵理解。对于第一种情况,对抽奖者而言中奖和不中奖的概率分别为0.01和0.99,这样可以算出来信息熵为0.08(以2为底,共2种情况,每种情况概率都为0.01和0.99)。而第二种情况信息熵为0.47(以2为底,共2种情况,每种情况概率都为0.1和0.9)。第三种情况信息熵为0(以2为底,共1种情况,概率都为1)。可以发现中奖或不中奖这两个概率相差越大,系统的信息熵越小。因为当某个事件占据绝对优势的时候,系统的混乱程度就下降了。以第一种情况为例,系统99%的概率都为一种事件,非常稳定。第三种情况其实是第一种情况的极限,只会有一种可能,完全稳定了。第二种情况事件概率差异小一些,更加混乱一些。进一步可以算得,对于只有两种可能事件的系统,这两个事件的概率相等都为1/2的时候,系统达到最大混乱状态,熵最大,为1对于三个可能事件的系统,同理当概率都为1/3的时候熵最大,约为1.58

从这个例子也可以得出两个结论:1.对于任何系统,可能的事件越多,系统越混乱,其熵越大;2.对于某个系统,当各可能事件的概率相等时达到最混乱状态,熵最大。用通俗易懂的语言解释就是各个事件“势均力敌”,打得很“焦灼”,这种状态是最混乱的。如果拿历史来举例,三国时代的熵是比较大的,因为各国势均力敌,相互混战;而秦朝的熵就是比较小的,因为全国大一统,其它势力基本不存在了,天下太平、稳定,也没有什么动乱了。

借上面两个例子可以好好理解一下信息熵的本质,表征的是系统的混乱程度,不确定性的大小。

2.信息熵在图像处理中的应用

(1)一维熵

在上面已经给出了信息熵的公式,对于影像而言,能随机变化的自然也就是影像中各像素的灰度了。 因此将信息熵的定义延伸到图像领域,那就是某个灰度级数在整个图像中出现的次数比上总像素个数,也就是定义中的Pi了。 对于8bit影像而言,灰度范围为0-255,共有256个。所以要计算一幅图像的信息熵其实很简单,对影像进行灰度直方图统计,然后求出各灰度级数的概率,最后套用信息熵公式即可。 采用Python代码实现如下。

# coding=utf-8
import cv2
import numpy as np


def calcEntropy(img):
    entropy = []

    hist = cv2.calcHist([img], [0], None, [256], [0, 255])
    total_pixel = img.shape[0] * img.shape[1]

    for item in hist:
        probability = item / total_pixel
        if probability == 0:
            en = 0
        else:
            en = -1 * probability * (np.log(probability) / np.log(2))
        entropy.append(en)

    sum_en = np.sum(entropy)
    return sum_en

代码中定义了calcEntropy()函数,传入的参数即是读取的8bit灰度影像。首先调用OpenCV的API统计了直方图(当然你可以自己用for循环写,但是效率肯定没有这个高),关于参数可参考这篇博客。 然后得到了总像素个数,对于直方图的每个灰度级数计算其出现概率,然后套用信息熵公式。 这里需要注意的是,如果出现概率为0的话需要单独写出来,不然程序在运行的时候会报错,因为log(0)等于负无穷。 最后,计算出各灰度级的熵之后,求和即可得到最终的结果。

利用上述函数对下图求信息熵(图片已转为灰度图)。 信息熵是5.9794。 信息熵是7.7776。 从信息熵的角度解释就是图二比图一在灰度层面更加混乱一些,也即可以认为图二的影像内容更加丰富一些。

(2)二维熵

上一部分的熵其实有点缺陷,缺点在于只从灰度层面统计了混乱度。但图像是二维的数据结构,所以应该充分利用空间信息来刻画图像内容的丰富程度,这样便有了二维熵的概念。二维信息熵的公式如下。

\[\left\{\begin{matrix} P_{i,j}=f(i,j)/W\cdot H\\ H=-\sum_{i=0}^{255}P_{i,j}logP_{i,j} \end{matrix}\right.\]

上式中W、H分别为影像的宽高,(i,j)为一个二元组,i表示某个滑动窗口内中心的灰度值,j为该窗口内除了中心像素的灰度均值。f(i,j)表示(i,j)这个二元组在整个图像中出现的次数。其余变量就跟一维相同了。 二维熵可以在图像所包含信息量的前提下,突出反映图像中像素位置的灰度信息和像素邻域内灰度分布的综合特征。

从公式上可以看出,二维熵相比于一维熵加入了对周围灰度信息的考虑,从而将空间信息纳入了计算。利用Python代码实现相比于一维熵稍微复杂一些。主要难点在滑动窗口部分。

# coding=utf-8
import cv2
import numpy as np
from collections import Counter
from multiprocessing import Pool


def calcIJ(img_patch):
    total_p = img_patch.shape[0] * img_patch.shape[1]
    if total_p % 2 != 0:
        tem = img_patch.flatten()
        center_p = tem[int(total_p / 2)]
        mean_p = (sum(tem) - center_p) / (total_p - 1)
        return (center_p, mean_p)
    else:
        print("modify patch size")


def calcEntropy2d(img, win_w=3, win_h=3, threadNum=6):
    height = img.shape[0]
    width = img.shape[1]

    ext_x = int(win_w / 2)
    ext_y = int(win_h / 2)

    # 考虑滑动窗口大小,对原图进行扩边,扩展部分灰度为0
    ext_h_part = np.zeros([height, ext_x], img.dtype)
    tem_img = np.hstack((ext_h_part, img, ext_h_part))
    ext_v_part = np.zeros([ext_y, tem_img.shape[1]], img.dtype)
    final_img = np.vstack((ext_v_part, tem_img, ext_v_part))

    new_width = final_img.shape[1]
    new_height = final_img.shape[0]

    # 依次获取每个滑动窗口的内容,将其暂存在list中,这里的写法不是最好的只是为了方便理解和便于并行。但在实际使用中没有太大必要重新建个list,都在一个循环里完成即可。
    patches = []
    for i in range(ext_x, new_width - ext_x):
        for j in range(ext_y, new_height - ext_y):
            patch = final_img[j - ext_y:j + ext_y + 1, i - ext_x:i + ext_x + 1]
            patches.append(patch)

    # 并行计算每个滑动窗口对应的i、j
    pool = Pool(processes=threadNum)
    IJ = pool.map(calcIJ, patches)
    pool.close()
    pool.join()
    # 串行计算速度稍慢,可以考虑使用并行计算加速
    # for item in patchs:
    #     IJ.append(calcIJ(item))

    # 循环遍历统计各二元组个数
    Fij = Counter(IJ).items()
    # 推荐使用Counter方法,如果自己写for循环效率会低很多
    # Fij = []
    # IJ_set = set(IJ)
    # for item in IJ_set:
    #     Fij.append((item, IJ.count(item)))

    # 计算各二元组出现的概率
    Pij = []
    for item in Fij:
        Pij.append(item[1] * 1.0 / (new_height * new_width))

    # 计算每个概率所对应的二维熵
    H_tem = []
    for item in Pij:
        h_tem = -item * (np.log(item) / np.log(2))
        H_tem.append(h_tem)

    # 对所有二维熵求和
    H = sum(H_tem)
    return H

对下面两幅图计算二维熵。 调低对比度和饱和度后,二维熵为8.3178。 原图二维熵为9.2895。

[2019-8-19更新]

之前有在Github反馈说二维熵的计算太慢了,所以又重新加速了下代码,相比于之前的MultiProcessing模块,这次使用的是Numba模块,相比于前者,加速效果非常好。关于这个模块的相关博客可以看这里。 主要修改的地方就是对calcIJ()函数用了Numba的Jit加速,部分代码如下。

@jit(nopython=True)
def calcIJ(img_patch):
    total_p = img_patch.shape[0] * img_patch.shape[1]
    if total_p % 2 != 0:
        center_p = img_patch[int(img_patch.shape[0] / 2), int(img_patch.shape[1] / 2)]
        mean_p = (np.sum(img_patch) - center_p) / (total_p - 1)
        return (center_p, mean_p)
    else:
        pass

完整代码见Github项目的entropy2dSpeedUp.py。加速后的代码运行示例img4.jpg3次取平均值,本机耗时为4.72s,使用MultiProcessing模块加速的平均耗时为33.91s。提升幅度达到86.1%。 此外利用更新后的代码计算我微单拍摄的原始影像,Github项目中的img5.jpg,大小为6000×4000,本机平均耗时为24.86s,也算在可接受范围内。若采用之前的代码可能就要等几分钟了。当然除了Python,也使用C++重新编写了代码,效率更高,代码还是在上面提到的这个Github项目里,C++文件夹下。

[更新结束]

最后,按照惯例将代码还有测试图片放到了Github上,点击查看,欢迎Star或Fork。

3.参考资料

  • [1]https://www.zhihu.com/question/22178202/answer/49929786
  • [2]https://blog.csdn.net/a6333230/article/details/81021922
  • [3]http://blog.sina.com.cn/s/blog_50e7455c01019l2k.html

本文作者原创,未经许可不得转载,谢谢配合

返回顶部