保姆级教程:用Python的SciPy库搞定超效率SBM模型(含非期望产出处理)
Python实战:超效率SBM模型的技术实现与商业洞察
当你在分析制造业碳排放效率时,是否遇到过这样的困境:明明A工厂和B工厂的DEA效率值都是1,却无法判断谁更优秀?或者在评估银行分支机构绩效时,传统模型无法处理"坏账率"这类非期望产出?这就是超效率SBM模型要解决的核心问题。
1. 环境配置与数据准备
在开始建模前,我们需要搭建一个稳定的Python环境。推荐使用Anaconda创建独立环境,避免依赖冲突:
conda create -n dea python=3.9 conda activate dea pip install numpy scipy pandas matplotlib典型的数据结构应该包含三类指标:
- 投入指标:如劳动力人数、资本投入、能源消耗
- 期望产出:如产值、利润
- 非期望产出:如污染物排放量、坏账率
用Pandas构建数据框架示例:
import pandas as pd data = { 'DMU': ['Factory_A', 'Factory_B', 'Factory_C'], 'Labor': [20, 25, 18], 'Capital': [500, 600, 450], 'Energy': [300, 350, 280], 'Output': [1000, 1200, 900], 'CO2': [200, 180, 210] } df = pd.DataFrame(data).set_index('DMU') inputs = df[['Labor', 'Capital', 'Energy']].values outputs = df[['Output', 'CO2']].values # CO2作为非期望产出注意:非期望产出在建模时需要特殊处理,通常有两种方式:
- 作为负向指标直接输入
- 在目标函数中设置方向系数
2. 超效率SBM模型核心算法
传统DEA模型的最大局限是效率值被限制在[0,1]区间,无法区分多个有效DMU之间的差异。超效率SBM通过松弛变量和前沿面调整解决了这个问题。
模型数学表达的核心思想:
min ρ = (1 - 1/m ∑(s_i^- /x_ik)) / (1 + 1/s ∑(s_r^+ /y_rk)) s.t. Xλ + s^- = x_k Yλ - s^+ = y_k λ ≥ 0, s^- ≥ 0, s^+ ≥ 0Python实现关键步骤:
from scipy.optimize import linprog import numpy as np def super_sbm(inputs, outputs, bad_outputs=None, crs=True): num_dmu, num_input = inputs.shape num_output = outputs.shape[1] efficiencies = [] for k in range(num_dmu): # 构建目标函数系数 c = np.zeros(num_input + num_output + num_dmu) c[:num_input] = 1/(num_input * inputs[k]) # 构建约束矩阵 A_eq = [] b_eq = [] # 投入约束 for i in range(num_input): row = np.zeros(num_input + num_output + num_dmu) row[i] = 1 row[num_input+num_output:] = -inputs[:,i] A_eq.append(row) b_eq.append(0) # 产出约束(处理非期望产出) for r in range(num_output): row = np.zeros(num_input + num_output + num_dmu) row[num_input+r] = 1 if bad_outputs and r in bad_outputs: row[num_input+num_output:] = outputs[:,r] else: row[num_input+num_output:] = -outputs[:,r] A_eq.append(row) b_eq.append(outputs[k,r] if (bad_outputs and r in bad_outputs) else 0) # 规模报酬约束 if crs: row = np.zeros(num_input + num_output + num_dmu) row[num_input+num_output:] = 1 A_eq.append(row) b_eq.append(1) # 求解线性规划 res = linprog(c, A_eq=A_eq, b_eq=b_eq, bounds=(0, None)) efficiency = 1/res.fun if res.success else np.nan efficiencies.append(efficiency) return efficiencies参数说明:
bad_outputs: 指定哪些产出指标是非期望产出(如[1]表示第二个产出)crs: 规模报酬不变(True)或可变(False)
3. 商业场景应用案例
3.1 制造业绿色效率评估
假设我们评估5家工厂的绿色生产效率,数据如下:
| 工厂 | 劳动力(人) | 资本(万元) | 能耗(吨标煤) | 产值(万元) | 废水(吨) |
|---|---|---|---|---|---|
| A | 50 | 200 | 300 | 1000 | 50 |
| B | 60 | 250 | 350 | 1200 | 45 |
| C | 40 | 180 | 280 | 900 | 55 |
实现代码:
inputs = np.array([[50,200,300], [60,250,350], [40,180,280]]) outputs = np.array([[1000,50], [1200,45], [900,55]]) # 废水作为非期望产出(第二个产出指标) eff_scores = super_sbm(inputs, outputs, bad_outputs=[1], crs=True) for dmu, score in zip(['A','B','C'], eff_scores): print(f"工厂{dmu}效率值:{score:.4f}")典型输出结果分析:
工厂A效率值:1.1254 工厂B效率值:1.0832 工厂C效率值:0.9765解读:
- 效率值>1表示超高效单位
- 工厂C效率值<1表示未达到前沿面
- 虽然工厂B产值最高,但A的综合效率更优
3.2 银行分支机构绩效评估
考虑三种产出指标:
- 期望产出:存款额、贷款额
- 非期望产出:不良贷款率
# 示例数据:员工数, 运营成本, 存款额, 贷款额, 不良率(%) bank_data = np.array([ [15, 80, 500, 300, 1.2], [20, 100, 600, 400, 1.5], [12, 70, 450, 250, 0.8] ]) inputs = bank_data[:,:2] outputs = bank_data[:,2:] # 不良率作为非期望产出(第三个产出指标) bank_eff = super_sbm(inputs, outputs, bad_outputs=[2], crs=False)4. 结果可视化与商业决策
效率值的可视化能更直观展示分析结果:
import matplotlib.pyplot as plt plt.figure(figsize=(10,6)) plt.bar(['A','B','C'], eff_scores, color=['green' if x>=1 else 'red' for x in eff_scores]) plt.axhline(1, color='gray', linestyle='--') plt.title('超效率SBM评估结果') plt.ylabel('效率值') plt.xlabel('决策单元') for i, v in enumerate(eff_scores): plt.text(i, v+0.02, f"{v:.3f}", ha='center') plt.show()商业决策建议:
- 标杆学习:识别效率值>1的超高效单位作为学习标杆
- 资源分配:对低效单位进行重点投入或流程优化
- 绩效激励:将效率值纳入KPI考核体系
- 战略调整:分析效率驱动因素,优化业务组合
5. 常见问题解决方案
Q1:如何处理量纲不一致的指标?
推荐标准化处理方法:
from sklearn.preprocessing import MinMaxScaler scaler = MinMaxScaler() scaled_inputs = scaler.fit_transform(inputs) scaled_outputs = scaler.fit_transform(outputs)Q2:模型结果不稳定怎么办?
- 检查数据是否有极端值
- 尝试不同的规模报酬假设(CRS/VRS)
- 增加DMU数量(一般不少于指标数的3倍)
Q3:如何解释效率值为负数?
这通常发生在:
- 数据输入错误(如将产出误作为投入)
- 非期望产出方向设置错误
- 线性规划无可行解
修正方法:
# 确保非期望产出方向正确 if bad_outputs: outputs[:,bad_outputs] = -outputs[:,bad_outputs]Q4:超效率值超过多少算异常?
经验阈值:
- 正常范围:0.5-2.0
- 警告范围:<0.3 或 >5.0
- 需要检查:>10.0
效率分析的实际价值在于相对比较而非绝对数值,关键是通过模型识别改进方向和标杆实践。
