R使用层次分析法进行综合指标等级划分

业务临时需要，需要确定多因素影响下的综合权重值，现使用层次分析法和拉格朗日多项式插值算法做简易值计算

1、建立层次分析结构模型，分析影响综合指标的各个因素，分层级，上层受下层影响，而同层各因素之间基本上相对独立（本实例只有两层）

2、连接oracle数据库，读取各个因素实际指标值

library(DBI)

library(ROracle)

drv=dbDriver('Oracle')

conn=dbConnect(drv,'AQTS_ZHZX','AQTS_ZHZX','192.160.0.0:1521/iims')

rs=dbSendQuery(conn,"with weather as

(

  select * from (

  select * from weather_fact_aqts w 

  where w.cjrq=(select max(cjrq) from weather_fact_aqts) 

  order by w.cjsd desc) where rownum=1

)

               

               select v.yfbm,v.dlbm,v.ldbm,weather.cjrq,weather.cjsd,v.wfshb,a.sgshb,

               weather.sbbh,weather.wd,weather.sd,weather.fl,weather.njd

               ,weather.fog_rank,weather.wind_rank from VIO_FACT_AQTS v

               join ACD_FACT_AQTS a on v.yfbm=a.sgyfbm and v.dlbm=a.lhbm and v.ldbm=a.glsbm

               join weather on 1=1") 

data=fetch(rs)

data=as.matrix(data)

3、影响综合因素的指标分两级：

#FIRST,A1:weather,A2:Flux,A3:acd,A4:vio

#SECOND1,B1:fog,B2:rain,B3:wind,B4:snow

#SECOND2,B2:FLUX1，FLUX2

4、构造第一层成对比较矩阵（采用成对比较法和1~9尺度）

#FIRST,A1:weather,A2:Flux,A3:acd,A4:vio

A=matrix(0,4,4)

diag(A)=1

A[1,2]=1/2;A[2,1]=2

A[1,3]=1/4;A[3,1]=4

A[1,4]=1/4;A[4,1]=4

A[2,4]=1/2;A[4,2]=2

A[2,3]=1/2;A[3,2]=2

A[3,4]=A[4,3]=1

5、构造第二层成对比较矩阵（构造针对上一层成对比较矩阵）

#SECOND,B1:fog,B2:rain,B3:wind,B4:snow

B=matrix(0,4,4)

diag(B)=1

B[1,2]=2;B[2,1]=1/2

B[1,3]=6;B[3,1]=1/6

B[1,4]=2;B[4,1]=1/2

B[2,3]=3;B[3,2]=1/3

B[2,4]=1;B[4,2]=1

B[3,4]=1/3;B[4,3]=3

#SECOND2,B2:FLUX1，FLUX2

C=matrix(1,2,2)

C[1,2]=3

C[2,1]=1/3

6、计算权向量并一致性校验

求权向量以及特征向量，进行一致性检验（因时间问题，这块没有加校验，后续修改）

7、求各层级因素权重，这里使用规范列平均法

Ac=colSums(A)

for (i in 1:nrow(A)) {

  A[,i]=A[,i]/Ac[i]

}

weigA=rowMeans(A)

#B

Bc=colSums(B)

for (i in 1:nrow(B)) {

  B[,i]=B[,i]/Bc[i]

}

weigB=rowMeans(B)

#C 1:flux,2:vehicle

Cc=colSums(C)

for (i in 1:nrow(C)) {

  C[,i]=C[,i]/Cc[i]

}

weigC=rowMeans(C)

########row 1:fog,2:rain,3:wind,4:snow,5:flux,6:vehicle,7:acd,8:vio

########new row 1:fog,2:rain,3:wind,4:snow,5:even,6:acd,7:vio

weit=matrix(0,8,4)

weit[1:4,1]=weigB

weit[5:6,2]=weigC

weit[7,3]=1

weit[8,4]=1

######

w=as.matrix(weigA)

p=weit

we=p%*%w

wenew=c(we[1:4],sum(we[5:6]),we[7],we[8])

8、构造各个因素的函数

#get the function of accident

acdx<-c(0,0.1,0.2,0.3,0.4)

acdy<-c(5,4,3,2,1)

acdf <- LagrangePolynomial(acdx,acdy)

acdv <- function(x,f){

  x=as.numeric(x)

  ranka=f(x)

  ranka[ranka>5]=5

  ranka[ranka<0]=0

  return (ranka)

}

#v=c(-0.1,0.2)

#vv=acdv(v,acdf)

#get the function of violation

viox<-c(0,0.1,0.2,0.3,0.4)

vioy<-c(5,4,3,2,1)

viof <- LagrangePolynomial(viox,vioy)

#get the function of flux flux

fluxfx<-c(0,0.15,0.3,0.45,0.6)

fluxfy<-c(5,4,3,2,1)

fluxff <- LagrangePolynomial(fluxfx,fluxfy)

#get the function of flux vehicel

fluxvx<-c(0,0.2,0.35,0.5,0.65)

fluxvy<-c(5,4,3,2,1)

fluxvf <- LagrangePolynomial(fluxvx,fluxvy)

#get the function of  fog rain snow

fogx<-c(0,50,100,200,500)

fogy<-c(1,2,3,4,5)

fogf <- LagrangePolynomial(fogx,fogy)

#get the function of wind

windx<-c(0,13.9,20.7,28.5,32.7)

windy<-c(5,4,3,2,1)

windf <- LagrangePolynomial(windx,windy)

9、对应各因素指标值与权重相乘求取综合指标级别，并写入数据库

rankc=cbind(data[,c(2:5)],as.matrix(acdv(data[,12],fogf))*wenew[1],as.matrix(acdv(data[,12],fogf))*wenew[2],

as.matrix(acdv(data[,11],fogf))*wenew[3],as.matrix(acdv(data[,12],fogf))*wenew[4],as.matrix(rep(3,nrow(data)))*wenew[5],

as.matrix(acdv(data[,7],acdf))*wenew[6],as.matrix(acdv(data[,6],viof))*wenew[7])

#comprehensive rank

alv=as.numeric(rankc[,5:11])

dim(alv)=c(nrow(rankc),7)

RANKOVER<-cbind(rankc[,1:4],ceiling(rowSums(alv)))

colnames(RANKOVER)<-c('DLBM','LDBM','CJRQ','CJSD','RANK')

RANKOVER=as.data.frame(RANKOVER)

#dbRemoveTable(conn, 'ZHZX_FACT_AQTS')

dbWriteTable(conn,'ZHZX_FACT_AQTS',RANKOVER, row.names = F, append = TRUE)

附：拉格朗日多项式插值算法

#LagrangePolynomial

LagrangePolynomial <- function(x,y) {

  len = length(x)

  if(len != length(y))

    stop("length not equal!")

  

  if(len < 2)

    stop("dim size must more than 1")

  

  #pretreat data abd alloc memery

  xx <- paste("(","a -",x,")")

  m <- c(rep(0,len))

  

  #combin express

  for(i in 1:len) {

    td <- 1

    tm <- "1"

    for(j in 1:len) {

      if(i != j) {

        td <- td*(x[i] - x[j])

        tm <- paste(tm,"*",xx[j])

      }

    }

    tm <- paste(tm,"/",td)

    m[i]<-tm #m[i] <- parse(text=tm)

  }

  

  #combin the exrpession

  m <- paste(m,"*",y)

  r <- paste(m,collapse="+")

  

  #combin the function

  fbody <- paste("{ return(",r,")}")

  f <- function(a) {}

  

  #fill the function's body

  body(f) <- parse(text=fbody)

  return(f)

}

（记录日常点滴）对算法以及R使用还不是很熟悉，后续完善！如有错误，希望指点>>>OK