代码收藏家技术教程 2023-12-13

使用Python实现AHP分析法的代码

AHP的简介

层次分析法（Analytic Hierarchy Process (AHP)）是的一种主观赋值评价方法，将与决策有关的元素分解成目标、准则、方案等多个层次，并在此基础上进行定性和定量分析，是一种系统、简便、灵活有效的决策方法。

简单来说，AHP常用于解决评价类问题，评价类问题例子：购买一个物品的时候会跟其他店的同款商品进行比较，最终确定一个商品。

AHP经常将一个决策（一个决定）分解为目标层，准侧层，方案层。例子如图所示。

而且AHP中涉及的矩阵都需要是一致矩阵。

AHP的局限

评价的决策层（方案层）和准则层不能过多，因为如果决策层和准则层过多，一致性检验难以实现（最大为n=15），而且判断矩阵和一致矩阵差异较大。从而难以通过一致性检验。

决策层的指标和数据已知，我们如何更好的利用已知的数据进行构造一致性矩阵。

AHP的基础知识

正互反矩阵成立条件

矩阵中的每个元素都要大于0

任意 $\text{[math]}$ ，满足这两个条件即为正互反矩阵。

一致矩阵成立条件

矩阵中的每个元素都要大于0

任意 $\text{[math]}$

任意 $\text{[math]}$

3. 一致性检验

AHP的的代码

在这里只需要修改准则层对于目标层的判断矩阵，以及对于方案层对于准则层矩阵的判断矩阵即可。其他都为现成的代码。可直接应用，但是需要按照对应进行修改方案层对于准则层矩阵的判断矩阵的个数。

全部代码

import numpy as np
import pandas as pd
import warnings


# 这里的意思是求权重之前我们需要进行一致性检验
# 参考：https://zhuanlan.zhihu.com/p/101505929
# 参考：https://blog.csdn.net/knighthood2001/article/details/127519604?spm=1001.2101.3001.6650.2&utm_medium=distribute.pc_relevant.none-task-blog-2%7Edefault%7ECTRLIST%7ERate-2-127519604-blog-98480769.pc_relevant_recovery_v2&depth_1-utm_source=distribute.pc_relevant.none-task-blog-2%7Edefault%7ECTRLIST%7ERate-2-127519604-blog-98480769.pc_relevant_recovery_v2&utm_relevant_index=3
# 一致性检验
def calculate_weight(data):
    RI = (0, 0, 0.52, 0.89, 1.12, 1.26, 1.36, 1.41, 1.46, 1.49, 1.52, 1.54, 1.56, 1.58, 1.59)
    # 转化为array类型的对象
    in_matrix = np.array(data)
    n, n2 = in_matrix.shape
    # 判断矩阵是否为方阵，而且矩阵的大小为n，n2
    if n != n2:
        print("不是一个方阵，所以不能进行接下来的步骤")
        return None
    for i in range(0, n):
        for j in range(0, n2):
            if np.abs(in_matrix[i, j] * in_matrix[j, i] - 1) > 1e-7:
                raise ValueError("不为正互反矩阵")
    eig_values, eig_vectors = np.linalg.eig(in_matrix)
    # eigvalues为特征向量，eigvectors为特征值构成的对角矩阵（而且其他位置都为0，对角元素为特征值）
    max_index = np.argmax(eig_values)
    # argmax为获取最大特征值的下标,而且这里是获取实部
    max_eig = eig_values[max_index].real
    # 这里max_eig是最大的特征值
    eig_ = eig_vectors[:, max_index].real
    eig_ = eig_ / eig_.sum()
    if n > 15:
        CR = None
        warnings.warn(("无法判断一致性"))
    else:
        CI = (max_eig - n) / (n - 1)
        if RI[n - 1] != 0:
            CR = CI / RI[n - 1]
        if CR < 0.1:
            print("矩阵的一致性可以被接受")
        else:
            print("矩阵的一致性不能被接受")
    return max_eig, CR, eig_


# 特征值法求权重
def calculate_feature_weight(matrix, n):
    # 特征值法主要是通过求出矩阵的最大特征值和对应的特征向量，然后对其特征向量进行归一化，最后获得权重
    eigValue, eigVectors = np.linalg.eig(matrix)
    max_index = np.argmax(eigValue)
    max_eig = eigValue[max_index].real
    eig_ = eigVectors[:, max_index].real
    # 返回的是特征向量，而且max_index为最大的特征值，在这里一般为n
    eig_ = eig_ / eig_.sum()
    # 这里返回的是特征向量
    return eig_


# 算术平均法求权重
def calculate_arithemtic_mean(matrix):
    n = len(matrix)
    matrix_sum = sum(matrix)
    normalA = matrix / matrix_sum  # 归一化处理
    average_weight = []
    for i in range(0, n):
        # 按照列求和
        temSum = sum(normalA[i])
        average_weight.append(temSum / n)
    return np.array(average_weight)


# 几何平均法求权重
def calculate_metric_mean(metrix):
    n = len(metrix)
    # 1表示按照行相乘，得到一个新的列向量,每行相乘获得一个列向量，所以用prod函数，
    vector = np.prod(metrix, 1)
    tem = pow(vector, 1 / n)
    # 开n次方
    # 归一化处理
    average_weight = tem / sum(tem)
    return average_weight


if __name__ == "__main__":
    # 准则重要性矩阵,即为我们在这里的时候要输入准则性矩阵
    criteria = np.array([[1, 2, 7, 5, 5],
                         [1 / 2, 1, 4, 3, 3],
                         [1 / 7, 1 / 4, 1, 1 / 2, 1 / 3],
                         [1 / 5, 1 / 3, 2, 1, 1],
                         [1 / 5, 1 / 3, 3, 1, 1]])

    # 对每个准则，方案优劣排序,即为方案层也会有一个一致矩阵，所以需要判断
    b1 = np.array([[1, 1 / 3, 1 / 8], [3, 1, 1 / 3], [8, 3, 1]])
    b2 = np.array([[1, 2, 5], [1 / 2, 1, 2], [1 / 5, 1 / 2, 1]])
    b3 = np.array([[1, 1, 3], [1, 1, 3], [1 / 3, 1 / 3, 1]])
    b4 = np.array([[1, 3, 4], [1 / 3, 1, 1], [1 / 4, 1, 1]])
    b5 = np.array([[1, 4, 1 / 2], [1 / 4, 1, 1 / 4], [2, 4, 1]])
    b = [b1, b2, b3, b4, b5]
    matrix_in = criteria
    max_eigen, CR, criteria_eigen = calculate_weight(matrix_in)
    print("准则层：最大特征值：{:.5f},CR={:<.5f},检验{}通过".format(max_eigen, CR, '' if CR < 0.1 else "不"))
    print("准则层权重为{}\n".format(criteria_eigen))
    max_eigen_list = []
    CR_list = []
    eigen_list = []
    for i in b:
        max_eigen, CR, eigen = calculate_weight(i)
        max_eigen_list.append(max_eigen)
        CR_list.append(CR)
        eigen_list.append(eigen)
    pd_print = pd.DataFrame(eigen_list, index=["准则" + str(i) for i in range(0, criteria.shape[0])],
                           columns=["方案" + str(i) for i in range(0, b[0].shape[0])])
    pd_print.loc[:, '最大特征值'] = max_eigen_list
    pd_print.loc[:, 'CR'] = CR_list
    pd_print.loc[:, '一致性检验'] = pd_print.loc[:, 'CR'] < 0.1
    print("方案层")
    print(pd_print)
    # 目标层
    # np.dot()函数为向量点积和矩阵乘法，即为AHP最后的目的是将准则层的的特征向量和方案层的特征向量进行矩阵乘法，而且最后是1*方案层的矩阵，
    # criteria_eigen的shape为（1,5），而且eight_list为（5,3）的矩阵
    # 而且reshape类似于转置矩阵的作用，所以使得原来为（5，）变成（1，5）
    object = np.dot(criteria_eigen.reshape(1, -1), np.array(eigen_list))
    print("\n目标层", object)
    print("最优选择方案{}".format(np.argmax(object)))

权重求法

求权重的方法有三种，分别为几何，算术，特征值求权重法。

几何法求权重

将矩阵按行相乘，求得一个列向量，之后再将向量进行开n次方，最后进行归一化处理（归一化处理每一个元素除以其所在列的和）

# 几何平均法求权重
def calculate_metric_mean(metrix):
    n = len(metrix)
    # 1表示按照行相乘，得到一个新的列向量,每行相乘获得一个列向量，所以用prod函数，
    vector = np.prod(metrix, 1)
    tem = pow(vector, 1 / n)
    # 开n次方
    # 归一化处理
    average_weight = tem / sum(tem)
    return average_weight

算术平均法求权重

将判断矩阵按照列归一化（每一个元素除以其所在列的和），将归一化的各列相加(按行求和)，将相加后得到的向量中每个元素除以n即可得到权重向量。而且也为列向量。

# 算术平均法求权重
def calculate_arithemtic_mean(matrix):
    n = len(matrix)
    matrix_sum = sum(matrix)
    normalA = matrix / matrix_sum  # 归一化处理
    average_weight = []
    for i in range(0, n):
        # 按照列求和
        temSum = sum(normalA[i])
        average_weight.append(temSum / n)
    return np.array(average_weight)

特征值求权重

求出矩阵A的最大特征值以及其对应的特征向量，对求出的特征向量进行归一化即可得到我们的权重。在这里我采用调用第三方库的方法直接求取。


# 特征值法求权重
def calculate_feature_weight(matrix, n):
    # 特征值法主要是通过求出矩阵的最大特征值和对应的特征向量，然后对其特征向量进行归一化，最后获得权重
    eigValue, eigVectors = np.linalg.eig(matrix)
    max_index = np.argmax(eigValue)
    max_eig = eigValue[max_index].real
    eig_ = eigVectors[:, max_index].real
    # 返回的是特征向量，而且max_index为最大的特征值，在这里一般为n
    eig_ = eig_ / eig_.sum()
    # 这里返回的是特征向量
    return eig_