Fisher最优分割法的结合应用

主成分分析与Fisher 最优分割法的结合应用
一. 主成分分析计算步骤
1.计算相关系数矩阵
⎥⎥
⎥
⎥⎥⎦⎤
⎢⎢⎢⎢
⎢⎣⎡=pp p p p p r r r r r r r r r R 2
1
22221
11211
在上式中，r ij （i ，j=1，2，…，p ）为原变量的xi 与xj 之间的相关系数，其计
算公式为
∑∑∑===----=
n
k n
k j kj
i ki
n
k j kj i ki
ij x x
x x
x x x x
r 1
1
2
2
1
)()
()
)((
因为R 是实对称矩阵（即r ij =r ji ），所以只需计算上三角元素或下三角元素即可。

2．计算特征值与特征向量
首先解特征方程0=-R I λ，通常用雅可比法（Jacobi ）求出特征值
),,2,1(p i i =λ，并使其按大小顺序排列，即0,21≥≥≥≥p
λλλ ；然后分别求出
对应于特征值i λ的特征向量),,2,1(p i e i =。

这里要求i e =1，即112
=∑=p
j ij e ，其中
ij e 表示向量i e 的第j 个分量。

3．计算主成分贡献率及累计贡献率 主成分i z 的贡献率为
),,2,1(1
p i p
k k
i
=∑=λ
λ
累计贡献率为
)
,,2,1(11p i p
k k
i
k k
=∑∑==λ
λ
一般取累计贡献率达85—95%的特征值m λλλ,,,21 所对应的第一、第二，…，第m （m ≤p ）个主成分。

4． 计算主成分载荷 其计算公式为
)
,,2,1,(),(p j i e x z p l ij i j i ij ===λ
得到各主成分的载荷以后，还可以进一步计算，得到各主成分的得分
⎥
⎥
⎥⎥⎦
⎤⎢⎢⎢⎢⎣⎡=nm n n m m z z z z z z z z z Z 2
1
22221
11211
二．Fisher 最优分割法的聚类步骤
1.定义类的直径
设某一类G 包含的样品有()()(){}()1,,...,i i j X X X j i +>，记为{},1,...,G i i j =+。

该
类的均值向量G X :
为
()1
1j
G t t i X X j i ==-+∑:
用(),D i j 表示这一类的直径，常用的直径有：
()()()'
,j
G G t t t i D i j X X X X =⎛⎫⎛⎫
=-- ⎪ ⎪⎝
⎭⎝⎭∑::
2．定义分类损失函数
用(),b n k 表示将n 个有序样品分为k 类的某一种分法，常记分发(),b n k 为：
{}{}{}11,1222,23,1,...,1,1,...,1,.................................1,...,,
k k k G i i i G i i i G i i n =+-=+-=+
其中分点为()12111...11k k k i i i n i i n ++=<<<<=-=+即。

定义上述分类法的损失函数为
()()11
,,1k
t t t L b n k D i i +==-⎡⎤⎣⎦∑
3.(),L b n k ⎡⎤⎣⎦的地推公式
费希尔算法最核心的部分是利用以下俩个地推公式：
()()(){}()()(){}
2,2min 1,1,,,min 1,1,j n k j n L b n D j D j n L b n k L P j k D j n ≤≤≤≤⎧=-+⎡⎤⎣
⎦⎪⎨
=--+⎡⎤⎡⎤⎪⎣
⎦⎣⎦⎩ 4．最优解的求法
若分类数()1k k n <<已知，求分类法(),P n k ，使它在损失函数意义下达最小，其求法如下：
首先找分点k j ，使地推公式达极小，即
()()(),1,1,k k L P n k L P j k D j n =--+⎡⎤⎡⎤⎣⎦⎣⎦
于是得第k 类{},1,...,k k k G i i n =+。

然后找1k j -，使它满足
()()()111,11,2,k k k k L P j k L P j k D j j ----=--+⎡⎤⎡⎤⎣⎦⎣⎦
得到地 k-1类{}11,11,...,1k k k k G i i j ---=+-。

类似的方法一次可以得到所有类
12,,...,k G G G ，这就是我们所求的最优解，即
(){}12,,,...,k P n k G G G =
总之，为了求最优解，主要是计算
(){},;1D i j i j n ≤<≤和(){},,1,L P i j i n i j n ≤≤≤≤⎡⎤⎣⎦
三．具体事例的分析
（中学生身体四项指标的主成分分析）在某中学随机抽取某年级30名学生，测量其身高()
X ，体重()X ，胸围()X 和坐高()X ，数据见下表：
1.主成分分析
首先，创建TXT文件，将数据纯入并保存，命名lyzbok.txt。

然后创建各个M文件（lyzstd.m总和和标准化矩阵；lyzfac.m计算相关系数矩阵，特征值和特征向量，特征根排序，贡献率，主成分数，主成分载荷等相关操作；lyzscore.m计算得分；lyzprint.m输出最终结果）。

最后，在Matlab窗口键入以下语句，并最终显示器结果。

lyzprint('lyzbook1.txt',30,4)
fid =
6
标准化结果如下:
v1 =
0.0570 0.0310 0.0577 0.0276
0.0158 0.0318 0.0133 0.0325
0.0277 0.0607 0.0277 0.0597
0.0301 0.0169 0.0300 0.0122
0.0536 0.0310 0.0566 0.0264
0.0131 0.0322 0.0178 0.0305
0.0274 0.0559 0.0311 0.0637
0.0293 0.0125 0.0300 0.0195
0.0617 0.0273 0.0550 0.0325
0.0189 0.0298 0.0133 0.0357
0.0297 0.0563 0.0266 0.0613
0.0331 0.0117 0.0292 0.0146
0.0574 0.0257 0.0607 0.0300
0.0139 0.0298 0.0178 0.0325
0.0258 0.0647 0.0281 0.0585
0.0304 0.0189 0.0330 0.0146
0.0613 0.0314 0.0592 0.0276
0.0173 0.0338 0.0167 0.0309
0.0308 0.0635 0.0296 0.0572
0.0331 0.0197 0.0323 0.0122
0.0547 0.0314 0.0573 0.0272 0.0119 0.0334 0.0159 0.0309 0.0254 0.0563 0.0277 0.0564 0.0293 0.0133 0.0311 0.0130 0.0590 0.0269 0.0558 0.0276 0.0166 0.0310 0.0144 0.0296 0.0293 0.0551 0.0277 0.0601 0.0320 0.0125 0.0296 0.0154 0.0578 0.0265 0.0596 0.0284 0.0166 0.0294 0.0148 0.0317 相关系数矩阵:
std =
1.0000 -0.1689 0.9831 -0.1795 -0.1689 1.0000 -0.1495 0.9644 0.9831 -0.1495 1.0000 -0.1734 -0.1795 0.9644 -0.1734 1.0000 特征向量(vec)及特征向量(val):
vec =
0.6725 0.2213 0.4889 -0.5097 0.2190 -0.6711 0.5123 0.4891 -0.6748 -0.2084 0.4969 -0.5042 -0.2107 0.6762 0.5016 0.4968 val =
0.0147 0 0 0 0 0.0376 0 0 0 0 1.6381 0 0 0 0 2.3096 特征根排序:
2.30963
1.6381
0.0375937
0.0146792
贡献率:
newrate =
0.5774 0.4095 0.0094 0.0037 主成分数：2
主成分载荷:
-0.7746 0.6258
0.7433 0.6557
-0.7663 0.6360
0.7550 0.6420
计算得分:
score =
-0.0445 0.1104 0.0659 18.0000 0.0257 0.0600 0.0857 10.0000 0.0475 0.1131 0.1606 2.0000 -0.0245 0.0568 0.0323 27.0000 -0.0419 0.1067 0.0648 20.0000 0.0231 0.0602 0.0833 12.0000 0.0446 0.1145 0.1591 4.0000 -0.0217 0.0581 0.0364 25.0000 -0.0451 0.1124 0.0673 16.0000 0.0243 0.0627 0.0870 8.0000 0.0448 0.1117 0.1565 5.0000
-0.0284 0.0564 0.0280 30.0000 -0.0492 0.1107 0.0615 21.0000 0.0222 0.0604 0.0826 13.0000 0.0507 0.1140 0.1647 1.0000 -0.0238 0.0618 0.0380 24.0000 -0.0487 0.1143 0.0656 19.0000 0.0222 0.0634 0.0856 11.0000 0.0439 0.1165 0.1604 3.0000 -0.0266 0.0620 0.0354 26.0000 -0.0425 0.1087 0.0663 17.0000 0.0266 0.0593 0.0859 9.0000 0.0435 0.1067 0.1502 7.0000 -0.0269 0.0552 0.0283 29.0000 -0.0476 0.1078 0.0602 23.0000 0.0215 0.0589 0.0804 15.0000 0.0424 0.1106 0.1530 6.0000 -0.0266 0.0569 0.0304 28.0000 -0.0493 0.1097 0.0605 22.0000
0.0215 0.0594 0.0809 14.0000
从数据中可以得出，第一主成分与第二主成分的贡献率之和高达98.69%，因此只需要用俩个主成分就能很好的概括这组数据，进而由俩个特征值所对应的特征向量可以写出第一和第二主成分：
****
11234
****
21234
0.67250.21900.67480.21070.22130.67110.20840.6762Z X X X X Z X X X X =+--=--+ 由上可知，第一和第二主成分都是标准化好的变量()*1,2,3,4i X i =的线性组合。

2．Fisher 最优分割法
在主成分分析的基础上可知第一主成分1Z 和第二主成分2Z 。

进而分别以第一主成分1Z 作为数据指标进行Fisher 最优分割法，及第二主成分2Z 作为数据指标进行Fisher 最优分割法。

以第一主成分1Z 作为数据指标进行Fisher 最优分割法：
clear
%输入标准化后的数据
0.0158 0.0277 0.0301 0.0536 0.0131 0.0274 0.0293 0.0617 0.0189 0.0297 0.0331 0.0574 0.0139 0.0258 0.0304 0.0613 0.0173 0.0308 0.0331 0.0547 0.0119 0.0254 0.0293 0.0590 0.0166 0.0293 0.0320 0.0578 0.0166 ];
b=[0.0310 0.0318 0.0607 0.0169 0.0310 0.0322 0.0559 0.0125
0.0273
0.0298 0.0563 0.0117 0.0257 0.0298 0.0647 0.0189
0.0314
0.0338 0.0635 0.0197 0.0314 0.0334 0.0563
0.0133 0.0269 0.0310 0.0551
0.0265
0.0294 ]; c=[ 0.0577 0.0133 0.0277 0.0300 0.0566 0.0178 0.0311 0.0300 0.0550
0.0133 0.0266 0.0292 0.0607 0.0178 0.0281
0.0330 0.0592 0.0167 0.0296 0.0323 0.0573 0.0159
0.0277 0.0311 0.0558 0.0144
0.0277
0.0296 0.0596
0.0148];
d=[ 0.0276
0.0325
0.0597
0.0122
0.0264
0.0305
0.0637
0.0195
0.0325
0.0357
0.0613
0.0146
0.0300
0.0325
0.0585
0.0146
0.0276
0.0309
0.0572
0.0122
0.0272
0.0309
0.0564
0.0130
0.0296
0.0601
0.0154
0.0284
0.0317];
y1=lyzzcf(a,b,c,d); %计算第一主成分
for i=1:30
X(i)=y1(i);
end
%计算直径D
D=zeros(30,30);
for i=1:29
for j=2:30
if(i<j)
sx=0;
d=0;
for t=i:j
sx=sx+X(t);
end
meanx=sx/(j-i+1);
for k=i:j
d=d+(X(k)-meanx).^2;
end
D(j,i)=d;
end
end
end
D;
%计算分类损失函数
LP=zeros(30,29);
J=zeros(30,29);
for k=2:29
for l=3:30
if(k<l)
min=1000;
if(k==2) %k=2时的最小损失函数
for j=2:l
temp=D(j-1,1)+D(l,j);
if(temp<min)
min=temp;
minj=j;
end
end
else %k>2时的最小损失函数
for j=k:l
temp=LP(j-1,k-1)+D(l,j); if(temp<min)
min=temp;
minj=j;
end
end
end
LP(l,k)=min;
J(l,k)=minj;
end
end
end
LP;
J;
%绘制分类损失函数随k变化趋势图
for i=1:28
y(i)=LP(30,i+1);
end
k=2:29;
plot(k,y);hold on
plot(k,y,'r*')
成5类就可以达到置信水平的要求。

以第二主成分
Z作为数据指标进行Fisher最优分割法：
2
其过程与第一主成分的过程基本相同，下面只给出分类损失函数随k变化趋势图。

成5类就也可以达到置信水平的要求。

3．综合主成分分析与Fisher最优分割法得出结论
以上讨论过了单从主成分分析考虑，该数据大致可分为四类，在此基础上又对第一和第二主成分分别作为数据指标进行了Fisher最优分割法，得出该数据在分成五类的时候同时达到了各个主成分的要求，故而将元数据分成五类是最佳的方法。

附录
1．主成分分析函数
lyzstd.m总和和标准化矩阵：
%lyzstd.m用总和和标准化矩阵
function std=lyzstd(vector)
lyzsum=sum(vector,1); %对列求和
[a,b]=size(vector); %矩阵大小，a为行数，b为列数
for i=1:a
for j=1:b
std(i,j)=vector(i,j)/lyzsum(j);
end
end
lyzfac.m计算相关系数矩阵，特征值和特征向量，特征根排序，贡献率，主成分数，主成分载荷等相关操作：
%lyzfac.m
function result=lyzfac(vector);
fprintf('相关系数矩阵:\n')
std=CORRCOEF(vector) %计算相关系数矩阵
fprintf('特征向量(vec)及特征值(val):\n')
[vec,val]=eig(std) %求特征值(val)及特征向量(vec)
newval=diag(val);
[y,i]=sort(newval); %对特征根进行排序，y为排序结果，i为索引
fprintf('特征根排序:\n')
for z=1:length(y)
newy(z)=y(length(y)+1-z);
end
fprintf('%g\n',newy)
rate=y/sum(y);
fprintf('\n贡献率:\n')
newrate=newy/sum(newy)
sumrate=0;
newi=[];
for k=length(y):-1:1
sumrate=sumrate+rate(k);
newi(length(y)+1-k)=i(k);
if sumrate>0.85 break;
end
end %记下累积贡献率大于85%的特征值的序号放入newi中
fprintf('主成分数：%g\n\n',length(newi));
fprintf('主成分载荷: \n')
for p=1:length(newi)
for q=1:length(y)
result(q,p)=sqrt(newval(newi(p)))*vec(q,newi(p));
end
end %计算载荷
disp(result)
lyzscore.m计算得分：
%lyzscore.m计算得分
function score=lyzscore(vector1,vector2);
sco=vector1*vector2;
lsum=sum(sco,2);
[newlsum,i]=sort(-1*lsum);
[newi,j]=sort(i);
fprintf('计算得分: \n')
score=[sco,lsum,j] %得分矩阵:sco为各主成分得分;lsum为综合得分;j 为排序结果
lyzprint.m输出最终结果：
%lyzprint.m
function print=lyzprint(filename,a,b); %filename为文本文件文件名;a为矩阵行数(样本数),b为矩阵列数(变量指标数)
fid=fopen(filename,'r')
vector=fscanf(fid,'%g',[a b]);
fprintf('标准化结果如下:\n')
v1=lyzstd(vector)
result=lyzfac(v1);
lyzscore(v1,result);
第二主成分得分对第一主成分得分散布图程序：
a=[-0.0445
0.0257
0.0475
-0.0245
-0.0419
0.0231
0.0446
-0.0217
-0.0451
0.0243
0.0448
-0.0284
-0.0492
0.0222
0.0507
-0.0238
-0.0487
0.0222
0.0439
-0.0266
-0.0425
0.0266
0.0435
-0.0269
-0.0476
0.0215
0.0424
-0.0266
-0.0493
0.0215 ];
b=[ 0.1104
0.0600
0.1131
0.0568
0.1067
0.0602
0.1145
0.0581
0.1124
0.0627
0.1117
0.0564
0.1107
0.0604
0.1140
0.0618
0.1143
0.0634
0.1165
0.0620
0.1087
0.0593
0.1067
0.0552
0.1078
0.0589
0.1106
0.0569
0.1097
0.0594];
plot(a,b,'*')
2．Fisher最优分割法
lyzzcf.m求第一主成分函数：function zcf=lyzzcf(a,b,c,d) n=size(a);
for i=1:n
zcf(i)=0.6725*a(i)+0.2190*b(i)-0.6748*c(i)-0.2107*d(i); end
lyzzcf1.m求第二主成分函数：
fuction zcf=lyzzcf1(a,b,c,d)
n=size(a);
for i=1:n
zcf(i)=0.2213*a(i)-0.6711*b(i)-0.2084*c(i)+0.6762*d(i); end。