简介DEA的MATLAB程序包内含DEA Solver是面向效率评价研究与实践的数据包络分析工具。它基于MATLAB环境支持多输入多输出系统的相对效率测算适合经管、运筹、管理工程等领域的学生、科研人员也适用于企业绩效评估等场景。压缩包共21个文件大小仅2.54MB主要包含11个xls结果表格、4个pdf文献/应用文档、2个m示例脚本、1个mat数据文件以及doc、txt和asv等配套说明与备份覆盖从模型定义、数据预处理、线性规划求解到结果输出的完整流程。已有442人学习下载。通过这套代码可以快速上手CCR、BCC等经典DEA模型理解效率得分、排序和有效单元的计算逻辑pdf与doc材料补充了方法背景和实际案例xls表单可直接查看输出结果m脚本易于修改扩展对课程教学、论文研究和企业效率诊断都很有价值。 去年帮学院一位老师处理商业银行效率评估数据她之前一直用 DEAP 和 MaxDEA 跑结果但后续的面板回归、画图、稳健性检验全在 MATLAB 里做。每次数据一改就得把 Excel 搬去 DEAP导出一份文本格式的结果再手工粘回 MATLAB一套流程下来大半天没了。当时我就琢磨与其来回倒腾不如直接用 MATLAB 自己写一个 DEA solver让数据清洗、效率计算、结果可视化全在一个环境里跑完。标题里的“内涵”我理解成“内置”或“包含”总之这套程序的核心就是自带一个可用的求解器。这个 solver 不是 Simulink 里的 Solver Configuration而是调用 MATLAB 自带的线性规划求解器 linprog不依赖任何第三方工具箱。DEAData Envelopment Analysis数据包络分析本质是求解一系列线性规划所以最顺理成章的做法就是把每个决策单元DMU的效率评估写成一次线性规划然后用 linprog 逐个解出来。这套方案我后来分享给几个学生和同事今天整理成一篇完整的博文。文章适合三类人看一是毕业论文或课题里要做效率评价但不想在 DEAP 和 MATLAB 之间反复搬运数据的人二是已经理解 DEA 原理、想自己动手实现一套完整模块的读者三是只想要一个“拿来就能跑”的工具跑完能出效率值、投影值和标杆报表的人。全文从模型原理讲到代码实现最后给出我实测过程中踩过的坑。1. 为什么我坚持在 MATLAB 里内置一个 DEA solver1.1 现成软件能跑但接不上我的分析流程市面上的 DEA 工具不少每样我都用过各有各的脾气。DEAP 是经典老牌程序免费、学术界认但输入输出全是文本文件批处理能力弱结果格式要靠手工整理。MaxDEA 功能很全中文支持好方便交互式操作但软件本身收费而且很难把每一步中间结果直接嵌入后续计量模型。Stata 的 dea 命令、R 的 Benchmarking 包其实也都能做但不是每个课题组都装了这些环境更多时候是“导师那边的数据在 MATLAB学生手头也是 MATLAB最后交作业还得是 MATLAB”。我用过一段时间之后做了个对比方案学习成本嵌入后续分析自定义程度费用DEAP低差低免费MaxDEA低差中收费R Benchmarking 包中中高免费MATLAB 自写 solver中好高已有授权表格不能表达的是工作流上的体感。自写 solver 后读 Excel、算效率、算投影、画图、导出结果整个流程不离开 MATLAB 一个环境。尤其当你有几十个模型、几百个 DMU或者要做窗口分析、Malmquist 指数时这个闭环优势非常明显。你只需要写一个循环把每期数据喂进去结果直接进 cell 数组或 table后续随你怎么玩。1.2 内置 solver 到底内置了什么很多人一听“自己写 solver”就发怵以为要重新实现单纯形法或内点法。其实完全不用。现代科学计算环境里“内置 solver”的正确姿势是把模型标准化然后调用成熟的优化内核。在 MATLAB 里这个内核就是 linprog。linprog 用的是一套经过几十年打磨的线性规划算法数值稳定性比自己写的代码可靠得多边界情况处理也完善。所以这套程序里“solver”只负责两件事把 DEA 模型的数学形式翻译成 linprog 要求的标准输入结构然后逐个 DMU 调用 linprog。真正的工作量不在求解而在模型构建和结果整理。这一点想通了整个实现难度就降下来一大半。另外注意这里的 solver 跟 Simulink 里的 Solver Configuration 完全是两个概念——Simulink 里那个是配置微分方程数值积分器的跟优化求解器不搭边。2. DEA 的数学内核CCR 和 BCC 到底在求什么2.1 从效率定义到 CCR 对偶模型DEA 的核心思想可以一句话概括不预设生产函数的具体形式让所有决策单元自己选一组对自己最有利的投入产出权重但全体 DMU 的效率上限被限制在 1 以内。每个 DMU 有 p 个投入和 q 个产出效率定义为加权产出除以加权投入。CCR 模型假设规模报酬不变数学上是一组分式规划。Charnes-Cooper 变换可以把分式规划转化成线性规划再取对偶得到包络形式也就是我们代码里真正要实现的版本min θs.t. Σ_j λ_j x_ij ≤ θ x_ik 对每个投入 iΣ_j λ_j y_rj ≥ y_rk 对每个产出 rλ_j ≥ 0这里的 θ 就是第 k 个 DMU 的效率值λ_j 是前沿面组合权重。通俗点说这套约束在问把其他所有 DMU 按 λ_j 比例组合能不能找到一个虚拟 DMU它的每个产出都不低于我而每个投入最多只有我的 θ 倍如果能那我就能压缩投入到原来的 θ 倍θ 越小说明我越没效率。几何含义也很直观把所有 DMU 的投入产出点画出来CCR 前沿面是过原点的一个锥面落在锥面上的 DMU 效率为 1落在锥面内部的效率小于 1。这就是“数据包络”名字的来源。2.2 BCC 只多了一个约束现实中很多场景并不满足规模报酬不变。比如银行合并之后效率可能不会线性增长医院病床数扩大一倍产出未必扩大一倍。这时用 BCC 模型更合适。BCC 和 CCR 唯一的区别是在约束里加了一个凸性条件Σ_j λ_j 1加了这个约束前沿面从“过原点的锥”变成了“凸包络面”允许规模报酬递增、不变、递减并存。两个模型配套使用还可以进一步分解规模效率CCR 测度的是综合技术效率BCC 测度的是纯技术效率两者相除得到规模效率。这也是为什么很多论文会同时跑两个模型。模型数学含义前沿面形状额外约束CCR规模报酬不变过原点的锥无BCC规模报酬可变凸包络面Σλ_j1如果你的研究对象是短期经营决策投入产出关系更贴近 VRS那就优先 BCC如果是长期宏观问题比如省份经济发展效率CCR 更常见。实际操作中两个都跑一遍也不费事毕竟只是多了一个等式约束。2.3 输入导向还是输出导向怎么选同样一个 CCR 或 BCC 模型还区分为输入导向和输出导向。输入导向回答的问题是“同样产出下投入最多能压缩多少”输出导向回答的是“同样投入下产出最多能扩多少倍”。代码里实现的是输入导向因为商业银行、医院这类场景中管理者对投入的控制力更强——压成本、控费用是每天都在做的事。选导向的原则其实很朴素如果你研究的问题里决策单元更容易控制投入就用输入导向如果决策单元更被动地接受投入、主动目标是扩大产出就用输出导向。两种导向算出来的 θ或 φ数值会不同但对最终排序的影响通常不大。严格的做法是把两种都报告出来看看稳健性。3. 程序架构solver 内置模块的划分与数据流3.1 文件结构与分工我习惯把代码拆成三个文件各自职责清晰调试时也能快速定位问题dea_main.m主脚本负责读取数据、设置参数、调用求解函数、调用报表函数dea_solver.m核心函数接收投入产出矩阵和模型选项内部调用 linprogdea_report.m结果处理生成效率表、投影值表、标杆表写入 Excel拆文件的好处是以后换数据不用动核心代码换模型不用动主脚本。如果只是临时跑一个小数据集把全部代码写进一个脚本也行但我不推荐因为你迟早要改参数拆开之后每个文件都是单点修改心理负担小得多。3.2 数据规范一张 Excel 表走天下数据组织用最朴素的格式一张 Excel 表前 p 列是投入指标后 q 列是产出指标一行一个 DMU不要有多余的合并单元格。读取代码很简单data readmatrix(bank_data.xlsx); % 如果第一行是表头改用 % data readmatrix(bank_data.xlsx, NumHeaderLines, 1); X data(:, 1:p); Y data(:, p1:end);p 和 q 在脚本里手动指定或者用变量名动态匹配。我建议手动指定因为大多数情况下指标数量在写脚本时已经确定了。注意检查数据里不能有 0 或负数经典 DEA 模型要求投入产出为正遇到缺失值最好整行删除或做填充不要让 NaN 混进 LP。3.3 主流程读取、校验、循环、汇总主脚本的逻辑非常直白读入数据拆成 X、Y校验行数、列数、非负性对 k 1..n构建第 k 个 DMU 的线性规划调用 linprog收集 θ 和 λ计算投影值生成报表这里最关键的决策是“逐 DMU 循环”。DEA 的问题结构决定了它不像最小二乘那样可以一次性矩阵化求解每个 DMU 的约束矩阵都不一样天然适合 for 循环。对几十到几百个 DMU 的小规模分析逐次调用 linprog 的速度完全够用一次跑几百个 DMU 也就几秒钟的事不需要考虑并行化。4. linprog 当 solver 用核心代码拆解与验证4.1 linprog 的调用约定linprog 的标准调用格式是x linprog(f, A, b, Aeq, beq, lb, ub, options)其中 f 是目标函数系数向量A、b 定义不等式约束 Ax ≤ bAeq、beq 定义等式约束 Aeqx beqlb、ub 是变量下界和上界。DEA 模型变量很简单每个 DMU 的线性规划变量只有 1n 个第一个是 θ后面 n 个是 λ_j。目标函数就是 min θ所以 f 的第一个元素为 1其余为 0。求解选项里我建议显式指定算法为 dual-simplex。这个算法对付 DEA 这种中等规模、约束结构稀疏的问题特别稳。老版本 MATLAB 可能没有这个选项那就用默认算法通常也没问题。4.2 构建第 k 个 DMU 的约束矩阵这是代码里最容易写错的地方我在草稿阶段也栽过跟头。关键是把矩阵维度和转置搞清楚。设 X 是 n×pY 是 n×q。变量向量是 v [θ; λ_1; ...; λ_n]长度 1n。投入约束Σ_j λ_j x_ij - θ x_ik ≤ 0对每个 i 写成一行。所以 A 矩阵的前 p 行里θ 的系数是 -x_ik即 -X(k,:)而 λ_j 的系数是第 i 个投入在所有 DMU 上的值也就是 X 的第 i 行等于 X 的对应行。产出约束-Σ_j λ_j y_rj ≤ -y_rk对每个 r 写成一行。这里 λ 的系数是 -Y 的第 r 行所以是 -Y。这一步最容易踩的坑就是忘记转置A(1:p, 2:end) 应该赋值 X而不是 X。如果不转置p 和 n 不相等时 MATLAB 直接报“矩阵维度不一致”相等时会给出完全错误的结果。我用一个很简单的随机数据验证过两种写法的差异转置是唯一正确的写法。4.3 完整可运行的 dea_solver.m下面这个函数可以直接抄走function [theta, lambda, exitflags] dea_solver(X, Y, model) % DEA_SOLVER 使用 linprog 求解 DEA 模型 % 输入: % X - n×p 投入矩阵n 个 DMUp 个投入指标 % Y - n×q 产出矩阵n 个 DMUq 个产出指标 % model - crs 或 vrs % 输出: % theta - n×1 效率值 % lambda - n×n 组合系数lambda(k,:) 是第 k 个 DMU 的参考组合 % exitflags - n×1 linprog 退出标志 n size(X, 1); p size(X, 2); q size(Y, 2); nVar 1 n; theta zeros(n, 1); lambda zeros(n, n); exitflags zeros(n, 1); options optimoptions(linprog, Display, off, Algorithm, dual-simplex); for k 1:n % 目标函数: min theta f zeros(nVar, 1); f(1) 1; % 不等式约束 A * v b A zeros(p q, nVar); b zeros(p q, 1); % 投入约束: sum_j lambda_j * X(j,i) - theta * X(k,i) 0 A(1:p, 1) -X(k, :); A(1:p, 2:end) X; b(1:p) 0; % 产出约束: -sum_j lambda_j * Y(j,r) -Y(k,r) A(p1:end, 2:end) -Y; b(p1:end) -Y(k, :); % 变量边界: theta 和 lambda 都非负无上界 lb zeros(nVar, 1); ub inf(nVar, 1); % VRS 模型追加凸性约束 Aeq []; beq []; if strcmp(model, vrs) Aeq zeros(1, nVar); Aeq(1, 2:end) 1; beq 1; end [x_opt, ~, exitflag] linprog(f, A, b, Aeq, beq, lb, ub, options); if exitflag 0 theta(k) x_opt(1); lambda(k, :) x_opt(2:end); else theta(k) NaN; end exitflags(k) exitflag; end end如果想改成输出导向本文还有配套的精品资源点击获取 MATLAB内置linprog实现DEA：从原理到代码的完整指南
简介DEA的MATLAB程序包内含DEA Solver是面向效率评价研究与实践的数据包络分析工具。它基于MATLAB环境支持多输入多输出系统的相对效率测算适合经管、运筹、管理工程等领域的学生、科研人员也适用于企业绩效评估等场景。压缩包共21个文件大小仅2.54MB主要包含11个xls结果表格、4个pdf文献/应用文档、2个m示例脚本、1个mat数据文件以及doc、txt和asv等配套说明与备份覆盖从模型定义、数据预处理、线性规划求解到结果输出的完整流程。已有442人学习下载。通过这套代码可以快速上手CCR、BCC等经典DEA模型理解效率得分、排序和有效单元的计算逻辑pdf与doc材料补充了方法背景和实际案例xls表单可直接查看输出结果m脚本易于修改扩展对课程教学、论文研究和企业效率诊断都很有价值。 去年帮学院一位老师处理商业银行效率评估数据她之前一直用 DEAP 和 MaxDEA 跑结果但后续的面板回归、画图、稳健性检验全在 MATLAB 里做。每次数据一改就得把 Excel 搬去 DEAP导出一份文本格式的结果再手工粘回 MATLAB一套流程下来大半天没了。当时我就琢磨与其来回倒腾不如直接用 MATLAB 自己写一个 DEA solver让数据清洗、效率计算、结果可视化全在一个环境里跑完。标题里的“内涵”我理解成“内置”或“包含”总之这套程序的核心就是自带一个可用的求解器。这个 solver 不是 Simulink 里的 Solver Configuration而是调用 MATLAB 自带的线性规划求解器 linprog不依赖任何第三方工具箱。DEAData Envelopment Analysis数据包络分析本质是求解一系列线性规划所以最顺理成章的做法就是把每个决策单元DMU的效率评估写成一次线性规划然后用 linprog 逐个解出来。这套方案我后来分享给几个学生和同事今天整理成一篇完整的博文。文章适合三类人看一是毕业论文或课题里要做效率评价但不想在 DEAP 和 MATLAB 之间反复搬运数据的人二是已经理解 DEA 原理、想自己动手实现一套完整模块的读者三是只想要一个“拿来就能跑”的工具跑完能出效率值、投影值和标杆报表的人。全文从模型原理讲到代码实现最后给出我实测过程中踩过的坑。1. 为什么我坚持在 MATLAB 里内置一个 DEA solver1.1 现成软件能跑但接不上我的分析流程市面上的 DEA 工具不少每样我都用过各有各的脾气。DEAP 是经典老牌程序免费、学术界认但输入输出全是文本文件批处理能力弱结果格式要靠手工整理。MaxDEA 功能很全中文支持好方便交互式操作但软件本身收费而且很难把每一步中间结果直接嵌入后续计量模型。Stata 的 dea 命令、R 的 Benchmarking 包其实也都能做但不是每个课题组都装了这些环境更多时候是“导师那边的数据在 MATLAB学生手头也是 MATLAB最后交作业还得是 MATLAB”。我用过一段时间之后做了个对比方案学习成本嵌入后续分析自定义程度费用DEAP低差低免费MaxDEA低差中收费R Benchmarking 包中中高免费MATLAB 自写 solver中好高已有授权表格不能表达的是工作流上的体感。自写 solver 后读 Excel、算效率、算投影、画图、导出结果整个流程不离开 MATLAB 一个环境。尤其当你有几十个模型、几百个 DMU或者要做窗口分析、Malmquist 指数时这个闭环优势非常明显。你只需要写一个循环把每期数据喂进去结果直接进 cell 数组或 table后续随你怎么玩。1.2 内置 solver 到底内置了什么很多人一听“自己写 solver”就发怵以为要重新实现单纯形法或内点法。其实完全不用。现代科学计算环境里“内置 solver”的正确姿势是把模型标准化然后调用成熟的优化内核。在 MATLAB 里这个内核就是 linprog。linprog 用的是一套经过几十年打磨的线性规划算法数值稳定性比自己写的代码可靠得多边界情况处理也完善。所以这套程序里“solver”只负责两件事把 DEA 模型的数学形式翻译成 linprog 要求的标准输入结构然后逐个 DMU 调用 linprog。真正的工作量不在求解而在模型构建和结果整理。这一点想通了整个实现难度就降下来一大半。另外注意这里的 solver 跟 Simulink 里的 Solver Configuration 完全是两个概念——Simulink 里那个是配置微分方程数值积分器的跟优化求解器不搭边。2. DEA 的数学内核CCR 和 BCC 到底在求什么2.1 从效率定义到 CCR 对偶模型DEA 的核心思想可以一句话概括不预设生产函数的具体形式让所有决策单元自己选一组对自己最有利的投入产出权重但全体 DMU 的效率上限被限制在 1 以内。每个 DMU 有 p 个投入和 q 个产出效率定义为加权产出除以加权投入。CCR 模型假设规模报酬不变数学上是一组分式规划。Charnes-Cooper 变换可以把分式规划转化成线性规划再取对偶得到包络形式也就是我们代码里真正要实现的版本min θs.t. Σ_j λ_j x_ij ≤ θ x_ik 对每个投入 iΣ_j λ_j y_rj ≥ y_rk 对每个产出 rλ_j ≥ 0这里的 θ 就是第 k 个 DMU 的效率值λ_j 是前沿面组合权重。通俗点说这套约束在问把其他所有 DMU 按 λ_j 比例组合能不能找到一个虚拟 DMU它的每个产出都不低于我而每个投入最多只有我的 θ 倍如果能那我就能压缩投入到原来的 θ 倍θ 越小说明我越没效率。几何含义也很直观把所有 DMU 的投入产出点画出来CCR 前沿面是过原点的一个锥面落在锥面上的 DMU 效率为 1落在锥面内部的效率小于 1。这就是“数据包络”名字的来源。2.2 BCC 只多了一个约束现实中很多场景并不满足规模报酬不变。比如银行合并之后效率可能不会线性增长医院病床数扩大一倍产出未必扩大一倍。这时用 BCC 模型更合适。BCC 和 CCR 唯一的区别是在约束里加了一个凸性条件Σ_j λ_j 1加了这个约束前沿面从“过原点的锥”变成了“凸包络面”允许规模报酬递增、不变、递减并存。两个模型配套使用还可以进一步分解规模效率CCR 测度的是综合技术效率BCC 测度的是纯技术效率两者相除得到规模效率。这也是为什么很多论文会同时跑两个模型。模型数学含义前沿面形状额外约束CCR规模报酬不变过原点的锥无BCC规模报酬可变凸包络面Σλ_j1如果你的研究对象是短期经营决策投入产出关系更贴近 VRS那就优先 BCC如果是长期宏观问题比如省份经济发展效率CCR 更常见。实际操作中两个都跑一遍也不费事毕竟只是多了一个等式约束。2.3 输入导向还是输出导向怎么选同样一个 CCR 或 BCC 模型还区分为输入导向和输出导向。输入导向回答的问题是“同样产出下投入最多能压缩多少”输出导向回答的是“同样投入下产出最多能扩多少倍”。代码里实现的是输入导向因为商业银行、医院这类场景中管理者对投入的控制力更强——压成本、控费用是每天都在做的事。选导向的原则其实很朴素如果你研究的问题里决策单元更容易控制投入就用输入导向如果决策单元更被动地接受投入、主动目标是扩大产出就用输出导向。两种导向算出来的 θ或 φ数值会不同但对最终排序的影响通常不大。严格的做法是把两种都报告出来看看稳健性。3. 程序架构solver 内置模块的划分与数据流3.1 文件结构与分工我习惯把代码拆成三个文件各自职责清晰调试时也能快速定位问题dea_main.m主脚本负责读取数据、设置参数、调用求解函数、调用报表函数dea_solver.m核心函数接收投入产出矩阵和模型选项内部调用 linprogdea_report.m结果处理生成效率表、投影值表、标杆表写入 Excel拆文件的好处是以后换数据不用动核心代码换模型不用动主脚本。如果只是临时跑一个小数据集把全部代码写进一个脚本也行但我不推荐因为你迟早要改参数拆开之后每个文件都是单点修改心理负担小得多。3.2 数据规范一张 Excel 表走天下数据组织用最朴素的格式一张 Excel 表前 p 列是投入指标后 q 列是产出指标一行一个 DMU不要有多余的合并单元格。读取代码很简单data readmatrix(bank_data.xlsx); % 如果第一行是表头改用 % data readmatrix(bank_data.xlsx, NumHeaderLines, 1); X data(:, 1:p); Y data(:, p1:end);p 和 q 在脚本里手动指定或者用变量名动态匹配。我建议手动指定因为大多数情况下指标数量在写脚本时已经确定了。注意检查数据里不能有 0 或负数经典 DEA 模型要求投入产出为正遇到缺失值最好整行删除或做填充不要让 NaN 混进 LP。3.3 主流程读取、校验、循环、汇总主脚本的逻辑非常直白读入数据拆成 X、Y校验行数、列数、非负性对 k 1..n构建第 k 个 DMU 的线性规划调用 linprog收集 θ 和 λ计算投影值生成报表这里最关键的决策是“逐 DMU 循环”。DEA 的问题结构决定了它不像最小二乘那样可以一次性矩阵化求解每个 DMU 的约束矩阵都不一样天然适合 for 循环。对几十到几百个 DMU 的小规模分析逐次调用 linprog 的速度完全够用一次跑几百个 DMU 也就几秒钟的事不需要考虑并行化。4. linprog 当 solver 用核心代码拆解与验证4.1 linprog 的调用约定linprog 的标准调用格式是x linprog(f, A, b, Aeq, beq, lb, ub, options)其中 f 是目标函数系数向量A、b 定义不等式约束 Ax ≤ bAeq、beq 定义等式约束 Aeqx beqlb、ub 是变量下界和上界。DEA 模型变量很简单每个 DMU 的线性规划变量只有 1n 个第一个是 θ后面 n 个是 λ_j。目标函数就是 min θ所以 f 的第一个元素为 1其余为 0。求解选项里我建议显式指定算法为 dual-simplex。这个算法对付 DEA 这种中等规模、约束结构稀疏的问题特别稳。老版本 MATLAB 可能没有这个选项那就用默认算法通常也没问题。4.2 构建第 k 个 DMU 的约束矩阵这是代码里最容易写错的地方我在草稿阶段也栽过跟头。关键是把矩阵维度和转置搞清楚。设 X 是 n×pY 是 n×q。变量向量是 v [θ; λ_1; ...; λ_n]长度 1n。投入约束Σ_j λ_j x_ij - θ x_ik ≤ 0对每个 i 写成一行。所以 A 矩阵的前 p 行里θ 的系数是 -x_ik即 -X(k,:)而 λ_j 的系数是第 i 个投入在所有 DMU 上的值也就是 X 的第 i 行等于 X 的对应行。产出约束-Σ_j λ_j y_rj ≤ -y_rk对每个 r 写成一行。这里 λ 的系数是 -Y 的第 r 行所以是 -Y。这一步最容易踩的坑就是忘记转置A(1:p, 2:end) 应该赋值 X而不是 X。如果不转置p 和 n 不相等时 MATLAB 直接报“矩阵维度不一致”相等时会给出完全错误的结果。我用一个很简单的随机数据验证过两种写法的差异转置是唯一正确的写法。4.3 完整可运行的 dea_solver.m下面这个函数可以直接抄走function [theta, lambda, exitflags] dea_solver(X, Y, model) % DEA_SOLVER 使用 linprog 求解 DEA 模型 % 输入: % X - n×p 投入矩阵n 个 DMUp 个投入指标 % Y - n×q 产出矩阵n 个 DMUq 个产出指标 % model - crs 或 vrs % 输出: % theta - n×1 效率值 % lambda - n×n 组合系数lambda(k,:) 是第 k 个 DMU 的参考组合 % exitflags - n×1 linprog 退出标志 n size(X, 1); p size(X, 2); q size(Y, 2); nVar 1 n; theta zeros(n, 1); lambda zeros(n, n); exitflags zeros(n, 1); options optimoptions(linprog, Display, off, Algorithm, dual-simplex); for k 1:n % 目标函数: min theta f zeros(nVar, 1); f(1) 1; % 不等式约束 A * v b A zeros(p q, nVar); b zeros(p q, 1); % 投入约束: sum_j lambda_j * X(j,i) - theta * X(k,i) 0 A(1:p, 1) -X(k, :); A(1:p, 2:end) X; b(1:p) 0; % 产出约束: -sum_j lambda_j * Y(j,r) -Y(k,r) A(p1:end, 2:end) -Y; b(p1:end) -Y(k, :); % 变量边界: theta 和 lambda 都非负无上界 lb zeros(nVar, 1); ub inf(nVar, 1); % VRS 模型追加凸性约束 Aeq []; beq []; if strcmp(model, vrs) Aeq zeros(1, nVar); Aeq(1, 2:end) 1; beq 1; end [x_opt, ~, exitflag] linprog(f, A, b, Aeq, beq, lb, ub, options); if exitflag 0 theta(k) x_opt(1); lambda(k, :) x_opt(2:end); else theta(k) NaN; end exitflags(k) exitflag; end end如果想改成输出导向本文还有配套的精品资源点击获取