Theory and Methodology Estimating and bootstrapping Malmquist indices 1 Leopold Simar a,2 , Paul W. Wilson b, * ,3 a Institut de Statistique and CORE, Universite Catholique de Louvain, Voie du Roman Pays 20, Louvain-la-Neuve, Belgium b Department of Economics, University of Texas at Austin, Austin, TX 78712-1173, USA Received 1 August 1996; accepted 1 October 1997 Abstract This paper develops a consistent bootstrap estimation procedure for obtaining confidence intervals for Malmquist indices of productivity and their decompositions. Although the exposition is in terms of input-oriented indices, the techniques can be trivially extended to the output orientation. The bootstrap methodology is an extension of earlier work described in Simar and Wilson (Simar, L., Wilson, P.W., 1998, Management Science). Some empirical examples are also given, using data on Swedish pharmacies. Ó 1999 Published by Elsevier Science B.V. All rights reserved. Keywords: DEA; Productivity; Resampling; Bootstrap; Malmquist indices 1. Introduction Fa ¨re et al. (1992) merge ideas on measurement of eciency from Farrell (1957) and on measure- ment of productivity from Caves et al. (1982) to develop a Malmquist index of productivity change. Caves et al. define their input-based Malmquist productivity index as the ratio of two input distance functions, while assuming no tech- nical ineciency in the sense of Farrell. Fare et al. extend the Caves et al. approach by dropping the assumption of no technical ineciency and devel- oping a Malmquist index of productivity that can be decomposed into indices describing changes in technology and eciency. We extend the Fare et al. approach by giving a statistical interpretation to their Malmquist productivity index and its components, and by presenting a bootstrap algo- rithm which may be used to estimate confidence intervals for the indices. This work will allow re- searchers to speak in terms of whether changes in productivity, eciency, or technology are signifi- cant in a statistical sense. In other words, our methods can be used to determine whether indicated changes in productivity, eciency, or European Journal of Operational Research 115 (1999) 459–471 * Corresponding author. Tel.: 512 471 3211; fax: 512 471 3510; e-mail: [email protected]. 1 Pontus Roos graciously provided the data used in the empirical examples. We alone, of course, are responsible for any remaining errors or omissions. 2 Research support from the contract ‘‘Projet d’Actions de Recherche Concert ees’’ (PARC No. 93/98–164) of the Belgian Government is gratefully acknowledged. 3 Research support from the Management Science Group, US Department of Veterans Aairs, is gratefully acknowledged. 0377-2217/99/$ – see front matter Ó 1999 Published by Elsevier Science B.V. All rights reserved. PII: S 0 3 7 7 - 2 2 1 7 ( 9 7 ) 0 0 4 5 0 - 5
Transcript
Theory and Methodology
Estimating and bootstrapping Malmquist indices 1
L�eopold Simar a,2, Paul W. Wilson b,*,3
a Institut de Statistique and CORE, Universit�e Catholique de Louvain, Voie du Roman Pays 20, Louvain-la-Neuve, Belgiumb Department of Economics, University of Texas at Austin, Austin, TX 78712-1173, USA
Received 1 August 1996; accepted 1 October 1997
Abstract
This paper develops a consistent bootstrap estimation procedure for obtaining con®dence intervals for Malmquist
indices of productivity and their decompositions. Although the exposition is in terms of input-oriented indices, the
techniques can be trivially extended to the output orientation. The bootstrap methodology is an extension of earlier
work described in Simar and Wilson (Simar, L., Wilson, P.W., 1998, Management Science). Some empirical examples
are also given, using data on Swedish pharmacies. Ó 1999 Published by Elsevier Science B.V. All rights reserved.
Keywords: DEA; Productivity; Resampling; Bootstrap; Malmquist indices
1. Introduction
FaÈre et al. (1992) merge ideas on measurementof e�ciency from Farrell (1957) and on measure-ment of productivity from Caves et al. (1982) todevelop a Malmquist index of productivitychange. Caves et al. de®ne their input-based
Malmquist productivity index as the ratio of twoinput distance functions, while assuming no tech-nical ine�ciency in the sense of Farrell. F�are et al.extend the Caves et al. approach by dropping theassumption of no technical ine�ciency and devel-oping a Malmquist index of productivity that canbe decomposed into indices describing changes intechnology and e�ciency. We extend the F�are etal. approach by giving a statistical interpretationto their Malmquist productivity index and itscomponents, and by presenting a bootstrap algo-rithm which may be used to estimate con®denceintervals for the indices. This work will allow re-searchers to speak in terms of whether changes inproductivity, e�ciency, or technology are signi®-cant in a statistical sense. In other words, ourmethods can be used to determine whetherindicated changes in productivity, e�ciency, or
European Journal of Operational Research 115 (1999) 459±471
3510; e-mail: [email protected] Pontus Roos graciously provided the data used in the
empirical examples. We alone, of course, are responsible for any
remaining errors or omissions.2 Research support from the contract ``Projet d'Actions de
Recherche Concert�ees'' (PARC No. 93/98±164) of the Belgian
Government is gratefully acknowledged.3 Research support from the Management Science Group,
US Department of Veterans A�airs, is gratefully acknowledged.
0377-2217/99/$ ± see front matter Ó 1999 Published by Elsevier Science B.V. All rights reserved.
PII: S 0 3 7 7 - 2 2 1 7 ( 9 7 ) 0 0 4 5 0 - 5
technology are real, or merely artifacts of the factthat we do not know the true production frontiersand must estimate them from ®nite samples.
The input-based Malmquist index of produc-tivity developed by F�are et al. measures produc-tivity change between times t1 and t2. This index, aswell as its component indices describing changes intechnology and e�ciency, consists of ratios of in-put distance functions (a more rigorous descrip-tion appears in Section 2). However, F�are et al. donot distinguish between the underlying true dis-tance functions and their estimates. For example,as a prelude to their Eq. (4) (page 88), they statethat ``the value of the distance function... is ob-tained as the solution to the linear programmingproblem...''. In fact, solving their linear program-ming problem yields an estimate of the distancefunction, not the distance function itself. F�are etal. are not alone in this regard; indeed, the litera-ture on nonparametric e�ciency measurement is®lled with such statements. Lovell (1993) andothers have labeled nonparametric, linear-pro-gramming based approaches to e�ciency mea-surement as deterministic, which seems to suggestthat these approaches have no statistical under-pinnings. Yet, if one views production data ashaving been generated from a distribution withbounded support over the true production set,then e�ciency, and changes in productivity, tech-nology, and e�ciency, are always measured rela-tive to estimates of underlying, true frontiers,conditional on observed data resulting from theunderlying (and unobserved) data-generatingprocess. Consequently, the estimates researchersare interested in involved uncertainty due to sam-pling variation. 4
Simar and Wilson (1998) develop a bootstrapprocedure which may be used to estimate con®-dence intervals for distance functions used tomeasure technical e�ciency, and demonstrate thatthe key to statistically consistent estimation ofthese con®dence intervals lies in the replication ofthe unobserved data-generating process. This pa-per extends those ideas to the case of Malmquistindices constructed from nonparametric distancefunction estimates using data from di�erent timeperiods.
In the next section, we de®ne the input-basedMalmquist productivity index and the distancefunctions from which it is constructed, and de-scribe how the productivity index can be decom-posed into indices of e�ciency change andtechnical shift. We also brie¯y discuss how thesemeasures can be estimated nonparametricallyusing linear programming techniques. While wefocus on input-based indices, one may triviallyextend our results to output-based measures bymerely modifying our notation. The bootstrapprocedure is presented in Section 3. In Section 4,we illustrate the bootstrap estimation using a panelof data on Swedish pharmacies previouslyexamined by F�are et al. Conclusions are given inSection 5.
2. Estimating Malmquist indices
To begin, consider ®rms which produce moutputs from n inputs. Let x 2 Rn
� and y 2 Rm�
denote input and output vectors, respectively. Theproduction possibilities set at time t is given by theclosed set
Pt � �x; y� j x can produce y at time tf g; �1�which may be described in terms of its sections
Xt�y� � x 2 Rn� j �x; y� 2 Pt
� ; �2�
i.e., its corresponding input requirement sets.Shephard (1970) discusses assumptions one mayreasonably make regarding Xt�x� (and hence Pt);typical assumptions, which we adopt, are (i) Xt�y�is convex for all y; t; (ii) all production requiresuse of some inputs, i.e., 0 62 Xt�y� if y P 0; y 6� 0;
4 The literature is also ®lled with references to the ``observed
best-practice frontier'', typically taken to be the boundary of
the convex, conical, or free-disposal hull of the observed data.
But, one may always ask what would happen if another
observation were obtained. Clearly, an additional observation
might push the observed best-practice frontier outward (but
never inward), although not beyond underlying boundary of
the true production set. Authors who use such terminology
have merely conditioned their analysis on an estimate of the
true production frontier. While possible, this is unnecessary.
460 L. Simar, P.W. Wilson / European Journal of Operational Research 115 (1999) 459±471
and (iii) both inputs and outputs are stronglydisposable, i.e.,
~x P x 2 Xt�y� ) ~x 2 Xt�y�and
~y P y) Xt�~y� � Xt�y�:Let subscript i; i � 1; . . . ;N ; denote a particu-
lar ®rm i; the N ®rms are each observed at (thesame) two points in time. The Shephard (1970)input distance function for ®rm i at time t1, rela-tive to the technology existing at time t2, is de®nedas
Dt1jt2i � sup h > 0 j xit1=h 2 Xt2�yit1�
� : �3�
The distance function Dt1jt2i gives a normalized
measure of distance from the ith ®rm's position inthe input/output space at time t1 to the boundaryof the production set at time t2 in the hyperplanewhere outputs remain constant. If t1 � t2, then wehave a measure of e�ciency relative to the con-temporaneous technology, and Dtjt
i P 1. If t1 6� t2,then Dt1jt2
i �<;�; >�1.FaÈre et al. (1992) write their Malmquist pro-
ductivity index as
Mi�t1; t2� � Dt2jt2i
Dt1jt1i
� Dt2jt1i
Dt2jt2i
� Dt1jt1i
Dt1jt2i
!�1=2�
; �4�
where t2 > t1. Values Mi�t1; t2� < 1 indicate im-provements in productivity between t1 and t2,while values Mi�t1; t2� > 1 indicate decreases inproductivity from time t1 to t2 (Mi�t1; t2� � 1would indicate no change in productivity). Theratio Dt2jt2
i =Dt1jt1i in Eq. (4) measures the change in
input technical e�ciency between periods t1 and t2,and de®nes an input-based index of e�ciencychange:
Ei�t1; t2� � Dt2jt2i
Dt1jt1i
: �5�
Values of Ei�t1; t2� less than (greater than) unityindicate improvements (decreases) in e�ciencybetween t1 and t2. Similarly, the remaining part ofthe right-hand side of Eq. (4) de®nes an input-based measure of technical change:
Ti�t1; t2� � Dt2jt1i
Dt2jt2i
� Dt1jt1i
Dt1jt2i
!�1=2�
: �6�
As with Mi�t1; t2� and Ei�t1; t2�, values of Ti�t1; t2�less than (greater than) unity indicate technicalprogress (regress) between times t1 and t2. 5
Unfortunately, the production set Pt is typi-cally unobserved; similarly, Xt�x� is also unob-served, as are the values of the distance functionswhich appear in the Malmquist index in Eq. (4)and its components in Eqs. (5) and (6). The indicesin Eqs. (4)±(6) represent true values which must beestimated. Substituting estimators for the corre-sponding true distance function values in Eqs. (4)±(6) yields estimatorsMi�t1; t2�, Ei�t1; t2�, andTi�t1; t2� of the productivity, e�ciency, and tech-nology change indices, respectively.
Estimation of the input distance functionscomprising Eqs. (4)±(6) requires estimation of Pt
and Xt�y�. Given a sample
S � �xit; yit� j i � 1; . . . ;N ; t � 1; 2f gof observations on N ®rms in 2 periods, there areseveral ways in which Pt may be estimated. Acommon approach is to estimate Pt by the conicalhull of the sample observations, which is tanta-mount to assuming constant returns to scale forthe production technology. The correspondingestimate of the input requirement set is
Xt�y� � x 2 Rn j y6Y tq; x P X tq; q 2 RN�
� ;
�7�where Y t � �y1t . . . yNt�, X t � �x1t . . . xNt�, with xit
and yit denoting �n� 1� and �m� 1� vectors ofobserved inputs and outputs, respectively, and q isan �N � 1� vector of intensity variables. This im-plies the distance function estimator
Dt1jt2i � sup k > 0 j xit1=k 2 Xt2�yit1�
n o; �8�
5 Note that output distance functions can be de®ned similar
to the input distance function in Eq. (3). Hence output-based
measures of productivity, e�ciency, and technical change can
be constructed by merely replacing the input distance functions
in Eqs. (4)±(6) with the corresponding output distance func-
tions.
L. Simar, P.W. Wilson / European Journal of Operational Research 115 (1999) 459±471 461
which may be computed by solving the linearprogram
Dt1jt2i
� �ÿ1
�min k j yit1 6Y t2 qi; kxit1 P X t2 qi; qi 2 RN
��
; �9�where t1�<;�; >�t2.
Alternatively, the production set Pt could beestimated by the convex hull of the sample obser-vations, which amounts to adding the constraint1*
qi � 1 in Eqs. (7)±(9), where 1*
is a �1� N� vectorof ones. This yields a distance function estimator~Dt1jt2
i which, when substituted into Eqs. (4)±(6),yields estimators ~Mi�t1; t2�, ~Ei�t1; t2�, and ~Ti�t1; t2�of the productivity, e�ciency, and technologychange indices, respectively. Other estimators ofPt are also possible; for instance, one could use thefree disposal hull of the sample observations sug-gested by Deprins et al. (1984). Use of the conicalhull as in FaÈre et al. (1992) implicitly assumes thatthe production technology exhibits constant re-turns to scale, allowing the Malmquist index inEq. (4) to be interpreted as an index of total factorproductivity. 6
3. Bootstrapping the Malmquist indices
The methodology for bootstrapping distancefunction estimators such as Eq. (9) presented inSimar and Wilson (1998) are easily adapted to thepresent case, except here the possible time-depen-
dence structure of the data must be taken intoaccount. As in our earlier work, we assume a data-generating process where ®rms randomly deviatefrom the underlying true frontier in a radial inputdirection. These random deviations from the con-temporaneous frontier at time t, measured by theShephard input distance function Dtjt
i in Eq. (3),are further assumed to result from ine�ciency. Ata given point in time, the marginal density of thesedeviations has support bounded on the left atunity.
Bootstrapping involves replicating this data-generating process, generating an appropriatelylarge number B of pseudosamples
S� � �x�it; y�it� j i � 1; . . . ;N ; t � 1; 2�
;
and applying the original estimators to thesepseudosamples. For each bootstrap replicationb � 1; . . . ;B, we use Eq. (9) to measure the dis-tance from each observation in the original sampleS to the frontiers estimated for either period fromthe pseudodata in S�. This is accomplished bysolving, for the case of the conical hull estimator,
Dt1jt2�i
� �ÿ1
�min k j yit1 6Y t2�qi; kxit1 P X t2�qi; qi 2 RN
��
;
�10�where Y t� � �y�1t . . . y�Nt� and X t� � �x�1t . . . x�Nt�. Fortwo time periods t1; t2, this yields bootstrap esti-mates fDt1jt1�
i �b�; Dt2jt2�i �b�; Dt1jt2�
i �b�; Dt2jt1�i �b�gB
b�1
for each ®rm i � 1; . . . ;N . 7 These estimates canthen be used to construct bootstrap estimatesM�
i �t1; t2��b�, E�i �t1; t2��b�, and T�i �t1; t2��b� (where
i � 1; . . . ;N and b � 1; . . . ;B) corresponding toEqs. (4)±(6), respectively, by replacing the truedistance function values in Eqs. (4)±(6) with theircorresponding bootstrap estimates.
Once these bootstrap values have been com-puted, we can correct for any ®nite-sample bias inthe original estimators of the distance functions or
6 Allowing for variable returns by using the convex hull or
free-disposal hull estimators does not guarantee solutions to
Eq. (8) for all observations when t1 6� t2, and in addition
prevents the Malmquist productivity index from being inter-
preted as an index of total factor productivity. However, if the
true technology contains regions with either increasing or
decreasing returns to scale, then the conical hull estimator will
fail to converge to Pt as N !1, and thus will give inconsistent
estimates. The convex and free-disposal hull estimators are
consistent regardless of whether returns to scale are constant or
variable in the sense that they converge to Pt as N !1, but
these estimators may converge more slowly than the conical
hull estimator if the true technology is one of constant returns
to scale everywhere. See Korostelev et al. (1995a) and
Korostelev et al. (1995b) for rates of convergence of the
various estimators of Pt.
7 Note that Dt1 jt1�i , Dt2 jt2�
i , and Dt2jt1�i can be computed by (i)
changing t2 to t1, (ii) changing t1 to t2, or (iii) reversing t1 and t2in Eq. (10), respectively. For purposes of the Malmquist index
and its components, we would typically order the time periods
so that t1 occurs before t2.
462 L. Simar, P.W. Wilson / European Journal of Operational Research 115 (1999) 459±471
the Malmquist indices, and con®dence intervals atthe desired level of signi®cance can be constructedusing the simple procedure outlined below. To il-lustrate, consider the set of bootstrap estimates forthe Malmquist index for ®rm i: fM�
i �t1; t2��b�gBb�1.
The e�ciency and technology change indices canbe analyzed similarly by merely changing M toeither E or T in the notation below. The bootstrapbias estimate for the original estimator Mi�t1; t2� is
dbiasB�Mi�t1; t2��
� Bÿ1XB
b�1
M�i �t1; t2��b� ÿMi�t1; t2�; �11�
which is the empirical bootstrap analog ofE�Mi�t1; t2�� ÿMi�t1; t2�. Therefore, a bias-cor-rected estimate of Mi�t1; t2� may be computed as
^Mi�t1; t2� � Mi�t1; t2� ÿdbiasB�Mi�t1; t2��
� 2Mi�t1; t2� ÿ Bÿ1XB
b�1
M�i �t1; t2��b�: �12�
Unfortunately, as Efron and Tibshirani (1993)note, the bias-corrected estimator in Eq. (12) mayhave higher mean-square error than the originalestimator. 8 The variance of the summation termon the right-hand side of the second line ofEq. (12) can be made arbitrarily small by in-
creasing B; yet, even if B!1, the bias-corrected
estimator^
Mi�t1; t2� will have variance equal to fourtimes that of the original estimator, Mi�t1; t2�, dueto the ®rst term on the second line of Eq. (12). Themean-square errors of the original and bias-cor-rected estimators can be compared by usingthe sample variance s2
�i of the bootstrap valuesM�
i �t1; t2��b�� B
b�1to estimate the variance of
Mi�t1; t2�. Then the estimated mean-square error
of^
Mi�t1; t2� is 4s2�i, and the estimated mean-square
error of Mi�t1; t2� is �s2�i � �dbiasB�Mi�t1; t2���2�. A bit
of algebra reveals that the bias-corrected estimatorwill likely have higher mean-square error than theoriginal estimator unless s2
�i <13�dbiasB�Mi�t1; t2���2.
In the case of individual distance function esti-mates, the corresponding bias may be substantial;but since Malmquist indices are de®ned as ratiosof distance functions, the overall bias of thesestatistics may be somewhat less than for individualdistance function estimates since the terms in boththe numerator and the denominator are biased inthe same direction. In any case, whether bias cor-rections such as Eq. (12) should be employed isalways an empirical question to be answered bythe data in any application.
To estimate con®dence intervals for theMalmquist index, note that the idea behind thebootstrap is to approximate the unknown distribu-tion of Mi�t1; t2� ÿMi�t1; t2�
ÿ �by the distribution
of M�i �t1; t2� ÿ Mi�t1; t2�
ÿ �conditioned on the origi-
nal data S. As noted above, our bootstrap proce-dure yields bootstrap values M�
i �t1; t2��b�� B
b�1,
which with the original estimate Mi�t1; t2�, can beused to obtain an empirical approximation to thesecond distribution.
If we knew the distribution of Mi�t1; t2�ÿÿ
Mi�t1; t2��, then it would be trivial to ®nd valuesaa, ba such that
Prob ÿ ba6 Mi�t1; t2� ÿMi�t1; t2�6 ÿ aa
� �� 1ÿ a �13�
for some small value of a, say 0.10 or 0.05. Sincewe do not know this distribution, we can use thebootstrap values to ®nd values a�a, b�a such that thestatement
Prob ÿ b�a6 M�i �t1; t2� ÿ Mi�t1; t2�6 ÿ a�ajS
� �� 1ÿ a �14�
is true with high probability. Mechanically, thisinvolves sorting the values
M�i �t1; t2��b� ÿ Mi�t1; t2�
� �; b � 1; . . . ;B
by algebraic value, deleting ��a=2� � 100�-percentof the elements at either end of this sorted array,and then setting ÿb�a and ÿa�a equal to the end-points of the resulting (sorted) array, with a�a6 b�a.When we say that Eq. (14) is ``true with highprobability'', we mean that this can be made so bymaking the number of bootstrap replications, B,
8 Note that we refer to^
Mi�t1; t2� as a bias-corrected, rather
than an unbiased, estimator, since Eq. (12) involves only a ®rst-
order correction of the bias in Mi�t1; t2�.
L. Simar, P.W. Wilson / European Journal of Operational Research 115 (1999) 459±471 463
large enough; as B!1, the probability thatEq. (14) is true approaches 1. 9
by substituting a�a and b�a for aa and ba in Eq. (13)and noting the conditioning in Eq. (15). Rear-ranging the terms in parentheses in Eq. (16) yieldsan estimated �1ÿ a�-percent con®dence interval
Mi�t1; t2� � a�a 6 Mi�t1; t2� 6 Mi�t1; t2� � b�a;
�17�and we say that the estimated Malmquist index issigni®cantly di�erent from unity (which wouldindicate no productivity change) if the interval inEq. (17) does not include unity. 10
The key to obtaining consistent bootstrap es-timates of the con®dence intervals lies in consis-tent replication of the data-generating process. Asdiscussed in Simar and Wilson (1998), resamplingfrom the empirical distribution of the data (i.e.,drawing with replacement from the set of originalobservations on inputs and outputs, or equiva-lently, from the set of original distance function
estimates) to construct the pseudosamples S� willlead to inconsistent bootstrap estimation of thecon®dence intervals. This results from the empir-ical distribution placing a positive probabilitymass at the boundary of the estimated productionset; this mass does not disappear as N !1, andso the empirical distribution provides an incon-sistent estimate of the underlying distribution ofine�ciencies measured by the input distancefunction.
Using a smooth bootstrap procedure as in Si-mar and Wilson (1998) overcomes this problemand yields consistent estimates. When bootstrap-ping distance function estimates from a singlecross-section of data, this may be accomplished byusing a univariate kernel estimator of the densityof the original distance function estimates, andthen drawing from this estimated density to con-struct the pseudosamples S� as in and Simar andWilson (1998). In the present case, however, wehave panel data, with the possibility of temporalcorrelation. For example, an ine�cient ®rm inperiod one may be more likely to be ine�cient inperiod two than a ®rm that is relatively more ef-®cient in period one. To preserve any temporalcorrelation present in the data, we use kernelmethods to estimate the joint density of�Dt1;t1
i ; Dt2;t2i �
� N
i�1.
The bivariate kernel density estimator with bi-variate kernel function K��� and bandwidth h isgiven by
f �z� � Nÿ1hÿ2XN
i�1
K�zÿ Z i�
h
� �; �18�
where z has dimension �1� 2� and Z i �� Dt1jt1
i Dt2jt2i � is the ith row of the �N � 2� matrix
containing the original data. 11 Note, however,that both Dt1;t1
i and Dt2;t2i are bounded from below
by unity. The density estimated from Eq. (18) canbe shown to be inconsistent and asymptotically
9 Of course, one could also choose ÿb�a and ÿa�a as the
100� s1th and 100� s2th percentiles of the sorted bootstrap
dard normal distribution function, and w�a� is the 100� ath
percentile of the standard normal distribution so that
U�w�a�� � a. This procedure involves a median-bias correction,
and was used in Simar and Wilson (1998), and is discussed in
Efron and Tibshirani (1993).10 In Simar and Wilson (1998), we constructed con®dence
intervals for distance functions using an estimate of bias
analogous to Eq. (11). The approach in this paper avoids
introducing the extra noise contained in this estimate. The bias
inherent in the distance function estimates is implicitly ac-
counted for here since we use the bootstrap values to construct
an empirical distribution of di�erences as in Eq. (15).
11 One might prefer to use di�erent bandwidths in each
direction; however, this is not necessary if the kernel function is
scaled by an estimate of the covariance matrix of the data as
discussed below.
464 L. Simar, P.W. Wilson / European Journal of Operational Research 115 (1999) 459±471
biased when the support of f is bounded, as in ourcase. 12
To overcome this problem, we adapt the uni-variate re¯ection method described by Silverman(1986) to our bivariate case. For the case of uni-variate data zif gN
i�1 bounded from below at unity,the re¯ection method involves using the univariatekernel density estimator to estimate the density ofthe original observations and their re¯ections
z0i� N
i�1about unity, where z0i � 2ÿ zi, 8 i �
1; . . . ;N . Truncating the resulting density estimate(based on 2N unbounded observations) on the leftat unity yields the desired density estimate for theunivariate case (see Simar and Wilson (1998) foran illustration). In the bivariate case, we proceedsimilarly, except that there are now two bound-aries in R2.
First, we form �N � 1� vectors
A � Dt1;t11 . . . Dt1;t1
N
h i0; �19�
and
B � Dt2;t21 . . . Dt2;t2
N
h i0: �20�
The values in A and B are bounded from below atunity. To re¯ect the distance function values aboutthe boundaries in two-dimensional space, we formthe �4N � 2� matrix represented in partitionedform by
D �
A B
2ÿ A B
2ÿ A 2ÿ B
A 2ÿ B
2666437775: �21�
The matrix D contains 4N pairs of values corre-sponding to the two time periods. 13 The temporalcorrelation of the original data A B� � is mea-sured by the estimated covariance matrix
R � r21 r12
r12 r22
24 35 �22�
of the columns of A B� �. 14 Note that R is alsonecessarily the estimated covariance matrix of there¯ected data 2ÿ A 2ÿ B� �. Moreover,
RR �r2
1 ÿr12
ÿr12 r22
24 35 �23�
must be the corresponding estimate of the covari-ance matrix of 2ÿ A B� � and A 2ÿ B� �.
Let Dj� denote the jth row of D. Then
g�z� � 1
4Nh2
X4N
j�1
Kjzÿ Dj�
h
� ��24�
is a kernel estimator of the density of the 4N re-¯ected data points represented by the rows of D,where z � z1 z2� �, and Kj��� is the bivariate nor-mal density function with shape R for j �1; . . . ;N ; 2N � 1; . . . ; 3N or shape RR for j � N �1; . . . ; 2N ; 3N � 1; . . . ; 4N . Then a consistent esti-mate of the density of the original data A B� �with bounded support is given by
12 The kernel function K��� in Eq. (18) must integrate to one.
At a point z, Eq. (18) estimates the underlying density as the
mean of N bivariate functions K��� centered on the original data
in Z; consequently, if the data are unbounded, then f �z�necessarily integrates to unity. However, when the data are
bounded as in our case, for data near the boundary, the
corresponding bivariate kernel functions extend over the
boundaries in one or both dimensions. Consequently, when
the density estimate obtained from Eq. (18) is integrated over
its bounded support, it will integrate to less than unity; the
estimated density will be too small near the boundaries.
Moreover, this problem remains as N !1.
13 To visualize the re¯ection, let ai; bi represent the ithelements of A and B, respectively �i � 1; . . . ;N�. Then for N®rms at times t1 and t2, we have N points �ai; bi� lying northeast
of the point �1; 1� in two-dimensional Euclidean space. These
points may be re¯ected by taking the ``mirror images'' about
vertical and horizontal lines passing through �1; 1�. This leaves
a boundary along the horizontal line passing through �1; 1� to
the left of this point, and along the vertical line passing through
�1; 1� below this point. Hence, an additional re¯ection is
required, obtained by taking the mirror image of the points to
the southeast of �1; 1� about the vertical line passing through
�1; 1�, or equivalently, by taking the mirror image of the points
to the northwest of �1; 1� about the horizontal line passing
through �1; 1�. Hence D in Eq. (21) has 4N rows.14 We use the sample covariance matrix to estimate the
covariance of A and B. Alternatively, one might use a robust
estimator of the covariance matrix, such as an M-estimator or
the minimum volume ellipsoid estimator proposed by Rous-
seeuw (1985).
L. Simar, P.W. Wilson / European Journal of Operational Research 115 (1999) 459±471 465
g��z� � 4g�z� for z1 P 1; z2 P 1;
0 otherwise:
(�25�
To generate the random deviates needed for thebootstrap, we do not have to actually estimate thedensity of the observations in D; rather, we use themethod suggested by Silverman (1986) and anal-ogous to that used for the univariate case in Simarand Wilson (1998). First, we randomly draw withreplacement N rows from D to form the �N � 2�matrix D� � �dij�, i � 1; . . . ;N , j � 1; 2 such thateach row of D has equal probability of selection.Let �d�j � Nÿ1
PNi�1 dij for j � 1; 2. Then compute
the �N � 2� matrix
C � �1� h2�ÿ1=2D� � h�� ÿ C
�d�1 0
0 �d�2
" # !
� C�d�1 0
0 �d�2
" #; �26�
where C is an �N � 2� matrix of ones, which givesan �N � 2� matrix of bivariate deviates from theestimated density of D, scaled to have the ®rst andsecond moment properties observed in the originalsample represented by A B� �. In addition, �� isan �N � 2� matrix containing N independentdraws from the kernel functions Kj��� in Eq. (24),with the ith row of �� representing (i) a draw froma normal density with shape R if D�i� was drawnfrom A B� � or 2ÿ A 2ÿ B� �; or (ii) a drawfrom a normal density with shape RR if D�i� wasdrawn from 2ÿ A B� � or A 2ÿ B� �.
Draws from a bivariate N�0; R� density can besimulated by generating independent, identicallydistributed pseudorandom N�0; 1� deviates �z1; z2�using the Box±Muller method (e.g., Press et al.,1986). The Cholesky decomposition of the �2� 2�matrix R yields the lower triangular matrix
L � `1 0
`2 `3
� �; �27�
where
LL0 � R; `1 � r1;
`2 � r12=r1; and `3 � r22 ÿ r2
12=r21
� �1=2
:
Then �`1z1; `2z1 � `3z2� � N�0; R�. Draws from aN�0; RR� density can be simulated similarly bycomputing �`1z1;ÿ`2z1 � `3z2� � N�0; RR�.
Finally, for each element cij of C, set
c�ij �cij; if cij P 1;
2ÿ cij; otherwise:
(�28�
The resulting �N � 2� matrix C� � �c�ij� consists oftwo column-vectors of simulated distance functionvalues. Pseudosamples S� are then constructed bysetting x�itj � c�ijxitj=Dtjjtj
i and y�itj � yitj for i �1; . . . ;N , j � 1; 2. 15
The only remaining issue is the choice of thebandwidth, h. Tapia and Thompson (1978), Sil-verman (1978, 1986), and HaÈrdle (1990) discussconsiderations relevant to the choice of h; in gen-eral, for a given sample size, larger values of hproduce more di�use (i.e., less e�cient) estimatesof the density, while very small values produceestimated densities with multiple modes. In theempirical examples which follow, we use Silver-man's (1986) suggestion for bivariate data by set-ting h � �4=5N�1=6
since we are using a bivariatenormal kernel scaled to have the same shape as thedata. 16
Finally, we note that in many applications, onemay wish to consider the evolution of e�ciency,technology, and productivity over more than twoperiods. Such was the case considered by FaÈre et al.(1992) and in our examples in the next section. The
15 Computing xitj=Dtj jtji scales the input vector back to the
ostensibly e�cient level indicated by the estimated frontier;
multiplying by c�ij simulates a random deviation away from this
frontier. If we were using output distance functions and the
output-based Malmquist index, we would retain the original
input vector in the pseudosample and generate a new output
vector.16 This bandwidth minimizes the approximate mean inte-
grated square error of the density estimate when the data are
bivariate normally distributed. Alternatively, one might use the
least-squares cross-validation procedure discussed by Silverman
(1986) to choose h. Since many of the elements in D will
typically equal one, the least cross-validation function may
su�er from degenerate behavior related to the discretization
problem described by Silverman. The results in Simar and
Wilson (1998) suggest, however, that the estimated con®dence
intervals are not very sensitive with respect to the choice of
bandwidth.
466 L. Simar, P.W. Wilson / European Journal of Operational Research 115 (1999) 459±471
Malmquist index (4) and its components inEqs. (5) and (6) are de®ned over two periods; inour bootstrap procedure, we account for the co-variance between successive periods in Eq. (22),but we ignore the covariance between periods thatare not adjacent to each other in time. This sug-gests that one might increase statistical e�ciencyby accounting for this information in our estima-tion procedure. Unfortunately, however, doing sowould require nonparametric estimation of a highdimensional density, rather than the bivariate es-timation we employ in Eqs. (23) and (24). It is wellknown that kernel methods, as with most non-parametric estimation methods, su�er from thecurse of dimensionality (see Silverman, 1986, for adiscussion of the magnitude of this problem for thecase of kernel density estimation). Therefore, anypossible gains in statistical e�ciency from consid-ering the covariance between distance functionestimates in nonsuccessive periods would almostcertainly be overwhelmed by the increase in mean-square error resulting from estimating a high-di-mensional density.
4. Empirical examples
FaÈre et al. (1992) describe annual data on 42Swedish pharmacies from 1980 to 1989, whichproduce four outputs from four inputs. F�are et al.assume constant returns to scale and estimatedistance functions using Eq. (10) to construct es-timates Ei�t1; t2�, Ti�t1; t2�, and Mi�t1; t2�, whichthey report in Tables 1±3 of their article. Main-taining the assumption of constant returns toscale, we applied the bootstrap methods outlinedin section three to obtain estimates of bias andvariance, and to test for signi®cant di�erencesfrom unity, while setting B � 2000. In comparingthe estimated biases and variances, we found thatthe bias correction in Eq. (12) for the e�ciency,technology, and productivity change indices would
increase mean-square error. Consequently, we do
not report the bias corrected estimates^Ei�t1; t2�,
^Ti�t1; t2�, and
^Mi�t1; t2�. Consistent with F�are
et al., we report the reciprocals of our originalestimates (which are identical to the values re-
ported by F�are et al.) in Tables 1±3 (respectively),so that numbers greater than unity denote progresswhile numbers less than unity denote regress. Inaddition, we use single asterisks (�) to indicatecases where the indices are signi®cantly di�erentfrom unity at the 0.10 level, and double asterisks���� to indicate cases where the indices are signi®-cantly di�erent from unity at the 0.05 level.
While examining changes in e�ciency, F�are etal. ®nd (page 96) ``®ve pharmacies (nos. 15, 32, 33,35, and 39) to be e�cient in all time periods'' asindicated by values of unity for e�ciency changebetween all successive pairs of years in Table 1 forthese ®ve pharmacies. Our bootstrap reveals anadditional six pharmacies (nos. 2, 14, 17, 20, 26,and 36) with statistically insigni®cant changes ine�ciency (at either the 0.10 or 0.05 levels), for allpairs of years. Three more pharmacies had sta-tistically insigni®cant changes in all but one pair ofyears (nos. 19, 31, and 38). In all, of the 237 es-timates of e�ciency change reported in Table 1that are not equal to unity, only 113 (47.7%) aresigni®cantly di�erent from unity at either 0.1 or0.05.
Turning to our results for the technical changeindex in Table 2, our bootstrap results generallysupport the statements made by F�are et al. withrespect to technical change. F�are et al. state that``between 1981 and 1982 almost all pharmaciesshowed technical progress'', and our results showthat our estimates of technical change are statis-tically signi®cant at the 0.05 level in all but fourinstances; of these four pharmacies, three showsigni®cant technical change at the 0.01 level. F�areet al. ®nd only one pharmacy (no. 33) whose ra-dially e�cient frontier point shows technical pro-gress in all periods, but our results indicate that thechanges for this pharmacy are not signi®cant infour periods.
Similarly, our results for the index of produc-tivity change in Table 3 generally support state-ments made by F�are et al. Where F�are et al. ®ndproductivity gains in 259 cases and productivitylosses in 119 cases, we ®nd signi®cant (at 0.10)gains in 217 cases, and signi®cant (at 0.10) losses in91 cases; 81.5 percent of the estimates shown inTable 3 are signi®cantly di�erent from unity at the0.10 level.
L. Simar, P.W. Wilson / European Journal of Operational Research 115 (1999) 459±471 467
As noted in Section 3, assuming constant re-turns to scale and using the conical hull of theobserved data to estimate the production set willyield statistically inconsistent distance functionestimates when the true technology has noncon-
stant returns to scale. The convex hull estimator ofthe production set, however, converges to the trueproduction set regardless of whether returns toscale are constant or otherwise. Therefore, lackinga formal test of returns to scale, using distance
Table 1
Changes in e�ciency (i.e., reciprocal of Ei�t1; t2�), 42 Swedish pharmacies. Numbers greater than one indicate improvements (constant
Note: Single asterisks (�) denote signi®cant di�erences from unity at 0.10; double asterisks (��) denote signi®cant di�erences from unity
at 0.05.
468 L. Simar, P.W. Wilson / European Journal of Operational Research 115 (1999) 459±471
function estimators based on the convex hull toconstruct the Malmquist index may be a saferapproach. As noted earlier, though, the Malmquistproductivity index loses its interpretation as anindex of total factor productivity when we do this.
5. Conclusions
Malmquist indices have been widely used inrecent years to examine changes in productivity,e�ciency, and technology not only within a variety
Table 2
Changes in technology (i.e., reciprocal of Ti�t1; t2�), 42 Swedish pharmacies. Numbers greater than one indicate improvements
Note: Single asterisks (�) denote signi®cant di�erences from unity at 0.10; double asterisks (��) denote signi®cant di�erences from unity
at 0.05.
L. Simar, P.W. Wilson / European Journal of Operational Research 115 (1999) 459±471 469
of industries, but across countries as well. In eachcase, researchers have provided point estimates,although clearly there must be uncertainty sur-rounding these estimates due to sampling varia-tion. Our methodology outlined in the preceding
sections provides a tractable approach for consis-tently estimating con®dence intervals. In addition,as illustrated in our empirical examples, ourbootstrap methodology provides a correction forthe inherent bias in nonparametric distance func-
Table 3
Changes in productivity (i.e., reciprocal of Mi�t1; t2�), 42 Swedish pharmacies. Numbers greater than one indicate improvements
Note: Single asterisks (�) denote signi®cant di�erences from unity at 0.10; double asterisks (��) denote signi®cant di�erences from unity
at 0.05.
470 L. Simar, P.W. Wilson / European Journal of Operational Research 115 (1999) 459±471
tion estimates (and hence in estimates of Ma-lmquist indices), as well as a method for checkingwhether the bias-correction will increase mean-square error.
Con®dence intervals such as those estimated inour empirical examples are essential in interpretingestimates of Malmquist indices. As with any esti-mator, it is not enough to know whether the Ma-lmquist index estimator indicates increases ordecreases in productivity, but whether the indi-cated changes are signi®cant in a statistical sense;i.e., whether the result indicates a real change inproductivity, or is an artifact of sampling noise.Our bootstrap procedure allows the researcher tomake these distinctions.
References
Caves, D.W., Christensen, L.R., Diewert, W.E., 1982. The
economic theory of index numbers and the measurement of
input, output and productivity. Econometrica 5, 1393±1414.
Deprins, D., Simar, L., Tulkens, H., 1984. Measuring labor
ine�ciency in post o�ces. In: Marchand, M., Pestieau, P.,
Tulkens, H. (Eds.), The Performance of Public Enterprises:
Concepts and Measurements. North-Holland, Amsterdam,
pp. 243±267.
Efron, B., Tibshirani, R.J., 1993. An Introduction to the
Bootstrap. Chapman and Hall, Inc., New York.
F�are, R., Grosskopf, S., Lindgren, B., Roos, P., 1992. Produc-
tivity changes in Swedish pharmacies 1980±1989: A non-
parametric approach. Journal of Productivity Analysis 3,
85±101.
Farrell, M.J., 1957. The measurement of productive e�ciency.
Journal of the Royal Statistical Society A 120, 253±281.
H�ardle, W., 1990. Applied Nonparametric Regression. Cam-
