模式识别实验二、贝叶斯分类参数估计matlab代码

u1=[1;2];%类条件分布概率参数u2=[5;3];d1=[1.50;01];d2=[30;03];n1=100;%训练样本数目n2=100;pw1=n1/(n1+n2);pw2=n2/(n1+n2);%先验概率x1=mvnrnd(u1,d1,n1);%训练样本x2=mvnrnd(u2,d2,n2);xt=[x1;x2];invsig1=inv(d1);invsig2=inv(d2);s1=abs(det(d1));s2=abs(det(d2));p1=mvnpdf(x1,u1',d1);p2=mvnpdf(x2,u2',d2);x1_test=mvnrnd(u1,d1,20);%测试样本x2_test=mvnrnd(u2,d2,20);x_test=[x1_test;x2_test];classify=zeros((n1+n2),1);%理论值fori=1:(n1+n2)g1(i)=(-0.5)*(xt(i,:)-u1')*invsig1*(xt(i,:)-u1')'-log(2*pi)-0.5*log(s1)+log(pw1);g2(i)=(-0.5)*(xt(i,:)-u2')*invsig2*(xt(i,:)-u2')'-log(2*pi)-0.5*log(s2)+log(pw2);if(g1(i)g2(i))classify(i)=1;elseif(g1(i)g2(i))classify(i)=2;endendne21=0;%应该是1类，错分到2类的数目ne12=0;%应该是2类，错分到1类的数目fori=1:n1if(classify(i)~=1)ne21=ne21+1;endendfori=n1+1:n1+n2if(classify(i)~=2)ne12=ne12+1;endendpe=(ne21+ne12)/(n1+n2)us1=mean(x1)%参数估计us2=mean(x2)si1=zeros(2);fori=1:n1si1=si1+(x1(i,:)-us1)'*(x1(i,:)-us1);endds1=si1/n1si2=zeros(2);fori=1:n2si2=si2+(x2(i,:)-us1)'*(x2(i,:)-us2);endds2=si2/n2invsigs1=inv(ds1);invsigs2=inv(ds2);ss1=abs(det(ds1));ss2=abs(det(ds2));x=xt(:,1);y=xt(:,2);[X,Y]=meshgrid(x,y);g=zeros((n1+n2),2);sclassify=zeros((n1+n2),1);%训练样本值fori=1:(n1+n2)g1(i)=(-0.5)*(xt(i,:)-us1)*invsigs1*(xt(i,:)-us1)'-log(2*pi)-0.5*log(ss1)+log(pw1);g2(i)=(-0.5)*(xt(i,:)-us2)*invsigs2*(xt(i,:)-us2)'-log(2*pi)-0.5*log(ss2)+log(pw2);if(g1(i)g2(i))sclassify(i)=1;elseif(g1(i)g2(i))sclassify(i)=2;endg(i,:)=[g1(i)g2(i)];endsne21=0;%应该是1类，错分到2类的数目sne12=0;%应该是2类，错分到1类的数目fori=1:n1if(sclassify(i)~=1)sne21=sne21+1;endendfori=n1+1:n1+n2if(sclassify(i)~=2)sne12=sne12+1;endendspe=(sne21+sne12)/(n1+n2)tclassify=zeros((n1+n2),1);%测试理论值fori=1:40g1(i)=(-0.5)*(x_test(i,:)-u1')*invsig1*(x_test(i,:)-u1')'-log(2*pi)-0.5*log(s1)+log(pw1);g2(i)=(-0.5)*(x_test(i,:)-u2')*invsig2*(x_test(i,:)-u2')'-log(2*pi)-0.5*log(s2)+log(pw2);if(g1(i)g2(i))tclassify(i)=1;elseif(g1(i)g2(i))tclassify(i)=2;endendtne21=0;%应该是1类，错分到2类的数目tne12=0;%应该是2类，错分到1类的数目fori=1:20if(classify(i)~=1)tne21=tne21+1;endendfori=21:40if(classify(i)~=2)tne12=tne12+1;endendtpe=(tne21+tne12)/40test_classify=zeros((n1+n2),1);%测试样本值fori=1:40g1(i)=(-0.5)*(x_test(i,:)-us1)*invsigs1*(x_test(i,:)-us1)'-log(2*pi)-0.5*log(ss1)+log(pw1);g2(i)=(-0.5)*(x_test(i,:)-us2)*invsigs2*(x_test(i,:)-us2)'-log(2*pi)-0.5*log(ss2)+log(pw2);if(g1(i)g2(i))test_classify(i)=1;elseif(g1(i)g2(i))test_classify(i)=2;endendteste21=0;%应该是1类，错分到2类的数目teste12=0;%应该是2类，错分到1类的数目fori=1:20if(test_classify(i)~=1)teste21=sne21+1;endendfori=21:40if(test_classify(i)~=2)sne12=sne12+1;endendtest_pe=(sne21+sne12)/40subplot(2,2,1)symsxxyy;gg=(-0.5)*[xx-u1(1),yy-u1(2)]*invsig1*[xx-u1(1),yy-u1(2)]'-0.5*log(s1)+log(pw1)-((-0.5)*[xx-u2(1),yy-u2(2)]*invsig2*[xx-u2(1),yy-u2(2)]'-0.5*log(s2)+log(pw2));ezplot(gg)holdonfori=1:n1+n2if(i=n1)if(classify(i)~=1)plot(x1(i,:)-1,'r*');elseplot(x1(i,:)-1,'ro');endelseif(classify(i)~=2)plot(x2(i-n1,1)+3,x2(i-n1,2),'b*');%为了让1类和2类在图像上区分开elseplot(x2(i-n1,:)+3,x2(i-n1,2),'bo');endendendtitle('理论决策面')subplot(2,2,2)symsxsys;ggs=(-0.5)*[xs-us1(1),ys-us1(2)]*invsigs1*[xs-us1(1),ys-us1(2)]'-0.5*log(ss1)+log(pw1)-((-0.5)*[xs-us2(1),ys-us2(2)]*invsigs2*[xs-us2(1),ys-us2(2)]'-0.5*log(ss2)+log(pw2));ezplot(ggs)holdonfori=1:n1+n2if(i=n1)if(sclassify(i)~=1)plot(x1(i,:)-1,'r*');elseplot(x1(i,:)-1,'ro');endelseif(sclassify(i)~=2)plot(x2(i-n1,1)+3,x2(i-n1,2),'b*');%为了让1类和2类在图像上区分开elseplot(x2(i-n1,1)+3,x2(i-n1,2),'bo');endendendtitle('样本决策面')subplot(2,2,3)symsxxyy;gg=(-0.5)*[xx-u1(1),yy-u1(2)]*invsig1*[xx-u1(1),yy-u1(2)]'-0.5*log(s1)+log(pw1)-((-0.5)*[xx-u2(1),yy-u2(2)]*invsig2*[xx-u2(1),yy-u2(2)]'-0.5*log(s2)+log(pw2));ezplot(gg)holdonfori=1:40if(i=20)if(tclassify(i)~=1)plot(x1_test(i,:)-1,'r*');elseplot(x1_test(i,:)-1,'ro');endelseif(tclassify(i)~=2)plot(x2_test(i-20,1)+3,x2_test(i-20,2),'b*');%为了让1类和2类在图像上区分开elseplot(x2_test(i-20,1)+3,x2_test(i-20,2),'bo');endendendtitle('测试理论决策面')subplot(2,2,4)symsxsys;ggs=(-0.5)*[xs-us1(1),ys-us1(2)]*invsigs1*[xs-us1(1),ys-us1(2)]'-0.5*log(ss1)+log(pw1)-((-0.5)*[xs-us2(1),ys-us2(2)]*invsigs2*[xs-us2(1),ys-us2(2)]'-0.5*log(ss2)+log(pw2));ezplot(ggs)holdonfori=1:40if(i=20)if(test_classify(i)~=1)plot(x1_test(i,:)-1,'r*');elseplot(x1_test(i,:)-1,'ro');endelseif(test_classify(i)~=2)plot(x2_test(i-20,1)+3,x2_test(i-20,2),'b*');%为了让1类和2类在图像上区分开elseplot(x2_test(i-20,1)+3,x2_test(i-20,2),'bo');endendendtitle('测试新决策面')