数据包络分析（DEA）详解

数据包络分析（DEA）详解

一、基本介绍
1.1原理
数据包络分析有多种模型，主要为：CCR模型，BBC模型、交叉模型、A&P模型。纵观该方法的各种模型，每一模型的具体数学推理过程基本一致，所得的标准线性规划求解公式也比较相似，它们之间的差异主要体现在每个模型所适用的条件有所不同。在国外，此方法较早地运用于银行医院、城市等方面效率的评价。近些年来，该方法在我国社会经济的许多领域也取得了不少应用成果。由于方法对评价对象的要求相对比较宽松，应用其评价相同类型DMU的相对有效性的优势地位，是其他方法难以取代的。

1.2CCR模型
CCR模型是指对DMU通过“投入一定数量的生产要素，并产出一定数量的产品”的经济系统来判断各个单元的相对合理性和有效性。从投入资源的角度来看，在当前产出的水准下，比较投入资源的使用情况，以此作为效益评价的依据，这种模式称为“投入导向模式”。

定义决策单元j的效率评价指数为：

对上式可以适当的取权系数v和u，使得，对第个决策单元进行效率评价，一般来说越大表明能够用相对较少的输入而取得相对较多的输出。故只需求得的最大值，即可探究在这n个DUM中相对来说是不是最优的。

1.3BCC模型
BCC模型是从产出的角度探讨效率，即在相同的投入水准下，比较产出资源的达成情况，这种模式称为“投入导向模式”。所得到的是“技术效益”，DEA=1称为“技术有效”，最优解是决策单元j的“技术效益”。

二、代码
2.1MATLAB代码

clear
clc
format long
data=[14.40 0.65 31.30 3621.00 0.00
16.90 0.72 32.20 3943.00 0.09
15.53 0.72 31.87 4086.67 0.07
15.40 0.76 32.23 4904.67 0.13
14.17 0.76 32.40 6311.67 0.37
13.33 0.69 30.77 8173.33 0.59
12.83 0.61 29.23 10236.00 0.51
13.00 0.63 28.20 12094.33 0.44
13.40 0.75 28.80 13603.33 0.58
14.00 0.84 29.10 14841.00 1.00]';
 
X=data([1:3],:);%X为输入变量
Y=data([4:5],:);%Y为输出变量
[m,n]=size(X);
s=size(Y,1);
A=[-X' Y'];%由于目标函数求最小，这里的-X就转化成了求最大
b=zeros(n,1);
LB=zeros(m+s,1);UB=[];
for i=1:n
   f=[zeros(1,m) -Y(:,i)'];
   Aeq=[X(:,i)',zeros(1,s)];
   beq=1;
   w(:,i)=linprog(f,A,b,Aeq,beq,LB,UB);%前3列为投入系数，后2列为产出系数
   E(i,i)=Y(:,i)'*w(m+1:m+s,i);%产出值*产出系数
end
theta=diag(E)';
fprintf('用DEA方法对此的相对评价结果为：\n');
disp(theta);

2.2Python代码

import gurobipy
import pandas as pd
 
# 分页显示数据, 设置为 False 不允许分页
pd.set_option('display.expand_frame_repr', False)
 
# 最多显示的列数, 设置为 None 显示全部列
pd.set_option('display.max_columns', None)
 
# 最多显示的行数, 设置为 None 显示全部行
pd.set_option('display.max_rows', None)
 
class DEA(object):
   def __init__(self, DMUs_Name, X, Y, AP=False):
      self.m1, self.m1_name, self.m2, self.m2_name, self.AP = X.shape[1], X.columns.tolist(), Y.shape[1], Y.columns.tolist(), AP
      self.DMUs, self.X, self.Y = gurobipy.multidict({DMU: [X.loc[DMU].tolist(), Y.loc[DMU].tolist()] for DMU in DMUs_Name})
      print(f'DEA(AP={AP}) MODEL RUNING...')
 
   def __CCR(self):
      for k in self.DMUs:
         MODEL = gurobipy.Model()
         OE, lambdas, s_negitive, s_positive = MODEL.addVar(), MODEL.addVars(self.DMUs),  MODEL.addVars(self.m1), MODEL.addVars(self.m2)
         MODEL.update()
         MODEL.setObjectiveN(OE, index=0, priority=1)
         MODEL.setObjectiveN(-(sum(s_negitive) + sum(s_positive)), index=1, priority=0)
         MODEL.addConstrs(gurobipy.quicksum(lambdas[i] * self.X[i][j] for i in self.DMUs if i != k or not self.AP) + s_negitive[j] == OE * self.X[k][j] for j in range(self.m1))
         MODEL.addConstrs(gurobipy.quicksum(lambdas[i] * self.Y[i][j] for i in self.DMUs if i != k or not self.AP) - s_positive[j] == self.Y[k][j] for j in range(self.m2))
         MODEL.setParam('OutputFlag', 0)
         MODEL.optimize()
         self.Result.at[k, ('效益分析', '综合技术效益(CCR)')] = MODEL.objVal
         self.Result.at[k, ('规模报酬分析', '有效性')] = '非 DEA 有效' if MODEL.objVal < 1 else 'DEA 弱有效' if s_negitive.sum().getValue() + s_positive.sum().getValue() else 'DEA 强有效'
         self.Result.at[k, ('规模报酬分析', '类型')] = '规模报酬固定' if lambdas.sum().getValue() == 1 else '规模报酬递增' if lambdas.sum().getValue() < 1 else '规模报酬递减'
         for m in range(self.m1):
            self.Result.at[k, ('差额变数分析', f'{self.m1_name[m]}')] = s_negitive[m].X
            self.Result.at[k, ('投入冗余率',  f'{self.m1_name[m]}')] = 'N/A' if self.X[k][m] == 0 else s_negitive[m].X / self.X[k][m]
         for m in range(self.m2):
            self.Result.at[k, ('差额变数分析', f'{self.m2_name[m]}')] = s_positive[m].X
            self.Result.at[k, ('产出不足率', f'{self.m2_name[m]}')] = 'N/A' if self.Y[k][m] == 0 else s_positive[m].X / self.Y[k][m]
      return self.Result
 
   def __BCC(self):
      for k in self.DMUs:
         MODEL = gurobipy.Model()
         TE, lambdas = MODEL.addVar(), MODEL.addVars(self.DMUs)
         MODEL.update()
         MODEL.setObjective(TE, sense=gurobipy.GRB.MINIMIZE)
         MODEL.addConstrs(gurobipy.quicksum(lambdas[i] * self.X[i][j] for i in self.DMUs if i != k or not self.AP) <= TE * self.X[k][j] for j in range(self.m1))
         MODEL.addConstrs(gurobipy.quicksum(lambdas[i] * self.Y[i][j] for i in self.DMUs if i != k or not self.AP) >= self.Y[k][j] for j in range(self.m2))
         MODEL.addConstr(gurobipy.quicksum(lambdas[i] for i in self.DMUs if i != k or not self.AP) == 1)
         MODEL.setParam('OutputFlag', 0)
         MODEL.optimize()
         self.Result.at[k, ('效益分析', '技术效益(BCC)')] = MODEL.objVal if MODEL.status == gurobipy.GRB.Status.OPTIMAL else 'N/A'
      return self.Result
 
   def dea(self):
      columns_Page = ['效益分析'] * 3 + ['规模报酬分析'] * 2 + ['差额变数分析'] * (self.m1 + self.m2) + ['投入冗余率'] * self.m1 + ['产出不足率'] * self.m2
      columns_Group = ['技术效益(BCC)', '规模效益(CCR/BCC)', '综合技术效益(CCR)','有效性', '类型'] + (self.m1_name + self.m2_name) * 2
      self.Result = pd.DataFrame(index=self.DMUs, columns=[columns_Page, columns_Group])
      self.__CCR()
      self.__BCC()
      self.Result.loc[:, ('效益分析', '规模效益(CCR/BCC)')] = self.Result.loc[:, ('效益分析', '综合技术效益(CCR)')] / self.Result.loc[:,('效益分析', '技术效益(BCC)')]
      return self.Result
 
   def analysis(self, file_name=None):
      Result = self.dea()
      file_name = 'DEA 数据包络分析报告.xlsx' if file_name is None else f'\\{file_name}.xlsx'
      Result.to_excel(file_name, 'DEA 数据包络分析报告')


————————————————
版权声明：本文为CSDN博主「饲养猿」的原创文章，遵循CC 4.0 BY-SA版权协议，转载请附上原文出处链接及本声明。
原文链接：https://blog.csdn.net/qq_48774513/article/details/120198871

Post Views: 7,632
Tags: , , ,