列向量互信息计算通用MATLAB代码

haidaowang

互信息的定义
正式地，两个离散随机变量 X 和 Y 的互信息可以定义为：
3 ?# m! n# `  t* a# w8 t
其中 p(x,y) 是 X 和 Y 的联合概率分布函数，而p(x)和p(y)分别是 X 和 Y 的边缘概率分布函数。

" h7 Q6 V8 _; A' d: A/ q, p其中 p(x,y) 当前是 X 和 Y 的联合概率密度函数，而p(x)和p(y)分别是 X 和 Y 的边缘概率密度函数。

: t7 }) P0 g# }* p7 U$ Q互信息量I(xi;yj)在联合概率空间P(XY)中的统计平均值。 平均互信息I(X;Y)克服了互信息量I(xi;yj)的随机性,成为一个确定的量。如果对数以 2 为基底，互信息的单位是bit。
0 O# }; @; ?* ~4 Z% V! [3 I
直观上，互信息度量 X 和 Y 共享的信息：它度量知道这两个变量其中一个，对另一个不确定度减少的程度。例如，如果 X 和 Y 相互独立，则知道 X 不对 Y 提供任何信息，反之亦然，所以它们的互信息为零。在另一个极端，如果 X 是 Y 的一个确定性函数，且 Y 也是 X 的一个确定性函数，那么传递的所有信息被 X 和 Y 共享：知道 X 决定 Y 的值，反之亦然。因此，在此情形互信息与 Y（或 X）单独包含的不确定度相同，称作 Y（或 X）的熵。而且，这个互信息与 X 的熵和 Y 的熵相同。（这种情形的一个非常特殊的情况是当 X 和 Y 为相同随机变量时。）
  R8 C0 v; a: s6 L2 z7 N" W
% l7 `0 W: s* G& k8 l此外，互信息是非负的（即 I(X;Y) ≥ 0; 见下文），而且是对称的（即 I(X;Y) = I(Y;X)）。

通用MATLAB代码
4 X- H' o( M7 A  _' y# p主函数main.m

) K1 E# T( h" |0 {! b/ @# w/ o- Rclc
u1 = rand(4,1);
; W% ~1 e; Z- Y1 w# ^% pu2 = [2;32;6666;5];
3 k# L$ i! E: t4 Qwind_size = size(u1,1);
: {4 p! }6 y+ R4 `# w) \( ~/ A2 I9 I; bmi = calmi(u1, u2, wind_size);

calmi.m
. _# ^6 m: X: V4 c" V, n4 E
%计算两列向量之间的互信息
%u1：输入计算的向量1
5 l4 ^+ l$ S6 d- G/ D%u2：输入计算的向量2
%wind_size：向量的长度
! ^1 O. A9 x* G$ O, _/ Yfunction mi = calmi(u1, u2, wind_size)
x = [u1, u2];
' h2 E. A  j! Y' V7 c" i$ Sn = wind_size;
. F( A( y' }- k' k1 {( V[xrow, xcol] = size(x);
bin = zeros(xrow,xcol);
pmf = zeros(n, 2);
4 q4 }+ [$ f1 \for i = 1:2
6 L: O9 o* w9 }8 X& a5 t    minx = min(x(:,i));
    maxx = max(x(:,i));
    binwidth = (maxx - minx) / n;
) }! M% ~+ x9 H# }( K6 V$ M    edges = minx + binwidth*(0:n);
    histcEdges = [-Inf edges(2:end-1) Inf];
    [occur,bin(:,i)] = histc(x(:,i),histcEdges,1); %通过直方图方式计算单个向量的直方图分布
; u" ~0 s  u8 x7 l: E' C    pmf(:,i) = occur(1:n)./xrow;
end
- B9 ^+ a7 |) N9 \. J7 W" B%计算u1和u2的联合概率密度
2 _& Y) d; x# |) yjointOccur = accumarray(bin,1,[n,n]);  %（xi，yi）两个数据同时落入n*n等分方格中的数量即为联合概率密度
2 Q1 B1 }% I, k% S( j, v' gjointPmf = jointOccur./xrow;
; }( q' {$ q0 i3 P# j/ }Hx = -(pmf(:,1))'*log2(pmf(:,1)+eps);
Hy = -(pmf(:,2))'*log2(pmf(:,2)+eps);
0 W# r0 a" p' ]% NHxy = -(jointPmf(: ))'*log2(jointPmf(: )+eps);
MI = Hx+Hy-Hxy;
. g) B& Z6 W& ^2 x/ w+ w# m; E0 Tmi = MI/sqrt(Hx*Hy);

& n3 x, j% @- O7 h6 [$ ?$ [' Q
  `, o; x( l7 `* n* v

NingW

列向量互信息计算通用MATLAB代码