*! version 1.0.3  30mar99
*
* syntax: outreg5 [varlist] using filename [, SE noLabel Bdec(int 3) Tdec(int 2)
*                noPAren noASter noCONs EForm CLAssic AT noNOBs noNI noOR noR2 
*                ADjr2 noNOTes COMma APpend REPLACE]
*

program define outreg5
   version 5.0
   * write formatted regression output to file

   local varlist "opt ex none"
   local using "req noprefix"
   local options "SE noLabel Bdec(int 3) Tdec(int 2) noPAren noASter noCONs EForm CLAssic AT(string) noNOBs noNI noR2 ADjr2 noNOTes COMma APpend REPLACE *"
   parse "`*'"

   tempname regobs n_k rsq adjrsq numi
   scalar `n_k'    = .
   scalar `rsq'    = .
   scalar `adjrsq' = .
   scalar `numi'   = .
   local   depvar  $S_E_depv
   local   robust  = "$S_E_vce"
   local   ivar    "."
   local   univar  = 1     /* univariate (single equation) regression */

   /* logistic reports coeffients in exponentiated form (odds ratios) */
   if "$S_E_cmd"=="logistic" {local eform "eform"}
   if ("`eform'"=="eform" | "$S_E_cmd"=="dprobit") {local cons "nocons"}
   if "`at'"!="" {local at = "at(`at')"}
   
   /* linear regression models */
   if ("$S_E_cmd"=="regress" | "$S_E_cmd"=="anova" | "$S_E_cmd"=="corc"    /*
      */| "$S_E_cmd"=="fit"  | "$S_E_cmd"=="glogit"| "$S_E_cmd"=="gprobit" /*
      */| "$S_E_cmd"=="regdw") {
      scalar `regobs' = _result(1)
      scalar `n_k'    = _result(5)
      scalar `rsq'    = _result(7)
      scalar `adjrsq' = _result(8)
   }

   /* panel data models */
   else if "$S_E_cmd"=="xtreg" | "$S_E_cmd"=="xtgee" | "$S_E_cmd"=="xtgls"  /*
        */ | "$S_E_cmd"=="xtpois" | "$S_E_cmd"=="xtprobit" {
      if ("$S_E_cmd"=="xtreg" &  ("$S_E_cmd2"=="xtreg_fe" | /*
         */ "$S_E_cmd2"=="xtreg_be" | "$S_E_cmd2"=="xtreg_re")) {
         scalar `regobs' = scalar(S_E_nobs)
         scalar `numi'   = scalar(S_E_n)
         local   ivar    $S_E_ivar
         if ("$S_E_cmd2"=="xtreg_fe" | "$S_E_cmd2"=="xtreg_be") {
            scalar `n_k' = $S_E_tdf
         }
      }
      else {
         scalar `regobs' = $S_E_nobs
         if "$S_E_cmd2"~="xtreg_ml" {
            local ivar $S_E_ivar
            if "$S_E_cmd"=="xtgls" {scalar `numi' = $S_E_ng}
            else {scalar `numi' = $S_E_ncl}
         }
      }
   }

   /* multi-equation models */
   else if "$S_E_cmd"=="sureg" | "$S_E_cmd"=="mvreg" | "$S_E_cmd"=="reg3" {
      local univar  = 0    /* multivariate (multiple equation) regression */
      local neq     $S_E_neq
      if "$S_E_cmd"=="reg3" {
         local eqlist  $S_E_eqns
         local r2 "nor2"  /* R2 not yet in compatible form in reg3 macros */
      }
      else {
         local eqlist  $S_E_elis
         scalar `n_k'  = $S_E_tdf
      }
      scalar `regobs' = $S_E_nobs
      local r2list  $S_E_r2
   }
   
   /* models reporting nobs etc. with S_E macros */
   else if "$S_E_cmd"=="areg" | "$S_E_cmd"=="bsqreg" | "$S_E_cmd"=="cnsreg"  /*
     */| "$S_E_cmd"=="dlogit2"|"$S_E_cmd"=="dmlogit2"| "$S_E_cmd"=="dprobit2"/*
     */| "$S_E_cmd"=="eivreg" | "$S_E_cmd"=="ereg"   | "$S_E_cmd"=="gnbreg"  /*
     */| "$S_E_cmd"=="heckman"| "$S_E_cmd"=="hlu"    | "$S_E_cmd"=="intreg"  /*
     */| "$S_E_cmd"=="nbreg"  | "$S_E_cmd"=="newey"  | "$S_E_cmd"=="nl"      /*
     */| "$S_E_cmd"=="poisson"| "$S_E_cmd"=="prais"  | "$S_E_cmd"=="qreg"    /*
     */| "$S_E_cmd"=="stereg" | "$S_E_cmd"=="stweib" |"$S_E_cmd"=="svylogit" /*
     */|"$S_E_cmd"=="svyprobt"| "$S_E_cmd"=="svyreg" | "$S_E_cmd"=="vwls"    /*
     */| "$S_E_cmd"=="weibull"| "$S_E_vce"=="Newey-West" {
      scalar `regobs' = $S_E_nobs
		capture scalar `n_k'    = $S_E_tdf
      if "$S_E_cmd"=="nl" {
         scalar `rsq'     = $S_E_rsq
         scalar `adjrsq'  = $S_E_rsqa
      }
      else if "$S_E_cmd"=="dmlogit2" {local univar  = 0}/* multivariate regression */
   }
   
   /* models using Stata internal [R]maximize */
   else if "$S_E_cmd"=="blogit" | "$S_E_cmd"=="bprobit"|"$S_E_cmd"=="clogit" /*
     */| "$S_E_cmd"=="cnreg"  | "$S_E_cmd"=="cox"     |"$S_E_cmd"=="dprobit" /*
     */| "$S_E_cmd"=="glm"    | "$S_E_cmd"=="logistic"|"$S_E_cmd"=="logit"   /*
     */| "$S_E_cmd"=="mlogit" | "$S_E_cmd"=="ologit"  |"$S_E_cmd"=="oprobit" /*
     */| "$S_E_cmd"=="probit" | "$S_E_cmd"=="rreg"    |"$S_E_cmd"=="stcox" /*
     */| "$S_E_cmd"=="tobit" {
      scalar `regobs' = _result(1)
      if "$S_E_cmd"=="mlogit" {local univar  = 0}/* multivariate regression */
   }
      
   quietly {
      preserve

      tempname coefcol varname
      tempname b vc
      mat `b' = get(_b)
      mat `vc' = get(VCE)
      if "`robust'"=="" {local robust "none"}
      
      if `univar' {
         coeftxt5 `varlist', `se' `label' bdec(`bdec') tdec(`tdec') `paren' `aster' `cons' `eform' `classic' `at' `nobs' `ni' `r2' `adjr2' `notes' `comma' `append' regobs(`regobs') n_k(`n_k') rsq(`rsq') adjrsq(`adjrsq') numi(`numi') ivar(`ivar') depvar(`depvar') robust(`robust') b(`b') vc(`vc') varname(`varname') coefcol(`coefcol')
   
         if "`append'"=="append" {
            appfile5 using `using', varname(`varname') coefcol(`coefcol') 
            outsheet v* `coefcol' using "`using'", nonames noquote `comma' replace
            drop v*
         }
         else {   
            noi outsheet `varname' `coefcol' using "`using'", nonames noquote `comma' `replace'
         }
      }
      else if ("$S_E_cmd"=="mlogit" | "$S_E_cmd"=="dmlogit2") {   
         /* multinomial logit */
         tempname b_cat vc_cat
         tab `depvar'
         local ncat = _result(2)-1 /* number of alternatives in depvar */
         local k = colsof(`vc')/`ncat'
         local oldepv `depvar'
         local cat = 1
         local catlist : roweq(`vc')
         local catword = 1
         while `cat' <= `ncat' {
            restore, preserve  /* to access variable labels after first equation */
            local catword = (`cat'-1)*`k' + 1 
            local catname: word `catword' of `catlist'
            if "`label'"=="" {local depvar "`catname'"}
            else {local depvar "`oldepv' `cat'"}
            if "$S_E_cmd"=="mlogit" {mat `b_cat'   = `b'[`cat',.]}
            else {mat `b_cat'   = `b'[.,"`catname':"]}
            mat `vc_cat'  = `vc'["`catname':","`catname':"]
   
            if `cat' == 1 & "`append'"~="append" {local apptmp ""}
            else {local apptmp "append"}
            coeftxt5 `varlist', `se' `label' bdec(`bdec') tdec(`tdec') `paren' `aster' `cons' `eform' `classic' `at' `nobs' `ni' `r2' `adjr2' `notes' `comma' `apptmp' regobs(`regobs') n_k(`n_k') rsq(`rsq') adjrsq(`adjrsq') numi(`numi') depvar(`depvar') robust(`robust') b(`b_cat') vc(`vc_cat') varname(`varname') coefcol(`coefcol')
            if `cat' == 1 & "`append'"~="append" {
               noi outsheet `varname' `coefcol' using "`using'", nonames noquote `comma' `replace'
            }
            else {
               appfile5 using `using', varname(`varname') coefcol(`coefcol') 
               outsheet v* `coefcol' using "`using'", nonames noquote `comma' replace
               drop v*
            }
            local cat = `cat' + 1
         }
      }
      else {   /* multiple equation regression */
         tempname b_eq vc_eq
         local eq = 1
         while `eq' <= `neq' {
            restore, preserve  /* to access variable labels after first equation */
            local eqname: word `eq' of `eqlist'
            eq ? `eqname'
            local depvar: word 1 of $S_1
            local r2str: word `eq' of `r2list'
            scalar `rsq' = real("`r2str'")
            mat `b_eq'   = `b'[.,"`eqname':"]
            mat `vc_eq'  = `vc'["`eqname':","`eqname':"]
   
            if `eq' == 1 & "`append'"~="append" {local apptmp ""}
            else {local apptmp "append"}
            coeftxt5 `varlist', `se' `label' bdec(`bdec') tdec(`tdec') `paren' `aster' `cons' `eform' `classic' `at' `nobs' `ni' `r2' `adjr2' `notes' `comma' `apptmp' regobs(`regobs') n_k(`n_k') rsq(`rsq') adjrsq(`adjrsq') numi(`numi') depvar(`depvar') robust(`robust') b(`b_eq') vc(`vc_eq') varname(`varname') coefcol(`coefcol')
            if `eq' == 1 & "`append'"~="append" {
               noi outsheet `varname' `coefcol' using "`using'", nonames noquote `comma' `replace'
            }
            else {
               appfile5 using `using', varname(`varname') coefcol(`coefcol') 
               outsheet v* `coefcol' using "`using'", nonames noquote `comma' replace
               drop v*
            }
            local eq = `eq' + 1
         }
      }
      
   }  /* for quietly */
end


program define coeftxt5
   version 5.0
   * convert vectors b and vc into coefficient and t-statistics text + variable name text

   local varlist "opt ex none"
   local options "SE noLabel BDEC(int 3) TDEC(int 2) noPAren noASter noCONs EForm classic at(string) noNOBs noNI noR2 ADjr2 noNOTes APpend regobs(string) n_k(string) rsq(string) adjrsq(string) numi(string) ivar(string) depvar(string) robust(string) b(string) vc(string) varname(string) coefcol(string)"
   parse "`*'"

   tempvar beta st_err tstat dotloc lead0 astrix mrgrow 
   tempfile bcopy 

    /* get rid of sigmas for xtreg_ml and _cut* for ologit in b and vc */
   if "$S_E_cmd"=="ologit" {
      local varnum = _result(3)  
      mat `b' = `b'[.,1..`varnum']
      mat `vc'  = `vc'[1..`varnum',1..`varnum']
      local cons "nocons"
   }
   if ("$S_E_cmd2"=="xtreg_ml") {  
      mat `b' = `b'[.,"`depvar':"]
      mat `vc'  = `vc'["`depvar':","`depvar':"]
   }
   else if "$S_E_cmd"=="dprobit" {
      tempname db dvc
      if "`at'"!="" {local at = "at(`at')"}
      dprbcal5 `varlist', `at' `classic' b(`b') vc(`vc') db(`db') 
   }
   
   mat `b' = `b''
   mat `vc' = vecdiag(`vc')
   mat `vc' = `vc''

   local bdecno = 10^(-`bdec')
   local tdecno = 10^(-`tdec')
   local zeros "00000000000000000000000"  
   local maxdec = max(`bdec', `tdec')
   local zeros = substr("`zeros'",1,`maxdec')
   
   if ("`label'"=="" & "`depvar'"!="" /*
     */ & "$S_E_cmd"!="mlogit" & "$S_E_cmd"!="dmlogit2") {
      local depname : var label `depvar'
      }   /* for mlogit, depvar is a category, not a variable */
   if "`depname'"=="" {local depname "`depvar'"}

   if (`numi'~=. & "`ni'"~="noni") {
      if ("`label'"=="" & "`ivar'"~=".") {local iname : var label `ivar'}
      if "`iname'"=="" {local iname "`ivar'"}
      if "`iname'"=="." {local iname "groups"}
   }

   /* in case varlist is specified:  */
   tempname arow newb newvc
   local vname : word 1 of `varlist'
   local i=1
   while "`vname'"~="" {
      local j = rownumb(`b',"`vname'")
      if `j'!=. {
         matrix `arow' = `b'[`j',.]
         matrix `newb' = `newb'\ `arow'
         matrix `arow' = `vc'[`j',.]
         matrix `newvc' = `newvc'\ `arow'
      }
      else if ("$S_E_cmd"!="sureg" & "$S_E_cmd"!="reg3") {
         di in red "`vname' not found in regression coefficients"
         exit 111
      }
      local i = `i'+1
      local vname : word `i' of `varlist'
   }

   if "`varlist'"~="" {
      if "`cons'"~="nocons" {     /* add the constant unless "nocons" is chosen */
         local j = rownumb(`b',"_cons")
         matrix `arow' = `b'[`j',.]
         matrix `newb' = `newb'\ `arow'
         matrix `arow' = `vc'[`j',.]
         matrix `newvc' = `newvc'\ `arow'
      }
      matrix `b' = `newb'
      matrix `vc' = `newvc'
   }
   else if "`cons'"=="nocons" {   /* delete the constant if "nocons" is chosen */
      local j_1 = rownumb(`b',"_cons")-1
      matrix `b' = `b'[1..`j_1',.]
      matrix `vc' = `vc'[1..`j_1',.]
   }

   local brows = rowsof(`b')
   gen double `beta' = matrix(`b'[_n, 1]) in 1/`brows'
   keep if `beta'~=.
   
   gen double `st_err' = matrix(`vc'[_n, 1]) in 1/`brows'
   replace `st_err' = sqrt(`st_err')
   gen double `tstat' = abs(`beta'/`st_err')
   if "`eform'"=="eform" {
      replace `beta' = exp(`beta')
      replace `st_err' = `beta'*`st_err'
   }
   if "$S_E_cmd"=="dprobit" {
      replace `beta' = matrix(`db'[_n,1]) in 1/`brows'
      replace `st_err' = matrix(`db'[_n,2]) in 1/`brows'
   }

   gen str31 `varname' = ""
   local bnames : rownames(`b')
   local bname : word 1 of `bnames'
   local i 1
   while "`bname'"~="" {
     if "`bname'"~="_cons" {
        if "`label'"=="" {local blabel : var label `bname'}
        else {local blabel ""}
     }
     else { 
        local blabel "Constant" 
     }
     if "`blabel'"=="" { 
        local blabel "`bname'" 
     }
     replace `varname' = "`blabel'" if `i'==_n
     local i = `i'+1
     local bname : word `i' of `bnames'
   }
     /* get rid of original data since labels accessed */
   keep `varname' `beta' `tstat' `st_err'  
   
   gen `mrgrow' = 2*_n + 1 + ("`append'"=="append")
   if "`aster'"~="noaster" {
      if `n_k'==. {
         scalar Tcrit1 = invnorm(0.995)
         scalar Tcrit5 = invnorm(0.975)
      }
      else {
         scalar Tcrit1 = invt(`n_k',0.99)
         scalar Tcrit5 = invt(`n_k',0.95)
      }
      gen str2 `astrix' = "*" if `tstat'>=Tcrit5
      replace `astrix' = "**" if `tstat'>=Tcrit1 
   }
   else {
      gen str2 `astrix' = ""
   }      
   if "`paren'"~="noparen" {
      local lparen "("
      local rparen ")"
   }
   else {
      local lparen ""
      local rparen ""
   }   
   if "`se'"~="se" {
      gen str12 `coefcol' = string(round(`tstat',`tdecno'))
      local t_bdec `tdec'
   }
   else {
      gen str12 `coefcol' = string(round(`st_err',`bdecno'))
      local t_bdec `bdec'
   }
   gen `dotloc' = index(`coefcol', ".")
   replace `coefcol' = substr(`coefcol'+"`zeros'", 1, `dotloc'+`t_bdec') if `dotloc' > 0
   replace `coefcol' = `coefcol'+"."+substr("`zeros'",1,`t_bdec') if `dotloc' == 0
   replace `coefcol' = "0" + `coefcol' if `dotloc'==1
   replace `coefcol' = "-0" + substr(`coefcol',2,length(`coefcol')) /* 
       */ if index(`coefcol', "-")==1 & `dotloc'==2
   replace `coefcol' = "`lparen'" + `coefcol' + "`rparen'" + `astrix'
   sort `mrgrow'
   save "`bcopy'", replace
   replace `mrgrow' = `mrgrow'-1
   replace `coefcol' = string(round(`beta',`bdecno'))
   replace `dotloc'  = index(`coefcol', ".")
   replace `coefcol' = substr(`coefcol'+"`zeros'", 1, `dotloc'+`bdec') if `dotloc' > 0
   replace `coefcol' = `coefcol'+"."+substr("`zeros'",1,`bdec') if `dotloc' == 0
   replace `coefcol' = "0" + `coefcol' if `dotloc'==1
   replace `coefcol' = "-0" + substr(`coefcol',2,length(`coefcol')) /* 
       */ if index(`coefcol', "-")==1 & `dotloc'==2
   sort `mrgrow'
   merge `mrgrow' using "`bcopy'"
   replace `varname' = " " if _merge==2
   drop `beta' `st_err' `tstat' `astrix' `dotloc' _merge
   sort `mrgrow'
   save "`bcopy'", replace

   drop `varname' `coefcol'
   local coefrow = 2*`brows'+1+("`append'"=="append")
   local totrows = `coefrow' + ("`nobs'"~="nonobs") + /*
    */(`numi'~=.) + ("`r2'"~="nor2"&`rsq'~=.) + /*
    */("`notes'"~="nonotes"&"`append'"~="append")*(1+("`aster'"~="noaster"))

   local minobs = max(1000,`brows')
   set obs `minobs'  /* obs must be set to >= 100 (depending on matsize?) */
   keep if _n <= `totrows'
   replace `mrgrow' = _n
   sort `mrgrow'
   merge `mrgrow' using "`bcopy'"
   sort `mrgrow'
   drop _merge
   
   replace `coefcol' = "`depname'" if _n==1 & ("`append'"~="append")
   replace `coefcol' = "`depname'" if _n==2 & ("`append'"=="append")
   if "`nobs'"~="nonobs" {
      local coefrow = `coefrow'+1
      replace `varname' = "Observations" if _n==`coefrow'
      replace `coefcol' = string(`regobs') if _n==`coefrow'
   }
   if (`numi'~=. & "`ni'"~="noni") {
      local coefrow = `coefrow'+1
      replace `varname' = "Number of `iname'" if _n==`coefrow'
      replace `coefcol' = string(`numi') if _n==`coefrow'
   }
   if "`r2'"~="nor2" & `rsq'~=. {
      local coefrow = `coefrow'+1
      if "`adjr2'"~="adjr2" {
         replace `varname' = "R-squared" if _n==`coefrow'
         replace `coefcol' = string(round(`rsq',`tdecno')) if _n==`coefrow'
            /* same decimal places for rsq as for tstats */
      }
      else {
         replace `varname' = "Adjusted R-squared" if _n==`coefrow'
         replace `coefcol' = string(round(`adjrsq',`tdecno')) if _n==`coefrow'
      }
      gen `dotloc'  = index(`coefcol', ".") if _n==`coefrow'
      replace `coefcol' = substr(`coefcol'+"`zeros'", 1, `dotloc'+`tdec') if _n==`coefrow' & `dotloc' > 0
      replace `coefcol' = `coefcol'+"."+substr("`zeros'",1,`tdec') if _n==`coefrow' & `dotloc' == 0
      replace `coefcol' = "0"+`coefcol' if _n==`coefrow' & `dotloc' == 1
   }
   if "`notes'"~="nonotes" & "`append'"~="append" {
      local coefrow = `coefrow'+1
      if "`se'"~="se" {
         if `n_k'~=. {
            local t_or_z "t"
         }
         else {
            local t_or_z "z"
         }
         if "`robust'"=="none" {
            replace `varname' = "Absolute value of `t_or_z'-statistics in parentheses" if _n==`coefrow'
         }
         else{
            replace `varname' = "`robust' `t_or_z'-statistics in parentheses" if _n==`coefrow'
         }
      }
      else {
         if "`robust'"=="none" {
            replace `varname' = "Standard errors in parentheses" if _n==`coefrow'
         }
         else{
            replace `varname' = "`robust' standard errors in parentheses" if _n==`coefrow'
         }
      }
      if "`aster'"~="noaster" {
         local coefrow = `coefrow'+1
         replace `varname' = "* significant at 5% level; ** significant at 1% level" if _n==`coefrow'
      }
    }

end  /* coeftxt5 */


program define appfile5
   version 5.0
   * append regression results to pre-existing file
   
   local using "req noprefix"
   local options "varname(string) coefcol(string)"
   parse "`*'"

   tempname vorder vartmp varsml
   tempfile tmpf1
   gen str80 `vartmp' = `varname'
   replace `vartmp' = "0" if _n==1
   replace `vartmp' = `vartmp'[_n-1]+"1" if trim(`vartmp')==""
   keep `varname' `coefcol' `vartmp'
   sort `vartmp'
   save `tmpf1', replace
   if index("`using'", ".")==0 {local using = "`using'"+".out"}
   insheet using "`using'", nonames clear
   gen str80 `vartmp' = v1
   replace `vartmp' = "01" if _n==1 & v2~="(1)"
   replace `vartmp' = "0" if _n==1 & v2=="(1)"
   replace `vartmp' = `vartmp'[_n-1]+"1" if `vartmp'==""
   gen `vorder' = _n+(1*(v2[1]~="(1)"))
   gen str10 `varsml' = `vartmp'
   count if (`varsml'=="R-squared" | `varsml'=="Adjusted R")
   local oldr2 = _result(1)
   count if (`varsml'=="Constant")
   local oldcon = _result(1)
   drop `varsml'
   sort `vartmp'

   merge `vartmp' using `tmpf1'
   gen str10 `varsml' = `vartmp'
   count if (`varsml'=="Constant" | `varsml'=="Constant1" /*
       */| `varsml'=="Observatio" | `varsml'=="R-squared"  /*
       */| `varsml'=="Adjusted R" | `varsml'=="Robust t-s" /*
       */| `varsml'=="Robust z-s" | `varsml'=="Absolute v" /*
       */| `varsml'=="Robust sta" | `varsml'=="Standard e" /*
       */| `varsml'=="* signific" | `varsml'=="Number of " /*
       */| `varsml'=="Newey-West" )
   local noncoef = _result(1)
   count if (`varsml'=="R-squared" | `varsml'=="Adjusted R")
   local newr2 = _result(1) - `oldr2'
   count if (`varsml'=="Constant")
   local newcon = _result(1) - `oldcon'
   count if _merge==2
   local newcoef = _result(1) - 2*`newcon' - `newr2'
   count if v2=="(1)"
   local endcoef = _N + 1*(_result(1)==0) - `noncoef' - `newcoef'
   replace `vorder' = `vorder'+_N-`endcoef'+2*`newcon' if `vorder'>`endcoef' & `vorder'~=.
   replace `vorder' = _n if `vorder'==.
   if `newcon' {  /* when constant in appended reg, but not prev reg */
      replace `vorder' = _N+1 if `vartmp'=="Constant" 
      replace `vorder' = _N+2 if `vartmp'=="Constant1"
   } 
   if `newr2' {  /* when R2 in appended reg, but not prev reg */
                 /* put after Constant & Constant1 (if any) */
      replace `vorder' = _N + 2*(`oldcon'|`newcon') + 0.5 /*
         */ if (`varsml'=="R-squared" | `varsml'=="Adjusted R")
   } 
   replace `vorder' = 1 if `vartmp'=="0"
   sort `vorder'
   replace v1 = `varname' if v1=="" & `varname'~=""
   drop `vorder' `varname' `vartmp' `varsml' _merge
   describe
   local numcol = _result(2)-1
   if (`numcol'==2) {
      replace v2 = "(1)" if _n==1
      replace `coefcol' = "(2)" if _n==1
   }
   else {
      replace `coefcol' = "(" + string(`numcol') + ")" if _n==1
   }
end  /* appfile5 */


program define dprbcal5
   version 5.0
   * calculate probability change for dprobit results 
   * rest of ado file adapted from dprobit.ado

   local varlist "opt ex none"
   local options "at(string) classic b(string) vc(string) db(string)"
   parse "`*'"

   tempname mns Vhat v xb fxb row

   if "`at'"!="" { 
      if colsof(`b')==colsof(matrix(`at')) { 
         if `at'[1,colsof(`b')] != 1 { 
            di in red "`at':  last (_cons) element not 1"
            exit 598
         }
         mat `mns' = `at'
      }
      else if colsof(`b')-1 == colsof(matrix(`at')) { 
         mat `mns' = (1)
         mat `mns' = `at',`mns' 
      }
      else { 
         di in red "`at':  conformability error" 
         exit 503
      }
   }
   else {mat `mns' = get(mns)}

   if "`classic'"!="" {
      local cons "_cons"
   }
   parse "$S_E_vl `cons'", parse(" ")

   Vhatdia5 `mns' `b' `vc' -> `Vhat'

   mat `xb' = `mns'*`b''
   scalar `xb' = `xb'[1,1]
   scalar `fxb' = normd(`xb')

   local lhs "`1'"
   mac shift

   local i 1
   local j 1
   local avar: word `j' of `varlist'
   local allvar = ("`avar'"=="")  /* i.e. `varlist' empty */
   while "``i''"!="" {
      if (`allvar' | "``i''"=="`avar'") {
         local isdum : word `i' of $S_E_dum
         if `isdum' & "`classic'"=="" {
            local c = /*
            */ normprob(`xb'+(1-`mns'[1,`i'])*`b'[1,`i']) - /*
            */ normprob(`xb'-(`mns'[1,`i'])*`b'[1,`i']) 
            vdummy5 `i' `mns' `b' `vc' -> `v'
         }
         else {
            local c = `b'[1,`i']*`fxb'
            scalar `v' = `Vhat'[1,`i']
         }
         local s = sqrt(`v')

         mat `row' = (`c',`s')
         mat `db' = `db' \ `row'

         if "``i''"=="`avar'" {
            local j=`j'+1
            local avar: word `j' of `varlist' 
         }
      }
      local i=`i'+1
   }
end


program define Vhatdia5 /* mns b V -> res */
   local xbar "`1'"
   local b "`2'"
   local V "`3'"
   local res "`5'"

   tempname tmp xb fxb dgdb

   mat `tmp' = `xbar'*`b''
   scalar `xb' = `tmp'[1,1]
   scalar `fxb' = normd(`xb')
   mat `dgdb' = `b'' * `xbar'
   mat `dgdb' = `xb'*`dgdb'
   local k = colsof(`V')
   mat `tmp' = I(`k')
   mat `dgdb' = `tmp' - `dgdb'
   mat `dgdb' = `fxb' * `dgdb'
   mat drop `tmp'

   mat `tmp' = `dgdb'*`V'
   mat `tmp' = `tmp'*`dgdb''
   mat `res' = vecdiag(`tmp')
end


program define vdummy5 /* i mns b V -> res */
   local i `1'
   local mns "`2'"
   local b "`3'"
   local V "`4'"
   local res "`6'"

   tempname x xb fxb Delp Vhat
   mat `x' = `mns'
   mat `x'[1,`i']=1
   mat `xb' = `x'*`b''
   scalar `fxb' = normd(`xb'[1,1])
   mat `Delp' = `fxb' * `x'

   mat `x'[1,`i']=0
   mat `xb' = `x'*`b''
   scalar `fxb' = normd(`xb'[1,1])
   mat `xb' = `fxb'*`x'
   mat `Delp' = `Delp' - `xb'

   mat `Vhat' = `Delp'*`V'
   mat `Vhat' = `Vhat'*`Delp''
   scalar `res' = `Vhat'[1,1]
end

exit