判別分析的MATLAB實現(xiàn)案例

%--------------------------------------------------------------------------%讀取examp10_01.xls中數(shù)據(jù)，進行距離鑒別%--------------------------------------------------------------------------%********************************讀取數(shù)據(jù)***********************************%讀取文獻examp10_01.xls的第1個工作表中C2:F51范圍的數(shù)據(jù)，即所有樣本數(shù)據(jù)，包括未判企業(yè)sample=xlsread('examp10_01.xls','','C2:F51');%讀取文獻examp10_01.xls的第1個工作表中C2:F47范圍的數(shù)據(jù)，即已知組別的樣本數(shù)據(jù)，training=xlsread('examp10_01.xls','','C2:F47');%讀取文獻examp10_01.xls的第1個工作表中B2:B47范圍的數(shù)據(jù)，即樣本的分組信息數(shù)據(jù)，group=xlsread('examp10_01.xls','','B2:B47');obs=[1:50]';%企業(yè)的編號%**********************************距離鑒別*********************************%距離鑒別，鑒別函數(shù)類型為mahalanobis，返回鑒別成果向量C和誤判概率err[C,err]=classify(sample,training,group,'mahalanobis');[obs,C]%查看鑒別成果err%查看誤判概率%--------------------------------------------------------------------------%加載fisheriris.mat中數(shù)據(jù)，進行貝葉斯鑒別%--------------------------------------------------------------------------%********************************加載數(shù)據(jù)***********************************loadfisheriris%把文獻fisheriris.mat中數(shù)據(jù)導(dǎo)入MATLAB工作空間%**********************************查看數(shù)據(jù)*********************************head0={'Obj','x1','x2','x3','x4','Class'};%設(shè)置表頭[head0;num2cell([[1:150]',meas]),species]%以元胞數(shù)組形式查看數(shù)據(jù)%*********************************貝葉斯鑒別********************************%用meas和species作為訓(xùn)練樣本，創(chuàng)立一種樸素貝葉斯分類器對象ObjBayesObjBayes=NaiveBayes.fit(meas,species);%運用所創(chuàng)立的樸素貝葉斯分類器對象對訓(xùn)練樣本進行鑒別，返回鑒別成果pre0，pre0也是字符串元胞向量pre0=ObjBayes.predict(meas);%運用confusionmat函數(shù)，并根據(jù)species和pre0創(chuàng)立混淆矩陣（包括總的分類信息的矩陣）[CLMat,order]=confusionmat(species,pre0);%以元胞數(shù)組形式查看混淆矩陣[[{'From/To'},order'];order,num2cell(CLMat)]%查看誤判樣品編號gindex1=grp2idx(pre0);%根據(jù)分組變量pre0生成一種索引向量gindex1gindex2=grp2idx(species);%根據(jù)分組變量species生成一種索引向量gindex2errid=find(gindex1~=gindex2)%通過對比兩個索引向量，返回誤判樣品的觀測序號向量%查看誤判樣品的誤判狀況head1={'Obj','From','To'};%設(shè)置表頭%用num2cell函數(shù)將誤判樣品的觀測序號向量errid轉(zhuǎn)為元胞向量，然后以元胞數(shù)組形式查看誤判成果[head1;num2cell(errid),species(errid),pre0(errid)]%對未知類別樣品進行鑒別%定義未判樣品觀測值矩陣xx=[5.8 2.7 1.8 0.735.6 3.1 3.8 1.86.1 2.5 4.7 1.16.1 2.6 5.7 1.95.1 3.1 6.5 0.625.8 3.7 3.9 0.135.7 2.7 1.1 0.126.4 3.2 2.4 1.66.7 3 1.9 1.16.8 3.5 7.9 1];%運用所創(chuàng)立的樸素貝葉斯分類器對象對未判樣品進行鑒別，返回鑒別成果pre1，pre1也是字符串元胞向量pre1=ObjBayes.predict(x)%--------------------------------------------------------------------------%加載fisheriris.mat中數(shù)據(jù)，進行Fisher鑒別%--------------------------------------------------------------------------%********************************加載數(shù)據(jù)***********************************loadfisheriris%把文獻fisheriris.mat中數(shù)據(jù)導(dǎo)入MATLAB工作空間%**********************************待判樣品*********************************%定義待判樣品觀測值矩陣xx=[5.8 2.7 1.8 0.735.6 3.1 3.8 1.86.1 2.5 4.7 1.16.1 2.6 5.7 1.95.1 3.1 6.5 0.625.8 3.7 3.9 0.135.7 2.7 1.1 0.126.4 3.2 2.4 1.66.7 3 1.9 1.16.8 3.5 7.9 1];%*********************************Fisher鑒別********************************%運用fisher函數(shù)進行鑒別，返回多種成果（見fisher函數(shù)的注釋）[outclass,TabCan,TabL,TabCon,TabM,TabG]=fisher(x,meas,species)%************************繪制兩個鑒別式得分的散點圖**************************%運用fisher函數(shù)進行鑒別，返回多種成果，其中ts為鑒別式得分[outclass,TabCan,TabL,TabCon,TabM,TabG,ts]=fisher(x,meas,species);%提取各類的鑒別式得分ts1=ts(ts(:,1)==1,:);%setosa類的鑒別式得分ts2=ts(ts(:,1)==2,:);%versicolor類的鑒別式得分ts3=ts(ts(:,1)==3,:);%virginica類的鑒別式得分plot(ts1(:,2),ts1(:,3),'ko')%setosa類的鑒別式得分的散點圖holdonplot(ts2(:,2),ts2(:,3),'k*')%versicolor類的鑒別式得分的散點圖plot(ts3(:,2),ts3(:,3),'kp')%virginica類的鑒別式得分的散點圖legend('setosa類','versicolor類','virginica類');%加標(biāo)注框xlabel('第一鑒別式得分');%給X軸加標(biāo)簽ylabel('第二鑒別式得分');%給Y軸加標(biāo)簽%************************只用一種鑒別式進行Fisher鑒別************************%令fisher函數(shù)的第4個輸入為0.5，就可以只用一種鑒別式進行鑒別[outclass,TabCan,TabL,TabCon,TabM,TabG]=fisher(x,meas,species,0.5)function[outclass,TabCan,TabL,TabCon,TabM,TabG,trainscore]=fisher(sampledata,training,group,contri)%FISHER鑒別分析.%class=fisher(sampledata,training,group)根據(jù)訓(xùn)練樣本training構(gòu)造鑒別式，%運用所有鑒別式看待判樣品sampledata進行鑒別.sampledata和training是具有相似%列數(shù)的矩陣，它們的每一行對應(yīng)一種觀測，每一列對應(yīng)一種變量.group是training對%應(yīng)的分組變量，它的每一種元素定義了training中對應(yīng)觀測所屬的類.group可以是一%個分類變量，數(shù)值向量，字符串?dāng)?shù)組或字符串元胞數(shù)組.training和group必須具有相%同的行數(shù).fisher函數(shù)把group中的NaN或空字符串作為缺失數(shù)據(jù)，從而忽視training%中對應(yīng)的觀測.class中的每個元素指定了sampledata中的對應(yīng)觀測所判歸的類,它和%group具有相似的數(shù)據(jù)類型.%%class=fisher(sampledata,training,group,contri)根據(jù)累積奉獻率不低于%contri，確定需要使用的鑒別式個數(shù)，默認(rèn)狀況下，使用所有鑒別式進行鑒別.contri%是一種在(0,1]區(qū)間內(nèi)取值的標(biāo)量，用來指定累積奉獻率的下限.%%[class,TabCan]=fisher(...)以表格形式返回所用鑒別式的系數(shù)向量，若contri%取值為1，則返回所有鑒別式的系數(shù)向量.TabCan是一種元胞數(shù)組，形如%'Variable''can1''can2'%'x1'[-0.2087][0.0065]%'x2'[-0.3862][0.5866]%'x3'[0.5540][-0.2526]%'x4'[0.7074][0.7695]%[class,TabCan,TabL]=fisher(...)以表格形式返回所有特性值，奉獻率，累積%奉獻率等.TabL是一種元胞數(shù)組，形如%'Eigenvalue''Difference''Proportion''Cumulative'%[32.1919][31.9065][0.9912][0.9912]%[0.2854][][0.0088][1]%%[class,TabCan,TabL,TabCon]=fisher(...)以表格形式返回混淆矩陣（包括總%的分類信息的矩陣）.TabCon是一種元胞數(shù)組，形如%'From/To''setosa''versicolor''virginica'%'setosa'[50][0][0]%'versicolor'[0][48][2]%'virginica'[0][1][49]%%[class,TabCan,TabL,TabCon,TabM]=fisher(...)以表格形式返回誤判矩陣.%TabM是一種元胞數(shù)組，形如%'Obj''From''To'%[71]'versicolor''virginica'%[84]'versicolor''virginica'%[134]'virginica''versicolor'%%[class,TabCan,TabL,TabCon,TabM,TabG]=fisher(...)將所用鑒別式作用%在各組的組均值上，得到組均值投影矩陣，以表格形式返回這個矩陣.TabG是一種元胞%數(shù)組，形如%'Group''can1''can2'%'setosa'[-1.3849][1.8636]%'versicolor'[0.9892][1.6081]%'virginica'[1.9852][1.9443]%[class,TabCan,TabL,TabCon,TabM,TabG,trainscore]=fisher(...)返回%訓(xùn)練樣品所對應(yīng)的鑒別式得分trainscore.trainscore的第一列為各訓(xùn)練樣品原本所%屬類的類序號，第i+1列為第i個鑒別式得分.%%Copyrightxiezhh.%$Revision:1.0.0.0$$Date:/10/0310:40:34$ifnargin<3error('錯誤：輸入?yún)?shù)太少，至少需要3個輸入.');end%根據(jù)分組變量生成索引向量gindex，組名元胞向量groups，組水平向量glevels[gindex,groups,glevels]=grp2idx(group);%忽視缺失數(shù)據(jù)nans=find(isnan(gindex));if~isempty(nans)training(nans,:)=[];gindex(nans)=[];endngroups=length(groups);gsize=hist(gindex,1:ngroups);nonemptygroups=find(gsize>0);nusedgroups=length(nonemptygroups);%判斷與否有空的組ifngroups>nusedgroupswarning('警告:有空的組.');end[n,d]=size(training);ifsize(gindex,1)~=nerror('錯誤:輸入?yún)?shù)大小不匹配，GROUP與TRAINING必須具有相似的行數(shù).');elseifisempty(sampledata)sampledata=zeros(0,d,class(sampledata));elseifsize(sampledata,2)~=derror('錯誤:輸入?yún)?shù)大小不匹配，SAMPLEDATA與TRAINING必須具有相似的列數(shù).');end%設(shè)置contri的默認(rèn)值為1，并限定contri在(0,1]內(nèi)取值ifnargin<4||isempty(contri)contri=1;endif~isscalar(contri)||contri>1||contri<=0error('錯誤:contri必須是一種在(0,1]內(nèi)取值的標(biāo)量.');endifany(gsize==1)error('錯誤:TRAINING中的每個組至少應(yīng)有兩個觀測.');end%計算各組的組均值gmeans=NaN(ngroups,d);fork=nonemptygroupsgmeans(k,:)=mean(training(gindex==k,:),1);end%計算總均值totalmean=mean(training,1);%計算組內(nèi)離差平方和矩陣E和組間離差平方和矩陣BE=zeros(d);B=E;fork=nonemptygroups%分別估計各組的組內(nèi)離差平方和矩陣.[Q,Rk]=qr(bsxfun(@minus,training(gindex==k,:),gmeans(k,:)),0);%各組的組內(nèi)離差平方和矩陣：AkHat=Rk'*Rk%判斷各組的組內(nèi)離差平方和矩陣的正定性s=svd(Rk);ifany(s<=max(gsize(k),d)*eps(max(s)))error('錯誤:TRAINING中各組的組內(nèi)離差平方和矩陣必須是正定矩陣.');endE=E+Rk'*Rk;%計算總的組內(nèi)離差平方和矩陣E%計算組間離差平方和矩陣BB=B+(gmeans(k,:)-totalmean)'*(gmeans(k,:)-totalmean)*gsize(k);end%求inv(E)*B的正特性值與對應(yīng)的特性向量EB=E\B;[V,D]=eig(EB);D=diag(D);[D,idD]=sort(D,'descend');%將特性值按降序排列V=V(:,idD);NumPosi=min(ngroups-1,d);%確定正特性值個數(shù)D=D(1:NumPosi,:);CumCont=cumsum(D/sum(D));%計算累積奉獻率%以表格形式返回所有特性值，奉獻率，累積奉獻率等.TabL是一種元胞數(shù)組head={'Eigenvalue','Difference','Proportion','Cumulative'};TabL=cell(NumPosi+1,4);TabL(1,:)=head;TabL(2:end,1)=num2cell(D);ifNumPosi==1TabL(2:end-1,2)={0};elseTabL(2:end-1,2)=num2cell(-diff(D));endTabL(2:end,3)=num2cell(D/sum(D));TabL(2:end,4)=num2cell(CumCont);%根據(jù)累積奉獻率的下限contri確定需要使用的鑒別式個數(shù)CumContGeConCumContGeCon=find(CumCont>=contri);CumContGeCon=CumContGeCon(1);V=V(:,1:CumContGeCon);%需要使用的鑒別式系數(shù)矩陣%以表格形式返回所用鑒別式的系數(shù)向量，若contri取值為1，%則返回所有鑒別式的系數(shù)向量.TabCan是一種元胞數(shù)組TabCan=cell(d+1,CumContGeCon+1);TabCan(1,1)={'Variable'};TabCan(2:end,1)=strcat('x',cellstr(num2str((1:d)')));TabCan(1,2:end)=strcat('can',cellstr(num2str((1:CumContGeCon)')));TabCan(2:end,2:end)=num2cell(V);%將訓(xùn)練樣品與待判樣品放在一起進行鑒別m=size(sampledata,1);gv=gmeans*V;stv=[sampledata;training]*V;nstv=size(stv,1);message='';outclass=NaN(nstv,1);fori=1:nstvobji=bsxfun(@minus,stv(i,:),gv);obji=sum(obji.