模式识别与智能信息处理实践实验报告

《模式识别与智能信息处理实践》实验报告书

专业：

学号：

姓名：

计算机科学与工程学院

2011.6.20

实验一聚类分析

一、实验目的

通过聚类分析实验，加深对聚类分析基本思想、方法的理解和掌握。

二、实验内容

了解动态、静态聚类算法的特点；

熟练掌握k-均值算法或层次聚类算法；

编写能对实际模式样本正确分类的算法程序。

掌握动态聚类算法的基本思想；

认识类别数、初始类心的选择对k-均值算法聚类结果的影响；

编写能对实际模式样本正确分类的k-均值算法程序。

三、方法手段

设类别数为k，选取k个初始聚类中心，按最小距离原则将各模式分配到k

类中的某一类，不断地计算类心和调整各模式的类别使每个模式特征矢量到其所

属类别中心的距离平方之和最小。

四、k-均值算法

（1）从D中随机取k个元素，作为k个簇的各自的中心。

（2）分别计算剩下的元素到k个簇中心的相异度，将这些元素分别划归到

相异度最低的簇。

（3）根据聚类结果，重新计算k个簇各自的中心，计算方法是取簇中所有

元素各自维度的算术平均数。

（4）将D中全部元素按照新的中心重新聚类。

（5）重复第4步，直到聚类结果不再变化。

五、k-均值程序运行结果

（1）改变初始类心，观察对聚类结果的影响

若选初始类心是[1 2 3]时的结果为其分为1类共39个，分为2类共61个，分为3类共50个，其中被分为第1类的样本为{51 53 78 101 103

104 105 106 108 109 110 111 112 113 116 117 118 119 121 123 125 126 129 130 131 132 133 135 136 137 138 140 141 142 144 145 146 148 149}，被分为第2

类的样本为{52 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 102 107 114 115 120 122 124 127 128 134 139 143 147 150}，被分为第3类的样本为{1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50}。

若选初始类心是[2 4 5]时其聚类结果为其分为1类共96个，分为2类共

22个，分为3类共个32，其中被分为第1类的样本为{51 52 53 54 55 56 57 59 60 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 95 96 97 98 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150}，被分为第2类的样本为{2 3 4 7 9 10 13 14 25 26 30 31 35 39 42 43 46 48 58 61 94 99}，被分为第3类的样

本为{1 5 6 8 11 12 15 16 17 18 19 20 21 22 23 24 27 28 29 32 33 34 36 37 38 40 41 44 45 47 49 50}。

可见，初始类心的选择对于K均值算法的影响较大。

（2）改变类别数k，比较其对类内距离平方和的大小的影响

若k=3且选初始类心是[1 2 3]时，其最终各类中心的距离为[1.1657 0.8230 4.0783]，若k=4且选初始类心是[1 2 3 4]时，其最终各类中心的距

离为[1.3049 0.4917 4.0783 1.3928]，可见，改变类别数k，其对类

内距离平方和的大小有直接的影响，且k越大, 其类内距离平方和距离越小。

六、实验总结

影响k-均值算法结果的因素有:初始类心的选择以及k值的选择，且k-均值算法适用于k值已知的情况，即类别号已知的情况。

七、Kmeans程序

function y=my_Kmeans(k,mid)

k=3;%类数

mid=[1 2 3]; %随便给三个聚类中心

%从文本文件读入数据放入X中

X=load('iris.txt','height','weight');

num=size(X,1);%获得X维数大小

for i=1:k

Z(i,:)=X(mid(i),:);%获取k个聚类中心的初始坐标

end

%计算新的聚类中心，K-均值算法的核心部分

temp=[];

while(~isequal(Z,temp)) %聚类中心是否变化,若不变化则停止循环

temp=Z;

class=cell(k,1);%初始化类样本class

value=cell(k,1);%初始化类样本的坐标value

for j=1:num

for t=1:k

D(t)=dist(X(j,:),Z(t,:)');%计算每个样本到类中心的距离

[minu,index]=min(D);%求出离聚类中心最小的一个样本

end

class{index}=cat(1,class{index},j);%将该样本归于一类

value{index}=cat(1,value{index},X(j,:));%存放该类样本的坐标

end

for i=1:k

Z(i,:)=mean(value{i});%计算k类样本的均值，更新聚类中心

end

end

celldisp(class);%显示Kmeans聚类结果

D %显示最终类间距离

实验二判别域代数界面方程法

一、实验目的

通过实验加深对判别域代数界面方程法基本思想、方法的认识，了解线性判别函数的分类能力，熟练掌握感知器算法，或H-K算法，或累积势函数分类法，以及它们各自的适用条件，加深对有监督训练、线性可分、解空间、权空间等概念的理解，编写能对实际模式样本正确分类的算法程序。

二、实验内容

编写能对实际模式样本正确分类的感知器算法，或H-K算法，或累积势函数分类法的算法程序，能对实际模式样本进行正确分类。

三、方法手段：

设已知类别模式特征矢量集为{x

1,x

2

,…,x

N

}，类别数为c。

感知器算法是基于一次准则函数的梯度下降法。从任意初始解矢量w

出发，

利用产生错分的样本对解矢量进行迭代校正：w

k =w

k-1

+ρx

i

，i=1,2,…,c，

k=0,1,2…，ρ=常数，从而解得线性判决函数d(x)=w’x的解矢量w.

四、感知器算法

1、算法思想

校正方法实际上是最优化技术中的负梯度下降法。该算法也是人工神经网络理论中的线性阈值神经元的学习算法。

2、算法原理步骤

设给定一个增广的训练模式集{}N x x x ,,,2

1，其中每个模式类别已知，它们分属1ω类和2ω类。

⑴ 置步数1=k ，令增量=ρ某正的常数，分别赋给增广权矢量初值)1(w

的

各分量较小的任意值。

⑵ 输入训练模式k x 。 ⑶ 计算判别函数值k x k w )('。 ⑷ 调整增广权矢量，规则是 ① 如果1ω∈k x 和0)(≤'k x k w ，则k x k w k w

ρ+=+)()1(；(偏小，故加) ② 如果2ω∈k x 和0)(≥'k x k w ，则k x k w k w ρ-=+)()1(；(偏大，故减) ③ 如果1ω∈k x 和0)(>'k x k w ，或2ω∈k

x 和0)(<'k x k w

，则)()1(k w k w

=+。（分类正确，不校正）

⑸ 令1+=k k 。如果N k ≤，返至⑵。如果N k >，检验判别函数x w

'对N x x x

,,,21是否都能正确分类。若是，结束；若不是，令1=k ，返至⑵。

如果训练模式已经符号规范化，即2ω∈k x 已乘以－1（包括增广分量1），

则校正权矢量的规则可统一为

???ρ+=+k x k w k w k w )()()1( )(0)()

(0)(错误分类若正确分类若≤'>'k k x k w x k w

在用全部模式训练完一轮后只要还有模式被判错，则需要进行第二轮迭代，用全部训练模式再训练一次，建立新的权矢量。如果对训练模式还有错分的，则

进行第三轮迭代，余类推，直至对所有训练模式均能正确分类为止，此时的w 即

为所求的权矢量。

五、感知器算法实验结果

（1）改变权矢量初值，观察对算法的影响

已知w1:X1=(0 0)’,X2=(1 0)’,X3=(1 1)’,X4=(1 2)’；

w2:X5=(1 3)’,X6=(2 3)’,X7=(3 2)’,X8=(3 3)’，

若步长因子ρ=1不变，权初始值取w=[0 0 0],则需迭代eth0=15步，若权初始值取w=[1 1 1]，则需迭代eth0=18步，可见，权初始值选取的不同，将直接导致算法的收敛速度不同。

（2）改变步长因子ρ，观察对算法的影响

若权初始值取w=[0 0 0]不变，取步长因子ρ=1，则需迭代eth0=15步，若取步长因子ρ=0.1，则需迭代eth0=16步，可见，步长因子ρ选取的不同，也将直接导致算法的收敛速度不同。

（3）算法的适用性

若已知w1:X1=(1 3)’,X2=(2 3)’,X3=(1 1)’,X4=(1 2)’；

w2:X5=(0 0)’,X6=(1 0)’,X7=(3 2)’,X8=(3 3)’，则导致该程序死循环，可见，感知器算法只适用两类样本线性可分的条件下（图1），若希望在两类样本线性不可分的条件下，构造一判决平面，则可考虑使用势函数。

图1 感知器算法分类结果

图2 线性不可分的情况

六、实验程序

function y=my_preception(W1,W2)

%感知器算法，对两类问题生成线性判别函数

%注意：前提是两类可线性判别，否则将死循环%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

W1=[0 0;1 0;1 1;1 2];%w1类中的样本

W2=[1 3;2 3;3 2;3 3];%w2类中的样本

%%%%%%初始化

w=[0 0 0]; %任取w初值

c=1; %任取校正增量系数c=1；

%%%%%%感知器算法部分%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%增广样本w1和w2乘(-1)

[m,n]=size(W1);

for i=1:m

W1(:,n+1)=1;%w1增广

end

[m,n]=size(W2);

for i=1:m

W2(:,n+1)=1;%w2增广

end

W2=-W2;%取反

%将增广向量转换成元组，便于处理

M=ones(1,m);

w1=mat2cell(W1,M,n+1);

w2=mat2cell(W2,M,n+1);

X=cat(1,w1,w2);%合并两类增广样本%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

eth0=0; %迭代次数

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%感知器算法核心部分

[m,n]=size(X);

temp1=[];temp2=1;%随意两个赋不等的值

while(~isequal(temp1,temp2)) %判断权值是否变换化，若不变化，则终止循环 eth0=eth0+1;

temp2=temp1;

clear temp1;%不清除此变量不行啊！！

for i=1:m

temp1{i,1}=w;

if(w*X{i}'<=0) %w乘X{i}的转置

w=w+c*X{i}; %若小于0更新权值

end

end

end %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

w %显示线性分类器权值

eth0 %显示迭代次数

实验三 Bayes 分类器设计

一、实验目的

通过实验，加深对统计判决与概率密度估计基本思想、方法的认识，了解影响Bayes 分类器性能的因素，掌握基于Bayes 决策理论的随机模式分类的原理和方法。

二、实验内容

设计Bayes 决策理论的随机模式分类器。

三、方法手段

Bayes 分类器的基本思想是依据类的概率、概密，按照某种准则使分类结果从统计上讲是最佳的。换言之，根据类的概率、概密将模式空间划分成若干个子空间，在此基础上形成模式分类的判决规则。准则函数不同，所导出的判决规则就不同，分类结果也不同。使用哪种准则或方法应根据具体问题来确定。

四、Bayes 算法

1.实验原理

多元正太分布的概率密度函数由下式定义

1

1

2211()exp ()()2(2)T d p X X X μμπ-??=--∑-????

∑ 由最小错误概率判决规则，可得采用如下的函数作为判别函数

()(|)(),

1,2,,i i i g x p X P i N ωω==

这里，()i P ω为类别i ω发生的先验概率，(|)i p X ω为类别i ω的类条件概率密度函数，而N 为类别数。

设类别i ω，i=1,2,……,N 的类条件概率密度函数(|)i p X ω，i=1,2,……,N 服从正态分布，即有(|)i p X ω~(,)i i N μ∑，那么上式就可以写为

112

2

()1()exp ()(),

1,2,,2(2)T i i d P g X X X i N ωμμπ-??=

--∑-=????

∑

由于对数函数为单调变化的函数，用上式右端取对数后得到的新的判别函数替代原来的判别函数()i g X 不会改变相应分类器的性能。因此，可取

111()()()ln ()ln ln(2)222

T i i i i i i d

g X X X P μμωπ-=--∑-+-∑-

显然，上式中的第二项与样本所属类别无关，将其从判别函数中消去，不会

改变分类结果。这样，判别函数()i g X 可简化为以下形式

111

()()()ln ()ln 22

T i i i i i i g X X X P μμω-=--∑-+-∑

2.实验步骤

1、求出两类样本的均值

11,2i

X i

X i N ω

μ∈

=

=∑

2、求每一类样本的协方差矩阵

1

1()()

,1,21i i

i N w w i jk

lj j lk k l i s x x j k N μμ==--=-∑

式中，l 代表样本在类中的序号，其中

lj x 代表i w 类的第l 个样本，第j 个特征值；

i

w j μ代表i w 类的i N 个样品第j 个特征的平均值

lk x 代表i w 类的第l 个样品，第k 个特征值；

i

w k μ代表i w 类的i N 个样品第k 个特征的平均值。

i w 类的协方差矩阵为

1112

2122i i

i

i

i ??∑∑∑= ?

∑∑?? 3、计算出每一类的协方差矩阵的逆矩阵1i -∑以及协方差矩阵的行列式i ∑ 4、求出每一类的先验概率

()/1,2i i P N N i ω≈=

5、将各个数值代入判别函数

111

()()()ln ()ln 22

T i i i i i i g X X X P μμω-=--∑-+-∑

判别边界为

12()()0g X g X -=

五、Bayes 分类器实验结果

已知（图1）数据a=[0.3760 0.0240 0.2440 -0.1740 0.0460 -0.3940 0.3760 0.7720 0.2660 0.5080 -0.4380 -0.0640 0.8160 0.5960 0.1120 0.3540 0.8380 -0.7680 0.4200 -0.7900]；其满足正态分布（图2）， 1. 最小错误率贝叶斯决策

图1 样本数据

图2 样本的类条件概率

根据最小错误率准侧，计算其后验条件概率（图3），通过程序运行出结果细胞分类结果为：

1 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 1 1 0，其中，0为判成正常细胞，1为判成异常细胞。

图3 后验条件概率

2. 最小风险贝叶斯决策

根据最小风险判别准侧，其损失函数赋值为r=[0 1000 0;200 0 0],则计算其条件风险概率（图4）通过程序运行出结果细胞分类结果为：

1 1 1 1 1 1 1 0 1 1 1 1 0 1 1 1 0 1 1 1，其中，0为判成正常细胞，1为判成异常细胞。

图4 条件风险概率

3.两类分类器结果不同原因分析

由最小错误率的贝叶斯判决和基于最小风险的贝叶斯判决得出图形中的分类结果可以看出，样本0.0240，0.2440等在前者中被分为“正常细胞”，在后者

被分为“异常细胞”，分类结果不同。因为在给予最小风险贝叶斯判决中，影响决策结果的因素多了损失r这一项，所以当结合最小风险贝叶斯决策表进行计算时，‘损失’起了主导作用，导致出现两者结果的不一致。

六、Bayes分类器程序代码

function y=my_bayes(n,a)

%%%%%%%%%%%%%最小错误率贝叶斯决策

% 构造实验数据

n=20; % 样本数

a=(round(100*rand(n,1))/100)*2.2-0.9;

% 样本数为n,特征数为1，数据在-0.9与1.3之间

figure

plot(1:n,a,'rx')

xlabel('样本数');

ylabel('生化化验值');

title('样本数据：生化化验值');

pause;

% 先验概率

P=[0.9 0.1];

% 作类条件概率密度曲线p(x|wi)

x=-0.9:0.01:1.3;

px(1,:)=(1/(sqrt(2*pi)*0.3))*exp(-0.5*(x/0.3).^2 );

px(2,:)=(1/(sqrt(2*pi)*0.1))*exp(-0.5*((x-1)/0.1 ).^2);

figure;

plot(x,px(1,:),'b',x,px(2,:),'r--');

xlabel('生化化验值');

ylabel('概率密度');

title('类条件概率密度曲线')

axis tight;

pause;

% 作后验概率曲线

for i=1:2

pwx(i,:)=px(i,:)*P(i)./(px(1,:)*P(1)+px(2,:)*P(2));

end

figure;

plot(x,pwx(1,:),'b',x,pwx(2,:),'r--');

xlabel('生化化验值');

ylabel('后验概率');

title('后验概率曲线')

axis tight;

pause;

% 计算给定生化化验值的类条件概率密度曲线

for j=1:n

s=a(j);

PXW1=spline(x,px(1,:),s);

PXW2=spline(x,px(2,:),s);

PXW=[PXW1,PXW2];

disp('样本')

s

%计算后验概率,判断输出

for i=1:2

Pwx(i)=PXW(i)*P(i)/(PXW(1)*P(1)+PXW(2)*P(2));

end

disp('后验概率 P(wi|x)=')

Pwx

plot(x,pwx(1,:),'b',x,pwx(2,:),'r--');

xlabel('生化化验值');

ylabel('后验概率');

title('后验概率曲线')

hold on

plot(s,Pwx(1),'or',s,Pwx(2),'ob');

axis tight;

hold off

if Pwx(1)>Pwx(2)

w(j,1)=s;

disp('正常人')

else

w(j,2)=s;

disp('感染病人')

end

pause;

end

disp('========================================') disp('正常人感染病人')

w

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%最小风险贝叶斯

n=20; % 样本数

figure

plot(1:n,a,'rx')

xlabel('样本数');

ylabel('生化化验值');

title('样本数据：生化化验值');

pause;

% 先验概率

P=[0.9 0.1];

% 作类条件概率密度曲线p(x|wi)

x=-0.9:0.01:1.3;

px(1,:)=(1/(sqrt(2*pi)*0.3))*exp(-0.5*(x/0.3).^2 );

px(2,:)=(1/(sqrt(2*pi)*0.1))*exp(-0.5*((x-1)/0.1 ).^2);

figure;

plot(x,px(1,:),'b',x,px(2,:),'r--');

xlabel('生化化验值');

ylabel('概率密度');

title('类条件概率密度曲线')

axis tight;

pause;

% 作后验概率曲线

for i=1:2

pwx(i,:)=px(i,:)*P(i)./(px(1,:)*P(1)+px(2,:)*P(2));

end

figure;

plot(x,pwx(1,:),'b',x,pwx(2,:),'r--');

xlabel('生化化验值');

ylabel('后验概率');

title('后验概率曲线')

axis tight;

pause;

% 损失函数

r=[0 1000 0;200 0 0];

%作条件风险曲线

for i=1:2

R(i,:)=r(i,1)*pwx(1,:)+r(i,2)*pwx(2,:);

end

figure;

plot(x,R(1,:),'b',x,R(2,:),'r--');

xlabel('生化化验值');

ylabel('条件风险');

title('条件风险曲线')

axis tight;

pause;

% 计算给定生化化验值的条件风险

for j=1:n

s=a(j);

PXW1=spline(x,px(1,:),s);

PXW2=spline(x,px(2,:),s);

PXW=[PXW1,PXW2];

disp('样本')

s

%计算后验概率

for i=1:2

Pwx(i)=PXW(i)*P(i)/(PXW(1)*P(1)+PXW(2)*P(2));

end

%计算条件风险,判断输出 for i=1:2

Rx(i,:)=r(i,1)*Pwx(1)+r(i,2)*Pwx(2);%%%%% end

disp('条件风险 R(ai|x)=')

Rx

plot(x,R(1,:),'b',x,R(2,:),'r--');

xlabel('生化化验值');

ylabel('条件风险');

title('条件风险曲线')

hold on

plot(s,Rx(1),'or',s,Rx(2),'ob');

axis tight;

hold off

if Rx(1)>Rx(2)

w(j,1)=s;

disp('正常人')

else

w(j,2)=s;

disp('感染病人')

end

pause;

end

disp('========================================') disp('正常人感染病人')

w

实验四特征提取与选择

一、实验目的

通过实验加深对特征提取与选择原理、方法的认识，编写基于离散K-L变换的特征提取程序，或分支定界特征挑选算法程序，掌握离散K-L变换原理、特点。编写基于离散K-L变换的特征提取程序。

二、实验内容

编写基于离散K-L变换的特征提取程序，主要使用主成分分析法（PCA）进行图像识别。

三、方法手段：

设n维特征矢量x=(x

1,x

2

,…,x

n

)’，其均值矢量m=E[x]，协方差矩阵

C=E[(x-m)(x-m)’]，设λ

1≥λ

2

≥…≥λ

n

是C的本征值，t

i

是λ

i

对应的本征矢量，

i=1,2,…,n. x经标准正交矩阵T’=(t

1,t

2

,…,t

n

)’正交变换后成为矢量

y=(y

1,y

2

,…,y

n

)’，即，y= T’x。这种正交变换称为离散K-L变换（间记为DKLT）。

DKLT的性质：

（1）变换后y各分量正交或不相关；

（2）变换后y各分量的方差更趋于不均匀，可分性强于x；

（3）最佳逼近性，即，若用y的前d

当n较大时，由DKLT的性质（2），可以用y的前d

特征脸方法是基于K-L变换的人脸识别方法，K-L变换是图像压缩的一种最优正交变换。高维的图像空间经过K-L变换后得到一组新的正交基，保留其中重要的正交基，由这些基可以张成低维线性空间。如果假设人脸在这些低维线性空间的投影具有可分性，就可以将这些投影用作识别的特征矢量，这就是特征脸方法用于人脸识别的基本思想。在人脸识别中，可以用离散K-L变换对人脸图像的原始空间进行转换，即构造人脸图像数据集的协方差矩阵，对之进行正交变换，求出协方差矩阵的特征向量，再依据特征值的大小对这些特征向量进行排序，每一个向量表示人脸图像中一个不同数量的变量，这些特征向量表示特征的一个集合，它们共同表示一个人脸图像。在人脸识别领域，人们常称这些特征向量为特征脸。每一个体人脸图像都可以确切地表示为一组特征脸的线性组合。这样我们首先通过有指导的训练（给定训练样本集已知分类）得到样本集在特征脸空间中的坐标。训练完成后，输入待辨识图像，求得在特征脸空间的坐标，采用最近邻法，就可以实现人脸识别。一般说K-L变换是均方差意义下的最佳变换，理论推导如下。

假设任何信号X 都可以被表示为：

∑∞

==1i i i u c X ，其中?

??≠==j i j

i u u j T

i

,0,1 (1.1)

K-L 变换是一种正交变换，即将一个向量X ，在某一种坐标系统中得描述，转换成用另一种基向量组成的坐标系表示。这组基向量是正交的，其中每个坐标基向量用u i 表示，c i 为其分量数值，变换的目的是如何将(1.1)表示的无限多维基向量系统改成有限d 维坐标系近似。如X 的近似值或者估计量表示为：

∑==d

i i i u c X 1

? (1.2)

要使向量X 的估计量误差小，即使引起的均方误差

[]

)?()?(X X X

X E T --=ε (1.3) 为最小，而K-L 变换就实现了这个目的，下面详细证明。

求某一个数据X 的相应的ci 值，可以通过X 与对应某个基u i 的点积来计算。由于不同的基之间是相互正交的，这个点积值就是c i 的值，即X

u c T i i =，如果我

们要求一组系数c i ,并将其表示成一个向量形式T c c C ,...),(21=,则我们可以

得：

UX X u u C T T =???

?

? ??=...21 (1.4)

则U 就是一个变换矩阵，其中每一行是某一个正交基向量的转置，由X 计算C 称为对X 的分解。反过来，如果我们希望用C 重构信号X ，则根据(1.2)，它是各个成分

之和。如果我们将重构信号表示为T d

x x x X ),...,,(?21=，则 ()C

U c c c u u u x x X T d d d =?

?????

?

??=?????

??= 21211,..,,? (1.5) 我们称∑∞

+==-1?d i i i u c X X 为残差，带进(1.3)中我们可以得到：

[

][][]

∑∑∑∞+=∞+=∞+==?=--=12

11)?()?(d i i d i d j T j j T i i T c E u c u c E X X X X E ε

另外：i T

j j T i j j j j T i

T i i u X u u c u c u

X u c ∑∑∞

=∞=====11，故：

[

][

]

∑∑

∑

∑

∞

+=∞

+=∞

+=∞

+=ψ===11112

][d i i

T i d i i T

T

i

d i i T

T

i

d i i

u u u XX E u u XX u E c E 欲使该均方误差为最小，就变成在保证正交变换的条件下，使ε最小的问题，设一函数： ，并对其进行求导可得：

∞+==-ψ,...,1,0)(d i u I i i λ (1.6)

可见向量u i 应是ψ矩阵的特征值i λ的特征向量，而此时截断误差为

∑∞

+==1d i i λε，如将按其大小顺序排列，即 ≥≥≥≥d λλλ21，则取前d 项特征值对应的特征向量组成的坐标系，可使向量的均方误差为最小。满足上述条件的变换就是K-L 变换。

五，K-L 变换人脸识别实验结果分析 1.数据准备

本实验数据是ORL 标准人脸数据库数据

1）给定样本集，在这里样本集的大小为400，一共是40人，每人10张图片。我选取每人5张图片一共200张作为训练集，训练集的大小是可选的，它的大小将直接影响到识别的正确率，在后面我对采用不同大小的训练集进行识别的正确率进行了统计。

2）确定采用何种产生矩阵，有总体散布矩阵{()()}T E x x μμ--和类间散布矩阵1

0()()()p T b i i i i S p ωμμμμ-==--∑，后面的例子程序中分别采用了这两种产生矩

阵，并进行了识别正确率的比较。

3）根据产生矩阵计算出特征脸，这里要用到奇异值分解，其基本原理就是通过计算较低维数矩阵的特征值与特征向量而间接求出较高维数矩阵的特征向量（特征脸），从而达到简化计算的目的。

4）分别求得训练集中各图像在特征脸空间中的坐标。

2.识别过程

首先求得待辨识图像相对于训练集平均脸的差值图像，然后求得该图像在特征脸空间中的坐标，最后采用最近邻法对图像进行归类。

3.实验结果

1）训练样本的过程中生成的一些特征脸（图1）。

∑

∑∞+=∞+=--ψ=11)1()(d i i T i i d i i T i i u u u u u g λ

图1 训练样本的过程中生成的一些特征脸2）识别的正确率统计

图2 基于两种产生矩阵的识别正确率统计

4.实验结论

通过对分别采用类间散布矩阵和总体散布矩阵为产生矩阵的K-L变换进行人脸识别的正确率的统计可以看出，应用总体散布矩阵作为产生矩阵的K-L变换进行人脸识别的正确率明显要高。我想主要是因为类间散布矩阵采用每个人的平均图像向量组成训练集，来生成特征脸空间，计算量较之采用总体散布矩阵作为产生矩阵的K-L变换计算量又进一步减少，但是它忽略了每个人的图像的多样性，造成识别率低。以总体散布矩阵作为产生矩阵的K-L变换进行人脸识别的计算量较之采用类间散布矩阵作为产生矩阵的K-L变换要大一些，但是它充分考虑了每个人图像的多样性，使得识别率相对较高。

六，实验程序

function y=mypca()

%%%%%%%%%%%%%%%%%%%%%%%%%PCA算法对人脸图像处理提取主成分程序

path = ['D:\PCA\']; % ['D:\PCA\ORL']

%读取图像

numimage=400; %40*10=400张人脸

imagepath=[path 'ORL\ORL001' '.bmp']; %第一张人脸文件的路径及文件名：D:\PCA\ORL\ORL001.bmp immatrix=imread(imagepath); % 读入第一张人脸文件，构成矩阵immatrix

[m,n]=size(immatrix); % 计算矩阵immatrix的行数m、列数n

DATA = uint8 (rand(m*n, numimage)); %随机生成m*n行、numimage列的矩阵，并取uint8

for i=1:numimage

s1=floor(i/100); % 取整,求第3位

tem=rem(i,100); % i除以100的余数，取后两位

s2=floor(tem/10); % 取第2位

s3=rem(tem,10); % 取第1位

imagepath=[path 'ORL\ORL' int2str(s1) int2str(s2) int2str(s3) '.bmp']; % 构成图像文件的路径即文件名

immatrix=imread(imagepath); % 读入每一张人脸文件，构成矩阵immatrix

imVector=reshape(immatrix,m*n,1); % 将矩阵immatrix转化为一个列向量，长度为m*n

DATA(:,i)=imVector; % 将列向量imVector依次加入到DATA矩阵的列中.DATA先随机生成过的end

clear i;clear j;

save DATA DATA; % 保存DATA

mn=mean(double(DATA'))'; % 计算DATA的行向量的均值

save mn mn; % 保存DATA的行向量的均值

%image substracted by mean of all train images

DATAzeromean=double(DATA)-repmat(mn,1,numimage);

save DATAzeromean DATAzeromean;

clear DATA;

L=DATAzeromean'*DATAzeromean;

[V,D]=eig(L);

enginvalue=diag(D);

[enginvalue,ix]=sort(enginvalue);%按升序排列矩阵元素

ix=flipud(ix);%从上到下翻转矩阵,即按降序

V=V(:,ix); %对V的特征向量位置调整

facespace=DATAzeromean*V; %脸空间

for t=1:numimage

facespace(:,t)=facespace(:,t)/norm(facespace(:,t));%Normalisation to unit length

end

subdim=200;

facespace=facespace(:,1:subdim);%选择子特征向量的协方差矩阵---facespace是400*200维，即前200个样本

projdata=facespace'*DATAzeromean;

save projdata projdata;

save facespace facespace;

datareconstruct=facespace*projdata;

fprintf('正在保存 Wakesplace中的图片数据\n');

save datareconstruct datareconstruct;

for i=1:numimage

imdata=datareconstruct(:,i);

imdata=reshape(imdata,m,n);

imwrite(imdata,[int2str(i) '.bmp'],'bmp');%得到重构图像1.bmp---400.bmp

end