[zz]求一维序列的信息熵（香浓熵）的matlab程序实例

对于一个二维信号，好比灰度图像，灰度值的范围是0-255，所以只要根据像素灰度值(0-255)出现的几率，就能够计算出信息熵。
    可是，对于一个一维信号，好比说心电信号，数据值的范围并非肯定的，不会是(0-255)这么肯定，若是进行域值变换，使其转换到一个整数范围的话，就会丢失数据，请高手指点，怎么计算。

好比数字信号是x(n),n=1~N
(1)先用Hist函数对x(n)的赋值范围进行分块，好比赋值范围在0~10的对应第  
    一块，10~20的第二块，以此类推。这以前须要对x(n)作一些归一化处理
(2)统计每一块的数据个数，并求出相应的几率
(3)用信息熵公式求解

以上求解方法得到的虽然是近似的信息熵，可是通常认为，这么作是没有问题的

求一维序列的信息熵的matlab程序代码以下：（已写成调用的函数形式）

测试程序：

fs=12000;
N=12000;
T=1/fs;
t=(0:N-1)*T;
ff=104;
sig=0.5*(1+sin(2*pi*ff*t)).*sin(2*pi*3000*t)+rand(1,length(t));

Hx=yyshang(sig,10)

 

%———————求一维离散序列信息熵matlab代码

function Hx=yyshang(y,duan)
%不以原信号为参考的时间域的信号熵
%输入：maxf:原信号的能量谱中能量最大的点
%y:待求信息熵的序列
%duan:待求信息熵的序列要被分块的块数
%Hx:y的信息熵
%duan=10;%将序列按duan数等分，若是duan=10,就将序列分为10等份
x_min=min(y);
x_max=max(y);
maxf(1)=abs(x_max-x_min);
maxf(2)=x_min;
duan_t=1.0/duan;
jiange=maxf(1)*duan_t;
% for i=1:10
% pnum(i)=length(find((y_p>=(i-1)*jiange)&(y_p<i*jiange)));
% end

pnum(1)=length(find(y<maxf(2)+jiange));
for i=2:duan-1
    pnum(i)=length(find((y>=maxf(2)+(i-1)*jiange)&(y<maxf(2)+i*jiange)));
end
pnum(duan)=length(find(y>=maxf(2)+(duan-1)*jiange));
%sum(pnum)
ppnum=pnum/sum(pnum);%每段出现的几率
%sum(ppnum)
Hx=0;
for i=1:duan
    if ppnum(i)==0
        Hi=0;
    else
        Hi=-ppnum(i)*log2(ppnum(i));
    end
    Hx=Hx+Hi;
end
end

%----------------

扩展阅读：

 

 

实验一：计算离散信源的熵

1、实验设备:

一、计算机

二、软件：Matlab

2、实验目的:

  一、熟悉离散信源的特色；

  二、学习仿真离散信源的方法

  三、学习离散信源平均信息量的计算方法

  四、熟悉 Matlab 编程；

3、实验内容:

   一、写出计算自信息量的Matlab 程序

   二、写出计算离散信源平均信息量的Matlab 程序。

   三、掌握二元离散信源的最大信息量与几率的关系。

   四、将程序在计算机上仿真实现，验证程序的正确性并完成习题。

4、实验报告要求

简要总结离散信源的特色及离散信源平均信息量的计算,写出习题的MATLAB实现语句。

 

信息论基础：

自信息的计算公式

              Matlab实现：I=log2(1/p) 或I=-log2(p)

熵（平均自信息）的计算公式

  

Matlab实现：HX=sum(-x.*log2(x))；或者h=h-x(i)*log2(x(i));

习题：

1. 甲地天气预报构成的信源空间为：

 

乙地信源空间为：

 

求此两个信源的熵。求各类天气的自信息量。

案：

运行程序：

p1=[1/2,1/4,1/8,1/8];%p1表明甲信源对应的几率

p2=[7/8,1/8];%p2表明乙信源对应的几率

H1=0.0;

H2=0.0;

I=[];

J=[];

for i=1:4

    H1=H1+p1(i)*log2(1/p1(i));

    I(i)=log2(1/p1(i));

end

disp('自信息量分别为：');

I

disp('H1信源熵为：');

H1

for j=1:2

    H2=H2+p2(j)*log2(1/p2(j));

     J(j)=log2(1/p2(j));

end

disp('自信息量分别为：');

J

disp('H2信源熵为：');

H2