基于R语言的综合评价（层次分析法AHP为例）

一.综合评价的基本概念

对一个事物的评价往往会涉及多个因素或者多个指标，评价是在多个因素相互作用下的一个综合判断。多指标综合评价方法具有以下的特点：包含若干个指标，分别说明被评价对象的不同方面，评价方法最终要对被评价对象作出一个整体性的评判，用一个总指标来说明被评价对象的一般水平。

那么评价指标的权重如何确定呢？有几种方法：德尔菲法确定权重，层次分析法确定权重。

德尔菲一般程序是先由各位专家单独对各个评价指标的重要程度进行评定，然后由综合评价人员对各位专家的评定结果进行综合，计算出平均数然后反馈给专家，如此反复进行几次，使各位专家的意见趋于一致，然后再确定出各评价指标的权数。

层次分析法(Analytic Hierarchy Process,AHP法）,核心问题是：权重构造。建立评价对象的综合评价指标体系，通过指标之间的两两比较确定出相对重要的程度，然后通过特征值法，最小二乘法，对数最小二乘法，上三角元素法等客观运算来确定各指标权数。

特征值法的运算步骤 

1.构造判断矩阵。通过指标之间的两两重要程度进行比较和分析判断，构造判断矩阵。层次分析法在对指标的相对重要程度进行测量的时候，引入九分位的相对重要的比例标度。即令A为判断矩阵，表示同一层次各个指标相对重要程度的判断值，它由若干位专家判定。则有A=(aij)m*m,矩阵A中的各元素aij表示横行指标Zj相对重要程度的两两比较值。

2.指标权数计算。层次分析法的信息判断矩阵是利用排序原理，求得各行的几何平均数，然后计算各评价指标的重要性权数，计算公式分别为：

ai=m√ai1*ai2*ai3*...aim=m√∏j=1aij

wi=ai/∑ai,i=1,2...,ai为几何平均数。

将各个评价指标的重要性权数用一个向量来表示，即W=(w1,w2,....wm)，该向量又称为判断矩阵的特征向量。

3.判断矩阵的一次性检验

与其他确定指标权重系数的方法相比，层次分析法的最大优点是通过一次性检验，保持和专家思想逻辑的一次性。

（1）计算判断矩阵的最大特征根：

       max=1/m∑(AW)i/Wi

（2）计算判断矩,阵一致性的指标

       CI=max-m/m-1

（3）计算判断矩阵一致性的比率

      CR =CI/RI≤0.1

综合评价方法较多，如综合评分法，综合指数法，秩和比法，层次分析法，TOPSIS法，模糊综合评价法，数据包络分析法等几种具有代表性的评价方法。我们今天主要说明两种方法：综合评分法以及层次分析法

1.综合评分法可以把各个指标得分直接相加得到一个最终总分，最后根据这个最终总分的高低来判定评价对象的优劣。这种方法的好处是对各个指标赋予同样的权重来同等对待，省去了确定指标权重的复杂步骤，但这同时也是一个缺点：不能区分各个指标的相对重要程度。因此经常把这个指标的相对重要程度不同给予不同的权重，然后用各个指标的得分乘以权重各个指标对各个不同方案的加权评分，每个方案的各指标加权评分除以权重所得到的商就是加权平均法，得分最多的方案就是最佳方案。

Si=∑WjZij=∑1/mZij=Zi

2.层次分析法

层次分析法的应用首先是确定各层次的目标体系。所谓分层，就是根据研究目标之间的内在联系和因果关系，逐步分解为层次的目标体系。

∑w i=1，∑wij=1

Zij为i个目标第j方案属性值向量的分量，利用公式：

Si=∑wijZij=Zi*Wi

对目标低层到高层进行计算，由此判断各个方案的优劣，Si值越大，这个方案越优。

 继续来看案例吧

区域自主创新能力的层次分析

数据：

构造函数msa.AHP

msa.AHP<-function(B){

 A=matrix(B,nrow=sqrt(length(B)),ncol=sqrt(length(B)),byrow=TRUE)

  print(A)

  m=ncol(A)

  ai=apply(A,1,prod)^(1/m)

  W=ai/sum(ai); 

  if(m>2){

    AW=A%*%W

    L_max=sum(AW/W)/m; 

    CI=(L_max-m)/(m-1); 

    RI=c(0,0,0.58,0.90,1.12,1.24,1.32,1.41,1.45,1.49,1.51)

    CR=CI/RI[m]

    cat('\n    L_max=',L_max,'\n')

    cat('    CI=',CI,'\n')

    cat('    CR=',CR,'\n')

    if(CR<=0.1) cat('  Consistency test is OK!\n\n')

    else cat('    Please adjust the judgment matrix!\n')

  }

  return(W)

}

A<-c(1,3,3,1,1/3,1,2,1/3,1/3,1/2,1,1/3,1,3,3,1)

judgeMatix<-matrix(A,ncol=4)

judgeMatix

权重计算

(Awet=msa.AHP(A))

B1=c(1,1,1,2,1,1,1,1,1,1,1/2,2,1,2,1,1,1,1/2,1,2,2,1/2,2,2,1,1,1,1,1,1/2,1,1,1,2,2,1,1,1/2,1,1/2,1,2,1,1/2,1/2,1,1/2,1/2,1)

(B1wet=msa.AHP(B1))

B2=c(1,1,1,2,1,1,1,1,1,1/2,2,1,1,1,1,1,1,2,1/2,2,1,1,1,1,1,1/2,1,1,1,2,1,1,1/2,1,1/2,1)

(B2wet=msa.AHP(B2))

B3=c(1,1,2,2,1,1,1,1,2,2,1/2,1,1,1/2,1/2,1/2,1/2,2,1,1/3,1,1/2,2,3,1)

(B3wet=msa.AHP(B3))

B4=c(1,1/2,1,1,1,1/2,2,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,2,1,1,1,1,1)

(B4wet=msa.AHP(B4))

结果如下：

> (Awet=msa.AHP(A))

          [,1] [,2] [,3]      [,4]

[1,] 1.0000000  3.0    3 1.0000000

[2,] 0.3333333  1.0    2 0.3333333

[3,] 0.3333333  0.5    1 0.3333333

[4,] 1.0000000  3.0    3 1.0000000

    L_max= 4.060434 

    CI= 0.02014454 

    CR= 0.02238282 

  Consistency test is OK!

[1] 0.3735942 0.1480936 0.1047180 0.3735942

> B1=c(1,1,1,2,1,1,1,1,1,1,1/2,2,1,2,1,1,1,1/2,1,2,2,1/2,2,2,1,1,1,1,1,1/2,1,1,1,2,2,1,1,1/2,1,1/2,1,2,1,1/2,1/2,1,1/2,1/2,1)

> (B1wet=msa.AHP(B1))

     [,1] [,2] [,3] [,4] [,5] [,6] [,7]

[1,]  1.0  1.0  1.0  2.0  1.0  1.0    1

[2,]  1.0  1.0  1.0  0.5  2.0  1.0    2

[3,]  1.0  1.0  1.0  0.5  1.0  2.0    2

[4,]  0.5  2.0  2.0  1.0  1.0  1.0    1

[5,]  1.0  0.5  1.0  1.0  1.0  2.0    2

[6,]  1.0  1.0  0.5  1.0  0.5  1.0    2

[7,]  1.0  0.5  0.5  1.0  0.5  0.5    1

    L_max= 7.485107 

    CI= 0.08085116 

    CR= 0.06125088 

  Consistency test is OK!

[1] 0.15552481 0.15552481 0.15552481 0.15552481 0.15552481

[6] 0.12758250 0.09479344

> B2=c(1,1,1,2,1,1,1,1,1,1/2,2,1,1,1,1,1,1,2,1/2,2,1,1,1,1,1,1/2,1,1,1,2,1,1,1/2,1,1/2,1)

> (B2wet=msa.AHP(B2))

     [,1] [,2] [,3] [,4] [,5] [,6]

[1,]  1.0  1.0  1.0  2.0  1.0    1

[2,]  1.0  1.0  1.0  0.5  2.0    1

[3,]  1.0  1.0  1.0  1.0  1.0    2

[4,]  0.5  2.0  1.0  1.0  1.0    1

[5,]  1.0  0.5  1.0  1.0  1.0    2

[6,]  1.0  1.0  0.5  1.0  0.5    1

    L_max= 6.329513 

    CI= 0.06590259 

    CR= 0.05314725 

  Consistency test is OK!

[1] 0.1858804 0.1656006 0.1858804 0.1656006 0.1656006 0.1314373

> B3=c(1,1,2,2,1,1,1,1,2,2,1/2,1,1,1/2,1/2,1/2,1/2,2,1,1/3,1,1/2,2,3,1)

> (B3wet=msa.AHP(B3))

     [,1] [,2] [,3] [,4]      [,5]

[1,]  1.0  1.0    2  2.0 1.0000000

[2,]  1.0  1.0    1  2.0 2.0000000

[3,]  0.5  1.0    1  0.5 0.5000000

[4,]  0.5  0.5    2  1.0 0.3333333

[5,]  1.0  0.5    2  3.0 1.0000000

    L_max= 5.325191 

    CI= 0.08129764 

    CR= 0.07258718 

  Consistency test is OK!

[1] 0.2516547 0.2516547 0.1258273 0.1332793 0.2375840

> B4=c(1,1/2,1,1,1,1/2,2,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,2,1,1,1,1,1)

> (B4wet=msa.AHP(B4))

     [,1] [,2] [,3] [,4] [,5] [,6]

[1,]    1  0.5    1    1    1  0.5

[2,]    2  1.0    1    1    1  1.0

[3,]    1  1.0    1    1    1  1.0

[4,]    1  1.0    1    1    1  1.0

[5,]    1  1.0    1    1    1  1.0

[6,]    2  1.0    1    1    1  1.0

    L_max= 6.080612 

    CI= 0.01612233 

    CR= 0.01300188 

  Consistency test is OK!

[1] 0.1314373 0.1858804 0.1656006 0.1656006 0.1656006 0.1858804

> 

由于指标的量纲不同，要进行无量钢化处理，规定效用值的范围是【0，100】对于正向指标Vij=(Vij-Vimin)/(Vimax-Vimin)*60+40,对于逆向指标Yij=(Vimax-Vij)/(Vimax-Vimin)*60+40，即Vij为正效应指标时，该指标值就越大，其效用值就越高，反之，同理。

构造函数：

fx<-function(x,a=40,b=60)(x-min(x)/max(x)-min(x))*b+a

Z=apply(citydata,2,fx)

S1=Z[,1:7]%*%B1wet

S2=Z[,8:13]%*%B2wet

S3=Z[,14:18]%*%B3wet

S4=Z[,19:24]%*%B4wet

S=cbind(S1,S2,S3,S4)%*%Awet

R=rank(-S)

data.frame(S1=S1,R1=rank(-S1),S2=S2,R2=rank(-S2),S3=S3,R3=rank(-S3),S4=S4,R4=rank(-S4),S=S,R=R)

B=barplot(S[,1],las=3)

结果如下：

S1 R1          S2 R2        S3 R3          S4 R4

1  12006.9293  2 44348.81419  1 28750.243 20  27668.5915  5

2  11969.7242  3 20780.45047  5 11984.601 29  19058.5719  8

3   1362.0888 22  5776.14175 19 39950.960 11   6156.9605 16

4   2090.1745 12  7310.02735 15 34876.252 16   1877.9753 24

5   1919.5381 14 11565.43826  8 37022.827 14   2295.4278 23

6   4721.7101  8 11209.67278  9 40329.049 10  13020.4897 10

7   1733.0044 16  7915.30907 14 17741.915 27   5494.1001 17

8   1794.5976 15  6886.34306 16 23640.209 24   3097.8206 20

9  22400.2120  1 42914.52212  2 40574.996  9  55529.0981  3

10  8640.0104  5 12336.41442  7 59630.548  2 129402.0112  2

11  6451.9606  6 22478.31510  4 40643.978  8  26682.2842  6

12  1498.9268 19    65.27659 30 31751.161 17   4547.2364 18

13  3170.7849  9 16474.67622  6 19838.445 26  19582.0302  7

14  1443.4736 20  1995.85466 28 20601.901 25   7416.0146 13

15  4947.6700  7  6434.25594 17 46995.443  5  45261.2230  4

16  1421.7331 21  2110.07774 27 45075.768  7   9407.5259 12

17  2139.6204 10  3442.88759 22 29935.337 19  10340.0860 11

18  1516.2422 18  2711.80908 25 30208.886 18   6446.9329 15

19  8718.6086  4 32155.82313  3 63529.246  1 170587.2931  1

20   628.1239 27  2940.36842 24 26614.641 22   2673.3329 22

21   557.7157 29 10980.71251 10  4435.302 30    459.8699 27

22  2071.6715 13  5857.04005 18 16130.619 28   3644.2968 19

23  1262.5608 24  2951.19351 23 53621.218  4  17552.2960  9

24   577.1137 28   856.18442 29 37896.852 12   2830.0411 21

25   542.7976 30  3471.24473 21 36117.588 15   1376.4672 25

26  2135.5043 11  8562.29210 12 27097.013 21   7130.1423 14

27   728.2675 26  2410.73505 26 37030.617 13    554.9793 26

28  1317.3158 23  5110.87264 20 46490.234  6    122.8069 30

29  1565.6284 17  7932.85125 13 57801.387  3    261.3034 28

30  1122.6515 25 10192.18622 11 24374.043 23    195.7004 29

           S  R

1  24400.989  5

2  15924.447  7

3   7848.068 14

4   6217.218 18

5   7164.411 16

6  12511.641 10

7   5730.109 21

8   5323.158 24

9  39718.237  3

10 59643.033  2

11 19963.813  6

12  5593.397 22

13 13017.556  9

14  5762.817 20

15 24631.892  4

16  9078.482 11

17  8306.984 12

18  6540.012 17

19 78402.371  1

20  4455.887 29

21  2470.794 30

22  4692.009 27

23 13081.281  8

24  5368.171 23

25  5013.257 25

26  7567.158 15

27  4714.200 26

28  6163.274 19

29  7910.182 13

30  4554.328 28

结论：从自主创新能力来看东部地区除河北，广西，海南三个省和自治区外，在全国范围内享有绝对优势。中西部的自主创新能力明显低于东部地区。