master

分支 (1)

管理

管理

master

gtfpch
/
gtfpch.ado

*! version 4.21, 6 Jan 2025
*! version 4.2, 12 Dec 2021
*! version 4.1, 21 Oct 2021
*! version 4.022, 18 Sep 2021
*! version 2.02, 20 June 2021
*! version 2.01, 20 Apr 2021
* use frame to store results
*! version 1.2
* 19 Jul 2020
* added hybrid Oriented
*! version 1.1
* 30 Apr 2020
* add biennial option

*! version 1.0
* By Kerry Du,Daoping Wang, Ning Zhang, 3 Dec 2019
**
*
capture program drop gtfpch
program define gtfpch, rclass prop(xt)
    version 16

   * recompile the mlib for different Stata version
	  _get_version gtfpch
	  _compile_mata, package(gtfpch) version(`package_version') verbose
    qui pwf
    local currentframe = r(currentframe)

    qui mata mata mlib index
    _xt, trequired
    local id=r(ivar)
    local time=r(tvar)

   * extract input and desriable output variables
    gettoken word 0 : 0, parse(" =:,")
    while `"`word'"' != ":" & `"`word'"' != "=" {
        if `"`word'"' == "," | `"`word'"'=="" {
                error 198
        }
		local invars `invars' `word'
		gettoken word 0 : 0, parse("=:,")
    }
    unab invars : `invars'

	gettoken word 0 : 0, parse(" =:,")
    while `"`word'"' != ":" & `"`word'"' != "=" {
        if `"`word'"' == "," | `"`word'"'=="" {
                error 198
        }
		local gopvars `gopvars' `word'
		gettoken word 0 : 0, parse(" =:,")
    }
    unab gopvars : `gopvars'

    syntax varlist [if] [in],    [dmu(varname) gx(varlist) gy(varlist) gb(varlist)  ///
	                              BIennial  SEQuential GLOBAL FGNZ RD  LUENberger ort(string) ///
								  WINdow(numlist intege max=1 >=1) SAVing(string)   Wmat(string) NONRadial ///
								  maxiter(numlist integer >0 max=1) tol(numlist max=1) frame(name) noPRINT NOCHeck]


   * check whether the new version is available
   if "`nocheck'"==""{ //default
      if "$gtfpchcheck"==""   gtfpchcheckupdate gtfpch
   }

    preserve

    sort `id' `time'
    if `"`frame'"'!=""{
       confirm new frame `frame'
    }

    marksample touse

    local bopvars `varlist'
    local ninp: word count `invars'
    local ngo:  word count `gopvars'
    local nbo:  word count `bopvars'
    local nvar: word count `invars' `gopvars' `bopvars'
    qui keep `invars' `gopvars' `bopvars' `touse' `id' `time' `dmu'  `gx' `gy' `gb'

    qui gen Row=_n
    label var Row "Row # in the original dataset"
    if "`nonradial'"!=""{
     local luenberger luenberger
	 *if "`ort'" =="" local ort h

     }

    if ("`ort'" == "") local ort = "out"
    else if ("`ort'" == "hybrid" | "`ort'" == "HYBRID" | "`ort'" == "h" | "`ort'" == "H") {
         local ort="h"
    }
    else {
        local ort = upper("`ort'")
        if ("`ort'" == "I" | "`ort'" == "IN" | "`ort'" == "INPUT") {
            local ort = "i"
        }
        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|h|hybrid) or nothing."
            exit 198
        }
    }

  if  "`luenberger'"=="" & "`ort'"=="h"{
      di as err "ort(hybrid) can not be used for Malmquist-Luenberger productivity index."
      exit 198
  }


  if "`gx'"!=""{
    local ngx: word count `gx'
    if `ngx'!=`ninp'{
        disp as error "# of input variables != # of variables specified in gx()."
        error 498
    }
    local gmat `gmat' `gx'
    local gmatname `gmatname' `gx'

  }
  else{

    if  `"`ort'"'!="out"{ // for ort(i) and ort(h)
      local invarscopy `invars'
      forv k=1/`ninp'{
        gettoken word invarscopy:invarscopy
        tempvar gx_`k'
        qui gen `gx_`k''=-`word'  // gx=-X for input orientied direction
        local gmat `gmat' `gx_`k''
        local gx `gx' `gx_`k''
        local gmatname `gmatname' -`word'
      }

    }
    else{
      forv k=1/`ninp'{
        tempvar gx_`k'
        qui gen `gx_`k''=0  // gx=0 for output orientied direction
        local gmat `gmat' `gx_`k''
        local gx `gx' `gx_`k''
        local gmatname `gmatname' 0
      }

    }


  }

  if "`gy'"!=""{
    local ngy: word count `gy'
    if `ngy'!=`ngo'{
        disp as error "# of desriable output variables != # of variables specified in gy()."
        error 498
    }
    local gmat `gmat' `gy'
    local gmatname `gmatname' `gy'

  }
  else{

    if  `"`ort'"'!="i"{ //for ort(out) and ort(h)
      local gopvarscopy `gopvars'
      forv k=1/`ngo'{
        gettoken word gopvarscopy:gopvarscopy
        tempvar gy_`k'
        qui gen `gy_`k''=`word'
        local gmat `gmat' `gy_`k''
        local gy `gy' `gy_`k''
        local gmatname `gmatname' `word'
      }

    }
    else{
      forv k=1/`ngo'{
        tempvar gy_`k'
        qui gen `gy_`k''=0
        local gmat `gmat' `gy_`k''
        local gy `gy' `gy_`k''
        local gmatname `gmatname' 0
      }
    }


  }

  if "`gb'"!=""{
    local ngb: word count `gb'
    if `ngb'!=`nbo'{
        disp as error "# of undesriable output variables != # of variables specified in gb()."
        error 498
    }
    local gmat `gmat' `gb'
    local gmatname `gmatname' `gb'

  }
  else{

    if `"`ort'"'!="i"{

    local bopvarscopy `bopvars'
    forv k=1/`nbo'{
      gettoken word bopvarscopy:bopvarscopy
      tempvar gb_`k'
      qui gen `gb_`k''=-`word'
      local gmat `gmat' `gb_`k''
      local gb `gb' `gb_`k''
      local gmatname `gmatname' -`word'
    }

    }
    else{
      forv k=1/`nbo'{
        tempvar gb_`k'
        qui gen `gb_`k''=0
        local gmat `gmat' `gb_`k''
        local gb `gb' `gb_`k''
        local gmatname `gmatname'  0
      }


    }


  }


    tempname weightvec wcheck
    if "`nonradial'"==""{
        if `"`wmat'"'!=""{
            disp as red "Warning: wmat() is only used when nonradial is specified."
        }

    }
    else{
        if `"`wmat'"'!=""{
            confirm matrix `wmat'
            local ncol=colsof(`wmat')
            if `ncol'!=`nvar'{
                dis as error `"# of column of `wmat' != # of input-output variables"'
                exit 498
            }

          forv i=1/`ncol'{
            local wival=`wmat'[1,`i']
            if `wival'<0{
              di as error `"The element of matrix `wmat' should not be negative."'
              exit 498
            }
          }

        //mat `weightvec'=`wmat'
        if "`ort'"=="out"{
            mat `weightvec'=`wmat'[1,1..`ninp']
            mata: st_numscalar("`wcheck'",sum(st_matrix("`weightvec'")))
            if `wcheck'!=0{
                disp as error "For output orientied Luenberger productivity indicator, the weight components for input reduction should be set to 0."
                exit 498
            }
        }
        else if "`ort'"=="i"{
            mat `weightvec'=`wmat'[1,(`ninp'+1)..(`ninp'+`nbo'+`ngo')]
            mata: st_numscalar("`wcheck'",sum(st_matrix("`weightvec'")))
            if `wcheck'!=0{
                disp as error "For input orientied Luenberger productivity indicator, the weight components for output reduction should be set to 0."
                exit 498
            }

        }
         mat `weightvec'=`wmat'
        }
        else{

           if "`ort'"=="out"{
              mat  `weightvec'=(J(1,`ninp',0),J(1,`ngo',1),J(1,`nbo',1))
           }
           else if "`ort'"=="i"{
              mat  `weightvec'=(J(1,`ninp',1),J(1,`ngo'+`nbo',0))
           }
           else{
              mat  `weightvec'=(J(1,`ninp',1),J(1,`ngo'+`nbo',1))
           }


        }

        forval i = 1/`=colsof(`weightvec')' {
             local wvalues `wvalues'  `=`weightvec'[1,`i']'
           }
        di
        disp " The weight vector is (`wvalues')"
        local rweightname (`wvalues')

    }


************************************** CRS case*************************************

    if "`nonradial'"==""{

        if "`luenberger'"==""{
            _malmqluen `invars'=`gopvars':`bopvars' if `touse', id(`id') time(`time') gx(`gx') ///
                               gy(`gy') gb(`gb')  ort(`ort')   ///
                              `biennial' `global' `sequential' window(`window') ///
                               maxiter(`maxiter') tol(`tol')
            local indexname "Malmquist-Luenberger Productivity Index"

        }
        else{

            _luen_ddf `invars'=`gopvars':`bopvars' if `touse', id(`id') time(`time') gx(`gx') ///
                               gy(`gy') gb(`gb')  ort(`ort')   ///
                              `biennial' `global' `sequential' window(`window') ///
                               maxiter(`maxiter') tol(`tol')
            local indexname "Luenberger Productivity Index (based on DDF)"
        }
    }
    else{

             _luen_nddf `invars'=`gopvars':`bopvars' if `touse', id(`id') time(`time') gx(`gx') ///
                               gy(`gy') gb(`gb')  ort(`ort')  wmat(`weightvec') ///
                              `biennial' `global' `sequential' window(`window') ///
                               maxiter(`maxiter') tol(`tol')
            local indexname "Luenberger Productivity Index (based on nonradial DDF)"

    }
            local resvars `r(rvars)'
*****************************************************************************************
  if "`rd'"==""&"`fgnz'"==""{

        di
        di " The directional vector is (`gmatname')"
        format `resvars' %9.4f
        qui keep if  `touse'
        qui cap bys `id' (`time'): gen Pdwise=`time'[_n-1]+"~"+`time' if _n>1
        qui cap bys `id' (`time'): gen Pdwise=string(`time'[_n-1])+"~"+string(`time') if _n>1
        label var Pdwise "Period wise"
        order Row `dmu' `id' Pdwise  `resvars'
        qui keep  Row `dmu' `id' Pdwise  `resvars'
        qui keep if !missing(Pdwise)

        if "`noprint'" == ""{
          disp _n(2) " Total Factor Productivity Change:`indexname'"
          disp "    (Row: Row # in the original data; Pdwise: periodwise)"
          list Row `dmu' `id'  Pdwise  `resvars', sep(0)
          di "Note: missing value indicates infeasible problem."
        }


        if `"`saving'"'!=""{
          save `saving'
          gettoken filenames saving:saving, parse(",")
          local filenames `filenames'.dta
          disp _n `"Estimated Results are saved in `filenames'."'

        }

        if `"`frame'"'!=""{
          disp _n `"Estimated Results are saved in `frame' frame."'
          frame copy `currentframe'  `frame'
        }

        return local file `filenames'
        return local frame `frame'
        *return local techtype `techtype'
        return local weight  `rweightname'
        return local gvec  `gmatname'
        restore

  }


  **************************** estiamte VRS and decompose scale effciency ************************
  else{


   foreach v of local resvars{

      rename `v' `v'_crs

    }

   if "`nonradial'"==""{

        if "`luenberger'"==""{
            _malmqluen `invars'=`gopvars':`bopvars' if `touse', id(`id') time(`time') gx(`gx') ///
                               gy(`gy') gb(`gb')  ort(`ort')  vrs ///
                              `biennial' `global' `sequential' window(`window') ///
                               maxiter(`maxiter') tol(`tol')
            //local indexname "Malmquist-Luenberger Productivity Index"

        }
        else{

            _luen_ddf `invars'=`gopvars':`bopvars' if `touse', id(`id') time(`time') gx(`gx') ///
                               gy(`gy') gb(`gb')  ort(`ort')  vrs  ///
                             `biennial' `global' `sequential' window(`window') ///
                               maxiter(`maxiter') tol(`tol')
            //local indexname "Luenberger Productivity Index (base on DDF)"
        }
    }
    else{

             _luen_nddf `invars'=`gopvars':`bopvars' if `touse', id(`id') time(`time') gx(`gx') ///
                               gy(`gy') gb(`gb')  ort(`ort')  wmat(`weightvec') ///
                              `biennial' `global' `sequential' window(`window') vrs ///
                               maxiter(`maxiter') tol(`tol')
            //local indexname "Luenberger Productivity Index (base on nonradial DDF)"

    }
    local resvars `r(rvars)'

    if "`luenberger'"==""{
      if "`rd'"!=""{

        qui gen SECH=TFPCH_crs/TFPCH
        label var SECH "Scale efficiecny change"
        qui replace TFPCH=TFPCH_crs
        local resvars `resvars' SECH

      }
      else{
        qui gen SECH=TECH_crs/TECH
        label var SECH "Scale efficiecny change"
        qui replace TFPCH=TFPCH_crs
		qui replace TECCH=TECCH_crs
        local resvars `resvars' SECH
      }

    }
    else{
      if "`rd'"!=""{

        //su TFPCH_crs TFPCH
        qui gen SECH=TFPCH_crs-TFPCH
        label var SECH "Scale efficiecny change"
        qui replace TFPCH=TFPCH_crs
        local resvars `resvars' SECH

      }
      else{
        //su TECH_crs TECH
        qui gen SECH=TECH_crs-TECH
        label var SECH "Scale efficiecny change"
        qui replace TFPCH=TFPCH_crs

		qui replace TECCH=TECCH_crs

        local resvars `resvars' SECH
      }

    }


    di
    di " The directional vector is (`gmatname')"
    format `resvars' %9.4f

    qui keep if `touse'
    qui cap bys `id' (`time'): gen Pdwise=`time'[_n-1]+"~"+`time' if _n>1
    qui cap bys `id' (`time'): gen Pdwise=string(`time'[_n-1])+"~"+string(`time') if _n>1
    label var Pdwise "Period wise"

    order Row `dmu' `id' Pdwise  `resvars'
    qui keep if !missing(Pdwise) & `touse'
    qui keep  Row `dmu' `id' Pdwise  `resvars'

    if "`noprint'"==""{
      disp _n(2) " Total Factor Productivity Change:`indexname'"
      disp "    (Row: Row # in the original data; Pdwise: periodwise)"

      list Row `dmu' `id'  Pdwise  `resvars' , sep(0)
      di "Note: missing value indicates infeasible problem."
    }

    if `"`saving'"'!=""{
      save `saving'
      gettoken filenames saving:saving, parse(",")
      local filenames `filenames'.dta
      disp _n `"Estimated Results are saved in `filenames'."'

    }
        if `"`frame'"'!=""{
          disp _n `"Estimated Results are saved in `frame' frame."'
          frame copy `currentframe'  `frame'
        }

        return local file `filenames'
        return local frame `frame'
        *return local techtype `techtype'
        return local weight  `rweightname'
        return local gvec  `gmatname'

        restore


  }


end


***********************************************************
/////////////////////subprograms///////////////////////////
***********************************************************
cap program drop _malmqluen
program define _malmqluen,rclass
    version 16

    gettoken word 0 : 0, parse("=:,")
    while `"`word'"' != ":" & `"`word'"' != "=" {
        if `"`word'"' == "," | `"`word'"'=="" {
                error 198
        }
    local invars `invars' `word'
    gettoken word 0 : 0, parse("=:,")
    }
    unab invars : `invars'

  gettoken word 0 : 0, parse("=:,")
    while `"`word'"' != ":" & `"`word'"' != "=" {
        if `"`word'"' == "," | `"`word'"'=="" {
                error 198
        }
    local gopvars `gopvars' `word'
    gettoken word 0 : 0, parse(" =:,")
    }
    unab gopvars : `gopvars'

    syntax varlist [if] [in], id(varname) time(varname)  [gx(varlist)    ///
                       BIennial  gy(varlist) gb(varlist) VRS ort(string) ///
                       GLOBAL SEQuential WINdow(numlist intege max=1 >=1) ///
                  maxiter(numlist integer >0 max=1) tol(numlist max=1 >0)   ]

  marksample touse
  local bopvars `varlist'

  local ninp: word count `invars'
  local ngo:  word count `gopvars'
  local nbo:  word count `bopvars'
  local nvar: word count `invars' `gopvars' `bopvars'

  confirm numeric var `invars' `gopvars' `bopvars'

  local techtype "contemporaneous"


   if "`global'"!=""{
     if "`sequential'"!=""{

       disp as error "global and sequential cannot be specified together."
       error 498

     }

     if "`window'"!=""{

       disp as error "global and window() cannot be specified together."
       error 498

     }

     if "`biennial'"!=""{

       disp as error "global and biennial cannot be specified together."
       error 498

     }

     local techtype "global"

  }


   if "`sequential'"!=""{

     if "`window'"!=""{

       disp as error "sequential and window() cannot be specified together."
       error 498

     }

     if "`biennial'"!=""{

       disp as error "sequential and biennial cannot be specified together."
       error 498

     }

     local techtype "sequential"

  }


     if "`window'"!=""{

       if "`biennial'"!=""{

         disp as error "biennial and window() cannot be specified together."
         error 498

       }

         local techtype "window"

     }

       if "`biennial'"!=""{

        local techtype "biennial"

       }


  local comvars: list invars & gopvars
  if !(`"`comvars'"'==""){
    disp as error "`comvars' should not be specified as input and desriable output simultaneously."
    error 498
  }

  local comvars: list invars & bopvars
  if !(`"`comvars'"'==""){
    disp as error "`comvars' should not be specified as input and undesriable output simultaneously."
    error 498
  }

  local comvars: list gopvars & bopvars
  if !(`"`comvars'"'==""){
    disp as error "`comvars' should not be specified as desriable and undesriable outputs simultaneously."
    error 498
  }


  if "`maxiter'"==""{
    local maxiter=-1
  }
  if "`tol'"==""{
    local tol=-1
  }


  tempvar period dmu

  qui egen `period'=group(`time')
  qui egen `dmu'=group(`id')


    qui su  `period'
    local tmax=r(max)

    tempvar flag temp DD D21 D12

    qui gen `DD'=.
    qui gen `D21'=.
    qui gen `D12'=.
    qui gen `temp'=.
    qui gen `flag'=0

    local gmat `gx' `gy' `gb'
    sort `period' `dmu'

  if `"`techtype'"'=="contemporaneous"{

      qui{
	  	*D11
        forv t=1/`tmax'{
            qui replace `flag'= (`period'==`t')
            _ddf  if `period'==`t' & `touse', rflag(`flag') gen(`DD') gv(`gmat') `vrs' ort(`ort') in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')
        }
		*D21
        local tt=`tmax'-1
        forv t=1/`tt'{
            qui replace `flag'=(`period'==`t'+1)
            _ddf  if `period'==`t' & `touse', rflag(`flag') gen(`D21') gv(`gmat') `vrs' ort(`ort') in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')
        }
		*D12
        forv t=2/`tmax'{
            qui replace `flag'=(`period'==`t'-1)
            _ddf  if `period'==`t' & `touse', rflag(`flag') gen(`D12') gv(`gmat') `vrs' ort(`ort') in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')
        }

    }


  }


    if `"`techtype'"'=="sequential"{

		*D11
        forv t=1/`tmax'{
            qui replace `flag'=(`period'<=`t')
            _ddf if `period'==`t' & `touse', rflag(`flag') gen(`DD') gv(`gmat') `vrs' ort(`ort') in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        }
		*D21
        local tt=`tmax'-1
        forv t=1/`tt'{
            qui replace `flag'=(`period'<=`t'+1)
            _ddf  if `period'==`t' & `touse', rflag(`flag') gen(`D21') gv(`gmat') `vrs' ort(`ort') in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')
        }
		*D12
        forv t=2/`tmax'{
            qui replace `flag'= (`period'<=`t'-1)
            _ddf if `period'==`t' & `touse', rflag(`flag') gen(`D12') gv(`gmat') `vrs' ort(`ort') in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')
        }


  }


     if `"`techtype'"'=="window"{

       local band=`window'

        forv t=1/`tmax'{
            qui replace `flag'=(`period'<=`t'+`band' & `period'>=`t'-`band')
            _ddf  if `period'==`t' & `touse', rflag(`flag') gen(`DD') gv(`gmat') `vrs' ort(`ort') in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')
        }
        local tt=`tmax'-1
        forv t=1/`tt'{
            qui replace `flag'= (`period'<=`t'+1+`band' &  `period'>=`t'-`band'+1)
            _ddf  if `period'==`t' & `touse', rflag(`flag') gen(`D21') gv(`gmat') `vrs' ort(`ort') in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')
        }

        forv t=2/`tmax'{
            qui replace `flag'=(`period'<=`t'-1+`band' & `period'>=`t'-1-`band')
            _ddf  if `period'==`t' & `touse', rflag(`flag') gen(`D12') gv(`gmat') `vrs' ort(`ort') in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')
        }


  }


  if `"`techtype'"'=="biennial" {

      qui{
        forv t=1/`tmax'{
            qui replace `flag'= (`period'==`t' | `period'==`t'+1)
            _ddf  if `period'==`t' & `touse', rflag(`flag') gen(`DD') gv(`gmat') `vrs' ort(`ort') in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')
        }

        forv t=2/`tmax'{
            qui replace `flag'=(`period'==`t'-1 | `period'==`t')
            _ddf  if `period'==`t' & `touse', rflag(`flag') gen(`D12') gv(`gmat') `vrs' ort(`ort') in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        }

    }


  }


  if `"`techtype'"'=="global"{

      qui replace `flag'=1
      _ddf  if `touse', rflag(`flag') gen(`temp') gv(`gmat') `vrs' ort(`ort') in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')
      qui bys `dmu' (`period'): gen TFPCH=`temp'/`temp'[_n-1]
      label var TFPCH "Total factor productivity change"
      cap drop `temp'

    sort `period' `dmu'
    forv t=1/`tmax'{
      qui replace `flag'=(`period'==`t')
      _ddf  if `touse' & `period'==`t', rflag(`flag') gen(`DD') gv(`gmat') `vrs' ort(`ort') in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')
    }

    qui bys `dmu' (`period'): gen TECH=`DD'/`DD'[_n-1]

	qui bys `dmu' (`period'): gen TECCH=TFPCH/TECH

    label var TECH  "Technical efficiency change"
	label var TECCH "Technological change"

	local resvars TFPCH TECH  TECCH


  }

  else if `"`techtype'"'=="biennial"{

    qui bys `dmu' (`period'): gen TFPCH=`D12'/`DD'[_n-1] if _n>1
    label var TFPCH "Total factor productivity change"
    cap drop `temp'

    sort `period' `dmu'
    forv t=1/`tmax'{
      qui replace `flag'=(`period'==`t')
      _ddf  if `touse' & `period'==`t', rflag(`flag') gen(`DD') gv(`gmat') `vrs' ort(`ort') in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

    }

    qui bys `dmu' (`period'): gen TECH=`DD'/`DD'[_n-1]
    qui bys `dmu' (`period'): gen TECCH=TFPCH/TECH

    label var TECH  "Technical efficiency change"
    label var TECCH "Techological change"
    local resvars TFPCH TECH  TECCH

  }

  else{

      //su `DD' `D12' `D21'
    qui {
      sort `dmu' `period'
      bys `dmu' (`period'): gen TECH=`DD'/`DD'[_n-1]
      bys `dmu' (`period'): gen TECCH=sqrt(`D12'/`DD'*`DD'[_n-1]/`D21'[_n-1])
      gen TFPCH= TECH*TECCH
      local resvars TFPCH TECH TECCH
        label var TFPCH "Total factor productivity change"
        label var TECH  "Technical efficiency change"
        label var TECCH "Techological change"
    }

  }

  return local rvars "`resvars'"


end

///////////////////////////////////////////////////////////////////////////
capture program drop _ddf
program define _ddf
    version 16

    syntax [if] [in], gen(string) INvars(varlist) OPvars(varlist) BADvars(varlist) [GVec(varlist) rflag(varname) ort(string) VRS maxiter(numlist) tol(numlist)]


    marksample touse
    markout `touse' `invars' `opvars' `badvars' `gvec'

    tempvar touse2
    local rflag=cond("`rflag'"=="","","if `rflag'")
    mark `touse2' `rflag' // `rflag' might be empty

    markout `touse2' `invars' `opvars' `badvars'

     local data `invars' `opvars' `badvars'
     local num1: word count `invars'
     local num2: word count `opvars'

    if "`vrs'"!=""{
      local rts=1
    }
    else{
      local rts=0

    }

    mata: _ddf("`data'",`num1',`num2',"`touse'", "`touse2'","`gvec'","`gen'",`rts',`maxiter',`tol')

	//di "`ort'"

    if "`ort'" =="out"{

      qui replace `gen'=1/(1+`gen') if `touse'
      }
    else{
      qui replace `gen'=(1-`gen') if `touse'
     }


end

//////////////////////////////////////////////////
cap program drop _luen_ddf
program define _luen_ddf,rclass
    version 16


    gettoken word 0 : 0, parse("=:,")
    while `"`word'"' != ":" & `"`word'"' != "=" {
        if `"`word'"' == "," | `"`word'"'=="" {
                error 198
        }
    local invars `invars' `word'
    gettoken word 0 : 0, parse("=:,")
    }
    unab invars : `invars'

  gettoken word 0 : 0, parse("=:,")
    while `"`word'"' != ":" & `"`word'"' != "=" {
        if `"`word'"' == "," | `"`word'"'=="" {
                error 198
        }
    local gopvars `gopvars' `word'
    gettoken word 0 : 0, parse(" =:,")
    }
    unab gopvars : `gopvars'

*********************************************************************************

    *local ninv: word count `invars'
    *local ngov: word count `gopvars'
    syntax varlist [if] [in], id(varname) time(varname)  [gx(varlist) ///
                   BIennial   gy(varlist) gb(varlist) VRS ort(string) ///
                             GLOBAL SEQuential WINdow(numlist intege max=1 >=1) ///
                 maxiter(numlist integer >0 max=1) tol(numlist max=1 >0)]

  marksample touse
  local bopvars `varlist'
  local ninp: word count `invars'
  local ngo: word count `gopvars'
  local nbo: word count `bopvars'
  local nvar: word count `invars' `gopvars' `bopvars'
  confirm numeric var `invars' `gopvars' `bopvars'

 local techtype "contemporaneous"


   if "`global'"!=""{
     if "`sequential'"!=""{

       disp as error "global and sequential cannot be specified together."
       error 498

     }

     if "`window'"!=""{

       disp as error "global and window() cannot be specified together."
       error 498

     }

     if "`biennial'"!=""{

       disp as error "global and biennial cannot be specified together."
       error 498

     }

     local techtype "global"

  }


   if "`sequential'"!=""{

     if "`window'"!=""{

       disp as error "sequential and window() cannot be specified together."
       error 498

     }

     if "`biennial'"!=""{

       disp as error "sequential and biennial cannot be specified together."
       error 498

     }

     local techtype "sequential"

  }


     if "`window'"!=""{

       if "`biennial'"!=""{

         disp as error "biennial and window() cannot be specified together."
         error 498

       }

         local techtype "window"

     }

       if "`biennial'"!=""{

        local techtype "biennial"

       }

  local comvars: list invars & gopvars
  if !(`"`comvars'"'==""){
    disp as error "`comvars' should not be specified as input and desriable output simultaneously."
    error 498
  }

  local comvars: list invars & bopvars
  if !(`"`comvars'"'==""){
    disp as error "`comvars' should not be specified as input and undesriable output simultaneously."
    error 498
  }

  local comvars: list gopvars & bopvars
  if !(`"`comvars'"'==""){
    disp as error "`comvars' should not be specified as desriable and undesriable outputs simultaneously."
    error 498
  }


  if "`maxiter'"==""{
    local maxiter=-1
  }
  if "`tol'"==""{
    local tol=-1
  }


    tempvar period dmu

  qui egen `period'=group(`time')
  qui egen `dmu'=group(`id')


    qui su  `period'
    local tmax=r(max)

    tempvar flag temp DD D21 D12

    qui gen `DD'=.
    qui gen `D21'=.
    qui gen `D12'=.
    qui gen `temp'=.
    local gmat `gx' `gy' `gb'

    qui gen `flag'=0

    sort `period' `dmu'

  if `"`techtype'"'=="contemporaneous"{

      qui{
        forv t=1/`tmax'{
            qui replace `flag'= (`period'==`t')
            _ddf_luen  if `period'==`t' & `touse', rflag(`flag') gen(`DD') gv(`gmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        }
        local tt=`tmax'-1
        forv t=1/`tt'{
            qui replace `flag'=(`period'==`t'+1)
            _ddf_luen  if `period'==`t' & `touse', rflag(`flag') gen(`D21') gv(`gmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        }

        forv t=2/`tmax'{
            qui replace `flag'=(`period'==`t'-1)
            _ddf_luen  if `period'==`t' & `touse', rflag(`flag') gen(`D12') gv(`gmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        }

    }


  }


  if `"`techtype'"'=="biennial"{

      qui{
        forv t=1/`tmax'{
            qui replace `flag'= (`period'==`t' | `period'==`t'+1)
            _ddf_luen  if `period'==`t' & `touse', rflag(`flag') gen(`DD') gv(`gmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        }

        forv t=2/`tmax'{
            qui replace `flag'=(`period'==`t'-1 | `period'==`t')
            _ddf_luen  if `period'==`t' & `touse', rflag(`flag') gen(`D12') gv(`gmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        }

    }


  }


    if `"`techtype'"'=="sequential"{


        forv t=1/`tmax'{
            qui replace `flag'=(`period'<=`t')
            _ddf_luen if `period'==`t' & `touse', rflag(`flag') gen(`DD') gv(`gmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')
        }
        local tt=`tmax'-1
        forv t=1/`tt'{
            qui replace `flag'=(`period'<=`t'+1)
            _ddf_luen  if `period'==`t' & `touse', rflag(`flag') gen(`D21') gv(`gmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')
        }

        forv t=2/`tmax'{
            qui replace `flag'= (`period'<=`t'-1)
            _ddf_luen if `period'==`t' & `touse', rflag(`flag') gen(`D12') gv(`gmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        }


  }


     if `"`techtype'"'=="window"{
    *local band=(`window'-1)/2
      local band=`window'

        forv t=1/`tmax'{
            qui replace `flag'=(`period'<=`t'+`band' & `period'>=`t'-`band')
            _ddf_luen  if `period'==`t' & `touse', rflag(`flag') gen(`DD') gv(`gmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        }
        local tt=`tmax'-1
        forv t=1/`tt'{
            qui replace `flag'= (`period'<=`t'+1+`band' &  `period'>=`t'-`band'+1)
            _ddf_luen  if `period'==`t' & `touse', rflag(`flag') gen(`D21') gv(`gmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        }

        forv t=2/`tmax'{
            qui replace `flag'=(`period'<=`t'-1+`band' & `period'>=`t'-1-`band')
            _ddf_luen  if `period'==`t' & `touse', rflag(`flag') gen(`D12') gv(`gmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        }


  }


  if `"`techtype'"'=="global"{

      qui replace `flag'=1
    _ddf_luen  if `touse', rflag(`flag') gen(`temp') gv(`gmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        qui bys `dmu' (`period'): gen TFPCH=`temp'[_n-1]-`temp'
    label var TFPCH "Total factor productivity change"
    cap drop `temp'

    sort `period' `dmu'
    forv t=1/`tmax'{
      qui replace `flag'=(`period'==`t')
      _ddf_luen  if `touse' & `period'==`t', rflag(`flag') gen(`DD') gv(`gmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

    }

    qui bys `dmu' (`period'): gen TECH=`DD'[_n-1]-`DD'
    *qui bys `dmu' (`period'): gen BPC=TFPCH-TECH
	qui bys `dmu' (`period'): gen TECCH=TFPCH-TECH

    label var TECH  "Technical efficiency change"
    *label var BPC "Best practice gap change"
	label var TECCH "Technological change"
    *local resvars TFPCH TECH  BPC
	local resvars TFPCH TECH  TECCH


  }
  else if `"`techtype'"'=="biennial"{

     qui bys `dmu' (`period'): gen TFPCH=`DD'[_n-1]-`D12' if _n>1
     label var TFPCH "Total factor productivity change"
     qui cap drop `temp'
     sort `period' `dmu'
    forv t=1/`tmax'{
      qui replace `flag'=(`period'==`t')
      _ddf_luen  if `touse' & `period'==`t', rflag(`flag') gen(`DD') gv(`gmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

    }

    qui bys `dmu' (`period'): gen TECH=`DD'[_n-1]-`DD'
    qui bys `dmu' (`period'): gen TECCH=TFPCH-TECH

    label var TECH  "Technical efficiency change"
    label var TECCH "Techological change"
    local resvars TFPCH TECH  TECCH


 }

  else{

      //su `DD' `D12' `D21'
    qui {
      sort `dmu' `period'
      bys `dmu' (`period'): gen TECH=`DD'[_n-1]-`DD'
      bys `dmu' (`period'): gen TECCH=0.5*(`DD'-`D12'+`D21'[_n-1]-`DD'[_n-1])
      gen TFPCH= TECH+TECCH
      local resvars TFPCH TECH TECCH
        label var TFPCH "Total factor productivity change"
        label var TECH  "Technical efficiency change"
        label var TECCH "Techological change"
    }

  }

  return local rvars "`resvars'"


end


///////////////////////////////////////////////////////////////
capture program drop _ddf_luen
program define _ddf_luen
    version 16

    syntax [if] [in], gen(string) INvars(varlist) OPvars(varlist) BADvars(varlist) [GVec(varlist) rflag(varname)  VRS maxiter(numlist) tol(numlist)]


    marksample touse
    markout `touse' `invars' `opvars' `badvars' `gvec'
    local rflag=cond("`rflag'"=="","","if `rflag'") //fix missing cond 18-9-2021
    tempvar touse2
    mark `touse2' `rflag' // `rflag' might be empty
    markout `touse2' `invars' `opvars' `badvars'


    local data `invars' `opvars' `badvars'
    local num1: word count `invars'
    local num2: word count `opvars'


    if "`vrs'"!=""{
      local rts=1
    }
    else{
      local rts=0

    }

    mata: _ddf("`data'",`num1',`num2',"`touse'", "`touse2'","`gvec'","`gen'",`rts',`maxiter',`tol')


end

////////////////////////////////////////////////////////////////////////
cap program drop _luen_nddf
program define _luen_nddf,rclass
    version 16


    gettoken word 0 : 0, parse("=:,")
    while `"`word'"' != ":" & `"`word'"' != "=" {
        if `"`word'"' == "," | `"`word'"'=="" {
                error 198
        }
    local invars `invars' `word'
    gettoken word 0 : 0, parse("=:,")
    }
    unab invars : `invars'

  gettoken word 0 : 0, parse("=:,")
    while `"`word'"' != ":" & `"`word'"' != "=" {
        if `"`word'"' == "," | `"`word'"'=="" {
                error 198
        }
    local gopvars `gopvars' `word'
    gettoken word 0 : 0, parse(" =:,")
    }
    unab gopvars : `gopvars'


    *local ninv: word count `invars'
    *local ngov: word count `gopvars'
    syntax varlist [if] [in], id(varname) time(varname) wmat(string)  [ gx(varlist) ///
                           BIennial  gy(varlist) gb(varlist) VRS ort(string) ///
                             GLOBAL SEQuential WINdow(numlist intege max=1 >=1) ///
                 maxiter(numlist integer >0 max=1) tol(numlist max=1 >0)]

  marksample touse
  local bopvars `varlist'

  local ninp: word count `invars'
  local ngo: word count `gopvars'
  local nbo: word count `bopvars'
  local nvar: word count `invars' `gopvars' `bopvars'
  confirm numeric var `invars' `gopvars' `bopvars'
 local techtype "contemporaneous"


   if "`global'"!=""{
     if "`sequential'"!=""{

       disp as error "global and sequential cannot be specified together."
       error 498

     }

     if "`window'"!=""{

       disp as error "global and window() cannot be specified together."
       error 498

     }

     if "`biennial'"!=""{

       disp as error "global and biennial cannot be specified together."
       error 498

     }

     local techtype "global"

  }


   if "`sequential'"!=""{

     if "`window'"!=""{

       disp as error "sequential and window() cannot be specified together."
       error 498

     }

     if "`biennial'"!=""{

       disp as error "sequential and biennial cannot be specified together."
       error 498

     }

     local techtype "sequential"

  }


     if "`window'"!=""{

       if "`biennial'"!=""{

         disp as error "biennial and window() cannot be specified together."
         error 498

       }

         local techtype "window"

     }

       if "`biennial'"!=""{

        local techtype "biennial"

       }

  local comvars: list invars & gopvars
  if !(`"`comvars'"'==""){
    disp as error "`comvars' should not be specified as input and desriable output simultaneously."
    error 498
  }

  local comvars: list invars & bopvars
  if !(`"`comvars'"'==""){
    disp as error "`comvars' should not be specified as input and undesriable output simultaneously."
    error 498
  }

  local comvars: list gopvars & bopvars
  if !(`"`comvars'"'==""){
    disp as error "`comvars' should not be specified as desriable and undesriable outputs simultaneously."
    error 498
  }

/*
  tempname weightvec wcheck
  if "`nonradial'"==""{
    if `"`wmat'"'!=""{
      disp as red "Warning: wmat() is only used when nonradial is specified."
    }

  }
  else{
    if `"`wmat'"'!=""{
      confirm matrix `wmat'
      local ncol=colsof(`wmat')
      if `ncol'!=`nvar'{
        dis as error `"# of column of `wmat' != # of input-output variables"'
        exit 498
      }

    //mat `weightvec'=`wmat'
    if "`ort'"=="out"{
      mat `weightvec'=`wmat'[1,1..`ninp']
      mata: st_numscalar("`wcheck'",sum(st_matrix("`weightvec'")))
      if `wcheck'!=0{
        disp as error "For output orientied Luenberger productivity index, the weight components for input reduction should be set to 0."
        exit 498
      }
    }
    else{
      mat `weightvec'=`wmat'[1,(`ninp'+1)..(`ninp'+`nbo'+`ngo')]
      mata: st_numscalar("`wcheck'",sum(st_matrix("`weightvec'")))
      if `wcheck'!=0{
        disp as error "For input orientied Luenberger productivity index, the weight components for output reduction should be set to 0."
        exit 498
      }

    }
     mat `weightvec'=`wmat'
    }
    else{
       if "`ort'"=="out"{
          mat  `weightvec'=(J(1,`ninp',0),J(1,`ngo',1)*(1/2/`ngo'),J(1,`nbo',1)*(1/2/`nbo'))
       }
       else{
          mat  `weightvec'=(J(1,`ninp',1/`ninp'),J(1,`ngo'+`nbo',0))
       }

    }

  }

*/

  if "`maxiter'"==""{
    local maxiter=-1
  }
  if "`tol'"==""{
    local tol=-1
  }


    tempvar period dmu

  qui egen `period'=group(`time')
  qui egen `dmu'=group(`id')


    qui su  `period'
    local tmax=r(max)

    tempvar flag temp DD D21 D12

    qui gen `DD'=.
    qui gen `D21'=.
    qui gen `D12'=.
    qui gen `temp'=.
    qui gen `flag'=0

    local gmat `gx' `gy' `gb'
    sort `period' `dmu'


  if `"`techtype'"'=="contemporaneous"{

      qui{
        forv t=1/`tmax'{
            qui replace `flag'= (`period'==`t')

            _nddf  if `period'==`t' & `touse', rflag(`flag') gen(`DD') gv(`gmat') wmat(`wmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        }
        local tt=`tmax'-1
        forv t=1/`tt'{
            qui replace `flag'=(`period'==`t'+1)
            _nddf  if `period'==`t' & `touse', rflag(`flag') gen(`D21') gv(`gmat') wmat(`wmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        }

        forv t=2/`tmax'{
            qui replace `flag'=(`period'==`t'-1)
            _nddf  if `period'==`t' & `touse', rflag(`flag') gen(`D12') gv(`gmat') wmat(`wmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        }

    }


  }


  if `"`techtype'"'=="biennial"{

      qui{
        forv t=1/`tmax'{
            qui replace `flag'= (`period'==`t' | `period'==`t'+1)

             _nddf  if `period'==`t' & `touse', rflag(`flag') gen(`DD') gv(`gmat') wmat(`wmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        }

        forv t=2/`tmax'{
            qui replace `flag'=(`period'==`t'-1 | `period'==`t')
            _nddf  if `period'==`t' & `touse', rflag(`flag') gen(`D12') gv(`gmat') wmat(`wmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        }

    }


  }


    if `"`techtype'"'=="sequential"{


        forv t=1/`tmax'{
            qui replace `flag'=(`period'<=`t')
            _nddf if `period'==`t' & `touse', rflag(`flag') gen(`DD') gv(`gmat') wmat(`wmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        }
        local tt=`tmax'-1
        forv t=1/`tt'{
            qui replace `flag'=(`period'<=`t'+1)
            _nddf  if `period'==`t' & `touse', rflag(`flag') gen(`D21') gv(`gmat') wmat(`wmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        }

        forv t=2/`tmax'{
            qui replace `flag'= (`period'<=`t'-1)
            _nddf if `period'==`t' & `touse', rflag(`flag') gen(`D12') gv(`gmat') wmat(`wmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        }


  }


     if `"`techtype'"'=="window"{
    local band=`window'

        forv t=1/`tmax'{
            qui replace `flag'=(`period'<=`t'+`band' & `period'>=`t'-`band')
            _nddf  if `period'==`t' & `touse', rflag(`flag') gen(`DD') gv(`gmat') wmat(`wmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        }
        local tt=`tmax'-1
        forv t=1/`tt'{
            qui replace `flag'= (`period'<=`t'+1+`band' &  `period'>=`t'-`band'+1)
            _nddf  if `period'==`t' & `touse', rflag(`flag') gen(`D21') gv(`gmat') wmat(`wmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        }

        forv t=2/`tmax'{
            qui replace `flag'=(`period'<=`t'-1+`band' & `period'>=`t'-1-`band')
            _nddf  if `period'==`t' & `touse', rflag(`flag') gen(`D12') gv(`gmat') wmat(`wmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        }


  }


  if `"`techtype'"'=="global"{

      qui replace `flag'=1
    _nddf  if `touse', rflag(`flag') gen(`temp') gv(`gmat') wmat(`wmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

        qui bys `dmu' (`period'): gen TFPCH=`temp'[_n-1]-`temp'
    label var TFPCH "Total factor productivity change"
    cap drop `temp'

    sort `period' `dmu'
    forv t=1/`tmax'{
      qui replace `flag'=(`period'==`t')
      _nddf  if `touse' & `period'==`t', rflag(`flag') gen(`DD') gv(`gmat') wmat(`wmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

    }

    qui bys `dmu' (`period'): gen TECH=`DD'[_n-1]-`DD'
    *qui bys `dmu' (`period'): gen BPC=TFPCH-TECH
	qui bys `dmu' (`period'): gen TECCH=TFPCH-TECH

    label var TECH  "Technical efficiency change"
    *label var BPC "Best practice gap change"
	label var TECCH "Technological Change"
    *local resvars TFPCH TECH  BPC
	local resvars TFPCH TECH  TECCH


  }
  else  if `"`techtype'"'=="biennial"{

    qui bys `dmu' (`period'): gen TFPCH=`DD'[_n-1]-`D12'
    label var TFPCH "Total factor productivity change"
    cap drop `temp'

    sort `period' `dmu'
    forv t=1/`tmax'{
      qui replace `flag'=(`period'==`t')
      _nddf  if `touse' & `period'==`t', rflag(`flag') gen(`DD') gv(`gmat') wmat(`wmat') `vrs'  in(`invars') op(`gopvars') bad(`bopvars') maxiter(`maxiter') tol(`tol')

    }

    qui bys `dmu' (`period'): gen TECH=`DD'[_n-1]-`DD'
    qui bys `dmu' (`period'): gen TECCH=TFPCH-TECH

    label var TECH  "Technical efficiency change"
    label var TECCH "Techological change"
    local resvars TFPCH TECH  TECCH


  }

  else{

      //su `DD' `D12' `D21'
    qui {
      sort `dmu' `period'
      bys `dmu' (`period'): gen TECH=`DD'[_n-1]-`DD'
      bys `dmu' (`period'): gen TECCH=0.5*(`DD'-`D12'+`D21'[_n-1]-`DD'[_n-1])
      gen TFPCH= TECH+TECCH
      local resvars TFPCH TECH TECCH
        label var TFPCH "Total factor productivity change"
        label var TECH  "Technical efficiency change"
        label var TECCH "Techological change"
    }

  }

  return local rvars "`resvars'"

end

//////////////////////////////////////////////////////////////
capture program drop _nddf
program define _nddf
    version 16

    syntax [if] [in], gen(string) INvars(varlist) OPvars(varlist) ///
	                  BADvars(varlist) wmat(string) GVec(varlist) ///
					  [rflag(varname)  VRS maxiter(numlist) tol(numlist)]


    marksample touse
    markout `touse' `invars' `opvars' `badvars' `gvec'

    tempvar touse2
     local rflag=cond("`rflag'"=="",""," if `rflag'")
    mark `touse2' `rflag' // `rflag' might be empty
    markout `touse2' `invars' `opvars' `badvars'


    local data `invars' `opvars' `badvars'
    local num1: word count `invars'
    local num2: word count `opvars'

    if "`vrs'"!=""{
      local rts=1
    }
    else{
      local rts=0

    }

    mata: _nddf("`data'",`num1',`num2',"`touse'", "`touse2'","`gvec'","`wmat'","`gen'",`rts',`maxiter',`tol')


end


*******************************************
///////////////////END//////////////////////
********************************************