找回密码
 注册
关于网站域名变更的通知
查看: 1385|回复: 2
打印 上一主题 下一主题

用MATLAB实现基于互信息的特征选择算法

[复制链接]

该用户从未签到

跳转到指定楼层
1#
发表于 2019-8-8 09:00 | 只看该作者 |只看大图 回帖奖励 |倒序浏览 |阅读模式

EDA365欢迎您登录!

您需要 登录 才可以下载或查看,没有帐号?注册

x
在概率论和信息论中,两个随机变量的互信息(Mutual Information,简称MI)或转移信息(transinformation)是变量间相互依赖性的量度。不同于相关系数,互信息并不局限于实值随机变量,它更加一般且决定着联合分布 p(X,Y) 和分解的边缘分布的乘积 p(X)p(Y) 的相似程度。互信息(Mutual Information)是度量两个事件集合之间的相关性(mutual dependence)。互信息最常用的单位是bit。
4 ?: _0 V3 {* g; ?; g5 P
* ?7 G+ M+ g4 Z8 M/ o( X+ d

" [9 O$ p- i" v, A7 P互信息的定义 9 p* e! u! q+ x( ~1 e5 z6 U
正式地,两个离散随机变量 X 和 Y 的互信息可以定义为:
, f* {8 V6 L# o6 a8 O+ m& Q& M" |% x; q0 H4 E4 z- Y
其中 p(x,y) 是 X 和 Y 的联合概率分布函数,而p(x)和p(y)分别是 X 和 Y 的边缘概率分布函数。
. G* p4 J5 @5 c  J, ?: w% y4 K" J/ k- a5 j1 Y5 f* y

1 _) l  P1 J% z7 y) c7 u3 q: k - K- U* T7 i* S( a0 ~+ @% D

" Y/ l' F- _  I9 D在连续随机变量的情形下,求和被替换成了二重定积分: 3 ~% l2 |+ J: J' z; z
' m: |; R4 Y, B, G

9 P) _/ s& o: ~
* m$ W+ G9 h; R# a ' h! [3 D! c  H: K
其中 p(x,y) 当前是 X 和 Y 的联合概率密度函数,而p(x)和p(y)分别是 X 和 Y 的边缘概率密度函数。. I/ ~5 b4 M3 @) a9 s6 S& D
# ^  @3 Z9 d$ a  ?" w
1 J6 |$ N# M5 E7 H
互信息量I(xi;yj)在联合概率空间P(XY)中的统计平均值。 平均互信息I(X;Y)克服了互信息量I(xi;yj)的随机性,成为一个确定的量。如果对数以 2 为基底,互信息的单位是bit。
9 M) `) {( x. D; V
+ J4 v+ R' x- `7 y' N
) t: N7 G: A4 t
直观上,互信息度量 X 和 Y 共享的信息:它度量知道这两个变量其中一个,对另一个不确定度减少的程度。例如,如果 X 和 Y 相互独立,则知道 X 不对 Y 提供任何信息,反之亦然,所以它们的互信息为零。在另一个极端,如果 X 是 Y 的一个确定性函数,且 Y 也是 X 的一个确定性函数,那么传递的所有信息被 X 和 Y 共享:知道 X 决定 Y 的值,反之亦然。因此,在此情形互信息与 Y(或 X)单独包含的不确定度相同,称作 Y(或 X)的熵。而且,这个互信息与 X 的熵和 Y 的熵相同。(这种情形的一个非常特殊的情况是当 X 和 Y 为相同随机变量时。)
. u( O9 B- W! t" V; q. Y( `2 E2 e, |+ b/ z# v
  a2 ^) @2 D/ g- i# O9 w! B
互信息是 X 和 Y 联合分布相对于假定 X 和 Y 独立情况下的联合分布之间的内在依赖性。于是互信息以下面方式度量依赖性:I(X; Y) = 0 当且仅当 X 和 Y 为独立随机变量。从一个方向很容易看出:当 X 和 Y 独立时,p(x,y) = p(x) p(y),因此:
1 j% z9 W! z8 |" o
8 |% \% v( w2 X4 l, p0 m+ W

( _# W; y! S* T$ K) d . o( r+ i( E. w! ?' t+ D

# s8 y- j5 [( H) q: p9 v此外,互信息是非负的(即 I(X;Y) ≥ 0; 见下文),而且是对称的(即 I(X;Y) = I(Y;X))。6 _! |6 j8 m5 R$ h8 C5 P

& }+ D5 A( x/ L
" W, N) C+ A3 u
互信息特征选择算法的步骤
( x4 Y4 R  B& G2 |! I  M①划分数据集
  H$ v" Z' O. a! ~②利用互信息对特征进行排序
( |! o- F% Q7 r③选择前n个特征利用SVM进行训练
/ q) k5 A& v5 [④在测试集上评价特征子集计算错误率 : N- Q) K4 i4 E% `
! L; p% A- B9 ~" j

. f. r0 L# ^, u' m& G7 d代码 ( u9 ?# N! c* s1 S
注意使用的数据集是dlbcl,大概五千多维,可以从UCI上下载,最终选择前100特征进行训练。5 @. `/ B  g$ A+ c: V9 I6 K

6 D- H) R0 J5 b# J$ k

5 ~2 a. a# N' ?7 _* a; m+ ?主函数代码:
+ w# `' [; |" V7 u+ |4 Y& w9 Q
3 S7 @5 Z( k& e+ a4 T2 n# |/ j

3 O+ u0 l& A/ E2 R* H
  • clear all
  • close all
  • clc;
  • [X_train,Y_train,X_test,Y_test] = divide_dlbcl();
  • Y_train(Y_train==0)=-1;
  • Y_test(Y_test==0)=-1;
  • % number of features
  • numF = size(X_train,2);
  • " l& ^+ b  ^. d% `
  • " `' M5 ^; C8 U5 O7 A$ o9 G

  • ; {9 D8 G) l# @; ]% L
  • [ ranking , w] = mutInfFS( X_train, Y_train, numF );
  • k = 100; % select the Top 2 features
  • svmStruct = svmtrain(X_train(:,ranking(1:k)),Y_train,'showplot',true);
  • C = svmclassify(svmStruct,X_test(:,ranking(1:k)),'showplot',true);
  • err_rate = sum(Y_test~= C)/size(X_test,1); % mis-classification rate
  • conMat = confusionmat(Y_test,C); % the confusion matrix
  • fprintf('\nAccuracy: %.2f%%, Error-Rate: %.2f \n',100*(1-err_rate),err_rate);  \% C  E5 a5 A; H- X1 U

6 a" O: [' E$ L$ y8 n& N! R# \0 c# P! t! W

( W! M3 ^7 ~8 Z: qmutInfFS.m0 |$ k0 Z2 R7 ?" h/ {4 b

* w, [# h& `$ b1 i

8 E. E% U5 a3 [
  • function [ rank , w] = mutInfFS( X,Y,numF )
  • rank = [];
  • for i = 1:size(X,2)
  •     rank = [rank; -muteinf(X(:,i),Y) i];
  • end;
  • rank = sortrows(rank,1);
  • w = rank(1:numF, 1);
  • rank = rank(1:numF, 2);

  • * e0 i$ y8 j2 a, Q* b$ T$ Z
  • end
    . v2 ]6 E% b( _6 h! P& \$ N, z$ D
# U+ X, w5 Y; S2 @/ `+ i

% V5 M4 n$ |; Z4 V9 a1 `% B
) t) @: N# G: B+ Z* B5 {/ [- ^muteinf.m
  C2 n- M+ @  S5 _- _5 C$ u( N/ f' L6 u

" x" S6 [$ H+ m$ {8 ^
  • function info = muteinf(A, Y)
  • n = size(A,1);%实例数量
  • Z = [A Y];%所有实例的维度值及标签
  • if(n/10 > 20)
  •     nbins = 20;
  • else
  •     nbins = max(floor(n/10),10);%设置区间的个数
  • end;
  • pA = hist(A, nbins);%min(A)到max(A)划分出nbins个区间出来,求每个区间的概率
  • pA = pA ./ n;%除以实例数量
  • , A# T' f" x& }2 F9 m7 N; j
  • i = find(pA == 0);
  • pA(i) = 0.00001;%不能使某一区间的概率为0
  • 1 C* `: ?6 z) |
  • od = size(Y,2);%一个维度
  • cl = od;
  • %下面是求实例不同标签的的概率值,也就是频率
  • if(od == 1)
  •     pY = [length(find(Y==+1)) length(find(Y==-1))] / n;
  •     cl = 2;
  • else
  •     pY = zeros(1,od);
  •     for i=1d
  •         pY(i) = length(find(Y==+1));
  •     end;
  •     pY = pY / n;
  • end;
  • p = zeros(cl,nbins);
  • rx = abs(max(A) - min(A)) / nbins;%每个区间长度
  • for i = 1:cl
  •     xl = min(A);%变量的下界
  •     for j = 1:nbins
  •         if(i == 2) && (od == 1)
  •             interval = (xl <= Z(:,1)) & (Z(:,2) == -1);
  •         else
  •             interval = (xl <= Z(:,1)) & (Z(:,i+1) == +1);
  •         end;
  •         if(j < nbins)
  •             interval = interval & (Z(:,1) < xl + rx);
  •         end;
  •         %find(interval)
  •         p(i,j) = length(find(interval));

  • ' U0 [3 q8 d/ Q6 W8 d
  •         if p(i,j) == 0 % hack!
  •             p(i,j) = 0.00001;
  •         end
  • . s' _+ O6 j. c
  •         xl = xl + rx;
  •     end;
  • end;
  • HA = -sum(pA .* log(pA));%计算当前维度的信息熵
  • HY = -sum(pY .* log(pY));%计算标签的信息熵
  • pA = repmat(pA,cl,1);
  • pY = repmat(pY',1,nbins);
  • p = p ./ n;
  • info = sum(sum(p .* log(p ./ (pA .* pY))));
  • info = 2 * info ./ (HA + HY);%计算互信息
    . ]9 f8 s0 Q  z. z! y. s
$ e7 o& J/ Y6 q) X" k  F2 ?' B
/ W4 R+ K0 j) E. r

! \9 H) F2 ^% K前100个特征的效果:; l: ?: d' K  w& l

' c9 u% \9 H  m9 m- x
* l3 J! w6 ^" F. Q
Accuracy: 86.36%, Error-Rate: 0.149 n( h5 O' h& g

8 v" c7 `& C: ^' Y6 i
7 D: U# }; b- M# Z$ z! j' W' H
选择前两个特征进行训练(压缩率接近100%,把上述代码中的K设为2即可)的二维图: & |+ ]4 o. H$ b% k

1 m% _3 n% |+ b) F8 b5 `4 U
( w$ b; D  {* @$ A& s4 O6 t1 }3 J" }, e/ F0 p

6 ]; X8 A4 \% _: q6 p. QAccuracy: 75.00%, Error-Rate: 0.255 E+ F, F0 Q8 ?
: N' Z& \- v0 a5 n

( T. b. A. _- a  _- t! X* a! y4 t/ l1 _+ W1 k) T
, M: ?( j8 C# |  W3 [9 O% v
您需要登录后才可以回帖 登录 | 注册

本版积分规则

关闭

推荐内容上一条 /1 下一条

EDA365公众号

关于我们|手机版|EDA365电子论坛网 ( 粤ICP备18020198号-1 )

GMT+8, 2025-10-5 19:35 , Processed in 0.140625 second(s), 26 queries , Gzip On.

深圳市墨知创新科技有限公司

地址:深圳市南山区科技生态园2栋A座805 电话:19926409050

快速回复 返回顶部 返回列表