用Matlab实现AHP的算法

.---可修编.1．MATLAB的基本容MATLAB（MATrixLABoratory，矩阵实验室的缩写）是一种特殊用途的计算机程序优化执行工程和科学计算。它开始为旨在执行矩阵数学程式的生活，但多年来它已发展成为一个灵活的计算系统基本上能够解决任何技术问题。MATLAB具有编程语言的基本特征，使用MATLAB也可以使用像BASIC、FORTRAN、C等传统编程语言一样，进行程序 设计 领导形象设计圆作业设计ao工艺污水处理厂设计附属工程施工组织设计清扫机器人结构设计 ，而且简单易学、编程效率高。正因为MATLAB的强大的功能，使得它在许多领域得到广泛应用。在科研与工程应用领域，MATLAB已被广泛地用于科学研究和解决各种具体的实际问题。许多科技工作者选用MATLAB做为计算工具，避免了繁琐的底层编程，从而可以把主要精力和时间花在科学研究和解决实际问题是上，提高了工作效率。1.1MATLAB矩阵矩阵是MATLAB的基本处理对象，因此根据本文所需，简单介绍所涉及MATLAB矩阵容。1.1.1MATLAB矩阵的建立1、直接输入法最简单的建立矩阵的方法是从键盘直接输入矩阵的元素。例如：A=[123;456;789]A=123456789也可以用回车键代替分号，按下列方式输入：A=[123456789]2、利用M文件建立矩阵比较大且复杂的矩阵，可以为它专门建立一个M文件，如同下例。利用M文件建立矩阵。启动有关编辑程序或MATLAB文本编辑器，并输入待建矩阵：MYMAT=[111,112,113,114,115,116,117,118,119;211,212,213,214,215,216,217,218,219];把输入的容以纯文本方式存盘（设文件名为mymatrix.m）。在MATLAB命令窗口中输入mymatrix，即运行该M文件，就会自动建立一个名为MYMAT的矩阵，可供以后使用。1.1.2矩阵的特征值与特征向量特征值和特征向量在科学研究和工程计算中都有非常广泛地应用。在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,0.5;1,1,0.25;0.5,0.25,2];[V,D]=eig(A)V=0.72120.44430.5315-0.68630.56210.4615-0.0937-0.69760.7103D=-0.01660001.48010002.5365求得的3个特征值是-0.0166、1.4801和2.5365，各特征值对应的特征向量为V的各列构成的向量。1.2MATLAB的M文件用MATLAB语言编写的程序称为M文件。M文件是由若干MATLAB命令组成在一起构成的，它可以完成某些操作，也可以实现某种算法。M文件可以根据调用方式的不同分为两类：命令文件（ScriptFile）和函数文件（FunctionFile）。它们的扩展名均为.m。函数文件由function语句引导，其基本结构为：function输出形参 表 关于同志近三年现实表现材料材料类招标技术评分表图表与交易pdf视力表打印pdf用图表说话 pdf =函数名（输入形参表）注释说明部分函数体语句我们通过举例说明如下：例2-2分别建立命令文件和函数文件，将求矩阵的一致性指标CI：CI=(λmax-n)/(n-1)程序1建立命令文件并以文件名CI.m存盘：max=input('pleaseinputmax:');n=input('pleaseinputn:');CI=(max-n)/(n-1)然后在MATLAB的命令窗口中输入CI即可。程序2建立函数文件CI.m。functionc=CI(max,n)c=(max-n)/(n-1)然后在MATLAB的命令窗口调用该函数文件。max=input('pleaseinputmax:');n=input('pleaseinputn:');c=CI(max,n)2．基于MATLAB的AHP实现2．1AHP的MATLAB的计算 流程 快递问题件怎么处理流程河南自建厂房流程下载关于规范招聘需求审批流程制作流程表下载邮件下载流程设计 框图根据层次 分析 定性数据统计分析pdf销售业绩分析模板建筑结构震害分析销售进度分析表京东商城竞争战略分析 法的一般步骤我们得到在MATLAB工具上实现的计算程序流程框图，如图2所示[16]输入准则层层数输入准则层对目标层判断矩阵元素判断矩阵归一化计算CR<0.1输入准则层第j个准则包含的方案数量n判断矩阵归一化计算输入方案层第j个准则层的判断矩阵元素bijCR<0.1另m=1;n=1输入方案层第m个准则层的相关矩阵元素计算权值矩阵权值的计算总CR计算开始结束CR<0.1图2以MATLAB实现的层次分析法的计算流程框图通过流程框图，层次分析的基本步骤如下：第一步：准则层对目标层的判断矩阵归一化且判断是否满足一致性；第二步：第一步满足时，将方案层对准则层的判断矩阵归一化并判断其一致性；第三步：当第一、二步满足时，求方案层的总排序权值与总CR并判断一致性。2．2平均随机一致性指标的MATLAB实现运用层次分析法决策者需要通过反复地解决决策问题，将同一层次的各元素与上一层次中某一准则的重要性进行比较，从而构造出两两判断比较矩阵A=(aij)nn(称为成对比较矩阵)。前面已经描述了九级标度法，此处运用其他描述，则这些成对比较矩阵应满足如下条件：(l)>0(2)·=l(3)=l按照事物逻辑要求，该矩阵还应具备一致性，即满足：·=前面已经给出由于客观事物的复杂性与决策者的认识的多样性，实际问题的成对比较矩阵不可能做到严格上的一致性，因而，借助平均随机一致性指标RI来相对判定其一致性程度。其中表1-4是已经计算好的1~15阶矩阵的RI值表，但未给出其实现过程，且各文献的RI值表不完全相同。究其原因除没有太大必要介绍外，真正去实现它却有如下三个难度：(1)随机两两判断矩阵中的元素要1~9和它们对应的倒数共17个整数与小数的均匀分布很难处理。(2)一般高级编程语言实现成对比较矩阵及相关计算，非常复杂，且占用存巨大，耗时多。(3)随机种子源不能控制。本文使用数学软件包MATLAB对其进行计算。其设计解决思路为：先用软件包中随机函数产生数1~17的均匀分布的n阶矩阵，然后在软件包中采用不同技巧将它转化为成对比较矩阵，最后用循环语句计算出RI值。结果如下表：表3-1计算的RI值阶数123456789RI000.54400.89801.13131.25151.34951.41901.4542以n=4为例，过程详见如下程序清单，其中随机成对比较矩阵的实现见相应注释部分。[15]MATLAB的程序M文件：functionri%计算RI值的命令文件n=4;ri=0;m=100;rand('seed',21)%控制随机发生器fori=1:ma=ceil(17*rand(n));%产生n阶l~17的随机阵a(find(a=8))=8.1;%消除0为分母b=1./(a-8);%产生一个辅助阵a(find(a>9))=b(find(a>9));%借助b，将9~17分别转化为~a(find(a=8.1))=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;%实现随机成对比较阵ak=size(a,1);%计算a的行维数ri=ri+(max(abs(eig(a)))-k)/(k-1);%计算100次RI值endri/m%计算平均RI值2．3AHP各环节的MATLAB实现以目标矩阵A=,准则层矩阵为P1=,P2=，P3=为例，运用MATLAB进行数据处理。2.3.1特征向量及其归一化的MATLAB实现MATLAB中求矩阵特征值和特征向量的函数是eig，其调用的格式为[V,D]=eig(A),其中，V为特征向量矩阵，D为特征值矩阵。层次分析法中需要求得是最大特征值及对应的归一化特征向量，而且考虑到eig函数在求得的特征值中可能会存在复数。因此，运用直接输入程序代码会产生一定的误差。在此需要对求得的V、D进行适当选择，定义一个M-filemaxeigvalvec.m来实现。function[maxeigval,w]=maxeigvalvec(A)%求最大特征值及对应的归一化特征向量%A为判断矩阵[eigvec,eigval]=eig(A);eigval=diag(eigval);%特征向量eigvalmag=imag(eigval);realind=find(eigvalmag=0.10disp([input('矩阵没通过一致性检验，请重新调整判断矩阵')]elsedisp([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);运行结果如下：矩阵通过一致性检验矩阵通过一致性检验矩阵通过一致性检验矩阵通过一致性检验2.3.3一致性检验及总排序的MATLAB实现通过层次单排序（权重）计算后，进行层次合成计算，在此本文定义tolsortvec.m函数计算层次总排序的权重并进行一致性检验。functiontw=tolsortvec(utw,dw,CIC,RIC)%求层次总排序权重并进行一致性检验%utw为上一层因素的总排序权重行向量%dw为下一层因素相对于上一层各因素的层次单排序权重矩阵%CIC为一致性指标列向量%RIC为随机一致性指标列向量tw=dw*utwCR=utw'*CIC/(utw'*RIC);ifCR>=0.10disp([input('层次总排序没通过一致性检验，请重新调整判断矩阵')]);elsedisp([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=0.42470.21230.06670.16460.02700.06980.0349层次总排序通过一致性检验其中tw是层次总排序结果。因此，根据数据建立如下的层次总排序表。表3-3层次总排序表（权重）APA层次P的总排序结果注0.63700.25830.1047P10.6667000.42470.3333000.2123P200.258300.066700.637000.164600.104700.0270P3000.66670.0698000.33330.0349注：按概率乘法，P层次总排序指标的权重值为N—P层次指标的权重与相应上一层次指标A—N层权重的积，且总排序权重值的和为1。2.3.4选择最优排序计算出层次总排序后，为了使决策者能迅速得出结果，本文对层次总排序进行最优排序。运用MATLAB键入如下指令：n=length(tw);fori=1:nt=max(tw);b(i)=t;[mn]=find(a==t);tw(n)=[];endb运行结果如下：b=0.42470.21230.16460.06980.06670.03490.0276利用MATLAB大大缩短了计算复杂矩阵的时间，为决策者节省了宝贵的时间，从而有更多的精力投入其他事务。3．基于MATLAB的AHP应用3．1挑选合适工作问题某毕业生选择工作，经双方恳谈，假设已有三个单位C1,C2,C3表示愿意录用他。该生对三个单位进行了解后，选取了一些中间指标进行考察，例如单位的研究课题，发展前途，待遇，同事情况，地理位置，单位名气等。根据层次分析法，试求该生工作优先排序（给出权值、计算程序），并给出最终选择决策。现以A、B、C表示选择工作的三个层次，建立如下结构模型：图3选择单位层次结构图根据成对比较法，得到相应判断矩阵如下表：表4-1A-B判断矩阵B1B2B3B4B5B6B1111411/2B2112411/2B311/21531/2B41/41/41/511/31/3B5111/3311B6222331表4-2B1~C判断矩阵11/41/241321/31表4-3B2~C判断矩阵C1C2C3C111/41/5C2411/2C3521表4-4B3~C判断矩阵C1C2C3C1131/3C21/311/7C3371表4-5B4~C判断矩阵C1C2C3C111/35C2317C31/51/71表4-6B5~C判断矩阵C1C2C3C1117C2117C31/71/71表4-7B6~C判断矩阵C1C2C3C1179C21/711C31/911现在在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);%随机一致性指标RIRI=[0,0,0.58,0.90,1.12,1.24,1.32,1.41,1.45,1.49,1.51];[Wa,LA]=eig(A)%求A的特征向量WA和特征根LAMaxn=input('pleaseinputlargesteigenvalue:');%输入最大特征根CIn=(Maxn-n)/(n-1);CRn=CIn/RI(n);%A的一致性比率CRnWA=Wa(:,1)/sum(Wa(:,1));%特征向量归一化ifCRn<0.10fprintf('A的CR%f通过一致性检验!\n',CRn);%控制文本格式elsefprintf('A的CR%f未通过一致性检验!\n',CRn);endfork=1:n%求B的特征向量WK和特征根LK[WB,LK]=eig(BS(1:3,(k-1)*m+1:(k-1)*m+3))Max(k)=input('pleaseinputlargesteigenvalue:');CIm(k)=(Max(k)-m)/(m-1);RIm(k)=RI(m);CRm(k)=CIm(k)/RIm(k);;%B的一致性比率CRmWK(:,k)=WB(:,1)/sum(WB(:,1));endfork=1:nifCRm(k)<0.10fprintf('B%d的CR%f通过一致性检验!\n',k,CRm(1,k));%控制文本格式elsefprintf('B%d的CR%f未通过一致性检验!\n',k,CRm(1,k));endenddisp('准则层对目标层权向量');disp(WA);disp('方案层对准则层权向量');disp(WK);E=WK*WAdisp('方案层组合权向量');disp(E);CI=CIm*WA;RI=RIm*WA;CR=CI/RI;%组合一致性比率CRifCR<0.10fprintf('组合一致性比率CR%f通过一致性检验!\n',CRn);elsefprintf('组合一致性比率CR%f未通过一致性检验!\n',CRn);end[MAX,CHOICE]=max(E);%最佳选择CHOICEMATLAB运行结果如下：Wa=-0.3396-0.1255-0.0563i-0.1255+0.0563i0.7354-0.1896+0.3838i-0.1896-0.3838i-0.4038-0.1884-0.5736i-0.1884+0.5736i-0.6464-0.4476-0.2693i-0.4476+0.2693i-0.42490.67240.67240.08340.3884-0.0605i0.3884+0.0605i-0.1063-0.0138+0.0429i-0.0138-0.0429i-0.0405-0.0592-0.0922i-0.0592+0.0922i-0.3298-0.1384+0.3417i-0.1384-0.3417i-0.12170.0035+0.1640i0.0035-0.1640i-0.6488-0.1467-0.0710i-0.1467+0.0710i0.13370.59200.5920LA=6.6178000000-0.1557+1.2808i000000-0.1557-1.2808i000000-0.0603000000-0.1230+0.5461i000000-0.1230-0.5461ipleaseinputlargesteigenvalue:6.6178A的CR0.099645通过一致性检验!WB=0.19990.1000+0.1731i0.1000-0.1731i0.9154-0.9154-0.91540.34930.1747-0.3025i0.1747+0.3025iLK=3.0183000-0.0091+0.2348i000-0.0091-0.2348ipleaseinputlargesteigenvalue:0.1999WB=0.14600.0730+0.1265i0.0730-0.1265i0.49940.2497-0.4325i0.2497+0.4325i0.8540-0.8540-0.8540LK=3.0246000-0.0123+0.2725i000-0.0123-0.2725ipleaseinputlargesteigenvalue:3.0246WB=0.3382-0.1691+0.2929i-0.1691-0.2929i0.1226-0.0613-0.1062i-0.0613+0.1062i0.93310.93310.9331LK=3.0070000-0.0035+0.1453i000-0.0035-0.1453ipleaseinputlargesteigenvalue:3.0070WB=0.3928-0.1964+0.3402i-0.1964-0.3402i0.91400.91400.91400.1013-0.0506-0.0877i-0.0506+0.0877iLK=3.0649000-0.0324+0.4448i000-0.0324-0.4448ipleaseinputlargesteigenvalue:3.0649WB=-0.7035-0.33920.6619-0.7035-0.9233-0.7495-0.10050.18040.0125LK=3.00000000.0000000-0.0000pleaseinputlargesteigenvalue:3.0000WB=-0.98440.98440.9844-0.1293-0.0647-0.1120i-0.0647+0.1120i-0.1189-0.0595+0.1030i-0.0595-0.1030iLK=3.0070000-0.0035+0.1453i000-0.0035-0.1453ipleaseinputlargesteigenvalue:3.0070B1的CR-2.413879通过一致性检验!B2的CR0.021207通过一致性检验!B3的CR0.006034通过一致性检验!B4的CR0.055948通过一致性检验!B5的CR0.000000通过一致性检验!B6的CR0.006034通过一致性检验!准则层对目标层权向量0.15070.17920.18860.04720.14640.2879方案层对准则层权向量0.13650.09740.24260.27900.46670.79860.62500.33310.08790.64910.46670.10490.23850.56950.66940.07190.06670.0965E=0.39520.29960.3052方案层组合权向量0.39520.29960.3052组合一致性比率CR0.099645通过一致性检验!CHOICE=12）、M—文件法：根据第三章节的层次分析法各环节的MATLAB实现定义如下三个函数：定义maxeigvalvec.m来实现最大特征根及对应的归一化特征向量。function[maxeigval,w]=maxeigvalvec(A)%求最大特征值及对应的归一化特征向量%A为判断矩阵[eigvec,eigval]=eig(A);eigval=diag(eigval);%特征向量eigvalmag=imag(eigval);realind=find(eigvalmag=0.10disp([input('矩阵没通过一致性检验，请重新调整判断矩阵')]elsedisp([input('矩阵通过一致性检验')]);end（3）定义tolsortvec.m函数计算层次总排序的权重并进行一致性检验。functiontw=tolsortvec(utw,dw,CIC,RIC)%求层次总排序权重并进行一致性检验%utw为上一层因素的总排序权重行向量%dw为下一层因素相对于上一层各因素的层次单排序权重矩阵%CIC为一致性指标列向量%RIC为随机一致性指标列向量tw=dw*utw;CR=utw'*CIC/(utw'*RIC);ifCR>=0.10disp([input('层次总排序没通过一致性检验，请重新调整判断矩阵')]);elsedisp([input('层次总排序通过一致性检验')]);end%main.m主程序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);%最佳选择CHOICEMATLAB运行结果如下：maxeigval=6.6178w=0.15070.17920.18860.04720.14640.2879maxeigval=3.0183w=0.13650.62500.2385maxeigval=3.0246w=0.09740.33310.5695maxeigval=3.0649w=0.27900.64910.0719maxeigval=3.0070w=0.24260.08790.6694maxeigval=3.0000w=0.46670.46670.0667maxeigval=3.0070w=0.79860.10490.0965矩阵通过一致性检验矩阵通过一致性检验矩阵通过一致性检验矩阵通过一致性检验矩阵通过一致性检验矩阵通过一致性检验矩阵通过一致性检验层次总排序通过一致性检验CHOICE=1由此分析得，通过层次分析法并且结合MATLAB，该生很容易便能作出了决策，选择工作1。而且容易得出如下表的层次总排序：表4-7层次总排序BCB1B2B3B4B5B6层次总排序0.15070.17920.18860.04720.14640.2879C10.13650.09740.24260.27900.46670.79860.3952C20.62500.33310.08790.64910.46670.10490.2996C30.23850.56950.66940.07190.06670.09650.3052从而得出结论：选择工作1.4．结束语利用层次分析法AHP，并结合被称为第四代计算机语言的MATLAB工具软件，来解决多目标、多准则的问题，大大缩减了决策者计算复杂矩阵的时间。目前很多领域都采用这种方法，如经济计划和管理、能源政策和分配、行为科学等。由此可见，层次分析法AHP结合MATLAB时很有发展前景的。但是层次分析法也有其局限性：（1）AHP方法中所用的指标体系需要有专家系统的支持，如果给出的指标不合理，得到的结果也就不准确；（2）AHP方法中经行多层比较的时候需要给出一致性比较，如果不满足一致性指标要求，则AHP方法就失去了作用。因此，需要更多的研究人员去探讨和解决。