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

Nov 3,2018   5121 words   19 min


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

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

1.熵的直观理解

如何理解信息熵呢?

(1)信息

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

(2)信息熵

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

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

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

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

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

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

通过以上例子可以发现一个事件其不确定性越高,熵越大,所包含的信息量越大。但需要注意的是,信息量越多,不见得对你越有用,例如下面这个例子。 例如有100个抽奖球,编号为1-100,而其中只有一个球是有奖的。抽一次获不获奖的不确定性是100。

如果你说,我也不知道几号是中奖球,等于我抽奖中奖的概率还是1/100。但这句话的信息量是非常大的,约为33.22。

如果你说,中奖球在60-70号之间,那么那么我中奖概率就由1/100变成了1/10,不中奖的概率为9/10。这句话的信息量约为0.469。

而如果你直接说,中奖球是66号,那么我中奖概率就是1,不中奖概率就是0。这句话的信息量是0。因为不存在不确定性了,系统的混乱程度降为0了。 这个例子可以好好理解一下信息熵的本质,表征的是系统的混乱程度,不确定性的大小。

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)二维熵

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

上式中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[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 = win_w / 2
    ext_y = 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[img_patch.shape[0] / 2, 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,也算在可接受范围内。若采用之前的代码可能就要等几分钟了。

[更新结束]

最后,按照惯例将代码还有测试图片放到了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

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

返回顶部