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.
65,113 characters · 20 sections · 107 citation commands
Modeling economies of scope in joint production: Convex regression of input distance function
\captionsetup[figure]{labelfont={bf},labelformat={default},labelsep=period,name={Fig.}} \captionsetup[table]{labelfont={bf},labelformat={default},labelsep=period,name={Table}}
\thispagestyle{empty}
\setcounter{page}{1} \setcounter{footnote}{0} \pagenumbering{arabic} \baselineskip 20pt
While the conventional microeconomic theory takes a single-product firm as a norm, there is growing empirical evidence that multiproduct firms are more common than usually assumed. Bernard2010 were among the first to argue that extensive reallocation of resources occurs through product switching by multiproduct firms: product switching is frequent, widespread, and influential in determining both aggregate and firm outcomes. A recent empirical study by Kuosmanen2023 suggest that approximately one-half of the Finnish manufacturing firms are multi-product firms, and approximately 20 percent of firms operate in multiple industry divisions.
Modeling the joint production of multiple outputs has proved a more vexing problem than most authors realize. In the deterministic setting where the input-output data contain no noise, the technology distance functions such as the input and output distance functions by Shephard1970 and the directional distance functions (DDF) by Chambers1996 (Chambers1996, Chambers1998b) provide theoretically sound representations of the technology, which can be conveniently estimated by linear programming techniques known as data envelopment analysis (DEA; Charnes1978, Charnes1978).\footnote{ The distance functions are reciprocals of Farrell's (Farrell1957) radial measures of technical efficiency. However, the distance functions are representations of technology, which can be used for analyzing marginal properties such as elasticities of substitution, transformation, or scale, even if one is not primarily interested in efficiency measurement. } By contrast, the situation becomes more difficult in the stochastic case, where data are perturbed by random noise. Most approaches known in the literature are based on hidden assumptions that are highly restrictive and unlikely to hold.
In mainstream economics, a commonly used approach to modeling multiproduct firms is to estimate a system of single-output production functions, dividing the firm's total inputs by outputs based on the revenue shares of outputs (Foster2008, Foster2008). However, such a proportional assignment of ratios of inputs to outputs seems completely arbitrary: a profit-maximizing multiproduct firm would demand inputs proportionate to the revenue shares of outputs only under very restrictive assumptions, which are unlikely to hold in the real world. De2016 propose an improved version of this approach, where they use the data of single-product firms to estimate production functions and then interpolate the estimated functions to model multiproduct firms. However, the linear interpolation approach assumes away possible synergies between multiple products: the resulting output isoquants are linear by construction. This seems counterintuitive due to economies of scope, in which physical synergies imply strictly convex output sets. Indeed, synergies are often seen as the main reason why firms find it beneficial to produce multiple products jointly rather than specialize in a single product (e.g., Panzar1981, Panzar1981). Therefore, assuming away the possible physical synergies from the outset seems restrictive.
In the econometric literature of stochastic frontier analysis (SFA), the modeling of joint production is usually based on Shephard's (Shephard1970) input and output distance functions (Lovell1994, Lovell1994; Coelli1999, Coelli1999; Kumbhakar2013, Kumbhakar2013). While the distance functions provide a general and theoretically sound framework for modeling joint production, the parametric functional forms usually employed in SFA tend to violate the theoretical properties of the distance function. For instance, it is easy to show that the Cobb-Douglas functional form for the distance function cannot support convex output sets at any parameter values, in fact, the output sets are not even bounded. The commonly used translog functional form is flexible enough to support convex output sets, but in that case, the output sets bend backward, violating the free disposability property. As a result, the translog distance functions cannot satisfy both convexity and free disposability axioms.
To our knowledge, the only theoretically sound approach to estimating technology distance functions of multiproduct firms is to resort to the shape-constrained nonparametric regression, following Kuosmanen2017a. Convex regression is a fully nonparametric method imposing shape constraints such as monotonicity, convexity, and homogeneity. The desired axiomatic properties of distance functions are thus guaranteed to hold. To disentangle technical inefficiency from random noise, Kuosmanen2017a employ the nonparametric kernel deconvolution method suggested by Hall2002. If the parametric distributional assumptions on inefficiency and/or noise are imposed, one can combine convex regression with the SFA techniques, as first outlined by Kuosmanen2012c, who refer to the blended application of nonparametric and parametric techniques as stochastic nonparametric envelopment of data (StoNED).
Thus far, the convex regression and StoNED estimation of distance functions have focused on either DDF by Kuosmanen2017a or a prior transformation on output variables by Schaefer2018. DDF is a general and theoretically valid representation of the technology, but its additive structure implies a somewhat limiting additive structure for the composite error term. In econometrics, the commonly used log-linear functional forms imply a multiplicative error term, where the size of the error is proportionate to the output. While a prior transformation on output variables is a meaningful attempt, the necessary orthogonality condition in Schaefer2018 can not be satisfied, and the approach is not theoretically sound, as we demonstrate below.
The purpose of this paper is to present a theoretically sound approach to model economies of scope in joint production using the convex regression approach. Our first methodological contribution is to develop a radial convex nonparametric least squares (CNLS) formulation to estimate the input distance function with the multiple input-multiple output specification. The orthogonality condition in radial CNLS will be stated and proven. We also extend the proposed radial CNLS approach to estimate output distance functions and the expected inefficiency.
Our second contribution is to compare the finite sample performance of radial CNLS versus the standard DEA and SFA estimators in terms of the input distance function estimation with multiple outputs in the controlled environment of the Monte Carlo study. The simulations reveal the superior performance of the radial CNLS estimator in the joint production setting involving multiple outputs.
The third contribution is to empirically apply the radial CNLS estimator to the Finnish electricity distribution firms' dataset. In contrast to the earlier DDF approach (Kuosmanen2017a, Kuosmanen2017a) or the CNLS with input requirement function approach (Kuosmanen2020d, Kuosmanen2020d), the proposed radial CNLS approach does not rely on any prior specification of the direction vector. To reduce the potential overfitting in shadow prices and increase the out-of-sample performance in energy regulation, we introduce the weight restriction constraints to radial CNLS. The results confirm that the input isoquants become more curved, implying better substitution possibilities. if we use two inputs rather than just project in the direction of one input as in the input requirement function. These results are highly relevant for the incentive regulation of electricity distribution firms in the real world.
The rest of this paper is organized as follows. Section (ref) introduces the theory on the input distance function and the regression model. Section (ref) presents a naive approach to estimating the input distance function and proposes the radial CNLS approach. The possible extensions are discussed in Section (ref). The Monte Carlo study and the empirical application of Finland's electricity distribution network regulation are demonstrated in Sections (ref) and (ref). Section (ref) concludes this paper with suggested avenues for future research. The formal proofs of Theorems and additional tables and figures are provided in the Online Supplement.
Consider a sample of firms where at least some are multiproduct firms, but a subset may specialize in just one or few products. We assume that each firm takes the output demands as given, and there is no noise in the observed outputs $\by \in \real_+^{s}$. The firm managers optimize some meaningful objective functions, e.g., minimizing cost. Note that cost minimization is equivalent to profit maximization when the output demands are taken as given. The observed inputs $\bx \in \real_+^{m}$ are subject to random perturbances $\varepsilon_i$, and the firm manager then adjusts their input demands according to
where $\bx^*$ is the optimal (cost minimizing) input vector.
Consider further the following production possibility set with the generated data $(\bx, \by)$.
where we assume that $T$ is a closed and nonempty set satisfying free disposability (monotonicity) and convexity axioms. See Fare1995 for further detailed discussions on the axiomatic production theory.
The Shephard's input distance function $D_i: \real_+^{m} \times \real_+^{s} \rightarrow \real_+^{1} \cup \{+ \infty\}$ is defined as Shephard1970:
where the construction of input distance function $D_i$ does not depend on behavioral hypotheses such as profit maximization or cost minimization.
Similar to Shephard1970, we consider $D_i$ as a general functional representation of multi-output technology: If $D_i(\bx, \by)=1$, then production plan $(\bx, \by)$ is technically feasible and efficient. If $D_i(\bx, \by)<1$, then $(\bx, \by)$ is technically feasible but inefficient. If $D_i(\bx, \by)>1$, then $(\bx, \by)$ is technically infeasible. Therefore, $(\bx, \by)$ for which $D_i(\bx,\by) \le 1$ are within the production possibility set and $(\bx, \by)$ for which $D_i(\bx, \by)=1$ are its efficient subset.
Linear homogeneity of $D_i$ is a fundamental property that does not relate to the properties of $T$, but arises from the definition of $D_i$. The linear homogeneity implies that
This property helps develop the following convex regression model of the input distance function. The value of the input distance function $D_i(\bx_i, \by_i)$ is then completely determined by the error term $\varepsilon_i$ (see Theorem (ref)).
This result is analogous to Proposition 2 in Kuosmanen2017a, adapted to the present setting of $D_i$. Note that even though the distance function is a deterministic representation of a deterministic technology, its value is a random variable when the observed data are perturbed by random noise and inefficiency.
To derive the regression equation, we can utilize the homogeneity property (ref) and normalize by $x_{1i}$ to obtain\footnote{ Note that any one of the input variables can be selected in the distance function normalization. }
where $\bx/x_1 = (1, x_2/x_1, \ldots, x_m/x_1)^\prime$. After taking the logs, we have $\ln x_{1i} = -\ln D_i(\bx_i/x_{1i}, \by_i) + \ln D_i(\bx_i, \by_i)$.
Using Theorem (ref), we substitute $\ln D_i(\bx_i, \by_i)$ with the error term $\varepsilon_i$ and obtain
A similar regression equation is commonly used in the SFA literature (e.g., Kumbhakar2013, Kumbhakar2013). The main difference is that we derive the regression equation without imposing any particular functional form for $D_i$. Note that the parametric functional forms usually employed in SFA tend to violate the theoretical properties of the distance function such as monotonicity and convexity of the output sets. In the SFA literature, $\varepsilon_i$ consists of asymmetric inefficiency $u \sim N^+(0, \sigma_u^2)$ and symmetric noise $v \sim N(0, \sigma_v^2)$. For simplicity, we here mainly focus on the conventional case of a symmetric error term, but the extension to stochastic frontier setting with an asymmetric composite error term is briefly examined in section (ref) below.
To estimate the regression equation (ref) developed in the previous section by convex regression, it would be tempting to formulate the CNLS problem as follows
The first two constraints of (ref) are used to characterize the regression equation (ref). The decision variable $\chi_i$ denotes the distance function $D_i(\bx_i/x_{1i}, \by_i)$. The second constraint can thus be interpreted as a multivariate linear regression among the normalized inputs, outputs, and distance function. The third constraint imposes the concavity on the distance function. The last set of constraints guarantees the monotonicity of the distance function.
Problem (ref) is henceforth referred to as a naive CNLS because the estimated CNLS residuals do not satisfy the usual sample orthogonality condition of the linear regression. In the parametric regression, the linear homogeneity of $D_i$ is satisfied by construction because the regression residuals $\hat{\varepsilon}$ satisfy the following sample orthogonality condition
which is the sample counterpart of the population orthogonality condition
This implies the equiproportionate radial orientation of $D_i$, but as noted above, the CNLS residuals do not generally satisfy these orthogonality conditions by construction. Therefore, we need to find a way of imposing the orthogonality for the input ratios $(\bx_i/x_{1i})$.
We note that Schaefer2018 propose an extended CNLS approach to estimate the ray production function in the multiple input-multiple output specification. This approach relies on the polar coordinate transformation to introduce the angle-related variables and generates a “multiple input-single output” model. While it is a meaningful attempt to model multiple outputs in the framework of CNLS, similar to problem (ref), the proposed formulation by Schaefer2018 does not impose the orthogonality conditions and therefore fails to satisfy the linear homogeneity property of the radial distance functions.
To gain intuition, consider first the special case of two inputs. We utilize the homogeneity property and normalize inputs by the Cobb-Douglas function $(x_{1i}^{1-d}x_{2i}^d)$, $0 < d < 1$, to obtain
where we can use any arbitrary parameter value $d$ (e.g., $d = 0.5$). We will later ensure that the specific choice of $d$ does not matter.
Taking logs and reorganizing, we now have the following regression equation
To estimate the equation (ref), we propose the following radial CNLS approach to model multiple input-multiple output setting
where problem (ref) is a semi-nonparametric partial linear formulation building on Johnson2011 (Johnson2011, Johnson2012a). Similar to problem (ref), the first two sets of constraints to represent the regression equation (ref), the last two sets of constraints impose the concavity and monotonicity on the input distance function, respectively.
In contrast to parameter $d$ that we simply postulate in equation (ref), $\delta$ is a decision variable to be estimated. Theorem (ref) shows that the semi-nonparametric partial linear formulation satisfies the orthogonality condition and hence $\delta$ can be a free decision variable in problem (ref).
Theorem (ref) shows that problem (ref) effectively enforces the equiproportionate radial orientation, in other words, the linear homogeneity, which is an essential property of $D_i$. We thus can consistently estimate the regression equation (ref) using the nonparametric radial CNLS formulation (ref).
Regarding the choice of parameter $d$, it is clearly possible or perhaps even likely that the coefficient $\delta$ obtained as the optimal solution to problem (ref) differs from our arbitrary specification of $d$. But it is unnecessary to guess the optimal value of $d$ beforehand. If we solve the problem (ref) for multiple different values of $d$, it is easy to verify that the parameter estimates of $\delta$ change, however, the residuals $\hat{\varepsilon}_i$ remain exactly the same irrespective of the choice of $d$. It turns out that $\delta$ is an over-identified parameter value that does not influence the optimal solution. More specifically, we need to consider $-\ln(\chi_i)+\delta(\ln x_{1i} - \ln x_{2i})$ as the estimate of $D_i$, not only the nonparametric part $-\ln(\chi_i)$. Note that the sum of two concave functions is itself a concave function. As the value of parameter $\delta$ changes, so does that of the nonparametric part, but their sum remains the same. This can be referred to as partial identification in the present context (see, e.g., Manski2009, Manski2009; Tamer2010, Tamer2010).
To find a unique estimate of $D_i$ in any given point $(\bx_i, \by_i)$, we can solve the following linear programming (LP) problem, directly analogous to the multiplier formulation of DEA
The adjustment by the CNLS residuals $\hat{\varepsilon}_i$ of problem (ref) effectively projects all observations to the estimated production frontier (in this case the average practice technology). The optimal solution to this LP problem provides the shadow prices $\bgamma$ and $\bbeta$ associated with the minimum extrapolation technology.
In this section, we consider the general case with more than two inputs, connect with the efficiency analysis, and estimate the output distance function.
To generalize the results of the previous section, consider the general case of $M$ inputs (i.e., $M > 2$). Let us normalize $D_i$ by the Cobb-Douglas function
In this case, the regression equation (ref) becomes
and the radial CNLS estimator can be stated as
This formulation is a more general CNLS estimation of the input distance function with multiple outputs. That is, problem (ref) is a special case of problem (ref), where $M=2$. Having solved the radial CNLS problem (ref), the same LP problem (ref) can be applied to obtain the unique estimate of $D_i$.
Given the residuals $\hat{\varepsilon}_i$ estimated by radial CNLS (ref), it is possible to estimate the expected inefficiency using the nonparametric kernel deconvolution or if one imposes further parametric distributional assumptions by quasi-likelihood or the method of moments. Following Kuosmanen2012c and Kuosmanen2017a, we in this section briefly summarize how to estimate the unconditional expected inefficiency $\E(u)$ and the conditional expected inefficiency $\E(u_i \mid \hat{\varepsilon}_i)$.
Nonparametric kernel deconvolution is a full nonparametric estimation of the expected inefficiency $\E(u)=\mu$. For the input distance function, the residuals $\hat{\varepsilon}_i$ are the consistent estimators of $\hat{\varepsilon}_i = \mu + (v_i - u_i)$. The density function of $\hat{\varepsilon}_i$ is defined as
where $K(\cdot)$ is a standard kernel function, $h$ is the bandwidth, $x$ is the projected data of $\hat{\varepsilon}$, and $n$ is the number of observations. The first derivative of the density function of the composite error term $(f_{\hat{\varepsilon}}^{\prime})$ is proportional to that of the inefficiency term ($f_u^{\prime}$) in the neighborhood of $\mu$ (Hall2002, Hall2002). Therefore, $\hat{\mu} = \arg \max_{x \in C}(\hat{f}_{\hat{\varepsilon}}^{\prime}(x))$ provides a nonparametric estimation of the expected inefficiency $\mu$, where $C$ is a closed interval in the right tail of $f_{\varepsilon}(\cdot)$.
Alternatively, the method of moments and the quasi-likelihood approaches for the residual decomposition build upon additional parametric distributional assumptions (e.g., the half-normal inefficiency and normal noise). For the method of moments approach, we first calculate the second and third central moments of the residuals $\hat{\varepsilon}_i$ by
Under the parametric distributional assumptions, we have $M_3 = (2/\pi)^{1/2}(1-4/\pi)\sigma_u^2$ and $M_2 = \sigma_v^2+\sigma_u^2(\pi-2)/\pi$. Given the estimated $\hat{M}_2$ and $\hat{M}_3$, we then calculate $\sigma_u$ and $\sigma_v$. Formally,
The quasi-likelihood estimation is an alternative to compute $\hat{\sigma}_u$ and $\hat{\sigma}_v$. We apply the standard maximum likelihood method to estimate the following likelihood function
where
The likelihood function $\ln L(\lambda)$ consists of a single parameter $\lambda$. $\Phi$ denotes the cumulative distribution function of normal distribution. After obtaining $\hat{\sigma}_u$ by these two parametric approaches, we can compute the expected inefficiency, $\mu = (2/\pi)^{1/2}\hat{\sigma}_u$.
Jondrow1982 propose a widely applied estimator for calculating the conditional expected inefficiency $\E[u_i \mid \hat{\varepsilon}_i]$. Specifically, for the input distance function, the conditional expected value of inefficiency $\E[u_i \mid \varepsilon_i]$ is formulated as
where $\mu_{*i}= -\hat{\varepsilon}_i \sigma_u^2/\sigma^2$, $\sigma_*^2 = \sigma_u^2\sigma_v^2/\sigma^2$, $\lambda = \sigma_u/\sigma_v$, and $\sigma^2 = \sigma_u^2 +\sigma_v^2$. $\phi$ and $\Phi$ are the standard normal density function and its cumulative distribution function, respectively. Note that the conditional mean inefficiency in Jondrow1982 could be further extended to conditional quantile inefficiency.
Consider now the general output distance function
where $D_i^O(\bx, \by)$ is homogeneous of degree 1 in $\by$, implying that
Given the case of $S$ outputs (i.e., $S \ge 2$) and the linear homogeneity, similar to equation (ref), the regression equation can be formulated as
We then estimate the output distance function using the following radial CNLS formulation
where the first two constraints characterize the regression model (ref). In contrast to problem (ref), the third set of constraints of problem (ref) impose the convexity of the output distance function. The last set of constraints also ensures the monotonicity of the output distance function.
The main challenge here is that the output isoquants of the nonparametric part $-\ln(\chi_i)$ are concave by construction, whereas those of the parametric part $\sum_{s=2}^{S}\delta_s(\ln y_{1i} - \ln y_{si})$ are convex by default. In contrast to the input distance function where both parts have the same curvature, in this case, there is no guarantee that the estimated output sets satisfy convexity.
If the data are well-behaved, the nonparametric part is flexible enough to offset the wrong curvature of the Cobb-Douglas part. In other words, the formulation (ref) could still satisfy the convexity of the output sets, however, this cannot be guaranteed always to hold. Note that possible violations of convexity of the output sets can always be fixed in the next step, where we apply the LP problem (ref) to estimate the output distance function. We leave a more detailed examination of the output distance function as an interesting avenue for future research.
In the Monte Carlo study, our main objective is to investigate the finite sample performance of the proposed radial CNLS approach and compare it with the naive CNLS approach and other commonly seen deterministic and stochastic frontier estimation approaches in estimating the input distance function.
We generate the two input-two output production datasets by using the following two sets of DGPs (Fare2010, Fare2010)
where the inputs $\bx$ in each DGP are independently and randomly drawn from the uniform distribution, $U[0, 1]$, and the error term $\varepsilon$ has three different specifications: $\varepsilon = v$ (i.e., only noise), $\varepsilon = -u$ (i.e., only inefficiency), and $\varepsilon = v - u$ (i.e., both noise and inefficiency). We generate noise $v \sim N(0, \sigma_v^2)$ and inefficiency $u \sim N^+(0, \sigma_u^2)$.
The first output $y_1$ in DGP I is randomly generated from the following two different gamma distributions
where $\alpha$ and $\beta$ are the shape and rate parameters of the gamma distribution. For DGP II, the first output $y_1$ is drawn from $U[e^{0.7}, e^{1.4}]$.
Given the different true parameter combinations ($\beta_k$, $k=0,\ldots,4$) in Table (ref), DGPs I and II consist of three polynomial and translog technologies, respectively. We thus consider 9 models with DGPs I and II in the experiments. Note that Model 3 in each DGP makes the true functions in terms of outputs more concave than the other two models, and Type-B in DGP I generates a more balanced outputs dataset (i.e., the value of generated $y_1$ is close to that of $y_2$).
In all experiments that follow, we resort to the Python/pyStoNED package (Dai2021b, Dai2021b) to solve the CNLS and DEA models and the Python/pySFA package to estimate the SFA models.\footnote{ A Python Package for Stochastic Frontier Analysis (pySFA): \url{https://github.com/gEAPA/pySFA}. } Note that the linear and nonlinear programming problems are solved by the off-the-shelf solvers MOSEK (9.3.11) and KNITRO (13.2), respectively. All experiments are run on Finland's high-performance computing cluster Puhti with Xeon @2.1 GHz processors, 2 CPUs, and 5 GB of RAM per task.
In experiment 1 we compare the performance of the proposed radial CNLS formulations (ref) and (ref) with the naive CNLS formulation (ref) in the absence of inefficiency (i.e., $\varepsilon = v$). We thus consider 135 scenarios with different numbers of observations, $n \in \{25, 50, 100, 200, 400\}$, and the variation of noise levels, $\sigma_v \in \{0.075, 0.15, 0.3\}$. To evaluate their finite sample performance, each scenario is run 500 times to calculate mean square error (MSE) and mean absolute deviation (MAD) (see the Online Supplement) in terms of input distance function estimation.
Figs. (ref) and (ref) depict the MSE results of radial CNLS and naive CNLS in DGP I with Type-A and DGP II. Radial CNLS outperforms naive CNLS in all scenarios considered, suggesting that radial CNLS is a more efficient method for estimating the input distance function. In contrast to arbitrarily choosing the input, normalizing the input distance function can satisfy the orthogonality conditions, which probably helps fit the true function with smaller expected errors. Several other interesting findings are also observed.
First, the curvature of the generated production function is likely to affect the performance of estimators. For instance, compared to the other two models in DGPs I and II, Model 3 has the smallest MSE values due to that the generated production function is more concave in terms of the outputs. Second, the MSE values increase as the noise becomes large. It is clearly evident from Figs. (ref) and (ref) that the estimation of the input distance function becomes worse when $\sigma_v$ increases. Third, compared to DGP I, DGP II has a smaller MSE difference between the two estimators. This is because the generated dataset by DGP II is more balanced (see Figure 1 in Fare2010, Fare2010), which dampens the disadvantage of the naive CNLS approach. Further, Figs.B1 and B2 in Appendix B demonstrate the MAD results of radial CNLS and naive CNLS in two sets of DGPs, showing similar findings as in the MSE comparison and supporting the main conclusion of experiment 1.
Given the potential influence of the curvature of the true output isoquant, we further compare two CNLS approaches based on DGP I with Type-B, where a more balanced dataset in terms of outputs is generated. As shown in Table (ref), the MSE values of naive CNLS are always larger than those of radial CNLS.
We proceed to compare the performance among radial CNLS, DEA, and SFA. Following Schaefer2018, we consider 180 scenarios in experiment 2 with different numbers of observations, $n \in \{25, 50, 100, 200, 400\}$, the noise-to-signal ratios, $\rho_{nts} \in \{0, 0.5, 1, 2\}$, where $\rho_{nts}=\sigma_v/\sigma_u$, and $\sigma_u=0.15$. Each scenario is also run 500 times to compute the MSE and MAD metrics. Note that the input-oriented DEA model with variable returns to scale and SFA with Cobb Douglas function (SFA-CD) and translog function (SFA-TL) specifications are considered in the benchmarking procedures.
Tables (ref) and B1 in Appendix B demonstrate the performance comparison of radial CNLS, DEA, and SFA with $n=400$. This experiment further confirms the superiority of the proposed radial CNLS approach in estimating input distance function compared to other available methods. Specifically, the radial CNLS approach has the lowest MSE and MAD values in all scenarios. Notably, in the case of the DGPs with inefficiency only, radial CNLS also works best compared to its counterparts. The MSE and MAD values of the SFA-CD approach are larger than those of others and increase much more sharply at the high noise-to-signal ratio of $\rho_{nts}=2$ (cf. Schaefer2018, Schaefer2018). This is because when estimating the Cobb-Douglas function using SFA, the considered DGPs are misspecified, especially in polynomial technologies. Note that the results show that radial CNLS and DEA are relatively robust to polynomial and translog technologies. Furthermore, as expected, the values of MSE and MAD increase as the noise increases (see, e.g., Henningsen2015, Henningsen2015; Schaefer2018, Schaefer2018; Ahn2023, Ahn2023).
The most significant real-world application of CNLS has thus far been incentive regulation of electricity distribution firms (Kuosmanen2012b, Kuosmanen2012b). This literature introduces the CNLS estimator that combines the virtues of DEA and SFA to the Finnish electricity distribution network regulation. But most existing applications of CNLS or other frontier estimation techniques focus on either variable cost or total cost and hence have demonstrable shortcomings (Kuosmanen2020d, Kuosmanen2020d). To consider variable cost and fixed cost simultaneously in the benchmark regulation, Kuosmanen2017a propose a new CNLS approach with DDF (i.e., DDF CNLS), and Kuosmanen2020d and Kuosmanen2022 further develop a CNLS approach to incorporate the input requirement function (i.e., IRF CNLS).
However, the DDF CNLS approach requires a prespecified direction vector, and the IRF CNLS approach needs to project in the direction of one input. The estimates by DDF CNLS or IRF CNLS are not invariant to the direction vector or the input choice. This motivates us to propose and apply the radial CNLS approach to energy regulation practice, where radial CNLS does not rely on any direction vector and is immune to the input selection.
Following Kuosmanen2020d, we consider the following the multiplicative cost frontier model with contextual variables
where $C$ is a non-decreasing convex conical cost function. $\bx$, $\by$, and $\bb$ are inputs, desirable outputs, and undesirable outputs, respectively. $\bz$ represents a set of contextual variables, and $\blambda$ denotes the associated marginal coefficients. $\varepsilon$ is a composite error term involving inefficiency and random effects. Accordingly, the input, output, and contextual variables of the cost frontier model (ref) are specified as follows. \setlist{nolistsep}
While radial CNLS (ref) is ready for estimating the multiple input-multiple output cost frontier model (ref), it might overfit training data and predictably perform poorly on testing data.\footnote{ Overfitting is a longstanding problem in convex regression and other general nonparametric regression, where the subgradients may become very large at the boundary of the convex hull of the design points (Liao2023, Liao2023). } Considering that the benchmark regulation generally involves future economic incentives, we set a lower bound and upper bound on each subgradient to control their magnitudes to increase the performance of radial CNLS (Kuosmanen2022, Kuosmanen2022). Other restrictions on fitted subgradients for reducing overfitting include the Lipchitz norm (e.g., Mazumder2019, Mazumder2019) and $L_2$ norm (e.g., Dai2023c, Dai2023c). One could even respecify the loss function to mitigate the effect of overfitting (e.g., Liao2023, Liao2023).
After taking into account these specifications and choosing $x_1$ as numeraire, we rephrase radial CNLS (ref) to estimate cost frontier model (ref) and have the following nonlinear programming problem
where the fourth set of constraints refers to the weight restrictions on subgradients, and $\mu_{it}$ denotes the marginal coefficient of the outages. The lower and upper bounds (e.g., $\bbeta_{lo}$ and $\bbeta_{up}$) in problem (ref) can be determined by either a data-driven approach (e.g., cross-validation) or decision-makers and/or stakeholders. A major practical advantage of the latter is that the weight restrictions can be easily communicated to decision-makers (e.g., the regulator) and stakeholders (e.g., the regulated firms and their customers) as they can literally see the weight restrictions themselves and comment if those are too loose, too restrictive, or just fine.
In practice, we first solve problem (ref) without weight restrictions (i.e., the fourth sets of constraints) to obtain the subgradients estimates (i.e., $\hat{\bbeta}$, $\hat{\bgamma}$, and $\hat{\mu}$) and then take 10% and 90% quartiles of each subgradient as its lower and upper bounds to reestimate problem (ref).
In this section we discuss what degree of difference in the estimates with and without weight restrictions and compare radial CNLS with its alternative. We thus apply the proposed radial CNLS approach to a panel of 77 Finnish electricity distribution firms in the years 2008--2020.\footnote{ The original dataset is applied to carry out the incentive regulation for Finnish electricity distribution networks (Kuosmanen2022, Kuosmanen2022), and its earlier version has been widely used in, e.g., Kuosmanen2012b, Kuosmanen2013, and Kuosmanen2020d. } The descriptive statistics for inputs, outputs, and contextual variables are summarized in Table B2 (Appendix B). See Kuosmanen2020d for a detailed introduction to these selected variables.
Table (ref) summarizes the results estimated by the basic radial CNLS model and weight-restricted radial CNLS model. Note that the estimated input distance function is a piece-wise linear function consisting of hyperplanes characterized by subgradients $\bbeta_i$, $\bgamma_i$, and $\mu_i$. When weight constraints are introduced to the radial CNLS model, the average estimates for the nonparametric part dramatically decrease, and the corresponding standard deviations also decline. This implies that a small part of firms has extremely estimated shadow prices (i.e., $\hat{\bbeta}_i$, $\hat{\bgamma}_i$, or $\hat{\mu}_i$), which can severely affect the fitted input distance function in a testing set and hence deteriorate the performance of economic incentives in energy regulation.
By comparing the 10%, 50%, and 90% quartile estimates, we observe that the changes in estimated shadow prices remain rather marginal, suggesting that most of the firms would not be highly affected by additional weight constraints. Furthermore, the estimated residual comparison also shows little impact of weight restrictions on the firm's relative performance in that most firms cluster around the 45-degree line (see Fig. (ref)). That is, while the applied weight-restricted CNLS model will affect the estimated input distance function, it can avoid extreme shadow prices and show its robustness in inefficiency estimation and endogeneity bias.
It is worth noting that the minimal estimated $\hat{\mu}_i$ is negative in the basic radial CNLS model but is positive in the weight-restricted radial CNLS model. Unsurprisingly, the marginal effect of outages at the firm level can be positive or negative in the context of the benchmark regulation (Kuosmanen2020d, Kuosmanen2020d). Recall that there is huge heterogeneity between Finnish electricity distribution firms in terms of variable and fixed costs as reflected by their standard deviation (see Table B2 in Appendix B). Firms with low outages can use a high level of quality and reliability of distribution as a competitive advantage in the yardstick competition. For these firms, the shadow price of outages is negative. Firms that face exceptionally severe weather shocks have higher operational costs due to large outages, deducing the positive shadow price of outages. This is because the standard compensations paid to customers for interruptions are also included in the variable cost. In practice, the threshold value for a positive outage effect is very high because the outage values must be exceptionally high compared to other firms such that the firm's operations appear competitive in the benchmark regulation between firms. Therefore, the impact of outages on the distribution firms is a “U-shaped” curve.
Furthermore, additional weight constraints also influence the estimated coefficient of the second contextual variable, energy loss. In the weight-restricted radial CNLS model, the value of the estimated coefficient of energy loss is 3.98, slightly lower than the value of 4.11 obtained from the basic radial CNLS model. Since the share of energy loss is on average higher in sparsely populated areas than in densely populated areas, the distorting effect of possible efficiency differences can be partially canceled out from the energy loss estimate with the help of weight restrictions.
Fig. (ref) demonstrates the input isoquants of the three largest distribution firms (Caruna, Elenia, and Helen) and the average level of all firms based on radial CNLS and IRF CNLS estimates. The horizontal axis is fixed cost (NKA), and the vertical axis is variable cost (KOPEX): both are rescaled in Million \officialeuro. Note that the shape of the input isoquant is determined by the output structure of the firm in the multiple output setting. Given the output structure of the three largest distribution firms, both two subfigures illustrate the considerable substitution possibilities between NKA and KOPEX. But the right figure shows that the average firm yields almost a Leontief-type input isoquant given the averaged output structure. Furthermore, it can be seen from Fig. (ref) that the input isoquants become more curved (i.e., more substitution possibilities) when we use two inputs in the radial CNLS model rather than projecting in the direction of one input in the IRF CNLS model.
Modeling of joint production is a more vexing problem than most authors realize. The axiomatic nonparametric approach is the only theoretically sound alternative when data are perturbed by random noise, and the parametric approach is likely to violate the theoretical properties of the distance function. Modeling radial input and output distance functions in CNLS/StoNED has also been an enigma for almost two decades. Finally, in this paper the problem has been solved in the case of the input distance function, but the output distance function partially remains.
This paper proposes a new radial CNLS approach to estimate the input distance function in the multiple input-multiple output specification. We demonstrate how to transform the input distance function correctly and show that the orthogonality conditions can be satisfied in the developed approach. We then discuss the possible extensions from the other three aspects: considering the general case with more than two inputs, connecting with the efficiency analysis, and estimating the output distance function.
Our simulations compare the finite sample performance of radial CNLS and other deterministic and stochastic frontier approaches regarding the input distance function estimation. Radial CNLS outperforms the naive CNLS, DEA, and SFA approaches in all scenarios considered, suggesting that radial CNLS is a more efficient method for estimating the input distance function.
The developed approach is further applied to the Finnish electricity distribution network regulation. To reduce the potential overfitting and increase the out-of-sample performance in energy regulation, we introduce weight restrictions on shadow prices to the framework of radial CNLS. In contrast to other traditional model specifications (e.g., DDF and IRF), radial CNLS does not rely on any direction vector and is immune to the input selection. The illustrated input isoquants suggest more substitution possibilities.
Future research could further examine joint production from the perspectives of nonparametric identification (see, e.g., Matzkin2007, Matzkin2007; Chiappori2015, Chiappori2015; Heckman2010, Heckman2010) and partial identification (see, e.g., Manski2009, Manski2009; Tamer2010, Tamer2010). Axioms of production theory (linear homogeneity) provide moment conditions (i.e., generalized method of moments), but the model remains inherently underidentified. Underidentification is conventionally seen as a failure in econometrics (Arellano2012, Arellano2012), but in many relevant applications, partial identification may be sufficient for the task at hand. One possible avenue of future research is to impose moment conditions supported by theory (linear homogeneity) and apply additional statistical criteria (e.g., maximum likelihood and least squares) to achieve (partial) identification.
The authors wish to acknowledge CSC – IT Center for Science, Finland, for computational resources. Sheng Dai gratefully acknowledges financial support from the OP Group Research Foundation [grant no. 20230008] and the Turku University Foundation [grant no. 081520].
\baselineskip 12pt
\baselineskip 20pt