2022-07-03-含有非期望产出的sbm模型python代码

传统的径向DEA模型无法考虑“松弛变量”对效率值的影响，也没有考虑同时使期望产出增加，非期望产出减少的技术变化，以此度量的效率值是不准确或有偏的，为了解决这一问题，Tone(2001)提出了基于投入产出松弛变量的环境效率评价模型，简称SBM模型，在此基础上，他进一步提出了SBM的拓展模型，从而实现了非期望产出条件下对环境效率的评价，含有非期望产出的SBM模型规划式如下图所示：

关于求解含有非期望产出的SBM模型的python代码如下：

from scipy.optimize import minimize
import numpy as np
import pandas as pd
import scipy.optimize as op
def sbmeff2(input_variable, desirable_output, undesirable_output, dmu,data,method = 'revised simplex'):
    """用于求解sbm模型
    
    Parameters:
    -----------
    input_variable:
        投入[v1,v2,v3,...] 
    desirable_output:
        期望产出[v1,v2,v3,...]
    undesirable_output:
        非期望产出[v1,v2,v3,...] 
    dmu:
        决策单元
    data:
        主数据
    method:
        求解方法.默认'revised simplex',可选'interior-point'

	Return:
	------
	res : DataFrame
		结果数据框[dmu	TE	slack...]
    """
 
    res = pd.DataFrame(columns = ['dmu','TE'], index = data.index)
    res['dmu'] = data[dmu]
    ## lambda有dmu个数个，S有变量个数个
    dmu_counts = data.shape[0]
     ## 投入个数
    m = len(input_variable)
    ## 期望产出个数
    s1 = len(desirable_output)
    ## 非期望产出个数
    s2 = len(undesirable_output)
    ## x[:dmu_counts] 为lambda
    ## x[dmu_counts:dmu_counts+1] 为 t
    ## x[dmu_counts+1 :dmu_counts + m + 1] 为投入slack
    ## x[dmu_counts+ 1 + m:dmu_counts + 1 + m + s1] 为期望产出slack
    ## x[dmu_counts + 1 + m + s1 :] 为非期望产出lack
    total = dmu_counts + m + s1 + s2 + 1
    cols = input_variable+desirable_output+ undesirable_output
    newcols = []
    for j in cols:
        newcols.append(j+'_slack')
        res[j+'_slack'] = np.nan
    for i in range(dmu_counts):
        
        ## 优化目标
        c = [0] * dmu_counts + [1] +  list(-1 / (m * data.loc[i, input_variable])) + [0] * (s1 + s2)
        
        ## 约束条件
        A_eq = [[0] * dmu_counts + [1] + [0] * m  + list(1/((s1 + s2) * data.loc[i, desirable_output])) + 
                                                    list(1/((s1 + s2) * data.loc[i, undesirable_output]))]
        
        ## 约束条件（1）
        for j1 in range(m):
            list1 = [0] * m
            list1[j1] = 1
            eq1 = list(data[input_variable[j1]]) + [-data.loc[i ,input_variable[j1]]] + list1 + [0] * (s1 + s2)
            A_eq.append(eq1)
        ## 约束条件（2）
        for j2 in range(s1):
            list2 = [0] * s1
            list2[j2] = -1
            eq2 = list(data[desirable_output[j2]]) + [-data.loc[i, desirable_output[j2]]] + [0] * m + list2 + [0] * s2
            A_eq.append(eq2)
        ## 约束条件（3）
        for j3 in range(s2):
            list3 = [0] * s2
            list3[j3] = 1
            eq3 = list(data[undesirable_output[j3]]) + [-data.loc[i, undesirable_output[j3]]] + [0] * (m + s1) + list3
            A_eq.append(eq3)         
        b_eq = [1] + [0] * (m + s1 + s2)               
        bounds = [(0, None)]*total
        ## 求解
        op1 = op.linprog(c = c,A_eq=A_eq,b_eq=b_eq,bounds=bounds,method = method)
        res.loc[i, 'TE'] = op1.fun
        res.loc[i, newcols] = op1.x[dmu_counts+1 :]
    return res

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

代码来自于CSDN的一位博主wonder1322，非常感谢他的分享。