先说个背景。我最早接触SBM模型的时候是在帮一家制造业公司做供应链效率评估。当时他们只有十几个决策单元DMU数据零散得不行指标口径还乱Excel里来回折腾了三天最后连“效率值怎么算出来的”都讲不清楚。后来我索性把整套流程迁到Python上pandas清洗、scipy求解、matplotlib可视化一步到位整个项目从接手到交付不到两周。这篇博文就是把那套流程完整拆开从数据清洗到SBM模型建模再到结果可视化每一步都给出可以照着跑的方案适合正在做效率评价、绩效评估、DEA分析但又不想只会用现成软件包的读者。1. 为什么是SBM模型而不是传统DEA1.1 径向模型查不准问题的根源DEA数据包络分析里最经典的是CCR和BCC两个模型它们属于径向模型。径向是什么意思就是假设所有投入或者所有产出都按同一比例缩放。比如某家银行的效率值是0.8径向模型的意思是这家银行的所有投入都等比例压到80%理论上就能跑到生产前沿面上。这个假设在实际业务里非常生硬。现实中的低效往往不是“全面低效”而是某一项投入严重冗余、另一项产出严重不足。比如一个物流仓库人力资源利用得很好但场地面积利用率极低。用BBC模型算下来效率值是0.85你根本不知道问题出在哪块场地还是哪块人力。径向模型没有给出松弛变量slack的信息相当于只告诉你“你离标杆差15%”但不告诉“你哪块比别人多花了钱”。1.2 SBM模型是怎么用松弛变量打准七寸的Tone在2001年提出的SBM模型Slacks-Based Measure基于松弛变量的效率测度解决了这个问题。SBM不再强求投入产出等比缩放而是直接把各项投入的冗余量s⁻和各项产出的不足量s⁺加权后放进了目标函数。效率值变成一个非径向的比值效率值 (1 − (1/m) × Σ(sᵢ⁻ / xᵢ₀)) / (1 (1/s) × Σ(sᵣ⁺ / yᵣ₀))其中m是投入个数s是产出个数xᵢ₀是第i项投入的实际值yᵣ₀是第r项产出的实际值。你可以把它理解成分子是“投入端平均节省比例”分母是“产出端平均增加比例”两者一除就能同时反映投入过多和产出不足两个维度的低效程度。这个设计带来的直接好处是效率值能被拆开看。一家公司效率值低你可以单独看是哪项投入的s⁻最大或者哪项产出的s⁺最大。这个信息在管理优化里非常重要——整改措施可以直接落到具体指标上。1.3 超效率SBM解决的痛点有效单元如何再排序经典SBM有个小尴尬效率值为1的DMU决策单元只能被判定为有效但多个有效DMU彼此之间谁更优分不出来。做效率排名的时候这是大问题——前几名全是1并列第一没法用。Tone在2002年顺势提出了Super-SBM超效率SBM。它的核心思路是计算某个有效DMU时把这个DMU从参考集里剔除让它去跟自己之外的样本对比。这样一来有效DMU的效率值可以突破1能继续排序。比如A和B效率都是1A剔除自己后效率值1.25B剔除自己后效率值1.08那A就排在B前面。这个思路跟“考试满分的人再多做附加题”是一个逻辑附加题分数越高说明原始成绩含金量越高。结合标题里的“超效率SBM模型”我这里会把经典SBM和超效率SBM一起讲因为超效率SBM依赖经典SBM的线性规划框架两者在代码实现上高度共享。2. 环境准备与数据清洗建模前最重要的工程环节2.1 搭建Python环境与核心依赖库如果你电脑上还没装Python直接从Python官网下载3.9以上版本装上就行安装时记得勾选“Add Python to PATH”。这一步不复杂但很多新手卡在环境变量上后面pip install统统报“不是内部或外部命令”基本都是这一步漏了。装好Python之后建议用虚拟环境隔离项目依赖实测下来最省心python -m venv sbm_env source sbm_env/bin/activate # Windows下是 sbm_env\Scripts\activate pip install pandas numpy scipy matplotlib openpyxl为什么要专门装这些库我逐个说明一下pandas负责数据读取、清洗、透视赛博表格numpy线性规划求解时处理矩阵和向量的基础scipy提供linprog线性规划求解器后面SBM和超效率SBM的核心数学模型就靠它解matplotlib输出效率分布、投影分析图openpyxl让pandas能读写Excel文件。你可能会问为什么不直接装一个DEA包比如dea或pyDEA我个人的经验是现成的DEA包在简单场景下确实方便但一遇到自定义约束比如非期望产出、面板数据、VRTS假设就很容易卡在包的限制上。自己用linprog搭一套SBM框架不过几十行代码后面想怎么改造就怎么改造。2.2 数据清洗的标准流程与细节SBM模型对数据质量非常敏感。数据里有缺失值、异常值或者量纲不在同一量级效率结果会变得很不稳定。我见过有人拿一份有一半缺失值的数据直接跑DEA跑出来40%的DMU效率值都是1结果完全没有区分度。建模前的数据清洗至少要完成下面四件事。第一步读入数据并做结构检查。项目开始的时候先用pandas把数据读进来确认每一列的类型和缺失情况import pandas as pd df pd.read_excel(sbm_data.xlsx, index_col0) print(df.info()) print(df.describe().T)这一步能看到什么df.info()告诉你每一列是不是数值类型有没有空值describe()告诉你每一列的均值、标准差、最小值、最大值方便第一时间发现“某项投入有的企业填了0”这种明显不合理的记录。第二步处理缺失值。SBM模型里任何一个DMU只要有一个指标缺失这个DMU就进不了线性规划求解。因为约束方程要求所有值非负且齐全。缺失值的处理方式有三种我按优先顺序排列如果缺失占比小于5%优先用同类DMU的中位数填充。中位数比均值稳健不容易被个别极端值带偏。如果缺失占比在5%到20%之间需要回溯原始业务口径去补确实补不了的建议删除该DMU。如果缺失占比超过20%这个指标本身就不可靠建议直接删除该指标别硬留。第三步剔除异常值。注意DEA模型里的“异常值”不能粗暴地用“超过3σ就删”这种统计口径判断。因为效率分析里“极好”和“极差”的样本恰恰会影响生产前沿面的形状。真正的异常值应该是指标口径错误比如某制造企业的“员工人数”填了50万人但营收只有5万元这明显是单位填错了或者多打了一个零。这类错误通过describe()看到最大最小值时基本能一眼识别处理方式就是修正或者删除但要记录删除原因保留审计轨迹。第四步检查指标方向性。SBM模型默认所有投入越小越好、所有产出越大越好。如果你的指标体系里有“不良率”“能耗”“投诉量”这类“越小越好”的产出指标就不能直接放进产出矩阵必须先做正向化处理。常用的正向化方式有两种取倒数不良率正向化 1 / 不良率。好处是简单但要注意不良率为0时会导致无穷大取最大值减当前值正向值 max(不良率) - 不良率。更稳但压缩了数据的变异程度。我在实际项目里常用第二种。因为取倒数会把分布变得极度偏斜影响线性规划求解的数值稳定性。2.3 量纲处理SBM模型到底要不要归一化先给结论SBM模型本身对量纲不敏感因为目标函数里sᵢ⁻ / xᵢ₀和sᵣ⁺ / yᵣ₀都是比值形式投入如果是“万元”分子分母同时带“万元”量纲抵消了。但这不代表你可以完全不做标准化。我踩过的坑是当一个指标的数值范围是0到1另一个指标的范围是几百万虽然数学上量纲能抵消但线性规划求解器在计算过程中会有数值精度问题。特别是样本量大的时候linprog默认的容差设置可能导致结果不稳定同一份数据不同求解顺序出来的效率值在第三位小数上有差异。所以我的习惯是建模前把投入指标和产出指标各自做min-max标准化统一压缩到0.1到1之间。这里为什么用0.1而不是0因为纯0值在比值计算中容易引起除零而SBM模型默认所有投入产出都为正数纯0会让线性规划退化。标准化公式很简单def minmax_positive(series): return (series - series.min()) / (series.max() - series.min()) * 0.9 0.1 for col in input_cols output_cols: df[col] minmax_positive(df[col])压缩到0.1到1之后所有指标数值保持在同一个数量级求解器跑起来非常稳也方便后面做投影分析对比。不过要专门提醒一句如果你需要把效率结果用于外部报告最好保留原始量纲做解释标准化后的数值只用于建模计算。3. SBM模型的线性规划化与Python实现3.1 把数学模型变成linprog认识的样子scipy.optimize.linprog只能解“最小化 cᵀx”形式的线性规划约束条件支持不等式Ax ≤ b和等式Aeq·x beq。SBM模型要硬套进去需要做一点数学改造。经典SBM的效率值求解原始形式是最小化 ρ (1 − (1/m) × Σ(sᵢ⁻ / xᵢ₀)) / (1 (1/s) × Σ(sᵣ⁺ / yᵣ₀))约束条件 Σⱼ λⱼ xᵢⱼ sᵢ⁻ xᵢ₀i 1,…,m Σⱼ λⱼ yᵣⱼ − sᵣ⁺ yᵣ₀r 1,…,s λⱼ ≥ 0, sᵢ⁻ ≥ 0, sᵣ⁺ ≥ 0这个目标函数是分式linprog不能直接处理。好在这个模型结构上是分式规划可以做一个Charnes-Cooper变换把分母消掉。思路不算复杂但很多人第一次上手时容易卡在变换细节上。考虑到大多数SBM应用的决策单元数量不超过50个我更推荐换一种更直白的思路直接对SBM的原始模型用“枚举松弛变量权重”的简化方法。具体来说对每个DMU₀求解效率值本质上是在找一组权重λ使得目标DMU对比前沿面时总松弛最小。经典SBM经过代数推导可以等价转化成下面这个线性规划这里给出我常用来直接编码的形式最小化 t − (1/m) × Σ(sᵢ⁻ / xᵢ₀) 约束条件 t (1/s) × Σ(sᵣ⁺ / yᵣ₀) 1 Σⱼ μⱼ xᵢⱼ sᵢ⁻ t · xᵢ₀ Σⱼ μⱼ yᵣⱼ − sᵣ⁺ t · yᵣ₀ μⱼ ≥ 0, sᵢ⁻ ≥ 0, sᵣ⁺ ≥ 0, t ≥ 0其中μⱼ t · λⱼ。这个变换的好处在于目标函数和约束全部线性化直接喂给linprog就行。其中变量向量的排列顺序是μⱼn个、sᵢ⁻m个、sᵣ⁺s个、t1个。如果你觉得这个变换稍微有点绕我提供另一个更省事的思路把所有DMU按效率值循环求解每次求解一个DMU时如果它的最优目标值非常接近1直接标记为有效对于小于1的DMU通过投影公式计算改进方向。实际业务分析中这个思路基本够用而且代码逻辑更容易读。3.2 核心代码经典SBM效率值计算下面这段代码是我项目里在用的精简版直接把上面的线性规划映射到linprog的接口上。你用的时候只需要替换X和Y两个矩阵X是投入矩阵行是DMU列是投入项Y是产出矩阵行是DMU列是产出项。import numpy as np from scipy.optimize import linprog def sbm_efficiency(X, Y, dm_index0): 计算第dm_index个DMU的SBM效率值 X: (n_samples, n_inputs) Y: (n_samples, n_outputs) 返回: 效率值, 投入松弛, 产出松弛 n X.shape[0] m X.shape[1] s Y.shape[1] x0 X[dm_index] y0 Y[dm_index] # 变量顺序: mu(长度n), sminus(长度m), splus(长度s), t(长度1) n_vars n m s 1 mu_idx slice(0, n) sminus_idx slice(n, n m) splus_idx slice(n m, n m s) t_idx n m s # 目标函数系数 c [0...0, -1/(m*x0[0]), ..., -1/(m*x0[m-1]), 0...0, 1] c np.zeros(n_vars) for i in range(m): c[n i] -1.0 / (m * x0[i]) c[t_idx] 1.0 # 等式约束 A_eq x b_eq A_eq [] b_eq [] # 约束1: t (1/s) * sum(splus / y0) 1 row np.zeros(n_vars) row[t_idx] 1.0 for r in range(s): row[n m r] 1.0 / (s * y0[r]) A_eq.append(row) b_eq.append(1.0) # 约束2: sum(mu_j * X[j]) sminus t * x0 for i in range(m): row np.zeros(n_vars) for j in range(n): row[mu_idx.start j] X[j, i] row[sminus_idx.start i] 1.0 row[t_idx] -x0[i] A_eq.append(row) b_eq.append(0.0) # 约束3: sum(mu_j * Y[j]) - splus t * y0 for r in range(s): row np.zeros(n_vars) for j in range(n): row[mu_idx.start j] Y[j, r] row[splus_idx.start r] -1.0 row[t_idx] -y0[r] A_eq.append(row) b_eq.append(0.0) # 变量边界: 全部 0 bounds [(0, None)] * n_vars result linprog(c, A_eqnp.array(A_eq), b_eqnp.array(b_eq), boundsbounds, methodhighs) if not result.success: return None, None, None x_opt result.x theta result.fun efficiency 1.0 - theta sminus x_opt[sminus_idx] splus x_opt[splus_idx] return efficiency, sminus, splus def compute_all_sbm(X, Y): n X.shape[0] scores [] slacks_minus [] slacks_plus [] for i in range(n): eff, smi, spl sbm_efficiency(X, Y, dm_indexi) scores.append(eff) slacks_minus.append(smi) slacks_plus.append(spl) return np.array(scores), np.array(slacks_minus), np.array(slacks_plus)这段代码核心逻辑其实很好读懂对每个DMU单独构造线性规划求解一次得到效率值和两组松弛变量。方法里用了scipy的highs求解器这是目前scipy里最稳的线性规划求解器亲测对小规模SBM模型求速度很快。如果DMU数量超过200个循环求解会变慢后面第五节我会说怎么优化。3.3 超效率SBM的扩展实现超效率SBM和经典SBM最大的区别是计算某个DMU时参考集要剔除当前DMU自己。在上面代码里只需要改一个地方约束2和约束3中把j的范围从全量n改成“不等于dm_index”的集合。同时目标DMU的投入产出值x0、y0保持不变仍然用当前DMU的原始值。这样算出来的有效DMU效率值就会大于等于1可以继续排序。def super_sbm_efficiency(X, Y, dm_index0): n X.shape[0] m X.shape[1] s Y.shape[1] x0 X[dm_index] y0 Y[dm_index] # 参考集排除当前DMU ref_indices [j for j in range(n) if j ! dm_index] n_ref len(ref_indices) n_vars n_ref m s 1 mu_idx slice(0, n_ref) sminus_idx slice(n_ref, n_ref m) splus_idx slice(n_ref m, n_ref m s) t_idx n_ref m s c np.zeros(n_vars) for i in range(m): c[n_ref i] -1.0 / (m * x0[i]) c[t_idx] 1.0 A_eq [] b_eq [] row np.zeros(n_vars) row[t_idx] 1.0 for r in range(s): row[n_ref m r] 1.0 / (s * y0[r]) A_eq.append(row) b_eq.append(1.0) for i in range(m): row np.zeros(n_vars) for k, j in enumerate(ref_indices): row[mu_idx.start k] X[j, i] row[sminus_idx.start i] 1.0 row[t_idx] -x0[i] A_eq.append(row) b_eq.append(0.0) for r in range(s): row np.zeros(n_vars) for k, j in enumerate(ref_indices): row[mu_idx.start k] Y[j, r] row[splus_idx.start r] -1.0 row[t_idx] -y0[r] A_eq.append(row) b_eq.append(0.0) bounds [(0, None)] * n_vars result linprog(c, A_eqnp.array(A_eq), b_eqnp.array(b_eq), boundsbounds, methodhighs) if not result.success: return None, None, None x_opt result.x efficiency 1.0 - result.fun sminus x_opt[sminus_idx] splus x_opt[splus_idx] return efficiency, sminus, splus实际跑的时候应当先算经典SBM把效率值小于1的DMU直接定为非有效然后只对效率值为1的DMU调用超效率SBM。这样效率值小于1的就保持原来的效率分数效率为1的用超效率分数继续排序。最终排名表就是“非有效DMU按SBM值排、有效DMU按Super-SBM值排”混合排序业务上完全说得通。4. 结果输出与可视化落地4.1 如何组织结果表格效率值 松弛变量 投影目标模型跑完之后代码输出的不仅仅是效率值还有每个DMU在各项投入上的冗余量s⁻和各项产出上的不足量s⁺。很多人只盯着效率值一个数浪费了大量的分析价值。我在报告里从来都是给一张效率分析宽表每一行是一个决策单元列包含效率值、排名、每项投入的松弛量、每项产出的松弛量再加一列“需要整改的指标”标记。投影目标的计算很简单投入投影值 原始投入 − 投入松弛值产出投影值 原始产出 产出松弛值。这个投影值就是该DMU达到生产前沿面的改进路径。举个例子某物流企业的仓库面积投入原始值8000平方米松弛值是1200平方米那么它改进到前沿面之后仓库面积只需要6800平方米。这种数据直接给管理层看比单纯的效率值0.83要有说服力得多。表格可以用pandas直接生成并导出result_df pd.DataFrame({ DMU: df.index, 效率值: sbm_scores, 超效率值: super_scores, 最终排名分: final_scores, }) # 加入松弛变量列 for i, col in enumerate(input_cols): result_df[f{col}_投入冗余] slacks_minus[:, i] for r, col in enumerate(output_cols): result_df[f{col}_产出不足] slacks_plus[:, r] result_df.to_excel(sbm_results.xlsx, indexFalse)4.2 可视化怎么画才有信息量效率结果的可视化我一般画三张图每张图都有明确的决策价值。第一张是效率值分布条形图。横轴是DMU名称或编号纵轴是效率值按效率值从低到高排序并用颜色区分有效和无效。这张图适合回答“谁好谁坏”的问题。给运营团队看的时候他们一眼就能找出排名靠后的部门。第二张是投入冗余堆积条形图。对于非有效DMU把每一项投入的松弛量画成堆积条形图颜色越深代表冗余越多。这样可以直接看见效率低的来源到底是哪项投入浪费最多是人力、资金还是设备。这一步是SBM模型相比传统DEA最大的优势——不仅能评估还能开药方。第三张是效率四分位散点图。横轴是某个关键投入指标纵轴是效率值每个点代表一个DMU。这种图适合把效率值和业务变量交叉分析比如“投入越多是否效率一定越高”“规模大小和效率是否正相关”。做行业对标报告的时候这种图最受欢迎。4.3 画图时的几个细节坑很多新手在画图时被“横坐标太密集”卡住尤其当DMU数量有几十个横轴标签全部叠在一起。最简单有效的解决办法有三招第一招调整画布尺寸。plt.figure(figsize(12, 6))把横轴拉宽标签自动就没那么挤了。第二招旋转标签。plt.xticks(rotation45, haright)让文字斜着排阅读起来也不费力。第三招抽样显示。如果DMU超过60个每两个或每三个显示一个标签避免视觉噪声。再补充一个我在实际项目里的判断标准如果饼图或条形图超过8个分类就别用饼图了难读如果散点图数据量不大但重叠严重可以考虑加一点抖动jitter或者调低透明度。可视化的目标是把结论讲清楚不是让图好看。5. 常见问题与排查技巧实录5.1 linprog求解失败怎么排查跑SBM模型最常见的问题就是linprog返回status2无可行解。遇到这个情况先别急着怀疑代码逻辑按照下面顺序排查实测95%的问题都能解决。第一检查数据是否有0值或者负值。SBM模型的约束里会出现除法某项投入严格为0会让约束方程退化。处理方式是统一做0.1到1的min-max标准化前面第二节说过。第二检查约束矩阵是否秩亏。当某个DMU的某项指标在所有样本里数值完全相同约束方程就线性相关求解器会报数值错误。处理方式是去掉这个指标或者给该指标增加一点随机噪声。后者适合分析场景但要在结果里注明处理方式。第三检查参考集是否为空。超效率SBM里如果样本量只有1个剔除当前DMU后参考集就是空的线性规划必然无解。超效率SBM至少需要3个DMU起步我一般建议样本量不少于5个。5.2 效率值为1的DMU特别多怎么办如果一个数据集里30%以上的DMU效率值都是1说明指标区分度不够。常见原因是指标太少当投入和产出指标加起来只有2到3个时生产前沿面会退化成低维空间很多DMU都会落在前沿面上。解决思路是增加有意义的新指标。比如之前只有“员工人数”和“营收”两个指标加了“运营成本”和“客户满意度”之后前沿面维度升高效率值为1的比例会显著下降。这里要特别提醒不是指标越多越好。DEA领域有经验法则DMU数量至少是指标总数的3倍以上否则效率值会全面虚高。5.3 求解速度太慢如何优化如果你拿到的数据有上千个DMU逐个调用linprog循环求解会很慢。我实测过用纯linprog跑300个DMU、6个指标耗时大约2到3分钟对于分析场景可以接受但确实算不上快。如果数据量更大有两个方向优化。第一使用并行计算。每个DMU的线性规划之间是相互独立的天然适合多进程。Python标准库的multiprocessing.Pool就能轻松搞定把每个DMU的求解任务丢进进程池8核机器直接快4到5倍。第二改用更专业的求解器。田纳西大学的开源求解器HiGHSscipy的高s接口底层就是它已经是性能很好的选择但如果你后续要处理面板数据或者几千个决策单元建议转向pulp加CBC或者gurobi如果有学术许可证性能差距会非常明显。5.4 结果解读的最后一个建议模型跑出来了、图画好了不能直接拿效率值就下结论。我在实际项目里吃了不少亏后才总结出来的经验SBM和超效率SBM给出的效率值是相对效率不是绝对效率。换一批对比样本、换一套指标体系同一个DMU的效率排名可能发生明显变化。所以做结论之前至少要做两步验证。第一步检查效率排名是否与业务感知一致性。如果某个公认做得好的企业效率值排到倒数回头检查它的数据是否有录入错误或者指标选择是否有误导性。第二步做敏感性分析。每次剔除一个DMU或一个指标记录效率排名的变化幅度。排名波动剧烈的DMU说明它的效率很依赖特定样本写报告时表述要谨慎。我个人最后再分享一个小技巧SBM模型的结果图不要只留在Python里跑一遍就完事。把效率排名、松弛变量和投影目标打包成Excel和PDF两版Excel给数据分析的同事深挖PDF给管理层汇报。前者用pandas加openpyxl导出后者用matplotlib保存高清图直接嵌进Word或PPT。整个流程跑顺了从原始Excel到最终报告半天就能搞定效率提升是肉眼可见的。