残缺区间群体决策模型（GDM）
0. 参考论文:
徐泽水： A consistency improving method in the analytic hierarchy process 1999年
刘芳：Acceptable consistency analysis of interval reciprocal comparison matrices 2009年
1. 主要思路（刘芳2012年）：
- 先把残缺区间矩阵\(U_k\)通过某种方法(LOP2)或者利用粒子群算法)进行补全成\(U_k^{’}\)，
- 检查补全以后的矩阵\(U_{k}^{'}\)是否具有一致性(\(U_{k}^{'}\)具有一致性 的充要条件 是 \(U_{k}^{'}\) 分解成单个的正互反判断矩阵B 和D，矩阵B和D的CR <= 0.1)，不具有则调整（按照徐泽水1999年提出的论文进行调整)
- 若\(U_{k}^{'}\)具有一致性，则进行群体决策模型构建，生成最终的区间判断矩阵U
- 通过最终的区间判断矩阵U生成区间权重，注意生成区间权重并没有归一化权重，根据区间权重生成可能度矩阵P.
2.主要函数构建：
consistency(A):求正互反判断矩阵的一致性指标，返回一个listem_get_w(A) :特征值求权重 —— 没有归一化权重gm_get_w(A)：几何平均求权重 — — 没有归一化权重get_w(B,D):分别获取B，D的权重（可以指定几何平均或者特征值求权重），然后组成区间权重向量（即小的在前，大的在后），这里返回的是一个矩阵，把每一个区间数看做矩阵的一行。fenjie(U ):把区间矩阵U分解成正互反判断矩阵B和Dadjust_w(A,lambda) :利用论文的方法进行调整，返回调整后符合一致性条件的一致性矩阵。degree_probability(a,b)函数计算两个区间数的可能度probability_matrix(w)给一个n*2 的区间数，求其可能度矩阵
介绍几个函数—-这几个函数都是利用粒子群算法把\(U_{k}\)变成\(U_{k}^{'}\)
以下是补全区间的关键元素。
index_matrix():此函数可以查找到区间正互反判断矩阵中那些元素是缺失的（缺失用0表示)，并生成一个矩阵，每一列都有其相关含义。hecheng(B,D):对应fenjie的逆向操作，把两个正互反矩阵进行合并成一个区间判断矩阵。`jianyan_yuesu(U):检验输入的区间判断矩阵U是否满足要求，即有没有输入上的出错。set_x(x,U):函数把x向量赋值给U中缺失的元素obj_fun(Ut):求补全好的区间矩阵的目标函数值随机初始化残缺元素
U_to_UU()函数是把区间矩阵U 重塑成(n * n ) 2形式的矩阵，以两列为单位，因为每两列是一个区间，第一列代表区间的下界，第二列代表区间的上界，故 形成2列矩阵，每一行代表一个区间数，UU_to_U()函数是通过找出UU中的缺失元素，然后赋予缺失元素一个随机值（ 该随机值保证了在1/9 ~ 9之间，且保证了UU中的下界与上界的关系，若某行只缺失一个数据，也能得到相应的保证），并把赋予好的完整矩阵变成区间判断矩阵的形式。U_runif()是U_to_UU()和UU_to_U()函数的融合，先利用U_to_UU()变成我们需要的格式，然后利用UU_to_U()产生随机值进行填充最后返回我们随机填充好的区间判断矩阵init_x()通过随机初始化U，找出解
GDM_PSO()主要的函数，通过利用粒子群算法把\(U_{k}\)变成\(U_{k}^{'}\)DGM_U()群体决策，把若干个完整的区间判断矩阵合并成一个区间判断矩阵。
rm(list = ls())
# 0。 一致性指标的求解
consistency = function(A){
lambda = Re(eigen(A)$values[1]) # 矩阵A的最大特征值
n = nrow(A)
RI = c(0,0,0.58,0.90,1.12,1.24,1.32,1.41,1.45);
CI = (lambda-n) / (n-1);
CR = CI / RI[n];
eig_w = eigen(A)$vectors[,1] / sum( eigen(A)$vectors[,1]);
return(list("eig_value"=lambda,"CI"=CI,"RI"=RI[n],"CR"=CR,'eig_w'=Re( eig_w )))
}
# 1. 特征值求权重
em_get_w = function(A){
n = nrow(A)
stopifnot(nrow(A) == ncol(A))
lambda = Re(eigen(A)$values[1]) # 矩阵A的最大特征值
n = nrow(A)
RI = c(0,0,0.58,0.90,1.12,1.24,1.32,1.41,1.45);
CI=(lambda-n)/(n-1);
CR=CI/RI[n];
eig_origin = eigen(A)$vectors[,1]
eig_w = eigen(A)$vectors[,1] # /sum(eigen(A)$vectors[,1]);
eig_w = Re(eig_w)
return(eig_w)
}
# 2. 几何平均求权重
gm_get_w =function(A){
n = nrow(A)
stopifnot(nrow(A) == ncol(A))
temp = apply(A, 1, function(x) prod(x)^(1/n) )
w = temp # /sum(temp)
return(w)
}
# 3. 合并权重
get_w = function(B,D,method = c('gm_get_w','em_get_w') ){
n = nrow(B)
stopifnot(n == nrow(D))
method <- match.arg(method)
f = get(method)
w_B = f(B)
w_D = f(D)
w_L = rep(0,length(w_B))
w_U = rep(0,length(w_B))
for(i in 1:length(w_B)){
w_L[i] = min(w_B[i],w_D[i])
w_U[i] = max(w_B[i],w_D[i])
}
w = matrix(c(w_L,w_U), ncol= 2 )
return(w)
}
# 4. 通过U进行分解，分解出B，D矩阵，
fenjie = function(U){
n = nrow(U)
stopifnot(ncol(U) == 2*n)
B = matrix(0,nrow = n,ncol = n)
D = matrix(0,nrow = n,ncol = n)
for(i in 1:n){
for(j in 1:n){
if(i<j){
B[i,j] = U[i,j*2]
D[i,j] = U[i,2*j-1]
}else if(i>j){
B[i,j] = U[i,2*j-1]
D[i,j] = U[i,j*2]
}else{
B[i,j] = U[i,2*j]
D[i,j] =U[i,2*j]
}
}
}
return(list('B'=B,'D'=D))
}
# 5. 徐泽水(1999年)的文章方法进行调整，直到满足一致性条件为准(CR <=0.1).
adjust_w <- function(A, lambda) {
k <- 0
n = nrow(A)
m = ncol(A)
stopifnot(n == m)
temp_CR <- consistency(A)$CR
temp_w <- consistency(A)$eig_w
while (temp_CR >= 0.1 && k < 1000) {
for (i in 1:n) {
for (j in 1:n) {
A[i, j] <- (A[i, j]^lambda) * (temp_w[i] / temp_w[j])^(1 - lambda)
}
}
temp_CR <- consistency(A)$CR
temp_w <- consistency(A)$eig_w
k <- k + 1
}
return(A)
}
# 6. degree_probability 函数计算两个区间数的可能度
degree_probability <- function(a, b) {
# 输入的a,b代表一个区间，即是一个二维向量，且小的在前面，大的元素在后面
stopifnot(length(a) == length(b), length(a) == 2, a[1] <= a[2], b[1] <= b[2])
temp <- 0
if (a[1] == a[2] && b[1] == b[2]) {
if (a[1] > b[1]) {
temp <- 1
} else if (a[1] == b[1]) {
temp <- 0.5
} else {
temp <- 0
}
} else if (a[1] == a[2] && b[1] != b[2]) {
if (a[1] > b[2]) {
temp <- 1
} else if (b[1] <= a[1] & a[1] <= b[2]) {
temp <- (a[1] - b[1]) / (b[2] - b[1])
} else {
temp <- 0
}
} else if (a[1] != a[2] && b[1] == b[2]) {
if (a[1] > b[1]) {
temp <- 1
} else if (a[1] <= b[1] & b[1] <= a[2]) {
temp <- (a[2] - b[1]) / (a[2] - a[1])
} else {
temp <- 0
}
} else if (a[1] != a[2] && b[1] != b[2]) {
if (a[1] < a[2] && a[2] <= b[1] && b[1] < b[2]) {
temp <- 0
} else if (a[1] <= b[1] && b[1] < a[2] && a[2] <= b[2]) {
s_t <- (a[2] - b[1]) * (a[2] - b[1]) * 0.5
s <- (b[2] - b[1]) * (a[2] - a[1])
temp <- s_t / s
} else if (a[1] <= b[1] && b[1] <= b[2] && b[2] <= a[2]) {
s_t <- ((a[2] - b[2]) + (a[2] - b[1])) * (b[2] - b[1]) * 0.5
s <- (b[2] - b[1]) * (a[2] - a[1])
temp <- s_t / s
} else if (b[1] < a[1] && a[1] < a[2] && a[2] < b[2]) {
# 可写等号b[1] <= a[1] && a[1]<a[2] &&a[2]<=b[2]
s_t <- ((a[1] - b[1]) + (a[2] - b[1])) * (a[2] - a[1]) * 0.5
s <- (b[2] - b[1]) * (a[2] - a[1])
temp <- s_t / s
} else if (b[1] < a[1] && a[1] < b[2] && b[2] < a[2]) {
s_tt <- (b[2] - a[1]) * (b[2] - a[1]) * 0.5
s <- (b[2] - b[1]) * (a[2] - a[1])
s_t <- s - s_tt
temp <- s_t / s
} else if (b[1] < b[2] && b[2] <= a[1] && a[1] < a[2]) {
temp <- 1
} else {
stop("运行出错,请检查")
}
} else {
stop("运行出错,请检查")
}
return(temp)
}
# 7. probability_matrix 给一个n*2 的区间数，求其可能度矩阵
probability_matrix <- function(Z) {
# probability_matrix函数输入一个n*2的矩阵，每一行代表输出各个方案的综合属性值得区间数
# 此函数输出各方案两两比较的可能度矩阵。
# degree_probability函数求两个区间数的可能度，
# a，b代表输入的区间数，输入这两个
P <- matrix(0, ncol = nrow(Z), nrow = nrow(Z))
for (i in 1:nrow(Z)) {
for (j in 1:nrow(Z)) {
P[i, j] <- degree_probability(Z[i, ], Z[j, ])
}
}
return(P)
}
# 8. index_matrix此函数可以查找到区间正互反判断矩阵中那些元素是缺失的（缺失用0表示)
index_matrix <- function(U) {
nnn <- which(U == 0, arr.ind = T)[, 1] # 行
mmm <- which(U == 0, arr.ind = T)[, 2] # 列
index <- matrix(0, nrow = length(nnn), ncol = 8)
# 矩阵index 每一行相当于一个变量
# 第一列代表 X索引，即变量X1,X2,...Xn
# 第二三列 代表变量Xi 所在U中的位置，用行和列表示
# 第四列，第五列，根据变量Xi在U的位置，
# 第6 列，判断Xi与Xj 是否为同一个位置，
index[, 1] <- 1:length(nnn)
index[, 2] <- nnn #
index[, 3] <- mmm
index
for (i in 1:nrow(index)) {
if (index[i, 3] %% 2 == 0) {
# 第3列的数字是偶数,检查U矩阵向前一列是否为0
if (U[index[i, 2], index[i, 3] - 1] == 0) {
index[i, 4] <- 0
index[i, 5] <- 0
} else {
index[i, 4] <- U[index[i, 2], index[i, 3] - 1]
index[i, 5] <- 0
}
} else {
if (U[index[i, 2], index[i, 3] + 1] == 0) {
index[i, 4] <- 0
index[i, 5] <- 0
} else {
index[i, 4] <- 0
index[i, 5] <- U[index[i, 2], index[i, 3] + 1]
}
}
}
for (i in 1:nrow(index)) {
if (index[i, 4] != 0 & index[i, 5] == 0) {
index[i, 5] <- 9
}
}
for (i in 1:nrow(index)) {
if (index[i, 5] != 0 & index[i, 4] == 0) {
index[i, 4] <- 1 / 9
}
}
index[, 6] <- ifelse((index[, 4] == index[, 5] & index[, 4] == 0), 1, 0)
jinjibiao <- c()
for (i in 1:nrow(index)) {
for (j in 1:nrow(index)) {
if (index[i, 6] == 1 & index[j, 6] == 1 & j != i) {
if (index[i, 2] == index[j, 2] & abs(index[i, 3] - index[j, 3]) == 1) {
if (!(index[i, 7] %in% jinjibiao)) {
jinjibiao <- c(jinjibiao, j)
index[i, 7] <- index[i, 1]
index[i, 8] <- index[j, 1]
}
}
}
}
}
for (i in 1:nrow(index)) {
if (index[i, 7] > index[i, 8]) {
t <- index[i, 8]
index[i, 8] <- index[i, 7]
index[i, 7] <- t
}
}
return(index)
}
# 9. 合成
hecheng <- function(B, D) {
n = nrow(B)
m = ncol(B)
stopifnot(n == m, n == nrow(D) )
UU <- matrix(0, nrow = n, ncol = n * 2)
for (i in 1:n) {
for (j in 1:(2*n) ) {
k <- (j + 1) %/% 2
if (i > k) {
if (j %% 2 == 0) {
UU[i, j] <- D[i, k]
} else {
UU[i, j] <- B[i, k]
}
} else if (i < k) {
if (j %% 2 == 0) {
UU[i, j] <- B[i, k]
} else {
UU[i, j] <- D[i, k]
}
} else {
UU[i, j] <- 1
}
}
}
return(UU)
}
####################################################################################
######################## 构建 目标函数 补全区间残缺矩阵 ##################
# 10. 检验U矩阵是否符合条件
jianyan_yuesu <- function(U) {
if (!all(U >= 1 / 9 & U <= 9)) {
stop("U矩阵不满足约束条件")
}
B <- fenjie(U)$B # 先对U进行分解，分解成B，D约束条件
D <- fenjie(U)$D
R <- matrix(0, nrow = nrow(D), ncol = ncol(D)) # 存储C与D比较的结果
for (i in 1:nrow(B)) {
for (j in 1:ncol(B)) {
if (i < j) {
R[i, j] <- (B[i, j] >= D[i, j])
} else if (i > j) {
R[i, j] <- (B[i, j] <= D[i, j])
} else {
R[i, j] <- all(B[i, j] == D[i, j], B[i, j] == 1)
}
}
}
if (!all(R == 1)) {
stop("U矩阵分解后的B,D矩阵不满足约束条件")
}
return(1) # 返回结果为1 即满足约束
}
# 11. 把x的值赋予U中的残缺元素
set_x <- function(x, U) {
U[which(U == 0)] <- x
return(U)
}
# 12，U补全后的矩阵记为Ut，计算补全后的目标函数值
obj_fun <- function(Ut) {
s_sum <- 0
B_sum <- 0
D_sum <- 0
B <- fenjie(Ut)$B
D <- fenjie(Ut)$D
n <- nrow(B)
if (jianyan_yuesu(Ut)) {
for (i in 1:n) {
for (j in 1:n) {
for (k in 1:n) {
B_sum <- B_sum + (log(B[i, j], 9) - log(B[i, k], 9) - log(B[k, j], 9))^2
D_sum <- D_sum + (log(D[i, j], 9) - log(D[i, k], 9) - log(D[k, j], 9))^2
}
}
}
s_sum <- s_sum + B_sum + D_sum
return(s_sum)
} else {
stop("!!!输入的矩阵不符合约束条件")
}
}
#########################################################
#######################################################
######################## 13. 初始化x #######
######### U_runif 通过U产生 符合要求的随机矩阵U2
# 13.1 重塑U ---把U变成UU（即(n*n )*2)形式的矩阵）
U_to_UU <- function(U) {
n <- nrow(U)
m <- ncol(U) / 2
stopifnot(m >= 4, m == nrow(U))
# 先重塑U,U的每两列进行列合并成UU，通过UU产生0所在位置的随机数，并返回0所在的行下标，
# 即矩阵B，D为0的下标
UU <- U[, 1:2]
for (j in 2:(ncol(U) / 2)) {
temp1 <- U[, c(2 * j - 1, 2 * j)]
UU <- rbind(UU, temp1)
}
return(UU)
}
# 13.2 UU 转变成U的形式，并赋予随机值
UU_to_U <- function(UU) {
n <- sqrt(nrow(UU))
if (!as.integer(n)) stop("n不是整数")
BD_index2 <- t(apply(UU, 1, function(x) {
temp <- c(0, 0)
if (x[1] == 0 && x[2] != 0) {
temp[1] <- runif(1, 1 / 9, x[2])
temp[2] <- x[2]
} else if (x[1] != 0 && x[2] == 0) {
temp[1] <- x[1]
temp[2] <- runif(1, x[1], 9)
} else if (x[1] == 0 && x[2] == 0) {
temp <- sort(runif(2, 1 / 9, 9))
} else {
temp <- x
}
return(temp)
}))
# 把BD_index2 转变为 U的形式
stopifnot(nrow(BD_index2) %% n == 0)
U2 <- BD_index2[1:n, ]
for (j in 2:(nrow(BD_index2) / n)) {
temp2 <- BD_index2[(n * j - n + 1):(n * j), ]
U2 <- cbind(U2, temp2)
}
return(U2)
}
# 13.3 直接通过U随机初始化U
U_runif <- function(U) {
UU <- U_to_UU(U)
# 对UU的每一行进行判断，产生出符合要求的随机数，并使BD_index2的形式与UU的形式一样，都为 (n*m) * 2矩阵
U2 <- UU_to_U(UU)
return(U2)
}
# 13.4 通过随机初始化U，找出解
init_x <- function(U) {
U2 <- U_runif(U)
return(U2[which(U == 0)])
}
###############################################################
## 14。通过粒子群算法求目标函数最小值，来补全残缺区间矩阵的元素
GDM_PSO <- function(U) {
Popsize <- 100 # 种群的大小
Dim <- sum(U == 0) # 粒子的维度
c1 <- 2 # 学习因子
c2 <- 2 # 学习因子
LB <- 1 / 9 # 下限
UB <- 9 # 上限
w_start <- 0.9
w_end <- 0.4
Vmax <- 5
Vmin <- -Vmax
MaxIter <- 3000
index <- index_matrix(U)
# 初始化粒子位置和速度
X <- matrix(0, nrow = Popsize, ncol = Dim)
for (i in 1:Popsize) {
X[i, ] <- init_x(U)
}
V <- matrix(runif(Popsize * Dim), nrow = Popsize, ncol = Dim)
# 设置当前位置为粒子的最好位置，并记录其最好值
PBest <- X
FPBest <- apply(X, 1, function(xx) obj_fun(set_x(xx, U)))
# FPBest
if (anyNA(FPBest)) {
stop("适应度函数出现NA 或NaN，请检查")
}
# 找出初始微粒群体的最好微粒
Best <- PBest[which.min(FPBest), ] # 最好的微粒
FGBest <- min(FPBest) # 以及最好的微粒对应的适应度值
FX <- c(FGBest, rep(0, MaxIter - 1)) # 用来记录每一代的最粒子对应的适应度值
iter <- 0 # 初始迭代次数
while (iter < MaxIter) {
# if(iter%%100 == 0) print(iter) # 显示循环进度条
iter <- iter + 1
# 更新权重的值
w_now <- (w_start - (w_start - w_end) / MaxIter * iter)
A <- matrix(rep(Best, each = Popsize), nrow = Popsize, ncol = Dim)
# 生成随机数
R1 <- matrix(runif(Popsize * Dim), nrow = Popsize, ncol = Dim)
R2 <- matrix(runif(Popsize * Dim), nrow = Popsize, ncol = Dim)
# 速度更新
V <- w_now * V + c1 * R1 * (PBest - X) + c2 * R2 * (A - X)
# 对进化后速度大于最大素的的微粒进行处理
changes <- V > Vmax
V[changes] <- Vmax
changes <- V < Vmin
V[changes] <- Vmin
# 微粒位置进行更新
X <- X + 1.0 * V
# 对进化后微粒位置 大于搜索空间的的微粒进行处理
changes <- X > UB
X[changes] <- UB
changes <- X < LB
X[changes] <- LB
# 并再次检验是否满足约束，如果不满足约束，则进行调整
for (ii in 1:nrow(index)) {
if (index[ii, 6] == 0) {
X[which(X[, ii] < index[ii, 4]), ii ] <- index[ii, 4]
X[which(X[, ii] > index[ii, 5]), ii] <- index[ii, 5]
} else {
X[which(X[, ii] < 1 / 9), ii ] <- 1 / 9
X[which(X[, ii] > 9), ii] <- 9
# 交换值 ，
if (ii %in% index[, 7] & index[ii, 7] != 0) {
kkkk <- ii
jjjj <- index[ii, 8]
Low <- ifelse(X[, kkkk] > X[, jjjj], X[, jjjj], X[, kkkk])
UP <- ifelse(X[, kkkk] > X[, jjjj], X[, kkkk], X[, jjjj])
X[, kkkk] <- Low
X[, jjjj] <- UP
}
}
}
# 重新计算新位置的适应度值
FPBest2 <- apply(X, 1, function(xx) obj_fun(set_x(xx, U)))
# 更新每个微粒最好的位置
P <- (FPBest2 < FPBest)
FPBest[P] <- FPBest2[P] # 适应度值更换
PBest[P, ] <- X[P, ] # 粒子位置更换
# 保存最好的粒子和适应度值
Best <- PBest[which.min(FPBest), ]
FGBest <- min(FPBest)
FX[iter] <- FGBest
}
# Best 为找到的最优解，把最优解赋值给区间型矩阵
Ut = set_x(Best,U) # Ut 即补全后的区间型矩阵
obj = obj_fun(Ut ) # 最优区间矩阵对于的目标函数值
return(list("best" = Best, "FGBest" = FGBest,"Ut" = Ut,"obj"= obj))
}
## 15. DGM_U()群体决策，把若干个完整的区间判断矩阵合并成一个区间判断矩阵。
DGM_U = function(...){
####### 第一部分是对参数的检查 #####################
# 输入的参数为U_{k}，U_{k}代表一个完整的区间乘积互反判断矩阵。
n = nrow(..1)
m = ncol(..1)
args <- list(...)
stopifnot( length(args) >=2, all( map_lgl(args,is.matrix)),2*n ==m )
# 检查所有的矩阵维度是否一样
library(purrr)
t = map(args,dim) %>% do.call(rbind, .)
stopifnot( nrow(unique(t)) == 1)
####### 第二部分 对参数的操作，达到想要的目的###########
## 1. 使用快捷的函数
k = 1 / length(args)
UB = reduce( map(args,function(x)fenjie(x)$B),`*`)^k
UD = reduce( map(args, function(x)fenjie(x)$D), `*`)^k
return( hecheng(UB,UD))
}3.测试
3.1 测试1
先利用一个区间残缺矩阵进行补全操作，然后计算其可能度矩阵P。
U = matrix(c(1,1,1,2,0,2,2,3,
1/2,1,1,1,0,0,2,5,
1/2,0,0,0,1,1,3,0,
1/3,1/2,1/5,1/2,0,1/3,1,1),nrow = 4,byrow = T)
U
#> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
#> [1,] 1.0000000 1.0 1.0 2.0 0 2.0000000 2 3
#> [2,] 0.5000000 1.0 1.0 1.0 0 0.0000000 2 5
#> [3,] 0.5000000 0.0 0.0 0.0 1 1.0000000 3 0
#> [4,] 0.3333333 0.5 0.2 0.5 0 0.3333333 1 1
best_U = GDM_PSO(U)
best_U
#> $best
#> [1] 1.5000000 0.7745967 1.5000000 0.6666667 0.6666667 0.3333333 1.2909944 3.0000000
#>
#> $FGBest
#> [1] 2.560765
#>
#> $Ut
#> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
#> [1,] 1.0000000 1.0 1.0000000 2.0 0.6666667 2.0000000 2 3
#> [2,] 0.5000000 1.0 1.0000000 1.0 0.6666667 1.2909944 2 5
#> [3,] 0.5000000 1.5 0.7745967 1.5 1.0000000 1.0000000 3 3
#> [4,] 0.3333333 0.5 0.2000000 0.5 0.3333333 0.3333333 1 1
#>
#> $obj
#> [1] 2.560765
library(magrittr)
library(purrr)
## # 矩阵U进行分解成B，D，并检验B，D的CR值
best_U$Ut %>%fenjie() %>% map(.,function(x)consistency(x)$CR)
#> $B
#> [1] 0.04881706
#>
#> $D
#> [1] 9.026064e-11
# 直接求出区间权重，以及根据权重求出区间可能度矩阵
( w = get_w(B= fenjie(best_U$Ut)$B ,D = fenjie(best_U$Ut)$D) ) #每一行对应第i个方案的区间权重
#> [,1] [,2]
#> [1,] 1.0745699 1.861210
#> [2,] 1.0745699 1.340343
#> [3,] 1.0382256 1.611855
#> [4,] 0.3860974 0.537285
( P = probability_matrix(w) ) # 可能度矩阵P
#> [,1] [,2] [,3] [,4]
#> [1,] 0.5000000 0.8310704 0.6801310 1.0
#> [2,] 0.1689296 0.5000000 0.2950183 1.0
#> [3,] 0.3198690 0.7049817 0.5000000 1.0
#> [4,] 0.0000000 0.0000000 0.0000000 0.53.2 测试2
多个区间矩阵的操作
U1 = matrix(c(1,1,2,4,0,0,0,0,
1/4,1/2,1,1,1,3,0,0,
0,0,1/3,1,1,1,1/2,1,
0,0,0,0,1,2,1,1),nrow = 4,byrow = T)
U2 = matrix(c(1,1,1,3,2,5,3,6,
1/3,1,1,1,3,0,1/2,4,
1/5,1/2,0,1/3,1,1,0,0,
1/6,1/3,1/4,2,0,0,1,1),nrow = 4,byrow = T)
U3 = matrix(c(1,1,0,0,1/3,2,2,4,
0,0,1,1,0,0,1/3,2,
1/2,3,0,0,1,1,1/2,3,
1/4,1/2,1/2,3,1/3,2,1,1),nrow = 4,byrow = T)
#########################################################################
############### 对U1进行操作 ###########
best_U1 = GDM_PSO(U1)
best_U1
#> $best
#> [1] 0.1111111 0.1111111 0.7177010 2.1337232 0.3849002 9.0000000 1.4055756 9.0000000 0.5114926 0.2640794 9.0000000 2.5980762
#>
#> $FGBest
#> [1] 2.824147
#>
#> $Ut
#> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
#> [1,] 1.0000000 1.000000 2.0000000 4 1.405576 9 0.5114926 9.000000
#> [2,] 0.2500000 0.500000 1.0000000 1 1.000000 3 0.2640794 2.598076
#> [3,] 0.1111111 0.717701 0.3333333 1 1.000000 1 0.5000000 1.000000
#> [4,] 0.1111111 2.133723 0.3849002 9 1.000000 2 1.0000000 1.000000
#>
#> $obj
#> [1] 2.824147
### 矩阵U进行分解成B，D，并检验B，D的CR值
best_U1$Ut %>%fenjie() %>% map(.,function(x)consistency(x)$CR)
#> $B
#> [1] 0.002877359
#>
#> $D
#> [1] 0.1582194#########################################################################
############### 对U2进行操作 ###########
best_U2 = GDM_PSO(U2)
best_U2
#> $best
#> [1] 0.3333333 0.7905694 3.0000000 2.0000000 0.5000000 1.2649111
#>
#> $FGBest
#> [1] 8.227604
#>
#> $Ut
#> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
#> [1,] 1.0000000 1.0000000 1.0000000 3.0000000 2.0000000 5 3.0 6.000000
#> [2,] 0.3333333 1.0000000 1.0000000 1.0000000 3.0000000 3 0.5 4.000000
#> [3,] 0.2000000 0.5000000 0.3333333 0.3333333 1.0000000 1 0.5 1.264911
#> [4,] 0.1666667 0.3333333 0.2500000 2.0000000 0.7905694 2 1.0 1.000000
#>
#> $obj
#> [1] 8.227604
## # 矩阵U进行分解成B，D，并检验B，D的CR值
best_U2$Ut %>%fenjie() %>% map(.,function(x)consistency(x)$CR)
#> $B
#> [1] 0.01941419
#>
#> $D
#> [1] 0.1401432#########################################################################
############### 对U3进行操作 ###########
best_U3 = GDM_PSO(U3)
best_U3
#> $best
#> [1] 0.4158431 0.4158454 2.4372290 1.3458265 2.4377710 3.1576567 0.3133446 0.7397628
#>
#> $FGBest
#> [1] 10.52634
#>
#> $Ut
#> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
#> [1,] 1.0000000 1.0000000 2.437229 2.437771 0.3333333 2.0000000 2.0000000 4
#> [2,] 0.4158431 0.4158454 1.000000 1.000000 0.3133446 0.7397628 0.3333333 2
#> [3,] 0.5000000 3.0000000 1.345827 3.157657 1.0000000 1.0000000 0.5000000 3
#> [4,] 0.2500000 0.5000000 0.500000 3.000000 0.3333333 2.0000000 1.0000000 1
#>
#> $obj
#> [1] 10.52634
## # 矩阵U进行分解成B，D，并检验B，D的CR值
best_U3$Ut %>%fenjie() %>% map(.,function(x)consistency(x)$CR)
#> $B
#> [1] 0.00605522
#>
#> $D
#> [1] 0.2106209发现U2 和U3 填充后没有满足一致性条件，于是需要调整
## U2调整后的矩阵
temp2 = best_U2$Ut %>% fenjie() %>% map2(.,list(0.6,0.88),adjust_w)
hecheng(temp2$B,temp2$D)
#> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
#> [1,] 1.0000000 1.0000000 1.0902384 3.0000000 2.2142197 5.000000 2.5464288 6.000000
#> [2,] 0.3333333 0.9172306 1.0000000 1.0000000 2.7801250 3.000000 0.5831896 4.000000
#> [3,] 0.2000000 0.4516264 0.3333333 0.3596961 1.0000000 1.000000 0.4911658 1.264911
#> [4,] 0.1666667 0.3927068 0.2500000 1.7147083 0.7905694 2.035972 1.0000000 1.000000
## U3调整后的矩阵
temp3 = best_U3$Ut %>% fenjie() %>% map2(.,list(0.6,0.88),adjust_w)
hecheng(temp3$B,temp3$D)
#> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
#> [1,] 1.0000000 1.0000000 2.575345 2.437771 0.4348642 2.0000000 1.5314698 4
#> [2,] 0.4158431 0.3918649 1.000000 1.000000 0.2970813 0.7397628 0.3218428 2
#> [3,] 0.5000000 2.2995685 1.345827 3.341806 1.0000000 1.0000000 0.6476240 3
#> [4,] 0.2500000 0.6529675 0.500000 3.107107 0.3333333 1.5441057 1.0000000 1 如果分解后的B和D没有达到一致性，则进行调整，由于这里用的是PSO算法进行残缺元素的填充，因此带有一定的随机值。于是我们用以下三个矩阵代代替补全后的矩阵(论文中给出的)
U_t1 = matrix(c(1,1,2,4,2,9,1,9,
1/4,1/2,1,1,1,3,0.5,2.25,
1/9,1/2,1/3,1,1,1,1/2,1,
1/9,1,0.4444,2,1,2,1,1),nrow = 4,byrow = T)
U_t2 = matrix(c(1,1,1.0815,3,2.0609,5,2.5209,6,
1/3,0.9246,1,1,2.5515,3,0.6414,4,
1/5,0.4852,1/3,0.3919,1,1,1.1237,1.3333,
1/6,0.3967,1/4,1.5587,0.75,0.8899,1,1),nrow = 4,byrow = T)
U_t3 = matrix(c(1,1,0.7256,2.1183,0.4837,2,1.2234,4,
0.4721,1.3782,1,1,0.6666,0.8008,0.4299,2,
1/2,2.0675,1.2488,1.5001,1,1,0.6449,3,
1/4,0.8174,1/2,2.3262,1/3,1.5507,1,1),nrow = 4,byrow = T)
( U = DGM_U(U_t1,U_t2,U_t3) ) # 多个区间矩阵融合为一个综合的区间矩阵
#> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
#> [1,] 1.0000000 1.0000000 1.1621208 2.9402860 1.2585998 4.481405 1.4555977 6.000000
#> [2,] 0.3401086 0.8604891 1.0000000 1.0000000 1.1936774 1.931622 0.5166013 2.620741
#> [3,] 0.2231443 0.7945333 0.5177063 0.8377192 1.0000000 1.000000 0.7129147 1.587388
#> [4,] 0.1666667 0.6870140 0.3815587 1.9355892 0.6299605 1.402705 1.0000000 1.000000
( w = get_w(fenjie(U)$B,fenjie(U)$D) )# 对综合区间矩阵进行求权重--注意：权重没归一化
#> [,1] [,2]
#> [1,] 1.2079392 2.9818703
#> [2,] 0.8534871 1.1454892
#> [3,] 0.6543922 0.8299693
#> [4,] 0.4473844 1.1686549
probability_matrix(w)#可能度矩阵
#> [,1] [,2] [,3] [,4]
#> [1,] 0.5 1.00000 1.0000000 1.0000000
#> [2,] 0.0 0.50000 1.0000000 0.7654600
#> [3,] 0.0 0.00000 0.5000000 0.4087181
#> [4,] 0.0 0.23454 0.5912819 0.50000004 .总结
如果区间矩阵是残缺的，那么用以下函数去补全：
U = matrix(c(1,1,1,2,0,2,2,3,
1/2,1,1,1,0,0,2,5,
1/2,0,0,0,1,1,3,0,
1/3,1/2,1/5,1/2,0,1/3,1,1),nrow = 4,byrow = T)
U
best_U = GDM_PSO(U)
best_U
## 补全以后要进行一致性检验：
best_U$Ut %>%fenjie() %>% map(.,function(x)consistency(x)$CR
## 一致性检验没通过，则用下面的方法进行调整
temp = best_U$Ut %>% fenjie() %>% map2(.,list(0.6,0.88),adjust_w)
U1 = hecheng(temp$B,temp$D) # 一定是一致性检验通过的区间判断矩阵若是多个完整的区间矩阵，要综合成一个矩阵，则可以直接用下面的方法：
# U_t1 ,U_t2,... ,是完整的区间判断矩阵，且符合一致性检验条件
( U = DGM_U(U_t1,U_t2,U_t3) ) # 多个区间矩阵融合为一个综合的区间矩阵，
( w = get_w(fenjie(U)$B,fenjie(U)$D) )# 对综合区间矩阵进行求权重--注意：权重没归一化
probability_matrix(w)#可能度矩阵