花了两个多星期,终于是把杨德平老师这本《经济预测与决策技术》学完了。。。有些分析我自己比较熟悉,就没有编写方法代码。。。明天起,开始学空间分析,争取能够实现经济学和地理学对空间分析的嵌入和整合。

本次发布的是“考虑非期望产出的DEA SBM模型”,这个模型的理论基础建立在Tone(2001)和钱争鸣(2013)两篇文章上,其逻辑在于基于DEA_CCR_Dual_Model估算出相关参数,然后结合Tone(2001)的算法,转化成DEA_SBM的相关参数。

代码的运行效果图如下:

原始代码自如下开始:

% 本代码由重庆科技学院范巧副教授于2018年7月28日23:50撰写;

% 本代码版权保护,不得用于科研或学习之外的其他用途;

% 小范学数量经济学之九:考虑非期望产出的DEA_SBM模型MATLAB模拟;

% 投入指标:劳动力、资本、综合资源投入水平;

% 期望产出指标:城市实际GDP;

% 非期望产出指标:工业废水排放量;工业SO2排放量;工业烟粉尘排放量;

% 决策单元个数:285个地级市;

% 参考文献一:Tone K. A slacks-based measure of efficiency in data envelopment analysis[J]. European Journal of Operational Research,2001,130(3):498-509.

% 参考文献二:钱争鸣,刘晓晨.中国绿色经济效率的区域差异与影响因素分析[J].中国人口·资源与环境,2013,23(07):104-109.

% 本代码基于软件MATLAB R2018a撰写。

% 录入投入指标、期望产出指标和非期望产出指标;

X0 = randn(m0,n);

Y_E0 = randn(m1,n);

Y_NE0 = randn(m2,n);

X = abs(X0); % 投入指标;

Y_E = abs(Y_E0); % 期望产出指标;

Y_NE = abs(Y_NE0); % 非期望产出指标;

J1 = size(X,2);

J2 = size(Y_E,2);

J3 = size(Y_NE,2);

% 检查决策单元的一致性;

if J1 == J2, J2 == J3;

('The decision making units are equals between inputs and outputs')

% 确定投入指出和产出指标的维数;

m = size(X,1);

Y = [Y_E;Y_NE];

s = size(Y,1);

s1 = size(Y_E,1);

% 确定线性规划的约束条件;

f = [zeros(1,J+m+s),1];

LB = zeros(J+m+s+1,1);

% 求解线性规划的最优解;

for j = 1:J;

Aeq = [X,eye(m),zeros(m,s),-X(:,j);

Y,zeros(s,m),-eye(s),zeros(s,1)];

beq = [zeros(m,1);Y(:,j)];

w = linprog(f,A,b,Aeq,beq,LB,UB);

% linprog(f,A,b) solves min f'*x such that A*x ≤ b.

% f, x, b, beq, lb, and ub are vectors, and A and Aeq are matrices.

% linprog defines a set of lower and upper bounds on the design variables x;

% the solution is always in the range lb ≤ x ≤ ub;

% Set Aeq = [] and beq = [] if no equalities exist.

% 确定DEA_CCR模型的相关参数;

Lamda_CCR = W(1:J,:);

s_minus_CCR = W(J+1:J+m,:);

s_plus_CCR = W(J+m+1:J+m+s,:);

Theta_CCR = W(J+m+s+1,:);

% 基于Tone (2001)文献的CCR向SBM转化的方法,确定DEA_SBM的相关参数;

Lamda_SBM = Lamda_CCR;

s_plus_SBM = s_plus_CCR;

Theta_CCR0 = 1./Theta_CCR;

x0 = Theta_CCR0.*(X*Lamda_CCR+s_minus_CCR);

y0 = Y*Lamda_CCR-s_plus_CCR;

s_minus_SBM = s_minus_CCR+(1-Theta_CCR).*x0;

% DEA_SBM_Dual_Model is efficient while:

% 1. the value of theta equals to 1;

% 2. Corresponding column elements in matrix s_minus and s_plus equal to 0.

% 3. Corresponding elements on the diag of Lamda equal to 1.

% 基于DEA_SBM模型,计算DEA效率评价值,包含全部产出、期望产出和非期望产出;

rho0 = Theta_CCR-mean(s_minus_CCR./x0);

rho1 = mean(s_plus_CCR./y0);

rho2 = mean(s_plus_CCR(1:m1,:)./y0(1:m1,:));

rho3 = mean(s_plus_CCR(m1+1:m1+m2,:)./y0(m1+1:m1+m2,:));

rho_total = rho0./(1+rho1); % 基于全部产出的DEA效率评价值;

rho_expected =  rho0./(1+rho2); % 基于期望产出的DEA效率评价值;

rho_non_expected = rho0./(1+rho3); % 基于非期望产出的DEA效率评价值;

[M10,N10] = sort(rho_total,'descend');

[M1,N1] = sort(N10);

[M20,N20] = sort(rho_expected,'descend');

[M2,N2] = sort(N20);

[M30,N30] = sort(rho_non_expected,'descend');

[M3,N3] =sort(N30);

plot(rho_total,'-p');

plot(rho_expected,'-o');

plot(rho_non_expected,'-*');

axis([1,J,-0.00001,1.00001])

xlabel('评价单元')

ylabel('效率评价值降序排序')

legend('基于全部产出的DEA效率评价值降序排序','基于期望产出的DEA效率评价值降序排序','基于非期望产出的DEA效率评价值降序排序')