Stata-DEA：数据包络分析一文读懂

👇 连享会 · 推文导航 | www.lianxh.cn

• 🍎 Stata：Stata基础 | Stata绘图 | Stata程序 | Stata新命令
• 📘 论文：数据处理 | 结果输出 | 论文写作 | 数据分享
• 💹 计量：回归分析 | 交乘项-调节  | IV-GMM | 时间序列 | 面板数据  | 空间计量 | Probit-Logit | 分位数回归
• ⛳ 专题：SFA-DEA | 生存分析  | 爬虫 | 机器学习 | 文本分析
• 🔃 因果：DID | RDD | 因果推断 | 合成控制法 | PSM-Matching
• 🔨 工具：工具软件 | Markdown | Python-R-Stata
• 🎧 课程：最新专题 |  计量专题 | 关于连享会

🍓 课程推荐：2024 机器学习与因果推断专题
主讲老师：司继春 (上海对外经贸大学) ；张宏亮（浙江大学）
课程时间：2024 年 11 月 9-10 日 ；16-17日
课程咨询：王老师 18903405450（微信）

 课程特色 · 2024机器学习与因果推断：

• 懂原理、会应用。本次课程邀请了两位老师合作讲授，目的在于最大限度地实现理论与应用的有机结合。为期四天的课程，分成两个部分：第一部分讲解常用的机器学习算法和适用条件，以及文本分析和大语言模型；第二部分通过精讲 4-6 篇发表于 Top 期刊的论文，帮助大家理解各类机器学习算法的应用场景，以及它们与传统因果推断方法的巧妙结合。
• 以 Top 期刊论文为范例。目前多数人的困惑是不清楚如何将传统因果推断方法与机器学习结合起来。事实上，即便是 MIT 和 Harvard 的大牛们也都在「摸着石头过河」。为此，通过论文精讲和复现来学习这部分内容或许是目前最有效的方式了。张宏亮老师此前在浙江大学按照这一模式教授了「因果推断和机器学习」课程，效果甚佳：学生们能够逐渐建立起研究设计的理念，并在构造识别策略时适当地嵌入机器学习方法。

作者： 谢雁翔 (南开大学)
邮箱： xyxmask1995@163.com

编者按： 本文部分参考如下内容，同时也特别感谢杜克锐教授和王道平博士对于部分命令的答疑。

Source:

• 数据包络分析 (DEA) 简明教程
• 厦门大学杜克锐主页
• Measuring technical efficiency and total factor productivity change with undesirable outputs in Stata

目录

• 1. 引言

• 2. DEA模型简介

• 2.1 技术效率的概念

• 2.2 径向距离模型

• 2.3 非径向距离模型

• 2.4 SBM 模型

• 2.5 Malmquist 指数分解

• 3. 相关命令的安装及解读

• 3.1 DEA命令基本操作

• 3.2 技术效率统计推断的系列命令基本操作

• 3.3 SBM 模型基本操作

• 3.4 方向距离函数与非径向方向距离函数基本操作

• 3.5 Malmquist-DEA 基本操作

• 4. 参考文献

• 5. 相关推文

温馨提示： 文中链接在微信中无法生效。请点击底部「阅读原文」。或直接长按/扫描如下二维码，直达原文：

1. 引言

数据包络分析 (DEA) 是一种数据驱动的非参数效率测度方法，并由美国著名运筹学家 A.Charnes、W.W.Cooper、E.Rhodes (1978) 首先提出。由于其不预设定具体函数形式和允许多种投入产出的优点，现已被学者们广泛用于评估决策单元的投入产出效率。

DEA 的基础是构建生产技术，在此基础上，基于生产单元到生产前沿的距离来估算效率。按照距离函数的不同，效率可分为径向效率和非径向效率。

Fare 等 (1994) 首先将 DEA 方法引入了 Malmquist-TFP 测算中。Chung 等 (1997) 提出方向距离函数 (DDF)，并在此基础上开发了 Malmquist-Luenberger 生产率指数，用来解决包含污染排放作为非期望产出情况下的全要素生产率测度问题。此后，考虑到松弛变量及非期望产出的影响，日本学者刀根熏（Tone）在 2001 年提出了 SBM-DEA (Slacks-Based Measure) 模型。

现有的 DEA 及其复杂拓展模型主要通过 Maxdea (价格高昂5000+) 完成，其他软件诸如 DEA-SOLVER、Mydea、以及 DEAP 均提供较为简单的测度模型。当然，我们也可以通过 Stata 实现！

对比现有 DEA 命令，dea 过于简单，本文将介绍效率测算的统计推断系列命令 teradial、tenonradial、teradialbc、nptestind、nptestrts。同时，由于上述命令仅能处理截面数据的静态效率计算，且不能考虑非期望产出，本文将进一步介绍由厦大杜克锐老师编写的用于测算静态径向/非径向技术效率的命令 teddf、sbmddf、ddfeff、nddfeff、以及用于测算全要素生产率的命令 gtfpch、malmq2。

2. DEA模型简介

2.1 技术效率的概念

在 DEA 分析中，技术效率是指一个生产单元 (DMU) 的生产水平达到该组数据最优技术水平的程度。技术效率可以从投入 (Debreu) 和产出 (Farrel) 两个角度来衡量，在投入既定的情况下，技术效率由产出最大化的程度来衡量。在产出既定的情况下，技术效率由投入最小化的程度来衡量。当然，在计算 TFP 的过程中，一般选择投入既定。

下面是一个单一投入产出的例子，来帮助理解技术效率的概念。

上表中， 反应为各个生产单元技术效率的高低， (标准化) 是将各单元的 除以其中的最大值，从而更好的比较这一数值。当涉及多个投入产出时，就会对各个投入与产出赋予一定的权重，然后分别加权，计算产出投入比，具体如下所示：

由上式可知，产出投入比为 。DEA 的数据驱动特征即通过数据本身获得权重，从而计算各个 DMU 的技术效率。由于 DEA 的数据驱动特征，要求投入产出数据不能存在缺失值。

2.2 径向距离模型

径向效率的度量命令 teradial 采用的是 Debreu–Farrell (Debreu，1951；Farrell，1957) 方法。具体如下：

• 假设有 个 ，对于 ，有 种投入，记为 ，有 种产出，记为 ；
• 然后，我们假设在技术条件 下，投入 产出为 ，数学表达为： are producible by ；
• 那么在技术条件 下，生产可能集表示为: 。投入的需求集表示为：。

以生产可能集为例，技术效率就表示为，某个给定数据点与生产可能集边界的距离。若以 DEA 模型来测量此种技术效率，则对于 个 ，每个 有 种投入， 种产出的数据集来说，Debreu–Farrell (Debreu，1951；Farrell，1957) 的以产出为导向的估计方法，可以通过下述线性规的方程式来表示：

其中， 是一个 的产出矩阵， 是一个 的投入矩阵。估计 是最小的包围面 (smallest convex free-disposal hull)。上述线性规划求解的是规模报酬不变 (CRS) 的技术效率。在其他关于规模报酬的假设下，只需改变 的约束，例如规模报酬可变 (VRS)，设置 即可。

命令 teddf 和 ddfeff 考虑了非期望产出，假设有 个决策单元 (DMU)，并且假设每个 使用投入向量 来联合生产期望产出向量 和非期望产出向量 ，则多产出的生产技术可以表示为：

根据 Chung 等 (1997) 的研究，在生产理论中，生产可能性集是一个有界集和闭集，且投入和期望产出具有强处置性。此外，为了让 表示环境技术，我们还需要增加弱可处置性（Weak disposability）和零结合假设（Nulljonitness）这两个额外的公理：

Chung 等 (1997) 首次将 DDF 应用于测度包含污染产出的情况，即用方向性产出距离函数来评估瑞典纸浆厂的绩效。在这里基本的方向性距离函数被定义为，寻求极大限度地增加期望产出，同时减少投入和非期望产出，即：

其中 为投入和产出应缩放的方向向量， 表示为无效率。根据研究的目的 (例如政策目标)，DDF 可使用不同大小的方向向量。结合之前的等式，我们就可以通过求解以下 DEA 模型计算出 的值，并且来测度每个 DMU 的环境无效率值。

2.3 非径向距离模型

命令 tenonradial 的非径向效率测量方法是 Russell (Färe 和 Lovell，1978；Färe 等，1994a) 法。对于以产出为导向的非径向量度定义为：

其线性规划方程式定义为：

命令 nptestrts 规模收益的两个假设检验，首先是：

如果假设 被拒绝，则可以进行下述假设检验：

也是就是说先查看此技术条件下是否是规模报酬不变的，若不是再看是否是 NIRS (NonDecreasing Returns to scale) 的。

命令 teddf 和 nddfeff 考虑了非期望产出，Zhou 等 (2012) 首先对考虑了非期望产出的 NDDF 给出了如下的正式定义：

其中， 为与投入和产出相关的标准化权重向量， 为方向向量， 为尺度因子向量。结合环境生产技术以及 NDDF 的定义，我们可以通过求解下面的 DEA 模型来计算 的值。

2.4 SBM 模型

SBM (Slacks-Based Measure) 效率测量方法是 DEA 效率测量方法中的非径向效率测度，其优点在于，直接度量多余的投入量与不足的产出量，投入与产出到生产前沿面的距离被称作松弛量 (Slacks)，来测算效率。具体原理如下：

• 假定 个决策单元 (DMU)，且每个决策单元有 3个向量，投人向量、期望产出向量和非期望产出向量，分别表示为 、、；
• 定义矩阵为 、、 分别为 、、。根据实际投入产出，假设 、 、 生产可能性集合为 ，即 种要素投入 所产生的期望和非期望产出的所有组合，可以定义为：

根据定义，加人非期望产出的 SBM-Undesirable，模型如下：

式中 、、 分别为第 个决策单元 DMU 的表示投入冗余量、正产出不足量及副产出超标量，、、 为其对应的向量， 为权重向量， 为目标函数且是严格递减的，且 。

• 当 时，即 、、 时，决策单元才是有效率的；
• 当 时，即 、、 三者中至少有一个不等于零时，决策单元是无效率的，存在着投人产出上改进的必要性。同时，由于该模型是一个非线性规划模型，可以根据 Charnes-Cooper 的转换方法将其转换成线性规划模型进行求解。

2.5 Malmquist 指数分解

Lovell 将构成 CRS 生产可能集的前沿技术称为基准技术，即为了计算 TFP 而定义的参照技术，并将构成 VRS 生产可能集的前沿技术称为最佳实践技术，即现实中存在的前沿技术。Malmquist 生产率指数应当定义在基准技术之上，基于 和 期参照技术的 Malmquist 生产率指数分别为：

因为基于 和 期参照技术定义的 Malmquist 生产率指数在经济含义上是对称的，按照 Fisher (1922) 理想指数思想，定义它们的几何平均为综合生产率指数:

对于 Malmquist 生产率指数本身，FGNZ 与 RD 之间并没有分歧，分歧在于对指数的分解上。FGNZ 的分解为：

RD 的分解为：

其中 、 和 分别表示技术效率变动、技术进步和规模报酬变动。显然 FGNZ 和 RD 在技术效率变动分解上是一致的，不同之处在于对技术进步和规模报酬变动的分解 FGNZ 本质错误在于技术进步定义，FGNZ 承认现实技术是 VRS 形式的，并据此计算了技术效率变动和规模报酬变动。但是在计算技术进步时却舍弃了现实的 VRS 技术，反而采用了假设的 CRS 技术。因此 FGNZ 分解中的技术进步并不是现实技术进步，而是参照技术的进步。由此也导致二者对规模报酬变动的不同界定。定义 为规模效率，规模报酬变动应当反映沿着同一生产前沿的规模效率变化回，即 RD 分解 ，而不是 分解 所表示的沿不同生产前沿的规模效率变化。

3. 相关命令的安装及解读

3.1 DEA命令基本操作

dea ivars = ovars [if] [in] [, options]

其中，ivars 是指投入变量，ovars 是指产出变量。options 主要包含以下内容：

• rts(crs|vrs|drs|nirs)：设定的不同的规模报酬，crs 为不变规模报酬，vrs 为可变规模报酬，drs 为递减规模报酬，nirs 为非递增规模报酬；
• ort(in|out)：设定模型类型，即 in 表示投入导向 DEA 模型，out 表示产出导向 DEA 模型，默认为投入导向 DEA 模型;
• stage(1|2)：设定识别效率松弛的公式，即 1 表示单阶段 DEA 模型，2 表示两阶段 DEA 模型，默认为后者；
• trace：将估计过程和结果显⽰在结果窗口，并保存到 dea.log 文件中；
• saving(filename)：将结果保存到名为 filename 的文件中。

具体示例代码如下：

use "coelli_table6.4.dta", clear 
dea i_x = o_q, rts(vrs) ort(in)   //投入导向规模报酬可变
dea i_x = o_q, rts(crs) ort(out) stage(2) saving(dea1_result)  //产出导向规模报酬不变
dea i_x = o_q, rts(vrs) ort(out) stage(2) saving(dea2_result)  //产出导向规模报酬可变

3.2 技术效率统计推断的系列命令基本操作

teradial 和 tenonradial

一次将效率测算的统计推断系列命令 teradial、tenonradial、teradialbc、nptestind、nptestrts 全部安装。

teradial    outputs = inputs [(ref_outputs = ref_inputs)] [if] [in] [, options]
tenonradial outputs = inputs [(ref_outputs = ref_inputs)] [if] [in] [, options]

outputs 为产出变量，inputs 为投入变量，ref_outputs 和ref_inputs 均是可选项，分别表示产出变量和投入变量的技术参考集，二者必须同时设定，若不设定则默认全部变量均为参考值。options 主要包含以下几项：

• rts (rtsassumption)：表示对规模报酬的假定，默认是 CRS，还可以设定 VRS 和 NIRS；
• base (basetype)：默认是产出导向 (output) 的计算，可设定投入导向；
• reference (varname)：定义变量的产出投入数据来构建技术参考边界。一般用全部变量，无需特别定义。如果 ref_outputs 和 ref_inputs 被定义，那么它们所涉及变量的数据点就被定义为技术参考边界；
• tename(newvar)：生成技术效率的新变量分别对应 Farrell 径向技术效率和 Russell 非径向技术效率；
• noprint：不报告估计过程。

具体示例如下：

set more off
use "ccr81.dta", clear 
teradial    y1 y2 y3 = x1 x2 x3 x4 x5, rts(crs)  base(output) tename(TErdCRSo) //产出导向，CCR模型
teradial    y1 y2 y3 = x1 x2 x3 x4 x5, rts(vrs)  base(output) tename(TErdVRSo) //产出导向，VRS径向，BCC模型
tenonradial y1 y2 y3 = x1 x2 x3 x4 x5, rts(crs)  base(output) tename(TEnrCRSo) //产出导向，CRS非径向
tenonradial y1 y2 y3 = x1 x2 x3 x4 x5, rts(nirs) base(output) tename(TEnrNRSo) //产出导向，NIRS非径向
tenonradial y1 y2 y3 = x1 x2 x3 x4 x5, rts(vrs)  base(output) tename(TEnrVRSo) //产出导向，VRS非径向
list TErdCRSo TErdVRSo TEnrCRSo TEnrNRSo TEnrVRSo in 1/7

      +-----------------------------------------------------------+
      |  TErdCRSo    TErdVRSo    TEnrCRSo    TEnrNRSo    TEnrVRSo |
      |-----------------------------------------------------------|
   1. | 1.0872573   1.0322941   1.1172095   1.0565404   1.0565404 |
   2. |  1.110133   1.1093141   1.3830886    1.277123    1.277123 |
   3. | 1.0790336   1.0684294   1.1705302   1.1165817   1.1165817 |
   4. | 1.1194336   1.1074129   1.4890863   1.4713013   1.4713013 |
   5. | 1.0758641           1   1.1967792   1.1967792           1 |
      |-----------------------------------------------------------|
   6. |  1.107752   1.1050748   1.3802144   1.3783777   1.3783777 |
   7. |  1.125782   1.1190873   1.5752878    1.547186    1.547186 |
      +-----------------------------------------------------------+

从上述结果可以看出，在相同规模报酬下，非径向技术效率大于等于径向技术效率。一般来说，当存在非零松弛时，径向测度方法要求投入产出同比例变动，会高估决策单元的效率。此处虽然结果相反，但径向并没有比非径向技术效率小太多。无论径向还是非径向，不同规模报酬下的技术效率表现为 CRS 大于等于 NIRS 大于等于 VRS。对于决策单元 1、2、3、4、7 而言，NIRS 和 VRS 下两种技术效率相等；对于决策单元 5 和 6，NIRS和 CRS 下两种技术效率相等。

nptestind

nptestind outputs = inputs [if] [in] [, options]

这个命令用来检验非参模型的独立性，且仅适用于径向技术效率模型的检验。基本语法与前述命令一致，只是在 options 里，有一个 alpha(#)，其设定了显著性水平，默认为 0.05。

通过 nptestind 判断效率与投入产出组合是否独立，选择相应的 bootstrap 进行检验。Bootstrap 方法的应用要依赖于特定的假设，若是基于产出导向的技术效率，则主要假设为该效率是否与产出组合独立；若是基于投入导向的技术效率，则主要假设为该效率是否与入组合独立。如果独立，则平滑同质性 bootstrap（smoothed homogeneous bootstrap) 技术可以应用；反之，则需要使用平滑异质性 bootstrap (heterogeneous bootstrap) 技术来提供有效的统计推断。

具体示例如下：

matrix testsindpv = J(2, 3, .)    //生成三行两列的矩阵均为空值
matrix colnames testsindpv = CRS NiRS VRS  //三列标题分别命名为CRS NIRS VRS
matrix rownames testsindpv = output-based input-based //两行分别命名为output input
nptestind y1 y2 y3 = x1 x2 x3 x4 x5, rts(crs)  base(output) reps(999) alpha(0.05)  //判断
matrix testsindpv[1,1] = e(pvalue) 
nptestind y1 y2 y3 = x1 x2 x3 x4 x5, rts(nirs) base(output) reps(999) alpha(0.05) noprint
matrix testsindpv[1,2] = e(pvalue) 
nptestind y1 y2 y3 = x1 x2 x3 x4 x5, rts(vrs)  base(output) reps(999) alpha(0.05) noprint
matrix testsindpv[1,3] = e(pvalue) 
nptestind y1 y2 y3 = x1 x2 x3 x4 x5, rts(crs)  base(input)  reps(999) alpha(0.05) noprint
matrix testsindpv[2,1] = e(pvalue)  
nptestind y1 y2 y3 = x1 x2 x3 x4 x5, rts(nirs) base(input)  reps(999) alpha(0.05) noprint
matrix testsindpv[2,2] = e(pvalue) 
nptestind y1 y2 y3 = x1 x2 x3 x4 x5, rts(vrs)  base(input)  reps(999) alpha(0.05) noprint
matrix testsindpv[2,3] = e(pvalue) 
matrix list testsindpv   //列出矩阵进行判断

testsindpv[2,3]
                    CRS       NiRS        VRS
output-based          .  .21821822  .04104104
input-based   .03203203    .004004  .25225225

可以看出，在基于产出导向的非参模型中，只有 VRS 技术的独立性检验 p 值是 0.041，即在该技术下，独立性假设在 5% 显著水平被拒绝。在基于投入导向的非参模型中，结果恰好相反，只有 VRS 技术下的独立性假设在 5% 水平没有被拒绝。由于我们以产出导向为例，所以按照前述理论，对于 VRS 技术下的非参模型，应该使用平滑异质性 bootstrap 方法。

teradialbc

teradialbc outputs = inputs [(ref_outputs = ref_inputs)] [if] [in] [, options]

teradialbc 是针对径向技术效率的统计推断，通过截面数据的计算可以得到技术效率的校正偏差、方差波动、以及确定上下界。命令语法与上述命令基本一致，这里重点介绍几个 options 选项：

• subsampling：要求使用分样本 bootstrap 技术来构建生产参考集，用于对比；
• smoothed：要求使于平滑 bootstrap 技术来构建生产参考集，默认选项；
• heterogeneous：要求使用异质性平滑 bootstrap 技术来构建生产参考集，若没有设定，为默认值；
• rep(#): 设定 bootstrap (可放回重复抽样构成新的生产集合估计) 重复次数，默认为 reps(999)，最小值为 reps(200)。如果使用校正偏差方法对置信区间进行估计，重复次数最好为最大值；
• tebc (newvar)：生成 1 个新变量 (误差校正效率)，该变量为校正偏差的径向技术效率；
• biassqvar (newvar)：生成 1 个新变量衡量方差，该变量为径向技术效率偏差平方与方差之比的三倍，该值越小方差越大；
• telower (newvar) ：生成 1 个新变量，该变量为径向技术效率估计值的下界；
• teupper (newvar)：生成 1 个新变量，该变量为径向技术效率估计值的上界。

具体示例如下：

. teradial   y1 y2 y3 = x1 x2 x3 x4 x5, rts(vrs) base(output)  ///
            tename(TErdVRSo1) . noprint   //估算出原始的径向技术效率

. teradialbc y1 y2 y3 = x1 x2 x3 x4 x5, rts(vrs) base(output)  ///
             heterogeneous reps(200) tebc(TErdVRSoBC2) biassqvar(TErdVRSoBC2bv) ///
             telower(TErdVRSoLB2) teupper(TErdVRSoUB2) noprint  //平滑异质性 bootstrap

. teradialbc y1 y2 y3 = x1 x2 x3 x4 x5, rts(vrs) base(output) ///
             subsampling reps(200) tebc(TErdVRSoBC3) biassqvar(TErdVRSoBC3bv)   ///
             telower(TErdVRSoLB3) teupper(TErdVRSoUB3) noprint  //分样本 bootstrap

. list TErdV* in 1/10 

    +-----------------------------------------------------------------------------------------------------------------------+
      |  TErdVRSo   TErdVRSo1   TErdVR~C2   TErdV~2bv   TErdV~LB2   TErdV~UB2   TErdVR~C3   TErdV~3bv   TErdV~LB3   TErdV~UB3 |
      |-----------------------------------------------------------------------------------------------------------------------|
   1. | 1.0322941   1.0322941   1.1809929   76.225244   1.1390678   1.2941726   1.1505643   3.6002157   1.0400416   1.4177618 |
   2. | 1.1093141   1.1093141   1.2230107   35.312512   1.1581699   1.3216823   1.1809096   5.1156057   1.1151088   1.2975177 |
   3. | 1.0684294   1.0684294   1.1576379   46.506792   1.1153948   1.2206362   1.1346832   6.3092316   1.0695672   1.2479747 |
   4. | 1.1074129   1.1074129    1.144445    21.58194   1.1257194    1.183656    1.136345   6.6563286   1.1098659    1.173818 |
   5. |         1           1           .           .           .           .           1   16.238698           1           1 |
      |-----------------------------------------------------------------------------------------------------------------------|
   6. | 1.1050748   1.1050748   1.3115154   18.700563   1.1962903   1.6766238   1.1853757   1.3266262   1.1103577   1.9223312 |
   7. | 1.1190873   1.1190873   1.1551476   23.303265   1.1324907     1.19003   1.1499285   7.1769983   1.1226502    1.180018 |
   8. | 1.1042363   1.1042363   1.3625102   45.553601   1.2274371   1.6945506   1.2628257   3.4194518   1.1068878   1.6552209 |
   9. | 1.1607212   1.1607212   1.2031425    23.99597   1.1763622   1.2399871   1.2134079   7.8238199   1.1622436    1.265083 |
  1.  | 1.0545148   1.0545148   1.1759478    67.57958   1.1288292   1.2655937   1.1454321   4.2965037   1.0573176   1.3261682 |
      +-----------------------------------------------------------------------------------------------------------------------+

从结果可以看出，对于分样本 bootstrap 方法而言，其方差波动值要远远小于平滑异质性 bootstrap 方法下的结果，这说明分样本 bootstrap 方法估计值的方差要更大，而且校正偏差的效率估计值的均方误差也要比原始效率的均方误差更大，即比起分样本 bootstrap 方法，平滑异质性 bootstrap 方法更适用于该样本。值得注意的是，第 5 个决策单元在平滑异质性 bootstrap 方法下没有结果，这是因为该观测值位于 bootstrap 构建的边界之内，造成线性规划无解。

nptestrts

nptestrts outputs = inputs [if] [in] [, options]

命令 nptestrts 用来检验非参模型的规模报酬且仅适用于径向技术效率模型的检验，重点介绍几个 options：

• base (basetype)：指定优化类型，基本类型可以是产出导向或投入导向；
• heterogeneous：异质性 bootstrap，默认为同质；
• reps (#)：bootstrap 抽样次数 replications，默认 999；
• alpha (#):确定置信区间，默认 0.5；
• testtwo：设定 test#2 被检验，否则仅检验 test#1；
• tecrsname(newvar)：在 CRS 的假设下，tecrsname(newvar) 使用技术效率的径向度量生成变量；
• tenirsname(newvar)：在 NIRS 的假设下，tenirsname(newvar) 使用技术效率的径向度量生成变量；
• tevrsname(newvar):在VRS的假设下，tevrsname(newvar) 使用技术效率的径向度量生成变量；
• sefficiency(newvar)：生成具有规模效率的变量；
• psefficient(newvar)：生成变量，其检验的 p 值表示如果统计上有效，则数据点的小数位数效率不等于 1；
• sefficient：生成 1 个新变量（有效），当决策单元在统计上规模有效时该变量等于 1，否则等于 0；
• nrsovervrs(newvar)：在 NIRS 和 VRS 下以技术效率的度量之比生成变量；
• pineffdrs(newvar)：生成带有测试的 p 值的变量，该数据的数据点由于 DRS 而导致规模无效；
• sineffdrs：生成 1 个新变量 (无效)，当决策单元由于规模报酬递减在统计上规模无效时，该变量等于 1，否则等于 0；
• noprint：禁止显示日志；
• nodots：抑制复制点的显示。

具体示例如下：

nptestrts y1 y2 y3 = x1 x2 x3 x4 x5, base(output) heterogeneous reps(200) ///
alpha(0.05) testtwo sefficient(SEffnt_het) sineffdrs(SiDRS_het)   //前文检验支持平滑异质性bootstrap
table SEffnt_het
table SiDRS_het

Radial (Debreu-Farrell) output-based measures of technical efficiency under assumption of CRS, NIRS, and VRS
technology are computed for the following data:

     Number of data points (K) = 70
     Number of outputs     (M) = 3
     Number of inputs      (N) = 5
     
  Reference set is formed by 70 data points, for which measures of technical efficiency are computed.
  ----------------------------------------------------------------------------------------
           Test #1
 
  Ho: mean(F_i^CRS)/mean(F_i^VRS) = 1
    and
  Ho: F_i^CRS/F_i^VRS = 1 for each of 70 data point(s)
 
  Bootstrapping reference set formed by 70 data points and computing radial (Debreu-Farrell) output-based measures of
  technical efficiency under assumption of CRS and VRS technology for each of 70 data points relative to the
  bootstrapped reference set
 
  Smoothed heterogeneous bootstrap (200 replications)
 ----+--- 1 ---+--- 2 ---+--- 3 ---+--- 4 ---+--- 5 
 ..................................................    50
 ..................................................   100
 ..................................................   150
 ..................................................   200
  p-value of the Ho that mean(F_i^CRS)/mean(F_i^VRS) = 1 (Ho that the global technology is CRS) = 1.0000:
 
  mean(hat{F_i^CRS})/mean(hat{F_i^VRS}) = 1.0164 is not statistically greater than 1 at the 5% significance level
 
  All data points are scale efficient
  ----------------------------------------------------------------------------------------
 
 . table SEffnt_het
 ----------------------------
 Indicator       |
 variable if     |
 statistically   |
 scale efficient |      Freq.
 ----------------+-----------
 scale efficient |         70
 ----------------------------
 
 . table SiDRS_het
 variable SiDRS_het not found
 r(111);

从上面的估计结果可以看出，使用平滑异质性 bootstrap 方法，由于 test#1 的 p 值为 1.0164，说明无法拒绝全域生产技术是 CRS 的原假设，进一步对 70 个决策单元的规模效率进行检验，发现所有的决策单元都是规模有效的，因此无需再进行 test#2 的检验。

上述变量 SEffnt_het 是表述决策单元规模有效的指示变量，由于所有决策单元均是规模有效的，因此频数为样本数 70；变量 SiDRS_het 是表述决策单元规模无效的指示变量，由于在上述样本估计中，不存在规模无效的样本，因此无法生成该变量，所以会显示该变量无法找到。

3.3 SBM 模型基本操作

teddf

安装方式：

*第一种安装方式
net install gtfpch, from("https://raw.githubusercontent.com/kerrydu/gtfpch/master/") force   //如果网络不好可以尝试第二种
*第二种安装方式：下载压缩文件并将其解压缩到计算机磁盘上
copy https://codeload.github.com/kerrydu/gtfpch/zip/master gtfpch-master.zip
unzipfile gtfpch-master.zip
net install gtfpch, from(`c(pwd)'/gtfpch-master) force
*第三种安装方式：如果上述方法均受网速限制，建议手动下载
viewe browse "https://github.com/kerrydu/gtfpch"  //打开网页下载
net install gtfpch, from(D:\Lianxh\TE_Lianxh\T2_DEA\gtfpch\gtfpch-master)force  //将下载的gtfpch-master.zip放入相应路径，更改路径安装，注：D:\Lianxh\TE_Lianxh\T2_DEA\gtfpch为更改路径

语法说明：

teddf inputvars = desirable_outputvars : undesirable_outputvars [if] [in] , dmu(varname) [options]

命令 teddf 必须使用面板数据，指定面板变量和时间变量。

• inputvars：投入变量；
• desirable_outputvars：期望产出；
• undesirable_outputvars：非期望产出；
• dmu(varname): 指定各决策单元变量 (DMU)，必填项；
• time(varname)：指定表示同时期的生产技术 (当年)，指定当年生产投入组合 convex combination 作为技术参考集。设定时期若不写这一项，则默认为全局生产技术 (global production technology)。
• gx(varlist)：指定投入调整的方向向量，gx() 中指定的变量顺序应与投入变量中的顺序相同，默认是 gx=(0,..,0)；
• gy(varlist)：指定期望产出调整的方向向量。默认 gy=Y，即期望产出观测值，gy() 中指定的变量顺序应与 desirable_outputvars 中的顺序相同；
• gb(varlist)：指定非期望产出调整的方向向量。默认 gb=-B，即非期望产出观测值，gb() 中指定的变量顺序应与 undesirable_outputvars 中的顺序相同；
• nonradial：指定使用非径向方向距离函数；
• wmat(name)：指定用于调整输入和输出变量的权重矩阵 wmat(name) 仅在指定非径向时才可以使用，默认值为 W =（1/n，...，1/n)；
• vrs:指定具有可变规模收益的生产技术。默认情况下，j 假定是规模报酬不变 (CRS) 的生产技术；
• sup：指定超效率 super-dea 度量，暂时处于调试阶段，还不能使用；
• tone:指定了 Tone(2004) 的 SBM 模型，计算松弛变量 (Slack) 分离分析单要素效率；
• window(#)：以 ＃ 周期带宽指定窗口生产技术，能够解决不可行解问题，一般窗口设置为 3；
• biennial：期指定两年一次的生产技术，能够解决不可行解问题；
• sequential：表示指定序列生产技术，技术不会退步假设，DEA 的模型设定根据自身要求，均为先验的；
• global：表示指定全局生产技术，解决不可行解问题；
• saving(filename[， replace])：以 filename.dta 格式储存结果；
• maxiter(#): 指定最大迭代次数，该迭代次数必须为大于 0 的整数。最大迭代次数的默认值为 16000；
• tol (real): 收敛标准容差，其值必须大于 0。若不写这项，则默认为 1e-8。

具体示例如下：

cls
*crs sbm
use example_ddf.dta,clear
qui teddf labor capital energy= gdp: co2, dmu(id) time(t) tone sav(crssbm,replace)   //SBM-Slack模型保存到crssbm文件，默认crs
qui merge 1:1 id t using crssbm
qui gen EE=1-(S_energy/energy)  //计算单要素投入，根据松弛程度计算能源效率
qui drop _merge 
save tee.dta,replace        //保存crs下计算的效率值
*vrs sbm
use example_ddf.dta,clear
qui teddf labor capital energy= gdp: co2, dmu(id) time(t) tone vrs sav(vrssbm,replace) 
qui merge 1:1 id t using vrssbm
qui gen EE2=1-(S_energy/energy)
qui drop _merge 
save pee.dta,replace        //保存vrs下计算的效率值
*合并数据
use tee.dta,clear
qui merge 1:1 id t using pee.dta,keepusing(EE2)
qui gen SEE=EE/EE2   //计算规模效率
list EE EE2 SEE in 1/10   

      +-----------------------------------+
      |        EE         EE2         SEE |
      |-----------------------------------|
   1. | .78440115   .89173482   .87963499 |
   2. | .78855801   .85038927   .92729064 |
   3. |         1           1           1 |
   4. | .87553169           1   .87553169 |
   5. | .95080883   .95581305   .99476444 |
      |-----------------------------------|
   6. | .95307049   .99590829    .9569862 |
   7. | .35996606   .79537908    .4525717 |
   8. | .35878288   .77911833   .46049858 |
   9. |         1           1           1 |
  1.  |         1           1           1 |
      +-----------------------------------+

sbmeff

安装方式：

*第一种安装方式
net install sbmeff, from("https://raw.githubusercontent.com/kerrydu/sbmeff/master/")
*第二种安装方式：下载压缩文件并将其解压缩到计算机磁盘上
copy https://codeload.github.com/kerrydu/sbmeff/zip/master sbmeff-master.zip
unzipfile sbmeff-master.zip
net install sbmeff, from(`c(pwd)'/sbmeff-master)

3.4 方向距离函数与非径向方向距离函数基本操作

teddf 和 ddfeff

teddf 和 ddfeff 均能够完成基于方向距离函数 (DDF) 效率的测算，dffeff 和 ndffeff 的安装方式如下：

*第一种安装方式
net install ddfeff, from("https://raw.githubusercontent.com/kerrydu/ddfeff/master/")
*第二种安装方式：下载压缩文件并将其解压缩到计算机磁盘上
copy https://codeload.github.com/kerrydu/ddfeff/zip/master ddfeff-master.zip
unzipfile ddfeff-master.zip
net install ddfeff, from(`c(pwd)'/ddfeff-master)   

由于 dffeff 与 teddf 语法结构一致，且 ddfeff 能够被 teddf 完全替代，此处仅通过 teddf 进行分析， ddfeff 的详细语法结构参见 help ddfeff。

具体示例代码如下：

cls
use example_ddf, clear
*利用teddf来测算：全域生产技术方向性距离函数的效率值
teddf labor capital energy= gdp: co2, dmu(id) sav(ddf_result,replace)   //全域生产技术
use ddf_result,clear
gen t = 10
replace t = 9 if mod(Row,2)       //因为是global生成时间t进一步区分匹配
merge 1:1 id t using example_ddf   
drop _merge 
gen ee=1/(1+Dval)   //计算环境技术效率ee，Dval即β值，表示距离生产前沿的距离，若为0零则在生产前沿上最有效，
list in 1/10   //等于1最有效，在生产前沿面上

*利用teddf来测算指定当期生产技术方向性距离函数的值
cls
use example_ddf, clear
qui teddf labor capital energy= gdp: co2, dmu(id) time(t) sav(ddf_result1,replace)
qui merge 1:1 id t using ddf_result1
drop _merge 

*改变投入产出方向向量g的大小，设定gy = 2Y，gb = -2B 用来计算 DDF 的值 gx = 0，结果如下：
cls
gen gdp_2 = 2*gdp
gen co2_2 = 2*co2
qui teddf labor capital energy= gdp: co2, dmu(id) time(t) gy(gdp_2) gb(co2_2) sav(ddf_result2,replace) //改变投入产出方向向量大小

teddf 和 nddfeff

teddf 和 nddfeff 均能够完成基于非径向方向距离函数 (NDDF) 效率的测算，nddfeff 与 teddf  语法结构一致，teddf 的nonradial选项即可完成对于非径向方向距离函数 (NDDF) 效率的测算，此处通过 nddfeff 进行分析，nddfeff的详细语法结构参见help nddfeff。

cls
use example_ddf, clear
*首先使用默认设定值计算出 NDDF 的值
qui nddfeff labor capital energy= gdp: co2, dmu(id) sav(nddf,replace)   //全域生产技术,可以得到每一个投入产出变量的β值，进而得到单要素效率，Dval为综合β值
use nddf
qui gen dval = (1/9)*(B_labor + B_capital + B_energy) + (1/3) * B_gdp + (1/3) * B_co2
qui drop  _Row 
qui order dval  Dval
list in 1/10
 
*设置序列生产技术（即：sequential production technology）来计算NDDF的值。seq设定不存在技术退步
cls
use example_ddf, clear
nddfeff labor capital energy= gdp: co2, dmu(id) time(t) seq sav(nddf_result,replace)  
 
*可持续发展效率stp真实值计算
use nddf_result,clear
gen stp=1-Dval    //1为最有效，Dval为综合β无效率值，减去其非效率项得到stp为综合可持续发展效率

 +----------------------------------------------------------------------------------------+
 |      dval        Dval   id     B_labor   B_capital    B_energy       B_gdp       B_co2 |
 |----------------------------------------------------------------------------------------|
 | .25892964   .25892964    1   .64426131   .10559285   .23856332   7.075e-12   .44731641 |
 | .24840085   .24840085    1   .62711309    .1311299   .21144199   2.713e-09   .42197421 |
 | .08336507   .08336507    2   7.913e-11   .34000855   .17484528   1.368e-12   .07847726 |
 | .05095057   .05095057    2   1.002e-12   .32490275   .12563158   1.116e-11   .00267359 |
 | .10600371   .10600371    3   .24085305   .17981374   .07702746   9.320e-13   .15211304 |
 |----------------------------------------------------------------------------------------|
 | .08122075   .08122075    3   .15788438   .20410305   .04692951   4.874e-16   .10735661 |
 | .43398236   .43398236    4   .62884844   6.367e-11   .54237384   .30964538   .60189427 |
 | .43087166   .43087166    4   .57959145   7.869e-14   .51791025   .34368102   .58310006 |
 | 1.085e-11   1.085e-11    5   8.740e-15   8.494e-13   8.500e-12   5.505e-13   2.889e-11 |
 | 1.408e-10   1.408e-10    5   7.071e-11   4.726e-11   1.351e-10   3.750e-12   3.342e-10 |
 +----------------------------------------------------------------------------------------+

可以看出利用 nddfeff 命令计算的结果中不仅包含了 NDDF 的值 ( 值反映距离技术前沿的远近)，也包括了投入方向、期望产出方向以及非期望产出方向下的 值。因为该命令默认下的权重矩阵设定为各 1/3，而该样本集有 3 个投入变量，则各投入变量的权重为 1/9。利用第一行数据结果来验证，则：，这个值可以作为综合绩效指数，因为其包括经济能源和碳排放综合非效率信息。

3.5 Malmquist-DEA 基本操作

由于命令 gtfpch 在 malmq2 的基础上考虑了非期望产出的影响，因此本文着重介绍gtfpch命令。安装过程如上文所述，不再赘述；其基本语法如下：

tfpch inputvars = desirable_outputvars : undesirable_outputvars [if] [in], [options]

• inputvars：投入变量；
• desirable_outputvars：期望产出；
• undesirable_outputvars：非期望产出；
• dmu(varname): 各决策单元变量 (DMU) 的指定，必填项；
• gx(varlist)：指定投入调整的方向向量，默认 gx=0，gx() 中指定的变量顺序应与 inputvars 中的顺序相同；
• gy(varlist)：指定期望产出调整的方向向量。默认 gy=Y，即期望产出观测值，gy() 中指定的变量顺序应与 desirable_outputvars 中的顺序相同；
• gb(varlist)：指定非期望产出调整的方向向量。默认 gb=-B，即非期望产出观测值，gb() 中指定的变量顺序应与undesirable_outputvars 中的顺序相同；
• sequential: 表示指定序列生产技术，技术不会退步、经济持续增长假设；
• biennial：期指定两年一次的生产技术，解决不可行解问题；
• global: 表示指定全局生产技术，解决不可行解问题；
• nonradial：指定使用非径向方向距离函数；
• wmat(name)：指定用于调整输入和输出变量的权重矩阵，wmat(name) 仅在指定非径向时才可以使用；
• luenberger：指定估计 Luenberger 生产率指数。默认为 Malmquist–Luenberger 生产率指数。当指定非径向时，将估算 Luenberger 生产率指数；
• ort(string)：默认为ort(out)，表示产出导向的生产率指数。ort(i) 表示投入导向的生产率指数。
• fgnz：指定 fgnz 方法分解 TFP，使用 Färe 等 (1994) 的方法分解 TFP 的变化；
• rd：指定 rd 方法分解 TFP，使用 Ray 和 Desli (1997)的方法分解 TFP 变化；
• saving(filename[， replace])：以 filename.dta 格式储存结果；
• maxiter(#): 指定最大迭代次数，该迭代次数必须为大于 0 的整数。最大迭代次数的默认值为 16000；
• tol(real): 收敛标准容差，其值必须大于 0。若不写这项，则默认为 1e-8；
• window(#)：指定带有 ＃ 个周期的窗口生产技术。用于解决不可行解，一般窗口为 3，参考 Charnes (1994)、Asmild (2004)。

MLPI 和 LPI

实证研究中，计算绿色生产率时，MLPI 常和径向方向距离函数 (Radial DDF) 一起使用，LPI 常和非径向方向距离函数 (Non-radial DDF) 一起用。以下是这两个指数的常用代码实现过程：

use example_ddf.dta, clear
xtset id t   //必须使用面板数据

*MLPI (radial DDF)  径向DDF，几何平均乘法
gtfpch labor capital energy= gdp: co2, sav(ddfr_result, replace)  //存在不可行解得问题
 
*使用窗口window解决不可行解的问题
gtfpch labor capital energy= gdp: co2,  window(3) sav(ddfr_result1, replace)
 
*LPI（non-radial DDF）非径向DDF，加法
qui gtfpch labor capital energy= gdp: co2, nonradial luenberger sav(ddfn_result, replace) 
qui gtfpch labor capital energy= gdp: co2, nonradial window(3) luenberger sav(ddfn_result1, replace) 

SMLPI

当 DMU 样本容量很小或者考虑技术的本质 (即技术不会出现退步)，基于技术非退步的假设，构造 1 个不仅包括当年，还包括当年之前的所有年份的序列环境技术集，使用 Sequential Malmquist–Luenberger Productivity Index (SMLPI)。

*SMLPI 的 Stata 命令是在 MLPI 命令基础上加上了sequential 这个 option，具体如下：
gtfpch labor capital energy= gdp: co2, sequential sav(ddfs_result, replace)   //存在无解

需要注意的是：sequential 与 window() 不能同时使用。

GMLPI

GMLPI 即 Global Malmquist–Luenberger Productivity Index，全局曼奎斯特罗恩伯格生产率指数。其能克服在使用 MLPI 时出现的不可行解的问题， DMU 的 A 点是在 MLPI 下出现不可行解的问题，因为 A 不能沿着方向向量的方向与 t+1 期的生产前沿面有交点，所以传统的解法不可行，要构造1个全局的或者说全域的环境生产技术集，使其把所有当期的环境技术都包括，构成新的单一的技术集。

*GMLPI 的 Stata 命令是在 MLPI 命令基础上加上了 global 这个 option，具体如下：
gtfpch labor capital energy= gdp: co2, global sav(ddfg_result, replace)  

BMLPI

BMLPI 即 Biennial Malmquist–Luenberger Productivity Index，两期曼奎斯特罗恩伯格生产率指数。由王兵等（2013）提出，和 GMLPI 类似，解决存在非期望产出测度绿色生产率的不可行性解问题。GMLPI 虽然可以解决不可行性解问题，但是在数据变化的时候，比如新增加 1 年，新时期加入数据集时，它整个前沿面都要重新构建，整个技术都会改变，该指数要被重新计算，计算量非常复杂，缺乏稳定性。BMLPI 是每两期作为一个生产技术来计算效率和生产率指数，因此增加了一个新的时间它也不用重复计算，前面算过的不变，只把后边新加入的算一次，并且它同时可以考虑技术退步的情况。

*BMLPI 的 Stata 命令是在 MLPI 命令基础上使用了biennial 这个 option，具体如下：
gtfpch labor capital energy= gdp: co2, biennial sav(ddfb_result, replace)

除此之外还有 MMLPI，即 Metafrontier Malmquist–Luenberger Productivity Index，共同前沿曼奎斯特罗恩伯格生产率指数。由 Oh (2010) 提出，它是基础的 MLPI 的延伸，用来解决组间异质性问题。由于每组之间的生产环境技术不同，不同的组的 DMUs 的绩效不能直接做比较，但是通过利用 MMLPI 就可以比较不同组的生产率的增长。此外，生产率的增长可以分解为 3 个部分，分别是技术效率变化、技术进步变化和技术差距变化。对于 MMLPI 的测度可以通过 GMLPI 分组计算。

4. 参考文献

• Asmild M, Paradi J C, Aggarwall V, et al. Combining DEA window analysis with the Malmquist index approach in a study of the Canadian banking industry[J]. Journal of Productivity Analysis, 2004, 21(1): 67-89. -Link-
• Banker R D, Charnes A, Cooper W W. Some models for estimating technical and scale inefficiencies in data envelopment analysis[J]. Management science, 1984, 30(9): 1078-1092. -Link-
• Chambers R G, Chung Y, Färe R. Benefit and distance functions[J]. Journal of economic theory, 1996, 70(2): 407-419. -Link-
• Charnes A, Cooper W W, Rhodes E. Measuring the efficiency of decision making units[J]. European journal of operational research, 1978, 2(6): 429-444. -Link-
• Charnes A, Cooper W W, Lewin A Y, et al. Extensions to DEA Models, Data Envelopment Analysis: Theory, Methodology, and Application [M]. Norwell, Massachusetts: Kluwer Academic Publishers, 1994.
• Chung Y H, Färe R, Grosskopf S. Productivity and undesirable outputs: a directional distance function approach[J]. journal of Environmental Management, 1997, 51(3): 229-240. -Link-
• Cobb C W, Douglas P H. A theory of production[J]. The American Economic Review, 1928, 18(1): 139-165. -Link-
• Debreu G. The coefficient of resource utilization[J]. Econometrica: Journal of the Econometric Society, 1951: 273-292. -Link-
• Fare R, Färe R, Fèare R, et al. Production frontiers[M]. Cambridge university press, 1994.
• Färe R, Lovell C A K. Measuring the technical efficiency of production[J]. Journal of Economic theory, 1978, 19(1): 150-162. -Link-
• Färe R, Primont D. Multi-output production and duality: theory and applications[M]. Springer Science & Business Media, 2012.
• Färe R, Grosskopf S, Pasurka Jr C A. Environmental production functions and environmental directional distance functions[J]. Energy, 2007, 32(7): 1055-1066. -Link-
• Farrell M J. The measurement of productive efficiency[J]. Journal of the Royal Statistical Society: Series A (General), 1957, 120(3): 253-281. -Link-
• Ray S C, Desli E. Productivity growth, technical progress, and efficiency change in industrialized countries: comment[J]. The American Economic Review, 1997, 87(5): 1033-1039. -Link-
• Tone K. A slacks-based measure of efficiency in data envelopment analysis[J]. European journal of operational research, 2001, 130(3): 498-509. -Link-
• Zhang and Choi Cooper, W. W., Seiford, L. M., & Tone, K. (2007). Data envelopment analysis: A comprehensive text with models, applications, references and DEA-solver software (2nd ed.). Boston: Kluwer Academic Publishers.

5. 相关推文

Note：产生如下推文列表的命令为：
  lianxh dea 效率, m
安装最新版 lianxh 命令：
  ssc install lianxh, replace

• 专题：其它
• Stata: 数据包络分析 (DEA) 简明教程

 

🍓 课程推荐：2024 机器学习与因果推断专题
主讲老师：司继春 (上海对外经贸大学) ；张宏亮（浙江大学）
课程时间：2024 年 11 月 9-10 日 ；16-17日
课程咨询：王老师 18903405450（微信）

 

尊敬的老师 / 亲爱的同学们：

连享会致力于不断优化和丰富课程内容，以确保每位学员都能获得最有价值的学习体验。为了更精准地满足您的学习需求，我们诚挚地邀请您参与到我们的课程规划中来。请您在下面的问卷中，分享您 感兴趣的学习主题或您希望深入了解的知识领域 。您的每一条建议都是我们宝贵的资源，将直接影响到我们课程的改进和创新。我们期待您的反馈，因为您的参与和支持是我们不断前进的动力。感谢您抽出宝贵时间，与我们共同塑造更加精彩的学习旅程！https://www.wjx.cn/vm/YgPfdsJ.aspx# 再次感谢大家宝贵的意见！

 

New！ Stata 搜索神器：lianxh 和 songbl  GIF 动图介绍
搜： 推文、数据分享、期刊论文、重现代码 ……
👉 安装：
  . ssc install lianxh
  . ssc install songbl
👉  使用：
  . lianxh DID 倍分法
  . songbl all

 

🍏 关于我们

• 连享会 ( www.lianxh.cn，推文列表) 由中山大学连玉君老师团队创办，定期分享实证分析经验。
• 直通车： 👉【百度一下： 连享会】即可直达连享会主页。亦可进一步添加 「知乎」,「b 站」,「面板数据」,「公开课」 等关键词细化搜索。