*! version 1.0.0 03MAY2014
capture program drop dea_ci
program define dea_ci, rclass
version 11.0
// syntax checking and validation-----------------------------------------------
// -----------------------------------------------------------------------------
// returns 1 if the first nonblank character of local macro `0' is a comma,
// or if `0' is empty.
if replay() {
dis as err "ivars = ovars required."
exit 198
}
// quick check
#del ;
syntax anything(name=inoutvars equalok) [if] [in] [using/]
,
envars(varlist) // environment variable list
[
reps1(integer 100) // perform # bootstrap replications; default is 100
reps2(integer 2000) // perform # replications; default is 2000
rts(string) // ignore case sensitive,[{CRS|CCR}|{BCC|VRS}|DRS|IRS]
ort(string) // ignore case sensitive,[{IN|INPUT}|{OUT|OUTPUT}]
stage(integer 1) // dea stage 1 or 2
trace // Whether or not to do the log
* // other dea options(varname is options)
];
#del cr
// get invars and outvars
local varnames `inoutvars'
gettoken word inoutvars : inoutvars, parse("=")
while ~("`word'" == "=") {
local invars `invars' `word'
gettoken word inoutvars : inoutvars, parse("=")
}
local outvars `inoutvars'
unab invars : `invars'
unab outvars : `outvars'
// default model - CRS(Constant Return to Scale)
if ("`rts'" == "") {
local rts = "CRS"
}
// default orientation - Input Oriented
if ("`ort'" == "") local ort = "IN"
else {
local ort = upper("`ort'")
if ("`ort'" == "I" | "`ort'" == "IN" | "`ort'" == "INPUT") {
local ort = "IN"
}
else if ("`ort'" == "O" | "`ort'" == "OUT" | "`ort'" == "OUTPUT") {
local ort = "OUT"
}
else {
di as err "option ort allows for case-insensitive " _c
di as err "(i|in|input|o|out|output) or nothing."
exit 198
}
}
if ("`trace'" == "trace") {
di as txt "invars:[`invars']"
di as txt "outvars:[`outvars']"
di as txt "envars:[`envars']"
di as txt "reps1:[`reps1']"
di as txt "rts:[`rts']"
di as txt "ort:[`ort']"
di as txt "options:[`options']"
}
if ("`using'" != "") use "`using'", clear
if (~(`c(N)' > 0 & `c(k)' > 0)) {
dis as err "dataset required!"
exit 198
}
// end of syntax checking and validation ---------------------------------------
// REP1 : loop count of stage1
// REP2 : loop count of stage2
set more off
preserve
if ("`if'" != "" | "`in'" != "") {
qui keep `in' `if' // filtering : keep in range [if exp]
}
_stage1, invars(`invars') outvars(`outvars') envars(`envars') ///
reps1(`reps1') ///
rts(`rts') ort(`ort') stage(`stage') `trace' `options'
return add
restore, preserve
_stage2
restore, preserve
set more on
end
********************************************************************************
* STAGE1 -
********************************************************************************
program define _stage1, rclass
syntax, invars(varlist) outvars(varlist) envars(varlist) reps1(integer) ///
rts(string) ort(string) stage(integer) [trace *]
// declare the tempname
tempname dearslt
tempname TM EM T
capture drop theta_h1 theta_hs theta_hh
// step1: calculate the theta_h with dea
qui dea `invars' = `outvars', rts(`rts') ort(`ort') stage(`stage') `options'
matrix `dearslt' = r(dearslt)
matrix `TM' = `dearslt'[1...,"theta"]
mata: _roundmat("`TM'", 1e-14) // round off
matrix colnames `TM' = theta_h
svmat `TM', names(col)
if ("`ort'" == "IN") {
qui replace theta_h = 1 / theta_h
}
gen float theta_hs = 0
gen float theta_hh = 0
if ("`trace'" == "trace") {
di _n(2) as txt "[theta_h values as dea result.]:"
list
}
// step2: truncated regression theta_h on environment variables
qui truncreg theta_h `envars', ll(1) noconstant
scalar sigma_h = e(sigma)
matrix `EM' = e(b)
if ("`trace'" == "trace") {
di _n(2) as txt "[ereturn list as truncreg result.]:"
ereturn list sigma
di _n(2) as txt "[matrix list e(b) as truncreg result.]:"
matrix list e(b)
}
forvalues i = 1/`=`reps1'' {
capture drop _epsilon _theta_s _theta_hi
// step3: 1st bootstrap stage, L1 times : computing dea
// step3.1: calculate the epsilon1 using sigma_h in step 2.
gen float _epsilon = sigma_h * invnormal(normal(1/sigma_h) ///
+ (1 - normal(1/sigma_h)) * uniform())
// step3.2: calculate theta_s: _theta_s = ZiB + ei
gen float _theta_s = _epsilon
foreach envar of varlist `envars' {
matrix `T' = `EM'[1,"eq1:`envar'"]
qui replace _theta_s = _theta_s + `envar' * `T'[1,1] // beta_h
}
if ("`trace'" == "trace") {
di _n(2) as txt "[`i'][_theta_s = ZiB + ei]:"
list
}
// step3.3: reset the dataset using theta_h in step1:and theta_s in step3.2:
if ("`ort'" == "IN") {
foreach invar of varlist `invars' {
qui replace `invar' = `invar' * (theta_h / _theta_s)
}
}
else { // if ("`ort'" == "OUT") {
foreach outvar of varlist `outvars' {
qui replace `outvar' = `outvar' * (theta_h / _theta_s)
}
}
if ("`trace'" == "trace") {
di _n(2) as txt "[`i'][reset the dataset using theta_h]:"
list
}
// step3.4: compute dea with changed dataset in step 3.3:
qui dea `invars' = `outvars', ///
rts(`rts') ort(`ort') stage(`stage') `options'
matrix `dearslt' = r(dearslt)
matrix `TM' = `dearslt'[1...,"theta"]
mata: _roundmat("`TM'", 1e-14) // round off
matrix colnames `TM' = _theta_hi
svmat `TM', names(col)
// step4: Loop-i times
qui replace theta_hs = theta_hs + _theta_hi
if ("`trace'" == "trace") {
di _n(2) as txt "[`i'][result of calculate the theta_hs]:"
list
}
} // end of forvalues
capture drop _epsilon _theta_s _theta_hi
// step3.5: calculate mean of theta_hs
qui replace theta_hs = theta_hs / `reps1'
// step4: compute bias-corrected efficiency estimates
qui replace theta_hh = 2 * theta_h - theta_hs
if ("`trace'" != "trace") {
di _n(2) as txt "[result of stage1]:"
list
}
end
********************************************************************************
* STAGE2 -
********************************************************************************
program define _stage2, rclass
list
end
// Start of the MATA Definition Area -------------------------------------------
version 10
mata:
mata set matastrict on
// if function name start with '_', call at the stata.
/**
* used to look for special values:
* @param mat - matrix name
* @param cmp - compare type(ge, gt, eq, lt, le: default is eq)
* @param value - comparing value.
*/
function _all (string scalar mat, string scalar cmp, real scalar value)
{
real matrix M;
real scalar bool;
M = st_matrix(mat)
if (cmp == "ge") {
bool = all (M :>= value);
}
else if (cmp == "gt") {
bool = all (M :> value);
}
else if (cmp == "lt") {
bool = all (M :< value);
}
else if (cmp == "le") {
bool = all (M :<= value);
}
else { // if (cmp == "eq") {
bool = all (M :== value);
}
st_numscalar("r(bool)", bool)
}
/**
* used to look for special values:
* @param mat - matrix name
* @param cmp - compare type(ge, gt, eq, lt, le: default is eq)
* @param value - comparing value.
*/
function _any (string scalar mat, string scalar cmp, real scalar value)
{
real matrix M;
real scalar bool;
M = st_matrix(mat)
if (cmp == "ge") {
bool = any (M :>= value);
}
else if (cmp == "gt") {
bool = any (M :> value);
}
else if (cmp == "lt") {
bool = any (M :< value);
}
else if (cmp == "le") {
bool = any (M :<= value);
}
else { // if (cmp == "eq") {
bool = any (M :== value);
}
st_numscalar("r(bool)", bool)
}
end
// End of the MATA Definition Area ---------------------------------------------