TOPSIS算法



title: TOPSIS算法
date: 2020-02-24 11:18:06
categories: 数学建模
tags: [评价模型, MATLAB]
mathjax: true


定义

? C.L.Hwang和K.Yoon于1981年首次提出TOPSIS (Technique for Order Preference by Similarity to an Ideal Solution),可翻译为逼近理 想解排序法,国内常简称为优劣解距离法。

? TOPSIS 法是一种常用的综合评价方法,能充分利用原始数据的 信息,其结果能精确地反映各评价方案之间的差距。

? 基本过程为先将原始数据矩阵统一指标类型(一般正向化处理) 得到正向化的矩阵,再对正向化的矩阵进行标准化处理以消除各指 标量纲的影响,并找到有限方案中的最优方案和最劣方案,然后分 别计算各评价对象与最优方案和最劣方案间的距离,获得各评价对 象与最优方案的相对接近程度,以此作为评价优劣的依据。该方法对数据分布及样本含量没有严格限制,数据计算简单易行。

例子

评价下表中20条河流的水质情况。 注:含氧量越高越好;PH值越接近7越好;细菌总数越少越好;植物性营养物量介于10‐20之间最佳,超 过20或低于10均不好。

步骤

第一步:将原始矩阵正向化

最常见的四种指标:

指标名称 指标特点 例子
极大型(效益型)指标 越大(多)越好 成绩、GDP增速、企业利润
极小型(成本型)指标 越小(少)越好 费用、坏品率、污染程度
中间型指标 越接近某个值越好 水质量评估时的PH值
区间型指标 落在某个区间最好 体温、水中植物性营养物量

所谓的将原始矩阵正向化,就是要将所有的指标类型统一转化为 极大型指标。(转换的函数形式可以不唯一 )

极小型指标->极大型指标

极小型指标转换为极大型指标的公式:

\[max-x \]

如果所有的元素均为正数,那么也可以使用 1/x

中间型指标->极大型指标

{\(x_{i}\)}是一组中间型指标序列,且最佳的数值为\(x_{best}\),那么正向化的公式如下:

\[M=\max \left\{\left|x_{i}-x_{\text {best}}\right|\right\}, \quad \tilde{x}_{i}=1-\frac{\left|x_{i}-x_{\text {best}}\right|}{M} \]

举个例子:

区间型指标->极大型指标

{\(x_{i}\)}是一组中间型指标序列,且最佳的区间为[a,b],那么正向化的公式如下:

\[M=\max \left\{a-\min \left\{x_{i}\right\}, \max \left\{x_{i}\right\}-b\right\}, \quad \tilde{x}_{i}=\left\{\begin{array}{ll} {1-\frac{a-x}{M}} & {, xb} \end{array}\right. \]

第二步:正向化矩阵标准化

为了消去不同指标量纲的影响, 需要对已经正向化的矩阵进行标准化处理。

假设有n个要评价的对象,m个评价指标(已经正向化了)构成的正向化矩阵如下:

\[A=\begin{bmatrix}a_{11}& a_{12}& \cdots & a_{1n}\\a_{21}& a_{22}& \cdots & a_{2n}\\\vdots & \vdots & \ddots & \vdots \\a_{m1}& a_{m2}& \cdots & a_{mn}\end{bmatrix} \]

那么,对其标准化的矩阵记为Z,Z中的每一个元素:

\[z_{ij}=a_{ij}/ \sqrt{\sum_{i=1}^{n}a_{ij}^{2} } \]

\[每一个元素/\sqrt{V其所在列的元素的平方和} \]

注意:标准化的方法有很多种,其主要目的就是去除量纲的影响,未来我们还可能见到更多 种的标准化方法,例如:(x‐x的均值)/x的标准差;具体选用哪一种标准化的方法在多数情况下 并没有很大的限制,这里我们采用的是前人的论文中用的比较多的一种标准化方法。

第三步:计算得分并归一化

计算评分的公式

\[\frac{x-min}{max-min} \]

解释原因

(1)比较的对象一般要远大于两个。(例如比较一个班级的成绩)
(2)比较的指标也往往不只是一个方面的,例如成绩、工时数、课 外竞赛得分等。
(3)有很多指标不存在理论上的最大值和最小值,例如衡量经济增 长水平的指标:GDP增速

更深本质

\[\frac{x-min}{max-min}=\frac{x-min}{(max-x)+(x-min)}\\ =\frac{x与最小值的距离}{x与最大值的距离+x与最小值的距离} \]

注意:要区别开归一化和标准化。归一化的计算步骤也可以 消去量纲的影响,但更多时候,我们进行归一化的目的是为 了让我们的结果更容易解释,或者说让我们对结果有一个更 加清晰直观的印象。例如将得分归一化后可限制在0‐1这个区 间,对于区间内的每一个得分,我们很容易的得到其所处的比例位置。

可以使用层次分析法给这m个评价指标确定权重:

\[\sum_{j=1}^{m}w_{j}=1 \]

层次分析法的主观性太强了,更推荐大家使用熵权法来进行客观赋值

代码

topsis.m

%%  第一步:把数据复制到工作区,并将这个矩阵命名为X
% (1)在工作区右键,点击新建(Ctrl+N),输入变量名称为X
% (2)在Excel中复制数据,再回到Excel中右键,点击粘贴Excel数据(Ctrl+Shift+V)
% (3)关掉这个窗口,点击X变量,右键另存为,保存为mat文件(下次就不用复制粘贴了,只需使用load命令即可加载数据)
% (4)注意,代码和数据要放在同一个目录下哦,且Matlab的当前文件夹也要是这个目录。
clear;clc
load data_water_quality.mat
%% 注意:如果提示: 错误使用 load,无法读取文件 'data_water_quality.mat'。没有此类文件或目录。
% 那么原因是因为你的Matlab的当前文件夹中不存在这个文件
% 可以使用cd函数修改Matlab的当前文件夹
% 比如说,我的代码和数据放在了: D:第2讲.TOPSIS法(优劣解距离法)\代码和例题数据
% 那么我就可以输入命令:
% cd 'D:第2讲.TOPSIS法(优劣解距离法)\代码和例题数据'
% 也可以看我更新的视频:“更新9_Topsis代码为什么运行失败_得分结果怎么可视化以及权重的确定如何更加准确”,里面有介绍

%%  第二步:判断是否需要正向化
[n,m] = size(X);
disp(['共有' num2str(n) '个评价对象, ' num2str(m) '个评价指标']) 
Judge = input(['这' num2str(m) '个指标是否需要经过正向化处理,需要请输入1 ,不需要输入0:  ']);

if Judge == 1
    Position = input('请输入需要正向化处理的指标所在的列,例如第2、3、6三列需要处理,那么你需要输入[2,3,6]: '); %[2,3,4]
    disp('请输入需要处理的这些列的指标类型(1:极小型, 2:中间型, 3:区间型) ')
    Type = input('例如:第2列是极小型,第3列是区间型,第6列是中间型,就输入[1,3,2]:  '); %[2,1,3]
    % 注意,Position和Type是两个同维度的行向量
    for i = 1 : size(Position,2)  %这里需要对这些列分别处理,因此我们需要知道一共要处理的次数,即循环的次数
        X(:,Position(i)) = Positivization(X(:,Position(i)),Type(i),Position(i));
    % Positivization是我们自己定义的函数,其作用是进行正向化,其一共接收三个参数
    % 第一个参数是要正向化处理的那一列向量 X(:,Position(i))   回顾上一讲的知识,X(:,n)表示取第n列的全部元素
    % 第二个参数是对应的这一列的指标类型(1:极小型, 2:中间型, 3:区间型)
    % 第三个参数是告诉函数我们正在处理的是原始矩阵中的哪一列
    % 该函数有一个返回值,它返回正向化之后的指标,我们可以将其直接赋值给我们原始要处理的那一列向量
    end
    disp('正向化后的矩阵 X =  ')
    disp(X)
end
%% 作业:在这里增加是否需要算加权
% 补充一个基础知识:m*n维的矩阵A 点乘 n维行向量B,等于这个A的每一行都点乘B
% (注意:2017以及之后版本的Matlab才支持,老版本Matlab会报错)
% % 假如原始数据为:
%   A=[1, 2, 3;
%        2, 4, 6] 
% % 权重矩阵为:
%   B=[ 0.2, 0.5 ,0.3 ] 
% % 加权后为:
%   C=A .* B
%     0.2000    1.0000    0.9000
%     0.4000    2.0000    1.8000
% 类似的,还有矩阵和向量的点除, 大家可以自己试试计算A ./ B
% 注意,矩阵和向量没有 .- 和 .+ 哦 ,大家可以试试,如果计算A.+B 和 A.-B会报什么错误。

%% 这里补充一个小插曲
% % 在上一讲层次分析法的代码中,我们可以优化以下的语句:
% % Sum_A = sum(A);
% % SUM_A = repmat(Sum_A,n,1);
% % Stand_A = A ./ SUM_A;
% % 事实上,我们把第三行换成:Stand_A = A ./ Sum_A; 也是可以的哦 
% % (再次强调,新版本的Matlab才能运行哦)

%% 让用户判断是否需要增加权重
disp("请输入是否需要增加权重向量,需要输入1,不需要输入0")
Judge = input('请输入是否需要增加权重: ');
if Judge == 1
    disp(['如果你有3个指标,你就需要输入3个权重,例如它们分别为0.25,0.25,0.5, 则你需要输入[0.25,0.25,0.5]']);
    weigh = input(['你需要输入' num2str(m) '个权数。' '请以行向量的形式输入这' num2str(m) '个权重: ']);
    OK = 0;  % 用来判断用户的输入格式是否正确
    while OK == 0 
        if abs(sum(weigh) - 1)<0.000001 && size(weigh,1) == 1 && size(weigh,2) == m   % 这里要注意浮点数的运算是不精准的。
             OK =1;
        else
            weigh = input('你输入的有误,请重新输入权重行向量: ');
        end
    end
else
    weigh = ones(1,m) ./ m ; %如果不需要加权重就默认权重都相同,即都为1/m
end


%% 第三步:对正向化后的矩阵进行标准化
Z = X ./ repmat(sum(X.*X) .^ 0.5, n, 1);
disp('标准化矩阵 Z = ')
disp(Z)

%% 第四步:计算与最大值的距离和最小值的距离,并算出得分
D_P = sum([(Z - repmat(max(Z),n,1)) .^ 2 ] .* repmat(weigh,n,1) ,2) .^ 0.5;   % D+ 与最大值的距离向量
D_N = sum([(Z - repmat(min(Z),n,1)) .^ 2 ] .* repmat(weigh,n,1) ,2) .^ 0.5;   % D- 与最小值的距离向量
S = D_N ./ (D_P+D_N);    % 未归一化的得分
disp('最后的得分为:')
stand_S = S / sum(S)
[sorted_S,index] = sort(stand_S ,'descend')

% A = magic(5)  % 幻方矩阵
% M = magic(n)返回由1到n^2的整数构成并且总行数和总列数相等的n×n矩阵。阶次n必须为大于或等于3的标量。
% sort(A)若A是向量不管是列还是行向量,默认都是对A进行升序排列。sort(A)是默认的升序,而sort(A,'descend')是降序排序。
% sort(A)若A是矩阵,默认对A的各列进行升序排列
% sort(A,dim)
% dim=1时等效sort(A)
% dim=2时表示对A中的各行元素升序排列
% A = [2,1,3,8]
% Matlab中给一维向量排序是使用sort函数:sort(A),排序是按升序进行的,其中A为待排序的向量;
% 若欲保留排列前的索引,则可用 [sA,index] = sort(A,'descend') ,排序后,sA是排序好的向量,index是向量sA中对A的索引。
% sA  =  8     3     2     1
% index =  4     3     1     2

Inter2Max.m

function [posit_x] = Inter2Max(x,a,b)
    r_x = size(x,1);  % row of x 
    M = max([a-min(x),max(x)-b]);
    posit_x = zeros(r_x,1);   %zeros函数用法: zeros(3)  zeros(3,1)  ones(3)
    % 初始化posit_x全为0  初始化的目的是节省处理时间
    for i = 1: r_x
        if x(i) < a
           posit_x(i) = 1-(a-x(i))/M;
        elseif x(i) > b
           posit_x(i) = 1-(x(i)-b)/M;
        else
           posit_x(i) = 1;
        end
    end
end

Mid2Max.m

function [posit_x] = Mid2Max(x,best)
    M = max(abs(x-best));
    posit_x = 1 - abs(x-best) / M;
end

Min2Max.m

function [posit_x] = Min2Max(x)
    posit_x = max(x) - x;
     %posit_x = 1 ./ x;    %如果x全部都大于0,也可以这样正向化
end

Positivization.m

% function [输出变量] = 函数名称(输入变量)  
% 函数的中间部分都是函数体
% 函数的最后要用end结尾
% 输出变量和输入变量可以有多个,用逗号隔开
% function [a,b,c]=test(d,e,f)
%     a=d+e;
%     b=e+f;
%     c=f+d;
% end
% 自定义的函数要单独放在一个m文件中,不可以直接放在主函数里面(和其他大多数语言不同)

function [posit_x] = Positivization(x,type,i)
% 输入变量有三个:
% x:需要正向化处理的指标对应的原始列向量
% type: 指标的类型(1:极小型, 2:中间型, 3:区间型)
% i: 正在处理的是原始矩阵中的哪一列
% 输出变量posit_x表示:正向化后的列向量
    if type == 1  %极小型
        disp(['第' num2str(i) '列是极小型,正在正向化'] )
        posit_x = Min2Max(x);  %调用Min2Max函数来正向化
        disp(['第' num2str(i) '列极小型正向化处理完成'] )
        disp('~~~~~~~~~~~~~~~~~~~~分界线~~~~~~~~~~~~~~~~~~~~')
    elseif type == 2  %中间型
        disp(['第' num2str(i) '列是中间型'] )
        best = input('请输入最佳的那一个值: ');
        posit_x = Mid2Max(x,best);
        disp(['第' num2str(i) '列中间型正向化处理完成'] )
        disp('~~~~~~~~~~~~~~~~~~~~分界线~~~~~~~~~~~~~~~~~~~~')
    elseif type == 3  %区间型
        disp(['第' num2str(i) '列是区间型'] )
        a = input('请输入区间的下界: ');
        b = input('请输入区间的上界: '); 
        posit_x = Inter2Max(x,a,b);
        disp(['第' num2str(i) '列区间型正向化处理完成'] )
        disp('~~~~~~~~~~~~~~~~~~~~分界线~~~~~~~~~~~~~~~~~~~~')
    else
        disp('没有这种类型的指标,请检查Type向量中是否有除了1、2、3之外的其他值')
    end
end

结果

正向化后的矩阵 X =  
    4.6900    0.7172    3.0000    1.0000
    2.0300    0.4069   35.0000    0.6940
    9.1100    0.5241    8.0000    0.9058
    8.6100    0.9655    8.0000    0.4443
    7.1300    0.6552    4.0000    0.6914
    2.3900    0.8414   16.0000    0.6007
    7.6900    0.8552   16.0000    0.6551
    9.3000    0.8690   27.0000         0
    5.4500    0.5724   49.0000    1.0000
    6.1900    0.8138   37.0000    0.7848
    7.9300    0.6345   45.0000    0.6992
    4.4000    0.8069   37.0000    0.5419
    7.4600    0.1448   31.0000    1.0000
    2.0100         0    7.0000    0.4546
    2.0400    0.5862   31.0000    1.0000
    7.7300    0.4069    2.0000    1.0000
    6.3500    0.6000   29.0000    0.1824
    8.2900    0.0276   15.0000    1.0000
    3.5400    0.8138         0    0.4088
    7.4400    0.4897   46.0000    0.2731

请输入是否需要增加权重向量,需要输入1,不需要输入0
请输入是否需要增加权重: 1
如果你有3个指标,你就需要输入3个权重,例如它们分别为0.25,0.25,0.5, 则你需要输入[0.25,0.25,0.5]
你需要输入4个权数。请以行向量的形式输入这4个权重: [0.25,0.25,0.25,0.25]
标准化矩阵 Z = 
    0.1622    0.2483    0.0245    0.3065
    0.0702    0.1408    0.2863    0.2127
    0.3150    0.1814    0.0655    0.2776
    0.2977    0.3342    0.0655    0.1361
    0.2466    0.2268    0.0327    0.2119
    0.0826    0.2912    0.1309    0.1841
    0.2659    0.2960    0.1309    0.2008
    0.3216    0.3008    0.2209         0
    0.1885    0.1981    0.4009    0.3065
    0.2141    0.2817    0.3027    0.2405
    0.2742    0.2196    0.3682    0.2143
    0.1522    0.2793    0.3027    0.1661
    0.2580    0.0501    0.2536    0.3065
    0.0695         0    0.0573    0.1393
    0.0705    0.2029    0.2536    0.3065
    0.2673    0.1408    0.0164    0.3065
    0.2196    0.2077    0.2373    0.0559
    0.2867    0.0095    0.1227    0.3065
    0.1224    0.2817         0    0.1253
    0.2573    0.1695    0.3763    0.0837

最后的得分为:

stand_S =

    0.0451
    0.0478
    0.0485
    0.0488
    0.0431
    0.0448
    0.0539
    0.0510
    0.0681
    0.0684
    0.0702
    0.0591
    0.0527
    0.0192
    0.0533
    0.0434
    0.0466
    0.0438
    0.0358
    0.0565


sorted_S =

    0.0702
    0.0684
    0.0681
    0.0591
    0.0565
    0.0539
    0.0533
    0.0527
    0.0510
    0.0488
    0.0485
    0.0478
    0.0466
    0.0451
    0.0448
    0.0438
    0.0434
    0.0431
    0.0358
    0.0192


index =

    11
    10
     9
    12
    20
     7
    15
    13
     8
     4
     3
     2
    17
     1
     6
    18
    16
     5
    19
    14