用Matlab实现AHP的算法

1．MATLAB的基本内容
MATLAB（MATrix LABoratory，矩阵实验室的缩写）是一种特殊用途的计算机程序优化执行工程和科学计算。

它开始为旨在执行矩阵数学程式的生活，但多年来它已发展成为一个灵活的计算系统基本上能够解决任何技术问题。

MATLAB具有编程语言的基本特征，使用MATLAB也可以使用像BASIC、FORTRAN、C等传统编程语言一样，进行程序设计，而且简单易学、编程效率高。

正因为MATLAB的强大的功能，使得它在许多领域得到广泛应用。

在科研与工程应用领域，MATLAB已被广泛地用于科学研究和解决各种具体的实际问题。

许多科技工作者选用MATLAB做为计算工具，避免了繁琐的底层编程，从而可以把主要精力和时间花在科学研究和解决实际问题是上，提高了工作效率。

MATLAB矩阵
矩阵是MATLAB的基本处理对象，因此根据本文所需，简单介绍所涉及MATLAB 矩阵内容。

MATLAB矩阵的建立
1、直接输入法
最简单的建立矩阵的方法是从键盘直接输入矩阵的元素。

例如：
A=[1 2 3;4 5 6;7 8 9]
A =
1 2 3
4 5 6
7 8 9
也可以用回车键代替分号，按下列方式输入：
A=[1 2 3
4 5 6
7 8 9 ]
2、利用M文件建立矩阵
比较大且复杂的矩阵，可以为它专门建立一个M文件，如同下例。

利用M文件建立矩阵。

启动有关编辑程序或MATLAB文本编辑器，并输入待建矩阵：
MYMAT=[ 111 ,112,113,114,115,116,117,118,119;
211,212,213,214,215,216,217,218,219];
把输入的内容以纯文本方式存盘（设文件名为）。

在MATLAB命令窗口中输入mymatrix，即运行该M文件，就会自动建立一个名为MYMAT的矩阵，可供以后使用。

矩阵的特征值与特征向量
特征值和特征向量在科学研究和工程计算中都有非常广泛地应用。

在MATLAB 中，计算矩阵A的特征值和特征向量的函数是eig(A),常用的调用格式有3种[14]：
E = eig( A ) ：求矩阵A的全部特征值，构成向量E。

[V,D]=eig(A)：求矩阵A的全部特征值，构成对角矩阵D，并求A得特征向量构成V的列向量。

[V,D]=eig(A,’nobablance’)：与第2种格式中先对A作相似变换后求矩阵A的特征值和特征向量，而格式3直接求矩阵A的特征值和特征向量。

例如：
A=[1,1,;1,1,;,,2];
[V,D]=eig(A)
V =
D =
0 0
0 0
0 0
求得的3个特征值是、和，各特征值对应的特征向量为V的各列构成的向量。

MATLAB的M文件
M文件可以根据调用方式的不同分为两类：命令文件（Script File）和函数文件（Function File）。

它们的扩展名均为.m。

函数文件由function语句引导，其基本结构为：
function 输出形参表=函数名（输入形参表）
注释说明部分
函数体语句
我们通过举例说明如下：
例2-2 分别建立命令文件和函数文件，将求矩阵的一致性指标CI：CI=(λmax-n)/(n-1)程序1 建立命令文件并以文件名存盘：
max=input('please input max:');
n=input('please input n:');
CI=(max-n)/(n-1)
然后在MATLAB的命令窗口中输入CI即可。

程序 2 建立函数文件。

function c=CI(max,n)
c=(max-n)/(n-1)
然后在MATLAB的命令窗口调用该函数文件。

max=input('please input max:');
n=input('please input n:');
c=CI(max,n)
2．基于MATLAB的AHP实现
2．1 AHP的MATLAB的计算流程框图
根据层次分析法的一般步骤我们得到在MATLAB工具上实现的计算程序流程框图，如图2所示[16]
图 2 以MATLAB实现的层次分析法的计算流程框图通过流程框图，层次分析的基本步骤如下：
第一步：准则层对目标层的判断矩阵归一化且判断是否满足一致性；
第二步：第一步满足时，将方案层对准则层的判断矩阵归一化并判断其一致性； 第三步：当第一、二步满足时，求方案层的总排序权值与总CR 并判断一致性。

2．2 平均随机一致性指标的MATLAB 实现
运用层次分析法决策者需要通过反复地解决决策问题，将同一层次的各元素与上一层次中某一准则的重要性进行比较，从而构造出两两判断比较矩阵A=(a ij )nn (称为成对比较矩阵)。

前面已经描述了九级标度法，此处运用其他描述，则这些成对比较矩阵应满足如下条件：
(l)ij a >0 (2) ij a ·ji a = l (3) ii a =l
按照事物逻辑要求，该矩阵还应具备一致性，即满足：
ij a ·jk a =ik a
前面已经给出由于客观事物的复杂性与决策者的认识的多样性，实际问题的成对比较矩阵不可能做到严格上的一致性，因而，借助平均随机一致性指标RI 来相对判定其一致性程度。

其中表1-4是已经计算好的1~15阶矩阵的RI 值表，但未给出其实现过程，且各文献的RI 值表不完全相同。

究其原因除没有太大必要介绍外，真正去实现它却有如下三个难度：
(1)随机两两判断矩阵中的元素要求是1~9和它们对应的倒数共17个整数与小数的均匀分布很难处理。

(2)一般高级编程语言实现成对比较矩阵及相关计算，非常复杂，且占用内存巨大，耗时多。

(3) 随机种子源不能控制。

本文使用数学软件包MATLAB 对其进行计算。

其设计解决思路为：先用软件包中随机函数产生数1~17的均匀分布的n 阶矩阵，然后在软件包中采用不同技巧将它转化为成对比较矩阵，最后用循环语句计算出RI 值。

结果如下表：
表3-1 计算的RI值
以n=4为例，过程详见如下程序清单，其中随机成对比较矩阵的实现见相应注释部分。

[15]
MATLAB的程序M文件：
function ri %计算RI值的命令文件
n=4; ri=0; m=100;
rand('seed',21) %控制随机发生器
for i=1:m
a=ceil (17*rand(n)); %产生n阶l~17的随机阵
a(find(a=8))=; %消除0为分母
b=1./(a-8); %产生一个辅助阵
a(find(a>9))=b(find(a>9)); %借助b，将9~17分别转化为1
2
~
1
9
a(find(a=)=8;
e=eye(n); %产生一个4阶单位阵c=1./a; %将a中每个元素换成相应倒数c=c'; %将c转置
c=tril(c,-1); %抽取c的下三角(不含主对角线) a=triu(a,1); %抽取a的上三角(不含主对角线) a=a+c+e; %实现随机成对比较阵a k=size(a,1); %计算a的行维数
ri=ri+(max(abs(eig(a)))-k)/(k-1); %计算100次RI值
end
ri/m %计算平均RI值
2．3 AHP各环节的MATLAB实现
以目标矩阵A=
135
1313
1511
⎛⎫
⎪
⎪
⎪
⎝⎭
,准则层矩阵为P1=
12
121
⎛⎫
⎪
⎝⎭
, P2=
135
113
15131
⎛⎫
⎪
⎪
⎪
⎝⎭
，
P3=
12
121
⎛⎫
⎪
⎝⎭
为例，运用MATLAB进行数据处理。

特征向量及其归一化的MATLAB实现
MATLAB中求矩阵特征值和特征向量的函数是eig，其调用的格式为[V,D]=eig(A),其中，V为特征向量矩阵，D为特征值矩阵。

层次分析法中需要求得是最大特征值及对应的归一化特征向量，而且考虑到eig函数在求得的特征值中可能会存在复数。

因此，运用直接输入程序代码会产生一定的误差。

在此需要对求得的V、D进行适当选择，定义一个M-file 来实现。

function[maxeigval,w]=maxeigvalvec(A) %求最大特征值及对应的归一化特征向量
%A为判断矩阵
[eigvec,eigval]=eig(A);
eigval=diag(eigval); %特征向量
eigvalmag=imag(eigval);
realind=find(eigvalmag<eps);
realeigval=eigval(realind); %实特征根
maxeigval=max(realeigval) %最大特征值
index=find(eigval==maxeigval);
vecinit=eigvec(:,index); %最大特征值对应的特征向量
w=vecinit./sum(vecinit) %特征向量归一化
在MATLAB中键入如下指令：
A= [1,3,5;1/3,1,3;1/5,1/3,1];
P1=[1,2;1/2,1];
P2= [1,3,5;1/3,1,3;1/5,1/3,1];
P3=[1,2;1/2,1];
[max(1),wA]= maxeigvalvec(A); [max(2)wP1]= maxeigvalvec(P1); [max(3),wP2]= maxeigvalvec(P2); [max(4),wP3]= maxeigvalvec(P3); MATLAB运行结果如下：maxeigval=
w =
maxeigval=
2
w=
maxeigval=
w =
maxeigval=
2
w=
一致性检验及单排序的MATLAB实现
由AHP的MATLAB的计算流程图知，必须对各层次间的判断矩阵进行层次单排序和一致性检验。

因此，定义函数来实现层次单排序的一致性检验。

function [RI,CI]=sglsortexamine(maxeigval,A)
%层次分析法单排序一致性检验
%maxeigval为最大特征值，A为判断矩阵
n=size(A,1);
RIT= [0,0,,,,,,,,,];
RI=RIT(n);
CI=(maxeigval-n)/(n-1);
CR=CI/RI;
if CR>=
disp([input('矩阵没通过一致性检验，请重新调整判断矩阵')]
else
disp([input('矩阵通过一致性检验')]);
end
在MATLAB中键入如下指令：
[RIA,CIA]= sglsortexamine(max(1),A);
[RIP1,CIP1]= sglsortexamine(max(2),P1);
[RIP2,CIP2]= sglsortexamine(max(3),P2);
[RIP3,CIP3]= sglsortexamine(max(4),P3);
运行结果如下：
矩阵通过一致性检验
矩阵通过一致性检验
矩阵通过一致性检验
矩阵通过一致性检验
一致性检验及总排序的MATLAB实现
通过层次单排序（权重）计算后，进行层次合成计算，在此本文定义函数计算层次总排序的权重并进行一致性检验。

function tw=tolsortvec(utw,dw,CIC,RIC)
% 求层次总排序权重并进行一致性检验
% utw为上一层因素的总排序权重行向量
% dw为下一层因素相对于上一层各因素的层次单排序权重矩阵
% CIC为一致性指标列向量
% RIC为随机一致性指标列向量
tw=dw*utw
CR=utw'*CIC/(utw'*RIC);
if CR>=
disp([input('层次总排序没通过一致性检验，请重新调整判断矩阵')]);
else
disp([input('层次总排序通过一致性检验')]);
end
在MATLAB中输入如下指令：
dw=zeros(7,3);
dw=(1:2,1)=wP1; dw=(3:5,2)=Wp2; dw=(6:7,3)=wP3;
CIC=[CIP1;CIP2;CIP3];
RIC=[RIP1;RIP2;RIP3];
tw= tolsortvec(wA,dw,CIC,RIC)’;
运行结果如下：
tw=
层次总排序通过一致性检验
其中tw 是层次总排序结果。

因此，根据数据建立如下的层次总排序表。

表3-3 层次总排序表（权重）
注：按概率乘法，P层次总排序指标的权重值为N—P层次指标的权重
与相应上一层次指标A—N层权重的积，且总排序权重值的和为1。

选择最优排序
计算出层次总排序后，为了使决策者能迅速得出结果，本文对层次总排序进行最优排序。

运用MATLAB键入如下指令：
n=length(tw);
for i=1:n
t=max(tw);
b(i)=t;
[m n]=find(a==t);
tw(n)=[];
end
b
运行结果如下：
b=
利用MATLAB大大缩短了计算复杂矩阵的时间，为决策者节省了宝贵的时间，从而有更多的精力投入其他事务。

3．基于MATLAB的AHP应用
3．1 挑选合适工作问题
某毕业生选择工作，经双方恳谈，假设已有三个单位C1,C2,C3表示愿意录用他。

该生对三个单位进行了解后，选取了一些中间指标进行考察，例如单位的研究课题，发展前途，待遇，同事情况，地理位置，单位名气等。

根据层次分析法，试求该生工作优先排序（给出权值、计算程序），并给出最终选择决策。

现以A、B、C表示选择工作的三个层次，建立如下结构模型：
图3 选择单位层次结构图
根据成对比较法，得到相应判断矩阵如下表：
表4-1 A-B判断矩阵
A B1B2B3B4B5B6
B1111411/2
B2112411/2
表4-2 B1~C判断矩阵
表4-3 B2~C判断矩阵
表4-4 B3~C判断矩阵
表4-5 B4~C判断矩阵
表4-6 B5~C判断矩阵
表4-7 B6~C判断矩阵
现在在MATLAB中分别用直接输入程序法和M—文件方法求解。

1）、直接输入代码法：
在MATLAB中输入如下程序：
A = [ 1,1,1,4,1,1/2;
1,1,2,4,1,1/2;
1,1/2,1,5,3,1/2;
1/4,1/4,1/5,1,1/3,1/3;
1,1,1/3,3,1,1;
2,2,2,3,3,1];
B1 = [1,1/4,1/2;4,1,3;2,1/3,1];
B2 = [1,1/4,1/5;4,1,1/2;5,2,1];
B3 = [1,3,1/3;1/3,1,1/7;3,7,1];
B4 = [1,1/3,5;3,1,7;1/5,1/7,1];
B5 = [1,1,7;1,1,7;1/7,1/7,1];
B6 = [1,7,9;1/7,1,1;1/9,1,1];
BS = [B1,B2,B3,B4,B5,B6];
m = length(B1);n = length(A); %随机一致性指标RI
RI = [0,0,,,,,,,,,];
[Wa,LA] = eig (A) %求A的特征向量WA和特征根LA Maxn=input('please input largest eigenvalue:'); %输入最大特征根
CIn = (Maxn- n) / (n - 1);
CRn = CIn / RI(n); %A的一致性比率CRn
WA=Wa(:,1)/sum(Wa(:,1)); %特征向量归一化
if CRn <
fprintf('A 的CR %f 通过一致性检验!\n',CRn); %控制文本格式
else
fprintf('A 的CR %f 未通过一致性检验!\n',CRn);
end
for k = 1:n %求B的特征向量WK和特征根LK [WB,LK] = eig( BS(1:3,(k-1)*m+1:(k-1)*m+3) )
Max(k)=input('please input largest eigenvalue:');
CIm(k) = (Max(k)- m) / (m - 1);
RIm(k) = RI(m);
CRm(k) = CIm(k) / RIm(k);; %B的一致性比率CRm
WK(:,k)= WB(:,1)/sum(WB(:,1));
end
for k = 1:n
if CRm(k)<
fprintf('B%d的CR %f 通过一致性检验!\n',k,CRm(1,k)); %控制文本格式else
fprintf('B%d的CR %f 未通过一致性检验!\n',k,CRm(1,k));
end
end
disp('准则层对目标层权向量'); disp(WA);
disp('方案层对准则层权向量'); disp(WK);
E = WK * WA
disp('方案层组合权向量');disp(E);
CI = CIm * WA;
RI = RIm * WA;
CR = CI / RI; %组合一致性比率CR if CR <
fprintf('组合一致性比率CR %f 通过一致性检验!\n',CRn);
else
fprintf('组合一致性比率CR %f 未通过一致性检验!\n',CRn);
end
[MAX,CHOICE] = max(E); %最佳选择
CHOICE
MATLAB运行结果如下：
Wa =
- + + -
- + - +
- +
+ - - +
+ - + -
- +
LA =
0 0 0 0 0
0 + 0 0 0 0
0 0 - 0 0 0
0 0 0 0 0
0 0 0 0 + 0
0 0 0 0 0 -
please input largest eigenvalue:
A 的CR 通过一致性检验!
WB =
+ -
- +
LK =
0 0
0 + 0
0 0 -
please input largest eigenvalue:
WB =
+ -
- +
LK =
0 0
0 + 0
0 0 - please input largest eigenvalue:
WB =
+ -
- +
LK =
0 0
0 + 0
0 0 -
please input largest eigenvalue:
WB =
+ -
- +
LK =
0 0
0 + 0
0 0 - please input largest eigenvalue:
WB =
LK =
0 0
0 0
0 0
please input largest eigenvalue:
WB =
- +
+ -
LK =
0 0
0 + 0
0 0 - please input largest eigenvalue:
B1的CR 通过一致性检验!
B2的CR 通过一致性检验!
B3的CR 通过一致性检验!
B4的CR 通过一致性检验!
B5的CR 通过一致性检验!
B6的CR 通过一致性检验!
准则层对目标层权向量
方案层对准则层权向量
E =
方案层组合权向量
组合一致性比率CR 通过一致性检验!
CHOICE =
1
2）、M—文件法：
根据第三章节的层次分析法各环节的MATLAB实现定义如下三个函数：
（1）定义来实现最大特征根及对应的归一化特征向量。

function[maxeigval,w]=maxeigvalvec(A) %求最大特征值及对应的归一化特征向量
%A为判断矩阵
[eigvec,eigval]=eig(A);
eigval=diag(eigval); %特征向量
eigvalmag=imag(eigval);
realind=find(eigvalmag<eps);
realeigval=eigval(realind); %实特征根
maxeigval=max(realeigval) %最大特征值
index=find(eigval==maxeigval);
vecinit=eigvec(:,index); %最大特征值对应的特征向量
w=vecinit./sum(vecinit) %特征向量归一化
（2）定义函数来实现层次单排序的一致性检验。

function [RI,CI]=sglsortexamine(maxeigval,A)
%层次分析法单排序一致性检验
%maxeigval为最大特征值，A为判断矩阵
n=size(A,1);
RIT= [0,0,,,,,,,,,];
RI=RIT(n);
CI=(maxeigval-n)/(n-1);
CR=CI/RI;
if CR>=
disp([input('矩阵没通过一致性检验，请重新调整判断矩阵')]
else
disp([input('矩阵通过一致性检验')]);
end
（3）定义函数计算层次总排序的权重并进行一致性检验。

function tw=tolsortvec(utw,dw,CIC,RIC)
% 求层次总排序权重并进行一致性检验
% utw为上一层因素的总排序权重行向量
% dw为下一层因素相对于上一层各因素的层次单排序权重矩阵
% CIC为一致性指标列向量
% RIC为随机一致性指标列向量
tw=dw*utw;
CR=utw'*CIC/(utw'*RIC);
if CR>=
disp([input('层次总排序没通过一致性检验，请重新调整判断矩阵')]); else
disp([input('层次总排序通过一致性检验')]); end
%主程序
clear;
A = [ 1,1,1,4,1,1/2;
1,1,2,4,1,1/2;
1,1/2,1,5,3,1/2;
1/4,1/4,1/5,1,1/3,1/3;
1,1,1/3,3,1,1;
2,2,2,3,3,1];
B1 = [1,1/4,1/2;4,1,3;2,1/3,1];
B2 = [1,1/4,1/5;4,1,1/2;5,2,1];
B3 = [1,3,1/3;1/3,1,1/7;3,7,1];
B4 = [1,1/3,5;3,1,7;1/5,1/7,1];
B5 = [1,1,7;1,1,7;1/7,1/7,1];
B6 = [1,7,9;1/7,1,1;1/9,1,1];
[max(1),wA]=maxeigvalvec(A);
[max(2),wB1]=maxeigvalvec(B1);
[max(3),wB2]=maxeigvalvec(B2);
[max(4),wB4]=maxeigvalvec(B4);
[max(5),wB3]=maxeigvalvec(B3);
[max(6),wB5]=maxeigvalvec(B5);
[max(7),wB6]=maxeigvalvec(B6);
[RIA,CIA]=sglsortexamine(max(1),A);
[RIB1,CIB1]=sglsortexamine(max(2),B1);
[RIB2,CIB2]=sglsortexamine(max(3),B2);
[RIB3,CIB3]=sglsortexamine(max(4),B3);
[RIB4,CIB4]=sglsortexamine(max(5),B4);
[RIB5,CIB5]=sglsortexamine(max(6),B5);
[RIB6,CIB6]=sglsortexamine(max(7),B6);
dw=zeros(18,6);
dw(1:3,1)=wB1;dw(4:6,2)=wB2;dw(7:9,3)=wB3;dw(10:12,4)=wB4; dw(13:15,5)=wB5;dw(16:18,6)=wB6;
CIC=[CIB1;CIB2;CIB3;CIB4;CIB5;CIB6];
RIC=[RIB1;RIB2;RIB3;RIB4;RIB5;RIB6];
tw=tolsortvec(wA,dw,CIC,RIC)';
n=length(A);
[MAX,CHOICE] = max(n); %最佳选择
CHOICE
MATLAB运行结果如下：
maxeigval =
w =
maxeigval =
w = maxeigval = w = maxeigval = w = maxeigval = w =
maxeigval =
w =
maxeigval =
w =
矩阵通过一致性检验
矩阵通过一致性检验
矩阵通过一致性检验
矩阵通过一致性检验
矩阵通过一致性检验
矩阵通过一致性检验
矩阵通过一致性检验
层次总排序通过一致性检验CHOICE =
1
由此分析得，通过层次分析法并且结合MATLAB，该生很容易便能作出了决策，选择工作1。

而且容易得出如下表的层次总排序：
表4-7 层次总排序
从而得出结论：选择工作1.
4．结束语
利用层次分析法AHP，并结合被称为第四代计算机语言的MATLAB工具软件，来解决多目标、多准则的问题，大大缩减了决策者计算复杂矩阵的时间。

目前很多领域都采用这种方法，如经济计划和管理、能源政策和分配、行为科学等。

由此可见，层次分析法AHP结合MATLAB时很有发展前景的。

但是层次分析法也有其局限性：（1）AHP方法中所用的指标体系需要有专家系统的支持，如果给出的指标不合理，得到的结果也就不准确；（2）AHP方法中经行多层比较的时候需要给出一致性比较，如果不满足一致性指标要求，则AHP方法就失去了作用。

因此，需要更多的研究人员去探讨和解决。