*! version 1.9.0 2012-07-29 fixes for o. in stata 11+
capture program drop fitstat
program define fitstat, rclass
local ecmdversion = c(version) // 190
version 9.0
local cmd `e(cmd)'
if ( "`cmd'"!="regress" & "`cmd'"!="logit" & "`cmd'"!="probit" ///
& "`cmd'"!="logistic" & "`cmd'"!="cloglog" & "`cmd'"!="ologit" ///
& "`cmd'"!="oprobit" & "`cmd'"!="gologit" & "`cmd'"!="mlogit" ///
& "`cmd'"!="clogit" & "`cmd'"!="tobit" & "`cmd'"!="cnreg" ///
& "`cmd'"!="intreg" & "`cmd'"!="poisson" & "`cmd'"!="nbreg" ///
& "`cmd'"!="zip" & "`cmd'"!="zinb" ///
& "`cmd'"!="ztp" & "`cmd'"!="ztnb" & "`cmd'"!="slogit" ///
& "`cmd'"!="mprobit" & "`cmd'"!="rologit" ///
) {
di in r "fitstat does not work with the last model estimated."
exit
}
tempname n_obs n_obs2 n_parm llF llR df x2 r2 b v_ystar v v_y v_x v_e
tempname factor r2mcf r2mcfa r2mz r2ml r2cu r2 r2ef r2adj aic aicn bicx
tempname bicp isbin isord isvare isystar islgt ispbt ssr dev lrx2 lrx2p
tempname depv2 maxtmp uniqgrp n_obs2 depvgrp ifin v_error bicdif fllR
tempname bnocon phat yphat ydev predpr predcat wrong hiprob fbicp fn_obs
tempname totcor totone maxrow cntr2 cntr2a chkobs d d2 prevval prevdf
tempname results erank Vb noomit
syntax [if] [in] ///
[, Saving(string) Using(string) SAVE Diff FORCE Bic NOSTATA ///
Debug(integer 0) ///
]
* 181 mat `b' = get(_b)
* 190: remove o. from beta matrix
mat `b' = get(_b)
mat `Vb' = e(V)
local colfn : colfullnames `b'
local coln : colnames `b'
local coleq : coleq `b'
_ms_omit_info `b' // remove columns of zeros from b
local cols = colsof(`b')
mat `noomit' = J(1,`cols',1) - r(omit)
forvalues i = 1/`cols' { // double clear for 0 std error
if `Vb'[`i',`i'] == 0 & `noomit'[1,`i']==1 {
mat `noomit'[1,`i']==0
}
}
mata: newbeta = ///
select(st_matrix(st_local("b")), (st_matrix(st_local("noomit"))))
mata: st_matrix(st_local("b"),newbeta)
// reassign column names to b
local i = 0
foreach var of local coln {
local ++i
local nmvar : word `i' of `coln'
local nmeq : word `i' of `coleq'
local nmfn : word `i' of `colfn'
_ms_parse_parts `var'
if "`nmvar'"=="`nmfn'" { // no eq nm
local nmeq ""
}
else { // eq nm
local nmeq "[`nmeq']"
}
if !`r(omit)' & `nmeq'_se[`nmvar']!=0 {
local coln2 `coln2' `nmfn'
}
}
version 11: matrix colnames `b' = `coln2'
local depv `e(depvar)'
if "`nostata'"!="nostata" { // stata bic and aic
tempname stataic statabic stataaic
qui estat ic
mat `stataic' = r(S)
scalar `statabic' = `stataic'[1,6]
scalar `stataaic' = `stataic'[1,5]
}
if ("`save'"!="") local saving 0
if ("`diff'"!="") local using 0
*-> get information from last model estimated
local wtis ""
local iswgt "no"
if "`e(wtype)'"!="" {
local iswgt "yes"
local wtis "[`e(wtype)'`e(wexp)']"
}
mat `v_e' = J(1,1,1)
sca `n_obs' = e(N)
sca `llF' = e(ll)
sca `llR' = e(ll_0)
sca `df' = e(df_m)
if "`cmd'"~="regress" & "`cmd'"~="tobit" & "`cmd'"~="cnreg" {
sca `x2' = e(chi2)
}
if "`cmd'"=="tobit" | "`cmd'"=="cnreg" {
sca `x2' = e(F)
}
if "`cmd'"=="regress" {
sca `ssr' = e(rss)
sca `r2' = e(r2)
sca `r2adj' = e(r2_a)
mat `v_e'[1,1] = (e(rmse)*e(rmse)*e(df_r)/e(N))
}
if "`cmd'"=="clogit" {
if "`e(wtype)'"=="iweight" |"`e(wtype)'"=="pweight" ///
| "`e(wtype)'"=="aweight" {
di in r "cannot use fitstat after clogit with `e(wtype)'"
exit
}
qui egen `uniqgrp' = rank(_n) if e(sample)==1, by(`e(group)')
qui replace `uniqgrp' = . if `uniqgrp' ~= 1
qui sum `uniqgrp' `wtis', meanonly
sca `n_obs' = r(sum_w)
}
*-> classify type of model
* binary?
local isbin "no"
if ("`cmd'"=="probit") local isbin "yes"
if ("`cmd'"=="logit") local isbin "yes"
if ("`cmd'"=="logistic") local isbin "yes"
if ("`cmd'"=="cloglog") local isbin "yes"
* ordinal?
local isord "no"
if ("`cmd'"=="oprobit") local isord "yes"
if ("`cmd'"=="ologit") local isord "yes"
* var(e) assumed?
local isvare "no"
if ("`cmd'"=="probit") local isvare "yes"
if ("`cmd'"=="logit") local isvare "yes"
if ("`cmd'"=="logistic") local isvare "yes"
if ("`cmd'"=="oprobit") local isvare "yes"
if ("`cmd'"=="ologit") local isvare "yes"
* is there a ystar?
local isystar "`isvare'"
if ("`cmd'"=="tobit") local isystar "yes"
if ("`cmd'"=="intreg") local isystar "yes"
if ("`cmd'"=="cnreg") local isystar "yes"
* logit model?
local islgt "no"
if ("`cmd'"=="logit") local islgt "yes"
if ("`cmd'"=="ologit") local islgt "yes"
if ("`cmd'"=="logistic") local islgt "yes"
* probit model?
local ispbt "no"
if ("`cmd'"=="probit") local ispbt "yes"
if ("`cmd'"=="oprobit") local ispbt "yes"
* zip/zinb?
local iszi "no"
if ("`cmd'"=="zip") local iszi "yes"
if ("`cmd'"=="zinb") local iszi "yes"
* provide pseudoR2 measures - jf 20Jun2005
local pseudo "yes"
if ("`cmd'"=="regress") local pseudo "no"
if ("`cmd'"=="rologit") local pseudo "no"
if ("`cmd'"=="slogit") local pseudo "no"
if ("`cmd'"=="mprobit") local pseudo "no"
* log-lik intercept only - jf 20Jun2005
local llio "yes"
if ("`cmd'"=="rologit") local llio "no"
if ("`cmd'"=="slogit") local llio "no"
if ("`cmd'"=="mprobit") local llio "no"
*-> compute numbers of variables and parameters
sca `n_parm' = colsof(`b')
local n_cat = e(k_cat)
* version 1.8.1 2007-09-18 - stata 10 changed returns for some commands
if c(stata_version) >= 10 {
local n_cat = e(k_out)
if ("`cmd'"=="ologit" | "`cmd'"=="oprobit") local n_cat = e(k_cat)
}
if "`n_cat'"=="." {
local n_cat = 2
}
local n_rhs = `n_parm' - `n_cat' + 1
if "`cmd'"=="mlogit" {
local n_rhs = (`n_parm'/(`n_cat'-1)) - 1
}
if "`cmd'"=="tobit" | "`cmd'"=="intreg" /// for var(e)
| "`cmd'"=="nbreg" | "`cmd'"=="cnreg" /// for alpha
| "`cmd'"=="zip" /// for 2nd intercept
{
local n_rhs = `n_rhs' - 1
}
if "`cmd'"=="zinb" {
local n_rhs = `n_rhs' - 2
} // 2nd int & alpha
local varlist : colnames(`b')
parse "`varlist'", parse (" ")
local rhsnam ""
if "`isord'"=="yes" {
local i 1
while "``i''"!= "" {
local j = substr("``i''",1,4)
* pre1.7.0 050216 if "`j'"!="_cut" {
* 050216 for Stata 9
if "`j'"!="_con" {
local rhsnam "`rhsnam' ``i''"
}
local i = `i' + 1
}
}
if "`isbin'"=="yes" | "`cmd'"=="mlogit" | "`cmd'"=="regress" /*
*/ | "`cmd'"=="tobit" | "`cmd'"=="intreg" | "`cmd'"=="cnreg" /*
*/ | "`cmd'"=="poisson" | "`cmd'"=="nbreg" {
local i 1
while "``i''"!= "" {
if "``i''"!="_cons" & "``i''"!="_se" {
local rhsnam "`rhsnam' ``i''"
}
local i = `i' + 1
}
}
*-> compute ll0 for zip and zinb (ll0 is for model with two intercepts)
if "`iszi'"=="yes" {
tempname res zill0
tempvar insamp
sca `df' = `n_parm' - 2
if "`cmd'"=="zinb" {
* remove alpha parameter
sca `df' = `n_parm' - 3
}
gen `insamp' = e(sample)
local inftyp "`e(inflate)'"
if "`inftyp'"=="logit" local inftyp ""
local doit "`e(cmd)' `e(depvar)' if `insamp' `wtis',inf(_cons) `inftyp'"
_estimates hold `res'
quietly `doit'
sca `llR' = e(ll) /* sic: not ll_0 */
_estimates unhold `res'
sca `x2' = -2*(`llR'-`llF')
} // iszi model
if "`bic'" == "" {
di _newline in g "Measures of Fit for " in y "`cmd'" in g " of " ///
in y "`depv'" _newline
}
if "`isbin'" == "yes" | "`isord'" == "yes" | "`cmd'" == "mlogit" ///
| "`cmd'" == "mprobit" {
* sum does not allow pweights or iweights
if "`e(wtype)'"=="pweight" | "`e(wtype)'"=="iweight" {
sca `cntr2a' = .
sca `cntr2' = .
if "`isbin'" == "yes" {
sca `r2ef' = .
}
if "`isbin'"=="yes" {
di _n in blu /*
*/ "(Efron's R2, Count R2, and Adj Count R2 not calculated if " /*
*/ in whi "`e(wtype)'" in blu " used)"
}
else {
di _n in blu /*
*/ "(Count R2 and Adj Count R2 not calculated if " /*
*/ in whi "`e(wtype)'" in blu " used)"
}
}
* can compute measures if aweight, fweight or no weight
else quietly {
egen `depvgrp' = group(`depv') if e(sample)==1
sca `maxrow' = 0
gen `hiprob' = 0 if e(sample)==1
gen `predcat' = 0 if e(sample)==1
local counter = 1
* BRM - Efron R2 and Count R2s
if "`isbin'" == "yes" {
* Efron R2
* convert 0/!0 to 0/1 for binary outcomes
gen `depv2' = (`depv'!=0) if e(sample)==1
predict `phat' if e(sample)==1
sum `depv2' `wtis' if e(sample)==1
local mn = r(mean)
gen `yphat' = (`depv2' - `phat')^2 if e(sample)==1
gen `ydev' = (`depv2' - `mn')^2 if e(sample)==1
sum `yphat' `wtis' if e(sample)==1
local num = r(mean)
sum `ydev' `wtis' if e(sample)==1
local den = r(mean)
scalar `r2ef' = 1 - (`num'/`den')
* Count R2 and Adjusted Count R2
replace `predcat' = round(`phat', 1) if e(sample)==1
gen `wrong' = (`predcat'~=`depv2') if e(sample)==1
sum `depv2' `wtis' if `depv2' == 1 & e(sample)==1
sca `totone' = r(N)
if `totone' > (`n_obs'-`totone') {
sca `maxrow' = `totone'
}
else {
sca `maxrow' = `n_obs' - `totone'
}
}
* ORM - Count R2s
if "`isord'" == "yes" {
forvalues counter = 1(1)`n_cat' {
* ctnam is argument for ologitp/oprobitp
tempvar ct`counter'
local ctnam "`ctnam' `ct`counter''"
* find frequency of the modal category (for adj count R2)
sum `depvgrp' if (`depvgrp'==`counter') & (e(sample)==1) `wtis'
if r(N) > `maxrow' {
sca `maxrow' = r(N)
}
}
predict `ctnam' if e(sample) == 1
forvalues counter = 1(1)`n_cat' {
replace `predcat' = `counter' ///
if (`ct`counter''>`hiprob') & (e(sample)==1)
replace `hiprob' = `ct`counter'' ///
if (`ct`counter''>`hiprob') & (e(sample)==1)
}
gen `wrong' = (`predcat'~=`depvgrp') if e(sample)==1
}
* MNLM & and MPROBIT - Count R2s
if "`cmd'" == "mlogit" ///
| "`cmd'" == "mprobit" {
while `counter' <= `n_cat' {
tempvar ct`counter'
sum `depv' `wtis' ///
if (`depvgrp'==`counter') & (e(sample)==1), meanonly
local outcm = r(mean)
predict `ct`counter'', outcome(`outcm')
sum `depv' if (`depv'==`outcm') & (e(sample)==1) `wtis'
if r(N) > `maxrow' {
sca `maxrow' = r(N)
}
replace `predcat' = `outcm' ///
if (`ct`counter''>`hiprob') & (e(sample)==1)
replace `hiprob' = `ct`counter'' ///
if (`ct`counter''>`hiprob') & (e(sample)==1)
local counter = `counter' + 1
}
gen `wrong' = (`predcat'~=`depv') if e(sample)==1
}
* compute count R2s
sum `wrong' if `wrong' == 0 & e(sample) == 1 `wtis'
sca `totcor' = r(N)
sca `cntr2' = `totcor'/`n_obs'
sca `cntr2a' = (`totcor'-`maxrow')/(`n_obs'-`maxrow')
}
} // is binary, ordinal or nominal
* CLOGIT - Count R2s
if "`cmd'"=="clogit" {
quietly {
tempname cntr2 clphat clhi clpick clpred clwrong clright /*
*/ posval motest1 motest2 mptest1
* generate predicted outcomes for each group
predict `clphat' if e(sample) == 1
egen `clhi' = max(`clphat') if e(sample) == 1, by(`e(group)')
gen `clpred' = 1 if `clhi'==`clphat' & e(sample)==1
replace `clpred' = 0 if `clhi'!=`clphat' & e(sample)==1
* check if multiple predicted 1's per group
egen `mptest1' = sum(`clpred') if e(sample), by(`e(group)')
* weight of multiple predicted 1's per group split among each 1
replace `clpred' = `clpred'/`mptest1' if `clpred'==1 & e(sample)==1
* calculate count R2
sum `clpred' `wtis' if `clpred'>0 & `clpred'~=. & `depv'==0 & e(sample)==1
sca `clwrong' = r(sum)
sum `clpred' `wtis' if `clpred'>0 & `clpred'~=. & `depv'!=0 & e(sample)==1
sca `clright' = r(sum)
sca `cntr2' = `clright' / (`clwrong'+`clright')
* warn user if multiple outcomes per group
gen `posval'=1 if `depv'!=0 & `depv'~=. & e(sample)==1
egen `motest1' = sum(`posval') if `posval'==1, by(`e(group)')
sum `motest1'
sca `motest2' = r(max)
if `motest2' > 1 {
noi di in blu "(Multiple outcomes per group for clogit: Count R2 may not be valid)"
}
} // quietly
} // clogit
*-> compute Var(y*) using cov_x
if "`isystar'"=="yes" {
sca `factor' = 1/(`n_obs'-1)
quietly mat accum `v_x' = `depv' `rhsnam' `wtis' if e(sample), ///
deviations noconstant
sca `v_y' = `factor' * `v_x'[1,1] // sum y/(n-1)
mat `v_x' = `factor' * `v_x'[2...,2...] // cov matrix of rhs vars
if "`e(cmd)'"=="intreg" {
mat `v_x' = `v_x'[2...,2...]
}
mat `bnocon' = `b'[1,1..`n_rhs'] // trim off _con
mat `v' = `bnocon' * `v_x' * `bnocon''
if "`ispbt'"=="yes" {
mat `v_e'[1,1] = 1
}
if "`islgt'"=="yes" {
mat `v_e'[1,1] = (_pi*_pi)/3
}
if "`e(cmd)'"=="intreg" {
mat `v_e'[1,1] = e(sigma)*e(sigma)
}
if "`e(cmd)'"=="tobit" | "`e(cmd)'"=="cnreg" {
mat `v_e'[1,1] = `b'[1,"sigma:_cons"]*`b'[1,"sigma:_cons"]
}
mat `v' = `v' + `v_e'
sca `v_ystar' = `v'[1,1]
}
*-> information measures
sca `aicn' = ((-2*`llF') + (2*`n_parm'))
sca `aic' = `aicn'/`n_obs'
sca `dev' = -2*`llF'
local devdf = `n_obs'-`n_parm'
sca `lrx2' = 2*(`llF'-`llR')
local lrx2df = e(df_m)
if "`iszi'"=="yes" {
local lrx2df = `df'
}
if "`cmd'"=="mprobit" {
scalar `lrx2' = e(chi2)
}
sca `lrx2p' = chiprob(`lrx2df',`lrx2')
if "`cmd'"=="slogit" {
scalar `lrx2' = e(chi2)
sca `lrx2p' = e(p)
}
sca `bicx' = `dev' - (`devdf'*ln(`n_obs'))
sca `bicp' = -`lrx2' + (e(df_m)*ln(`n_obs'))
capture return scalar dev_df = `devdf'
capture return scalar lrx2_df = `lrx2df'
if "`iszi'"=="yes" {
* df = # rhs in both equations
sca `bicp' = -`lrx2' + (`df'*ln(e(N)))
}
* these equal usual R2 in regress
if "`pseudo'"=="yes" {
sca `r2mcf' = 1 - (`llF'/`llR')
sca `r2mcfa' = 1 - ((`llF'-`n_parm')/`llR')
sca `n_obs2' = 2/`n_obs'
sca `r2ml' = 1 - exp(2*(`llR'-`llF')/`n_obs')
sca `r2cu' = (`r2ml')/(1-exp(2*`llR'/`n_obs'))
}
if "`isystar'"=="yes" {
sca `r2mz' = (`v_ystar' - `v_e'[1,1])/`v_ystar'
sca `v_error' = `v_e'[1,1]
}
*-> output
* define output macros with format where / is parse character:
* output_text/tempname_of_scalar/abbrev_for_r()
* /df?/width/half(1, 2 or ? for either)
local stat2 "Log-Lik Intercept Only/llR/ll_0/no/11/1"
local stat3 "Log-Lik Full Model/llF/ll/no/11/2"
local stat4 "D/dev/dev/yes/11/1"
local stat6 "LR/lrx2/lrx2/yes/11/2"
local stat8 "Prob > LR/lrx2p/lrx2_p/no/9/2"
local stat9 "McFadden's R2/r2mcf/r2_mf/no/9/1"
local stat10 "McFadden's Adj R2/r2mcfa/r2_mfadj/no/9/2"
local stat11 "ML (Cox-Snell) R2/r2ml/r2_ml/no/9/1"
local stat12 "Cragg-Uhler(Nagelkerke) R2/r2cu/r2_cu/no/7/2"
local stat13 "McKelvey & Zavoina's R2/r2mz/r2_mz/no/9/?"
local stat14 "Efron's R2/r2ef/r2_ef/no/9/?"
local stat15 "Variance of y*/v_ystar/v_ystar/no/11/1"
local stat16 "Variance of error/v_error/v_error/no/11/2"
local stat17 "Count R2/cntr2/r2_ct/no/9/1"
local stat18 "Adj Count R2/cntr2a/r2_ctadj/no/9/2"
local stat19 "R2/r2/r2/no/9/1"
local stat20 "Adjusted R2/r2adj/r2_adj/no/9/2"
local stat21 "AIC/aic/aic/no/11/1"
local stat22 "AIC*n/aicn/aic_n/no/11/2"
local stat23 "BIC/bicx/bic/no/11/1"
local stat24 "BIC'/bicp/bic_p/no/11/2"
local stat25 "BIC used by Stata/statabic/statabic/no/11/1"
local stat26 "AIC used by Stata/stataaic/stataaic/no/11/2"
if "`llio'"=="no" {
local stat2 ""
local stat3 "Log-Lik Full Model/llF/ll/no/11/1"
local stat4 "D/dev/dev/yes/11/2"
}
if "`cmd'"=="slogit" local stat24 ""
if "`cmd'"=="mprobit" | "`cmd'"=="slogit" {
local stat6 "Wald X2/lrx2/lrx2/yes/11/1"
}
if "`cmd'"=="mprobit" | "`cmd'"=="slogit" {
local stat8 "Prob > X2/lrx2p/lrx2_p/no/9/2"
}
*-> check models for using and saving
if "`saving'" != "" {
if (length("`saving'")>4) {
di in red "saving() name must be < 5 characters long"
exit 198
}
}
if "`using'"!= "" {
if (length("`using'")>4) {
di in red "using() name must be < 5 characters long"
exit 198
}
qui capture mat list fs_`using'
if _rc == 111 {
di in red "Incorrect using() option: no indices saved as `using'"
exit 111
}
* if the using option has been set, prepare comparison macros
sca `fn_obs' = fs_`using'[1,1]
if `fn_obs' == -9999 {
sca `fn_obs' = .
}
* estimate command saved as row name
local fcmd : rownames(fs_`using')
local fcmd : word 1 of `fcmd'
* test if current and using are same type of model
local same "no"
if ("`fcmd'"=="logit"&"`e(cmd)'"=="logistic") | ///
("`fcmd'"=="logistic"&"`e(cmd)'"=="logit") {
local same "yes"
}
if ("`fcmd'"!="`e(cmd)'" & "`same'"=="no") & "`force'"=="" {
di in r "Current model estimated by `e(cmd)'" ///
in r ", but saved model estimated by `fcmd'."
di in r "To make the comparisons anyway, use the " ///
"force option."
exit
}
if ("`fcmd'"!="`e(cmd)'" & "`same'"=="no") {
di in y "Warning: Current model estimated by " /*
*/ "`e(cmd)', but saved model estimated by `fcmd'"
di
}
if `n_obs' != `fn_obs' & "`force'"=="force" {
di in y "Warning: N's do not match."
di
}
di in g _col(30) %9s "Current" _col(48) " Saved" ///
_col(65) "Difference"
di in g "Model:" _col(30) %9s in y "`e(cmd)'" ///
%9s _col(48) in y "`fcmd'"
di in g "N:" _col(30) %9.0f in y `n_obs' ///
%9.0f _col(48) in y `fn_obs' ///
%9.0f _col(66) in y (`n_obs' - `fn_obs')
* test if N's match
if `n_obs' != `fn_obs' & "`force'"=="" {
di _newline in r "N's do not match. To make the comparisons, use the " ///
in y "force" in r " option."
exit
}
* test of ll0 match
sca `fllR' = fs_`using'[1,2]
if `fllR'==-9999 {
sca `fllR' = .
}
if abs(`llR'-`fllR') > .0001 & (`llR'-`fllR') ~= . & "`force'"=="" {
di _newline in r /*
*/ "Log-likelihoods for intercept model differ for current and saved models."
di in r "To make the comparisons, use the " in y "force" in r " option."
exit
}
} /* if "`using'"!="" */
capture return scalar n_parm = `n_parm'
capture return scalar n_rhs = `n_rhs'
capture return scalar N = `n_obs'
*-> create matrices and macros containing results
local nstats = 24 // number of stats reported in output matrix
if "`nostata'"!="stata" local nstats = `nstats' + 2
mat `results' = J(1, `nstats'+2, .)
mat `results'[1,1] = `n_obs'
mat `results'[1,`nstats'+1] = `n_parm'
mat `results'[1,`nstats'+2] = `n_rhs'
local cnames "N" // start bulding column names for output matrix
forvalues i = 2(1)`nstats' {
* read individual stats
tokenize "`stat`i''", parse(/)
* outtxt - name of statistic
local outtxt`i' "`1'"
* scalnm - name of scalar holding value of statistic
local scalnm`i' "`3'"
* abnm - name of statistic for return value
local abnm`i' "`5'"
* df - whether statistic has df
local df`i' "`7'"
* outwid - width of output - TO DO: try to delete this
local outwid`i' "`9'"
* outhalf - need to be in first or second half
local outhalf`i' "`11'"
* cnames - macro of column names for output
local cnames "`cnames' `abnm`i''"
if "`df`i''"=="yes" {
local cnames "`cnames' `abnm`i''_df"
}
* see if scalar exists, if not, no output
capture return scalar `abnm`i'' = ``scalnm`i'''
if _rc != 0 {
if "`outtxt`i''" != "" {
mat `results'[1, `i'] = .x
if "`df`i''"=="yes" {
mat `results'[1,`i'+1] = .x
}
}
}
if _rc == 0 {
* add to results matrix
mat `results'[1, `i'] = ``scalnm`i'''
if "`df`i''"=="yes" {
mat `results'[1,`i'+1] = ``scalnm`i''df'
}
} // if _rc == 0
if "`df`i''"=="yes" {
local i = `i' + 1
}
} // forvalues i
*-> loop through to print results if using option not used
if "`using'"=="" {
* arrange numbers in rows and columns for printing
local trow = 1
local tcol = 1
forvalues i = 2(1)`nstats' {
* figure out if result should be in output
if `results'[1,`i']!=.x & "`outtxt`i''" != "" & ///
("`bic'"=="" | `i' >= 21) {
* figure out if need to advance to next row of table
if "`outhalf`i''" == "1" & `tcol' == 3 {
local ++trow
local tcol = 1
}
if "`outhalf`i''" == "2" & `tcol' == 1 {
local tcol = 3
}
local op`trow'_`tcol' = "`outtxt`i''"
if "`df`i''" == "yes" {
local temp = ``scalnm`i''df'
local op`trow'_`tcol' = "`op`trow'_`tcol''(`temp')"
}
local op`trow'_`tcol' = "`op`trow'_`tcol'':"
local ++tcol
local op`trow'_`tcol' = ``scalnm`i'''
* take care of output width
local fmt`trow'_`tcol' = "%`outwid`i''.3f"
// move to next column
local ++tcol
if `tcol' == 5 { // if needed, move to next row
local tcol = 1
local ++trow
}
}
} // forvalues i
* print rows and columns
local numrows = `trow'
if `tcol' == 1 {
local numrows = `trow' - 1
}
forvalues i = 1(1)`numrows' {
di in g "`op`i'_1'" _col(28) in y %11.3f `op`i'_2' ///
_col(42) in g "`op`i'_3'" _col(69) in y %11.3f `op`i'_4'
}
} // using option not used
*-> loop through to print results if using option has been used
if "`using'"!="" {
local trow = 1
forvalues i = 2(1)`nstats' {
* figure out if result should be in output
if `results'[1,`i']!=.x & ("`bic'"=="" | `i' >= 21) {
* calculate differences
sca `prevval' = fs_`using'[1,`i']
tempname difvals
sca `difvals' = ``scalnm`i'''-`prevval'
* for LR prob give LR prob -- computed during LR loop
if "`outtxt`i''"=="Prob > LR" {
sca `difvals' = `LRchip'
}
* print results
if "`df`i''" == "no" {
di in g "`outtxt`i''" ///
_col(28) in y %11.3f ``scalnm`i''' ///
_col(46) in y %11.3f `prevval' ///
_col(64) in y %11.3f `difvals'
}
if "`df`i''" == "yes" {
local tmp_curdf = ``scalnm`i''df'
sca `prevdf' = fs_`using'[1,(`i'+ 1)]
local tmp_prevdf = `prevdf'
local tmp_difdf = abs(`tmp_curdf' - `tmp_prevdf')
* for D remove the negative sign from the difference
if "`outtxt`i''"=="D" {
sca `difvals' = abs(`difvals')
}
* for LR statistic compute LR test prob
if "`outtxt`i''"=="LR" | "`outtxt`i''"=="Wald LR" {
tempname LRdif LRdf LRchip
sca `difvals' = abs(``scalnm`i''' - `prevval')
sca `LRdf' = abs(``scalnm`i''df'-`prevdf')
sca `LRchip' = chiprob(`LRdf', `difvals')
}
di in g "`outtxt`i''" ///
_col(28) in y %11.3f ``scalnm`i''' in g "(`tmp_curdf')" ///
_col(46) in y %11.3f `prevval' in g "(`tmp_prevdf')" ///
_col(64) in y %11.3f `difvals' in g "(`tmp_difdf')"
}
}
} // forvalues i
}
*-> BIC comparisons
if "`using'" != "" {
if abs(`llR'-`fllR') < .0001 { /* check that LL0s are the same */
sca `fbicp' = fs_`using'[1, 24]
sca `bicdif' = abs(`bicp'-`fbicp')
if `bicdif'~=. {
local sup "very strong"
if `bicdif'<= .0000000001 {
local sup "no"
}
else if `bicdif' <=2 {
local sup "weak"
}
else if `bicdif' <=6 {
local sup "positive"
}
else if `bicdif' <=10 {
local sup "strong"
}
local modpref "current"
if `bicp' > `fbicp' {
local modpref "saved"
}
if `bicdif'< .0000000001 & `bicdif'>-.0000000001 /*
*/ {
local modpref "either"
}
di _newline in g "Difference of " %8.3f in y `bicdif' in g /*
*/ " in BIC' provides " in y "`sup'" /*
*/ in g " support for " in y "`modpref'" in g " model."
}
} /* if `llR' == `fllR' */
} /* if "`using'" != "" */
if "`bic'"=="" & "`using'" != "" {
di _n in g "Note: p-value for difference in LR is only valid if models are nested."
}
local cnames "`cnames' n_parm n_rhs"
mat rownames `results' = `e(cmd)'
mat colnames `results' = `cnames'
if "`saving'" != "" {
capture mat drop fs_`saving'
mat fs_`saving' = `results'
di _n in blu "(Indices saved in matrix fs_`saving')"
}
quietly `cmd'
end
exit
* version 1.8.1 2007-09-18 returns for some commands changed in Stata 10
* version 1.8.0 19Jun2006 df correction in lrtest for zinb
* version 1.7.9 26Oct2005 fix slogit and mprobit output for LR to Wald X2
* version 1.7.8b 20Jun2005 small fix to new output routine
* version 1.7.8 20Jun2005 new output routine
* version 1.7.7 17Jun2005 streamlined/stata9 improvements to code
* version 1.7.6 15Jun2005 added rologit
* version 1.7.5 10Jun2005 clogit, tobit, cnreg, Cragg-Uhler (Nagelkerke) R2 for dif, ztp, ztnb
* version 1.7.4 10Jun2005 align Variance of y*
* version 1.7.2 13Apr2005
* version 1.7.1 28Feb2005 mprobit; fixed bic'; added count R2
* version 1.7.0 16Feb2005 slogit mprobit ztp, ztnb, format on others
* version 1.6.9 12Feb2005 add stata bic and aic
* version 1.6.8 07Feb2005 add ztp and ztnb
* version 1.6.7 22Jan2005 change LRX2