Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.
33,507 characters · 17 sections · 33 citation commands
Stochastic Frontier meets Breakdown Frontier
\nonstopmode
\address{Universidad ORT Uruguay}
{ Keywords: Stochastic Frontier Analysis, Technical inefficiency, Breakdown Frontier, Partial Identification.
JEL subject classification: C18, C2, C23, C51, D24 }
In this paper, we aim at responding the question, what are the maximum relaxations of the imposed assumptions on the unobservables which still sustain our conclusion in terms of efficiency of production units, given the observed data? Since its first introduction by aigner1977formulation and meeusen1977efficiency, stochastic frontier analysis has been applied to study the productivity and efficiency of production units in various economic sectors, such as banking ferrier1990measuring,adams1999semiparametric,kumbhakar2005measuring,malikov2016cost, healthcare zuckerman1994measuring,rosko2001cost,greene2004distinguishing,mutter2013investigating,comans2020cost and agriculture Bravo2007,Bravo2017, Trestini, Qushim2013,Ozden,Otieno,Gatti,Nwigwe, Qushim2018,Martinez,xin2025environmental,Lanfranco2013,Garcia2019,Garcia2022,Aguirre2024a,Aguirre2024b. The use of stochastic frontier models is still present in the current empirical practice, and significant extensions have been made to the methodology since the seminal works from aigner1977formulation and meeusen1977efficiency to account for more specific contexts, exploit more information, and to obtain more robust results. To name a few,schmidt1984production,greene2005fixed, colombi2014closed and kumbhakar2014technical contributed in the realm of panel data model extensions, including the estimation of the persistent and transient components of productive inefficiency. cornwell1990production,wang2002one,caudill1995frontier extended the analysis for models with determinants of inefficiency. simar2017nonparametric,wang2024flexible,centorrino2024nonparametric,zheng2024robust expanded the methodology for semi-parametric and non-parametric stochastic frontier models (SFM). See nguyen2022efficiency for a thorough review as well as software implementation details. Although these significant advances had made conclusions related to (in)efficiency become more robust to miss-specification, no sensitivity analysis per-se has been introduced in the literature. For sensitivity analysis masten2020inference we refer to a multidimensional approach that evaluates the robustness of empirical conclusions to simultaneous relaxations of a set of identifying assumptions. Typical approaches start with assumptions to identify a parameter, however, this approach begins with a specific conclusion and determines the weakest set of assumptions required for that conclusion to hold given the observed data. This is formalized in a breakdown frontier, which maps the tradeoffs between relaxing different assumptions while still maintaining the validity of the desired conclusion. This is relevant, firstly, because all the previous mentioned papers require certain identifying assumptions (no matter how slack they are) to provide a conclusion about (in)efficiency. Secondly, this analysis determines the threshold at which a straightforward model remains robust, tolerating sufficient deviation from baseline assumptions. Consequently, researchers may retain simpler models, maximizing both interpretability and implementation efficiency. Following masten2020inference we bridge this gap introducing a Breakdown frontier analysis in this context. The breakdown frontier generalizes the concept of an “identification breakdown point” introduced by horowitz1995identification, a one dimensional breakdown frontier. As noted by masten2020inference the breakdown frontier approach requires six main steps: (a) specifying a parameter of interest (typically in our context the inefficiency for a particular production unit), (b) specifying a set of baseline assumptions (for example truncated normality of the inefficiency term), (c) defining a class of assumptions indexed by a sensitivity parameter, which delivers a nested sequence of identified sets with the baseline assumptions obtained at one extreme and the no assumptions bounds obtained at the other, (d) characterizing identified sets for the parameter of interest as a function of the sensitivity parameter, (e) using those identified sets to define the breakdown frontier for a conclusion of interest, and (f) developing estimation and inference procedures for that frontier based on its characterization. We will follow these steps with particular focus on the assumptions from the seminal model from aigner1977formulation but this approach can be extended with proper modifications to the other models mentioned in this introduction above.
Let $y_i$ be the logarithm of an output of interest. Then, the classic model for stochastic frontiers states:
Where $f(\theta_x,x_i)$ is a known function of the inputs $x_i$ up to parameters $\theta_x \in \mathbb{R}^p$, $v_i$ is an error term, and $e_i$ is the unobservable component associated with technical (in)efficiency, where $e_i\geq 0$. It is common, in the stochastic frontier literature, to try to identify:
If one is interested in the estimates of individual (in)efficiency of a production unit, this measure captures it.
It is also common in the literature to impose additional structure on $v_i$ and $e_i$ to identify this parameter of interest. A typical set of assumptions is normality of $v_i$ and some form of truncated distribution on $e_i$, potentially with mean restrictions and heterosckedasticity. The intuition behind the sensitivity analysis is the following. Assume that there exist two unknown constants $b$ and $c$ that, respectively, measure the distance between the assumed distribution of $v$ and $e$ to the underlying or true distribution of these error terms. Then, in the context of the model presented above, suppose that the average value of technical inefficiency for a particular $u$ was estimated at $0.5$. The objective of the sensitivity analysis is to find the combinations of $c$ and $b$, which actually imply violations of the assumptions, such that $0.5$ is still a robust result. This is the so-called robust region (RR), in other words, a measure of how far the true distribution has to be from the imposed assumption in order to invalidate the estimated conclusion. The objective of the sensitivity analysis is to find, conditional on the estimated $0.5$ and a given value of $c$ that can be conceived as a tolerable violation of the efficiency measure, which is the value of the unknown constant $b$. This is the breakdown frontier, that is graphically depicted in Figure (ref).
The blue region are the set of combinations of $c$ and $b$, i.e. violations of the distributional assumptions, for which the conclusion of average (in)efficiency equal to some level holds. The yellow regions are the combinations for which we have not sufficient evidence in favor of the conclusions, given the imposed assumptions. The negative slope of the frontier implies that there is a tradeoff in reducing the violation of one of the assumptions, either on the error term ($v$) or in the efficiency term ($e$) that implies a larger violation of the other assumption. Or what is the same, if a relaxation of one of the assumptions is permitted, this comes at the benefit of a lower violation of the other assumption. The nonlinearity of the decreasing curve means that for a given relaxation of one assumption, say on $b$, it is required a reduction in the violation of the other ($c$) for the conclusions to hold, that varies depending on the value of $b$. Its convexity implies that a relaxation of $b$ when we allowed a small violation of that assumption permits a sizable reduction in the violation of the other assumption. And conversely, the same relaxation of $b$ when it is already broadly relaxed only allows for a small gain in the reduction of the violation of $c$.
Our starting point is to assume a relaxation of these distributional assumptions as a function of a parameter that measures the distance between these typically imposed baseline assumptions and the no-assumptions (and no identification of the parameter of interest) case.
Note that the case that $c=0,b=0$ would be the classic stochastic frontier model, which in the case that $v_i$ follows a Normal distribution and $e_i$ a truncated normal, we are under the case of aigner1977formulation. Moving away from $c=0$ and $b=0$ implies relaxing the structure of the distributional assumptions of the model.
In this section, we present how the identified set for the parameter of interest looks like under the relaxation of the classic assumptions ($b=0,c=0$). This result is collected in the following lemma.
Where $u=v-e$. Suppose that the average value of inefficiency for a particular $u$ from the model with $c,b=0$ was $0.5$, then this can be used to find the combinations of $c,b$ (of violations of the assumptions) such that $0.5$ is a robust result. This is the robust region, which would be a measure of how far away has to be the true values from the imposed assumption in order to reject our observed conclusion. This yields a measure of sensitivity.
Suppose that $e_i,v_i$ are independent of each other and also assume that $e_i\sim TN(\mu,\sigma^2_e)$ and $v_i\sim N(0,\sigma^2_v)$. Then:
Define then for the previously mentioned densities, $$ \int_{0}^{\infty } exp^{-e} f_{v_i|e_i}(u+e|e;\theta_v)de \equiv \Delta_1(u;\sigma_v), \int_{0}^{\infty } exp^{-e} f_{v_i|e_i}(u+e|e;\theta_v)f_{e_i}(e,\theta_e)de \equiv \Delta_1(u;\sigma_v,\sigma_e,\mu),$$ $$\int_{0}^{\infty } exp^{-e}f_{e_i}(e,\theta_e)de \equiv \Delta_1(u;\sigma_e,\mu), \int_{0}^{\infty } f_{v_i|e_i}(u+e|e;\theta_v)de \equiv \Delta_2(u; \sigma_v),$$ $$\int_{0}^{\infty } f_{v_i|e_i}(u+e|e;\theta_v)f_{e_i}(e,\theta_e)de \equiv \Delta_2(u;\sigma_v,\sigma_e,\mu).$$ Where their closed form expressions can be found as a result of Lemma (ref).
By an application of equation (ref):
We now formally define the breakdown frontier in this context. We also introduce the robust region, the area below the breakdown frontier. Suppose we run the usual Stochastic Frontier Model (with $c=0,b=0$) and we obtain that $E[exp^{-e_i}|v_i-e_i=u]=e_o$, we then for the breakdown frontier analysis begin with the conclusion that $E[exp^{-e_i}|v_i-e_i=u]\geq e_o$ for a fixed $e_o \in [0,1]$ (since the inefficiency is between 0 and 1). Relative to the baseline assumptions, what are the weakest assumptions that allow us to draw this conclusion, given the observed distribution of the data? Specifically, since larger values of $c$ and $b$ correspond to weaker assumptions, what are the largest values of $c$ and $b$ such that we can still definitively conclude that $E[exp^{-e_i}|v_i-e_i=u]\geq e_o$? This can be answered in two steps following masten2020inference. First, it requires gathering all values of $c$ and $b$ such that the conclusion holds. This is called the robust region. Since the lower bound of the identified set by Lemma (ref) and Equation (ref) is: $$ E[exp^{-e_i}|v_i-e_i=u] \geq \frac{\Delta_1(u; \sigma_v,\sigma_e, \mu)-c \Delta_1(u; \sigma_v)-b\Delta_1(u;\sigma_e, \mu)+bc}{\Delta_2(u; \sigma_v,\sigma_e, \mu)+c \Delta_2(u; \sigma_v)+b} $$ And since $c$ and $b$ can take any number between $0$ and plus infinity, we have that in our context the robust region is:
Note that when $c=0,b=0$, $e_o=\frac{\Delta_1(u; \sigma_v,\sigma_e, \mu)}{\Delta_2(u; \sigma_v,\sigma_e, \mu)}$ and since the previous function of $c,b$ is increasing in both of its arguments $c=0,b=0$ is included in the robust region and thus is not empty. The breakdown frontier is the set of points on the boundary of the robust region. Namely:
In order to express this in a simple graphical manner, we would like to see what are the values of $c$ for any value of $b$ that is in the frontier. This would be an equivalent way of characterizing the previous set. Thus, since the value of $e_o$ is fixed, solving for $b$:
This allows us to get the following analytical expression for the $BF(e_o)$:
In this section we study estimation and inference on the breakdown frontier defined above. The breakdown frontier is a known functional of the $\Delta$ objects which are known or estimable using classic maximum likelihood estimation as it is the case in aigner1977formulation model. Hence we propose plug-in estimators of the breakdown frontier. We can get $\sqrt{n}$-consistency and asymptotic distributional results using population smoothing procedures\footnote{ We can alternatively use the delta method for directionally differentiable functionals fang2019inference and implement a bootstrap procedure to construct asymptotically valid lower confidence bands for the breakdown frontier following the numerical delta method from hong2018numerical. Masten2019 denotes that there is no clear dominance between the numerical delta method and the population smoothing.} that recovers an approximation of this maximum that is not only directionally differentiable but fully differentiable, and therefore, standard asymptotic normality is obtained, and thus, the classic Nonparametric Bootstrap can be performed fang2019inference. Note that the main elements of the sensitivity frontier are the elements that can be estimated fron the "Null" model, a.k.a the model with $c=0,b=0$. Thus, we need to estimate $\theta=(\theta_x,\mu,\sigma_v,\sigma_e) \in \mathbb{R}^{p+3}$. This can be done via optimizing the following log-likelihood (for an $i.i.d$ sample of size $n$):
Where we obtain estimates for $\theta=(\theta_x,\mu,\sigma_v,\sigma_e)$ namely $\widehat{\theta}_n=(\widehat{\theta}_{xn},\widehat{\mu}_n,\widehat{\sigma}_{vn},\widehat{\sigma}_{en})$. Note that $ f(y_i,x_i;\theta)$ is the likelihood of observation $i$. By standard asymptotic results,
Where $ \mathcal{H}_{\theta}=-E[\frac{\partial ^2}{\partial \theta \partial \theta'} lnf(y_i,x_i;\theta)]$ or in more detail,
The variance covariance can be estimated consistently with $ \widehat{\mathcal{H}}_{\theta n}=-\frac{1}{n}\sum_{i=1}^n [\frac{\partial ^2}{\partial \widehat{\theta}_n \partial\widehat{\theta}_n'} lnf(y_i,x_i;\widehat{\theta}_n)]$. Note that the $\Delta$´s are all continuously differentiable functions of $\theta$, thus all of them have an asymptotically Gaussian behavior. Concretely, if $$\Delta=(\Delta_1(u;\sigma_v),\Delta_2(u;\sigma_v), \Delta_2(u;\sigma_v,\sigma_e,\mu), \Delta_1(u;\sigma_v,\sigma_e,\mu), \Delta_1(u;\sigma_e,\mu)) \in \mathbb{R}_+^{5}$$ And $$\widehat{\Delta}_n=(\Delta_1(u;\widehat{\sigma}_{vn}),\Delta_2(u;\widehat{\sigma}_{vn}), \Delta_2(u;\widehat{\sigma}_{vn},\widehat{\sigma}_{en},\widehat{\mu}_n), \Delta_1(u;\widehat{\sigma}_{vn},\widehat{\sigma}_{en},\widehat{\mu}_n), \Delta_1(u;\widehat{\sigma}_{en},\widehat{\mu}_n))$$ Then,
Where,
Furthermore, as long as $c-\Delta_1(u;\sigma_e,\mu)-e_o \neq 0$ it will then, by another application of the delta method, also be true that:
Where,
The next step is to obtain an asymptotic distribution for an estimator of $BF(e_o)$, but this is a non-differentiable function and thus the standard results break down. Here is where population smoothing enters. Population smoothing consists of using the smooth approximation of the maximum, that will allow to obtain a differentiable approximation to which we can apply standard bootstrap results. This approach is advantageous because we can use a non-parametric bootstrap to construct confidence intervals for the breakdown frontier. Consequently, we avoid estimating complex, convoluted variance components and bypass the need to handle a non-pivotal statistics.
Where as $\rho \rightarrow \infty$, $\frac{\sum_{k=1}^{K}m_{k}e^{\rho m_{k}}}{\sum_{k=1}^{K}e^{\rho m_{k}}} \rightarrow \max\{m_1,..m_K\}$. Thus, in our context, we focus on:
This is a differentiable function and thus, asymptotically normal point-wise in $c$ (or uniform under extra regularity conditions), then confidence intervals can be constructed by relying on the now asymptotic normality of:
In appendix (ref) we collect the influence function of our estimator which can then be used if of interest to construct a de-biased estimator of the form of:\footnote{See hines2022demystifying,fisher2021visually.} $$ \frac{\widehat{b}(c,e_o)_n e^{\rho \widehat{b}(c,e_o)_n}}{1+e^{\rho \widehat{b}(c,e_o)_n}} +\frac{1}{n}\sum_{i}^n \widehat{\psi}_{\frac{b(c,e_o) e^{\rho b(c,e_o)}}{1+e^{\rho b(c,e_o)}}}(x_i,y_i).$$ Where $\widehat{\psi}$ is the estimator of the influence function of the estimator.
If we pick $\rho$ to be a function of the data that grows at a rate slower than the rate our estimators converge, we can in a similar spirit as Bontempsetal2024 still obtain a consistent statistic for the maximum\footnote{ In fact we can construct a pivotal test statistic to invert and obtain confidence regions this is a result from Bontempsetal2024. Up to \today \quad the results from Bontempsetal2024 are not publicly available but we had knowledge of them from a presentation in IAAE-2024.} Thus, we can claim:
In this section, we provide a empirical illustration of the tool proposed in the previous sections, following the data used in nguyen2022efficiency, consisting of output and input information about rice producers in the Philippines.\footnote{Downloaded from the built in package \url{https://frontier.r-forge.r-project.org/}. Details of the data can be found in \url{https://www.rdocumentation.org/packages/frontier/versions/1.1-8/topics/riceProdPhil}.} The dataset includes the information about 43 rice producers in Tarlac, Philippines from 1990 to 1997. Out of the original dataset, we use one output $y_i$ (freshly threshed rice in tones) and three inputs $x_i$, the area planted (in hectares), labor used (man-days of family and hired labor), and fertilizer used (active ingredients in kilograms).
We run the SPF model with the output specification as $y_i=f(\theta_x,x_i)+v_i-e_i$, where the error term $v_i$ is assumed $v_i\sim N(0,\sigma^2_v)$, and the technical (in)efficiency $e_i$ is $e_i\sim TN(\mu,\sigma^2_e)$. It is also assumed that they are independent of each other. With the estimated frontier and for an arbitrary value of the composite error term $u_i=v_i-e_i$, say $u_i=0$, i.e. at the sample average, we compute the estimated efficiency as $E[exp^{-e_i}|v_i-e_i]=e_0=\frac{\Delta_1(u; \sigma_v,\sigma_e, \mu)}{\Delta_2(u; \sigma_v,\sigma_e, \mu)}=0.8872$. Then, allowing $b$ and $c$ to vary between zero and plus infinity, we can compute the breakdown frontier following equation (ref), which uses the estimated values of our $\sqrt{n}-$consistent plug-in estimators $\widehat{\Delta}_n$. Or equivalently, using equation (ref) and given $e_0=0.8872$, we find the values of $b(c,e_0)$ for each proposed value in the range of relaxations of $c$. The resulting breakdown frontier is presented in Figure (ref).
The downward slopping curve shows the tradeoffs between the relaxations of the distributional assumption of the efficiency $e$ and of the error term $v$, that maintain the conclusions ($e_0=0.8872$) to hold. The area under the curve is the robust region given in equation (ref), which is interpreted as all the combinations of violations of the imposed assumptions for which the mentioned conclusion hold. It is of particular interest to note that violations of the assumption on $e_i$ refers to a miss-specification of the distribution of the true (in)efficiency, however, violations on the assumption of $v_i$ includes both a miss-specification of the distribution but also a violation of the independence between the inefficiency and the error term.
The second part of the discussion in section 2.5 boils down to constructing confidence intervals for the breakdown frontier, allowed by the population smoothing procedures and performed with a nonparametric bootstrapping technique. In particular, we obtain confidence intervals for the breakdown frontier characterized by equation (ref), which required the use of population smoothing to sort the non-differentiability of this function. The confidence intervals are then constructed relying on the asymptotically normal distribution of $softBF(e_0)$ given in equation (ref), and by means of bootstrapping.
Heterosckedasticity is easily incorporated since the null model can be estimated with different assumptions on the variances and the mean of the errors.
The previous analysis focused on the distributional assumptions, but the breakdown frontier can be extended to incorporate relaxations of input exogeneity. A way of doing this would be to incorporate bounds on the following distances $|E(e_i|x)-E(e_i)|, |V(e_i|x)-V(e_i)|$.
The analysis has been exemplified in the case $\Delta$ can be analytically derived, but there will be cases in which expressions like $\int_{0}^{\infty} exp^{-e}f_{v_i|e_i}(u+e|e;\theta_v) f_{e_i}(e;\theta_e)de$ will not have a closed form solution (or the closed form solutions are not easy to implement). In such cases, we could leverage approximation results and obtain estimates for $\Delta$. For example if we obtain an estimate for $\theta_v,\theta_e$, $\widehat{\theta}_v,\widehat{\theta}_e$ we could generate via simulations a sample of size $S$ of realizations from $f_{v_i,e_i}(v,e;\widehat{\theta}_v, \widehat{\theta}_e)=f_{v_i|e_i}(v|e;\widehat{\theta}_v)f_{e_i}(e;\widehat{\theta}_e)$ and estimate $\int_{0}^{\infty} exp^{-e}f_{v_i|e_i}(u+e|e;\theta_v) f_{e_i}(e;\theta_e)de$ with $\sum_{s=1}^{S}exp^{-e_s}f_{v_i|e_i}(u+e_s|e_s;\widehat{\theta}_v) f_{e_i}(e_s;\widehat{\theta}_e)$. Under extra regularity conditions on the simulated sample we could (properly adjusting the variances) obtain similar asymptotic results.
The empirical application focuses mainly in the cross-sectional setting. Nevertheless, our procedure can be easily adapted to the panel data case in two directions. Firstly, the cross sectional analysis can be implemented year by year and thus construct a set of breakdown frontiers for each year. Secondly, in the presence of fixed effects, where the equation becomes $ y_{i,t}=f(\theta_x,x_{i,t})+v_{i,t}-e_{i,t}+c_i$ our model can be implemented on first differences $\Delta y_{i,t,t+1}=\Delta f(\theta_x,x_{i,t,t+1})+\Delta v_{i,t,t+1}-\Delta e_{i,t,t+1}$ in the spirit of chen2014consistent,belotti2012consistent. Future research can extend this methods to the cases of persistent and transitory shocks from colombi2014closed and kumbhakar2014technical.
In this paper we study sensitivity analysis of Stochastic Frontier Models. We elaborate relaxations of the baseline assumptions in the Stochastic Frontier Models and characterize the identified set under this relaxations. Furthermore, we derive the breakdown frontier for a relevant parameter of interest, the average inefficiency of a production unit. We show an application of the procedures on a well known dataset, and make the code available for the interested practitioner. For a particular average level of inefficiency we find a downward slopping curve which exhibits the tradeoffs between the relaxations of the distributional assumption of the efficiency $e$ and of the error term $v$, that maintain the conclusions. The nonlinearity of the decreasing curve means that for a given relaxation of one assumption, say on $b$, it is required a reduction in the violation of the other ($c$) for the conclusions to hold, that varies depending on the value of $b$. Its convexity implies that a relaxation of $b$ when we allowed a small violation of that assumption permits a sizable reduction in the violation of the other assumption. And conversely, the same relaxation of $b$ when it is already broadly relaxed only allows for a small gain in the reduction of the violation of $c$.