diff --git a/domin.ado b/domin.ado index b9abe04..0d08fe5 100644 --- a/domin.ado +++ b/domin.ado @@ -1,4 +1,4 @@ -*! domin version 3.5.0 8/14/2023 Joseph N. Luchman +*! domin version 3.5.1 1/12/2024 Joseph N. Luchman // version information at end of file **# Pre-program definition @@ -33,7 +33,7 @@ if replay() { **# Program definition and argument checks syntax varlist(min = 1 ts) [in] [if] [aw pw iw fw] , /// [Reg(string) Fitstat(string) Sets(string) All(varlist fv ts) /// - noCONditional noCOMplete EPSilon CONSmodel REVerse noESAMPleok mi miopt(string)] // remove direct mi support - make it into a wrapper + noCONditional noCOMplete EPSilon CONSmodel REVerse noESAMPleok mi miopt(string)] if strlen("`epsilon'") { local esampleok "esampleok" @@ -64,23 +64,21 @@ else if strlen("`epsilon'") & strlen("`e(title)'") local title "Epsilon-based `r else local title "Custom user analysis" matrix `domwgts' = r(domwgts) -ereturn post `domwgts' [`weight'`exp'], depname(`dv') obs(`=`r(N)'') esample(`touse') +ereturn post `domwgts' [`weight'`exp'], depname(`dv') obs(`=r(N)') esample(`touse') **# Ereturn scalars -ereturn scalar N = `=r(N)' - matrix `domwgts' = r(domwgts)*J(colsof(r(domwgts)), 1, 1) if strlen("`epsilon'") ereturn scalar fitstat_o = `=`domwgts'[1,1]' -else ereturn scalar fitstat_o = `=`domwgts'[1,1] + r(allfs) + r(consfs)' +else ereturn scalar fitstat_o = `=r(fullfs) + r(allfs) + r(consfs)' if `:list sizeof all' ereturn scalar fitstat_a = `=r(allfs)' if strlen("`consmodel'") ereturn scalar fitstat_c = `=r(consfs)' -if strlen("`setcnt'") { +if strlen("`setcount'") { - ereturn hidden scalar setcnt = `setcnt' + ereturn hidden scalar setcount = `setcount' - forvalues x = 1/`setcnt' { + forvalues x = 1/`setcount' { fvunab set`x': `set`x'' @@ -90,10 +88,10 @@ if strlen("`setcnt'") { } -else ereturn hidden scalar setcnt = 0 +else ereturn hidden scalar setcount = 0 **# Ereturn macros -ereturn hidden local dtitle "`title'" +ereturn hidden local disp_title "`title'" ereturn hidden local reverse "`reverse'" @@ -141,10 +139,10 @@ if !strlen("`conditional'") { ereturn matrix cdldom = `cdldom' } -matrix `ranks' = r(ranks) +matrix `ranks' = r(ranks) ereturn matrix ranking = `ranks' -matrix `stdzd' = r(stdzd) +matrix `stdzd' = r(stdzd) ereturn matrix std = `stdzd' Display @@ -174,12 +172,12 @@ tokenize `diivs' local dv = abbrev("`e(depvar)'", 10) -display _newline "{txt}General dominance statistics: `e(dtitle)'" _newline /// +display _newline "{txt}General dominance statistics: `e(disp_title)'" _newline /// "{txt}Number of obs{col 27}={res}{col 40}" %12.0f e(N) display "{txt}Overall Fit Statistic{col 27}={res}{col 36}" %16.4f e(fitstat_o) -if !missing(e(fitstat_a)) display "{txt}All Subsets Fit Stat.{col 27}={res}{col 36}" %16.4f e(fitstat_a) +if !missing(e(fitstat_a)) display "{txt}All Sub-models Fit Stat.{col 27}={res}{col 36}" %16.4f e(fitstat_a) if !missing(e(fitstat_c)) display "{txt}Constant-only Fit Stat.{col 27}={res}{col 36}" %16.4f e(fitstat_c) @@ -222,7 +220,7 @@ if `cpttest' { if e(estimate) == "dominance" & `=`cpttest'*`cdltest'' { - display _newline "{res}Strongest dominance designations" _newline + display _newline "{txt}Strongest dominance designations" _newline tempname bestdom cdl gen decision @@ -261,7 +259,7 @@ if e(estimate) == "dominance" & `=`cpttest'*`cdltest'' { local names `:colnames e(b)' - mata: display((select(vec(tokens(st_local("names"))':+((st_matrix("`bestdom'"):==1):*" completely dominates "):+tokens(st_local("names")))', /// + mata: display(("{txt}", select(vec(tokens(st_local("names"))':+((st_matrix("`bestdom'"):==1):*" completely dominates "):+tokens(st_local("names")))', /// regexm(vec(tokens(st_local("names"))':+((st_matrix("`bestdom'"):==1):*" completely dominates "):+tokens(st_local("names")))', /// "completely dominates")) , /// select(vec(tokens(st_local("names"))':+((st_matrix("`bestdom'"):==2):*" conditionally dominates "):+tokens(st_local("names")))', /// @@ -275,9 +273,9 @@ if e(estimate) == "dominance" & `=`cpttest'*`cdltest'' { } -if `=e(setcnt)' { +if `=e(setcount)' { - forvalues x = 1/`=e(setcnt)' { + forvalues x = 1/`=e(setcount)' { display "{txt}Variables in set`x': `e(set`x')'" @@ -285,7 +283,7 @@ if `=e(setcnt)' { } -if strlen("`e(all)'") display "{txt}Variables included in all subsets: `e(all)'" +if strlen("`e(all)'") display "{txt}Variables included in all sub-models: `e(all)'" end @@ -381,17 +379,17 @@ void domin_2mata( /*set processing*/ if ( strlen(sets) ) { - /*setup parsing sets - bind on parentheses*/ + /*setup parsing sets - bind on parentheses*/ parser = tokeninit(" ", "", "()") tokenset(parser, sets) /*get all sets and remove parentehses*/ iv_sets = tokengetall(parser) iv_sets = substr(iv_sets, 2, strlen(iv_sets):-2) - st_local("setcnt", strofreal( length(iv_sets) ) ) + st_local("setcount", strofreal( length(iv_sets) )) for (set=1; set<=length(iv_sets); set++) { - st_local("set"+strofreal(set), iv_sets[set] ) + st_local("set"+strofreal(set), iv_sets[set]) rc = _stata("fvunab waste: " + iv_sets[set], 1) if ( rc ) { @@ -457,7 +455,13 @@ void domin_2mata( stata("generate byte " + marks[1] + " = e(sample)", 1) - stata("count if " + marks[1], 1) + if (strmatch(weight, "fweight*") | strmatch(weight, "iweight*")) { + + stata("summarize " + marks[1] + " if " + marks[1] + + " [" + weight + "]", 1) + + } + else stata("count if " + marks[1], 1) obs = st_numscalar("r(N)") @@ -475,7 +479,13 @@ void domin_2mata( stata("markout " + marks[1] + " " + marks[2] + " " + invtokens(ivs) + " " + all, 1) - stata("count if " + marks[1], 1) + if (strmatch(weight, "fweight*") | strmatch(weight, "iweight*")) { + + stata("summarize " + marks[1] + " if " + marks[1] + + " [" + weight + "]", 1) + + } + else stata("count if " + marks[1], 1) obs = st_numscalar("r(N)") @@ -570,7 +580,7 @@ void domin_2mata( } - cons_fitstat = 0 + cons_fitstat = 0 if ( strlen(consmodel) ) { @@ -677,8 +687,9 @@ void domin_2mata( st_local("touse", marks[1]) st_local("dv", dv) - st_numscalar("r(allfs)", model_specs.get("all_fitstat")) - st_numscalar("r(consfs)", model_specs.get("cons_fitstat")) + st_numscalar("r(allfs)", all_fitstat) + st_numscalar("r(consfs)", cons_fitstat) + st_numscalar("r(fullfs)", full_fitstat) } @@ -943,8 +954,13 @@ end - updates to dominance.mata (to 0.1.0) - updated 'epsilon' method to use st_view() instead of st_data() - 'epsilon' method accepts weights - ** + // 3.5.1 - January 12, 2024 + - fixed references to 'all subsets' - should be 'all sub-models' + - fixed error in docuementation of -rforest- Example 11; needed -noesampleok- + - removed redundant e(N) setting + - fixed N computation with -fweight- and -iweight- --- + future domin ** planned ** - depreciated by -domin2-'s bmaregress interface; eventually will be defunct and enveloped by its interface */ diff --git a/domin.sthlp b/domin.sthlp index b3c0819..56d88d0 100644 --- a/domin.sthlp +++ b/domin.sthlp @@ -1,5 +1,5 @@ {smcl} -{* *! version 3.5.0 August 14, 2023 J. N. Luchman}{...} +{* *! version 3.5.1 January 12, 2024 J. N. Luchman}{...} {cmd:help domin} {title:Title} @@ -469,7 +469,7 @@ when overall model fit statistics are used that decrease with better fit (e.g., {phang} {stata end} {p_end} {phang} {stata rforest price mpg headroom weight, type(reg)} {p_end} {phang} {stata matrix list e(importance)} {p_end} -{phang} {stata domin price mpg headroom weight, reg(fitdom, reg_fd(rforest, type(reg)) postestimation(myRFr2) fitstat_fd(e(r2))) fitstat(e(fitstat))} {p_end} +{phang} {stata domin price mpg headroom weight, reg(fitdom, reg_fd(rforest, type(reg)) postestimation(myRFr2) fitstat_fd(e(r2))) fitstat(e(fitstat)) noesampleok} {p_end} {marker remark}{...} @@ -554,7 +554,7 @@ when they do not follow the traditional {it:depvar indepvars} format. As long a be expressed in some way that can be evaluated in {it:depvar indepvars} format, any analysis could be dominance analyzed. -{pstd}Any program used as a wrapper by {cmd:domin} must accept an {help if} statement in its {help syntax}. +{pstd}Any program used as a wrapper by {cmd:domin} must accept an {help if} statement, a comma, and at least one (possibly optional) option argument in its {help syntax}. It is recommended that wrapper programs parse the inputs as a {it:varlist} as well (see Example #9a). {pstd}A third wrapper program, {cmd:fitdom}, takes inspiration from the @@ -608,7 +608,7 @@ Please cite as: {title:Author} {p 4}Joseph N. Luchman{p_end} -{p 4}Principal Scientist{p_end} +{p 4}Research Fellow{p_end} {p 4}Fors Marsh{p_end} {p 4}Arlington, VA{p_end} {p 4}jluchman@forsmarsh.com{p_end} diff --git a/domme.ado b/domme.ado index 206af26..4c96150 100644 --- a/domme.ado +++ b/domme.ado @@ -1,12 +1,14 @@ -*! domme version 1.1.2 8/14/2023 Joseph N. Luchman +*! domme version 1.2.0 1/12/2024 Joseph N. Luchman +// version information at end of file -capture include dominance0.mata, adopath // include Mata-based dominance code +**# Pre-program definition +quietly include dominance.mata, adopath -program define domme, eclass // ~ history and version information at end of file ~ +program define domme, eclass - version 15.1 + version 15 - if replay() & !strlen("`0'") { //replay results - error if "by"; domme allows nothing in "anything" - allow it to replay only when there is nothing in options + if replay() & !strlen("`0'") { if ("`e(cmd)'" != "domme") error 301 @@ -17,1037 +19,796 @@ program define domme, eclass // ~ history and version informa exit 0 } - - /*define syntax*/ - syntax [anything(id="equation names" equalok)] [in] [if] /// - [aw pw iw fw], Reg(string) Fitstat(string) [Sets(string) /// - noCOMplete noCONditional /// - REVerse all(string) EXTraconstr(numlist) /// - ROPts(string) ADDConstr(numlist)] //addconstr() and extraconstr() undocumented - for use in possible extensions - - /*exit conditions*/ - capture which lmoremata.mlib //is -moremata- present? - - - if _rc { //if -domin- cannot be found, tell user to install it. - - display "{err}Module {cmd:moremata} not found. Install " /// - "{cmd:moremata} here {stata ssc install moremata}." - - exit 198 - - } - - /*general set up*/ - mata: model_specs = domme_specs() //initiate instance of domme_specs() structure - tempname ranks gendom stzd_gendom cdldom cptdom //temporary matrices for results - - local two "`anything'" //rename the input equation-to-independent variable mapping in "anything": done to generalize and simplify loop below - - local drop_exit = 0 //local macro indicating that -domme- should drop constraints made so they do not persist due to error/when something goes wrong before the end of command's successful execution - - //**process "anything" and obtain individual constraints**// - while strlen("`two'") { //process the equation-to-independent variable mapping in "anything" if something is present... - - gettoken one two: two, bind //parse the equation-to-independent variable mapping to bind all parenthetical statements together and pull out first parenthetical statement/equation - - if !regexm("`one'", "=") { //exit if there is no equal sign to make an equation - - display "{err}Equation {cmd:`one'} is missing a {cmd:=} to " /// - "distinguish equation and independent variable names." - - local drop_exit = 1 //indicate that constraints will be dropped - - continue, break //stop the -while- loop... - - } - - local one = regexr("`one'", "[/(]", "") //remove left paren from equation statement - - local one = regexr("`one'", "[/)]", "") //remove right paren from equation statement - - gettoken dv ivlist: one, parse("=") //further parse the focal equation to separate out dependent from independent variables - - if ( `: list sizeof dv' != 1 ) | regexm("`dv'", "=") { //multiple dependent variables/equations or no dependent variable where one should be... exit - - display "{err}Invalid equation name specified for {cmd:(`dv'`ivlist')}." - - local drop_exit = 1 //indicate that constraints will be dropped - - continue, break //stop the -while- loop... - - } - - local ivlist = regexr("`ivlist'", "=", "") //remove the equal sign from the independent variable list - - if ( `: list sizeof ivlist' == 0 ) { //empty independent variable list... exit as an empty list is an error - - display "{err}Empty set of independent variables specified for " /// - "equation {cmd:`dv'}." - - local drop_exit = 1 //indicate that constraints will be dropped +**# Program definition and argument checks +syntax [anything(id="equation names" equalok)] [in] [if] [aw pw iw fw], /// + Reg(string) Fitstat(string) /// + [Sets(string) All(string) /// + noCOMplete noCONditional REVerse /// + ROPts(string)] - continue, break //stop the -while- loop... +**# Wrapper process for -domme- processing in Mata +mata: domme_2mata("`reg'", "`fitstat'", "`sets'", "`all'", /// + "`conditional'", "`complete'", "`ropts'", "`reverse'", /// + "`weight'`exp'", "`in' `if'", "`anything'") - } - - local dv = trim("`dv'") //remove white spaces in dependent variable which can affect putting the string into equations - - foreach iv of local ivlist { //loop over the independent variables in the first equation... - - if substr("`dv'", 1, 1) == "~" { //method for parsing (G)SEM parameters with Stata defaults - - local dv "/" //(g)sem's DV is a forward slash - - local iv = subinstr("`iv'", "+", " ", .) //remove all plus signs required for covariance - - if `: list sizeof iv' == 3 /// //formatting if this is a covariance statement (situation where there would be 3 items in the parsed "independent variable") - local iv "`: word 1 of `iv''(`: word 2 of `iv'',`: word 3 of `iv'')" - - else local iv "`: word 1 of `iv''(`: word 2 of `iv'')" //formatting if this is a standard parameter - - } - - capture constraint free //find an unused constraint number - capture to ensure they're not all used - - if !_rc { //if an unused constraint can be found... - - local constr `r(free)' //use the free constraint by putting into local macro - - constraint `constr' _b[`dv':`iv'] = 0 //establish this constraint as one that domme will use - - local add "`dv':`iv'" //set up the equation, independent variable/parameter estimate combination for display - - local ivs = "`ivs' `add'" //record the parameter estimate for display +**# Ereturn processing +tempname domwgts cdldom cptdom stdzd ranks - local constrs "`constrs' `constr'" //add the constraint associated with the parameter estimate to the list of constraints - - mata: model_specs.param_list.put(st_local("add"), st_local("constr")) // add the parameters associated with the name to the parameter mapping - - } - - else { //if no unused constraints remain... +if !strlen("`conditional'") matrix `cdldom' = r(cdldom) +if !strlen("`complete'") matrix `cptdom' = r(cptdom) - display "{err}{cmd:domme} cannot make any more constraints as the " /// - "{help constraint dir} is full (see {help constraint drop})." +if strlen("`e(title)'") local title "`e(title)'" +else if !strlen("`e(title)'") & strlen("`e(cmd)'") local title "`e(cmd)'" +else local title "Custom user analysis" - local drop_exit = 1 //indicate that constraints will be dropped +matrix `domwgts' = r(domwgts) +ereturn post `domwgts' [`weight'`exp'] , obs(`=r(N)') esample(`touse') - continue, break //stop the -while- loop... +**# Ereturn scalars +ereturn scalar fitstat_o = `=r(fullfs) + r(allfs) + r(consfs)' - } +if `=r(allfs)' != 0 ereturn scalar fitstat_a = r(allfs) +if `=r(consfs)' != 0 ereturn scalar fitstat_c = r(consfs) - if `drop_exit' continue, break //stop the within equation -forvalues- loop... - } +if strlen("`setcount'") { - if `drop_exit' continue, break //stop the overall -while- loop... + ereturn hidden scalar setcount = `setcount' + + forvalues set = 1/`setcount' { + + ereturn local set`set' = trim("`set`set''") } - //**process and obtain sets**// - local two "`sets'" //rename the set-based equation-to-independent variable mapping: done to generalize and simplify loop below +} - local setcount = 0 //used if there are sets to start a count and number them +else ereturn hidden scalar setcount = 0 - while strlen("`two'") & !`drop_exit' { //if nothing's wrong so far and there are sets... +**# Ereturn macros +ereturn hidden local disp_title "`title'" - gettoken one two: two, bind //parse the sets of equation-to-independent variable mappings in parentheses +ereturn hidden local reverse "`reverse'" - local 2 = trim("`one'") //rename the set to generalize and simplify loop below +if strlen("allset") ereturn local all = "`allset'" - if ( substr("`2'", 1, 1) != "[" ) | /// - ( substr("`2'", strlen("`2'"), strlen("`2'")) != "]" ) { //if a left paren does not begin and a right paren does not end the first set... +if strlen("`ropts'") ereturn local ropts `"`ropts'"' - display "{err}Set {cmd:`2'} not bound by brackets, {cmd:[ ]}." +ereturn local reg "`reg'" - local drop_exit = 1 //indicate that constraints will be dropped +ereturn local fitstat "`fitstat'" - continue, break //stop the -while- loop... +ereturn local cmd "domme" - } - - local setcount = `setcount' + 1 //increment setcount in the case that the first (and subsequent) set is valid - - while strlen("`2'") & !`drop_exit' { //mirrors the processing of the equation-to-independent variable mapping for "anything" if they exist within a set... - - local 2 = subinstr("`2'", "[", "", 1) //remove left bracket +ereturn local title `"Dominance analysis for multiple equations"' - local 2 = trim(subinstr("`2'", "]", "", 1)) //remove right bracket +ereturn local cmdline `"domme `0'"' - gettoken 1 2: 2, bind //parse the equation-to-independent variable mapping to bind all brackets together and pull out first equation +**# Ereturn matrices +if !strlen("`conditional'") ereturn matrix cdldom `cdldom' - if ( substr("`1'", 1, 1) != "(" ) | /// - ( substr("`1'", strlen("`1'"), strlen("`1'")) != ")" ) { //if a left bracket does not begin and a right bracket does not end the first equation... +if !strlen("`complete'") ereturn matrix cptdom `cptdom' - display "{err}Equation {cmd:`1'} in set {cmd:`one'} not bound by " /// - "parentheses; {cmd:( )}." +matrix `ranks' = r(ranks) +ereturn matrix ranking `ranks' - local drop_exit = 1 //indicate that constraints will be dropped +matrix `stdzd' = r(stdzd) +ereturn matrix std `stdzd' - continue, break //stop the -while- loop... +Display - } - - if !regexm("`1'", "=") { //exit if there is no equal sign to make an equation +end - display "{err}Equation {cmd:`1'} in set {cmd:`one'} is missing a " /// - "{cmd:=} to distinguish equation and independent variable names." - local drop_exit = 1 //indicate that constraints will be dropped +/*Display program*/ +program define Display - continue, break //stop the -while- loop... +version 15 - } +/*set up*/ +tempname gendom stzd_gendom ranks - local 1 = subinstr("`1'", "(", "", 1) //remove left paren +matrix `gendom' = e(b) +matrix `stzd_gendom' = e(std) +matrix `ranks' = e(ranking) - local 1 = trim(subinstr("`1'", ")", "", 1)) //remove right paren +local diivs: colnames e(b) +local eqivs: coleq e(b) - gettoken dv ivlist: 1, parse("=") //further parse the focal equation to separate out dependent from independent variables +mata: st_local("cdltest", strofreal(cols(st_matrix("e(cdldom)")))) +mata: st_local("cpttest", strofreal(cols(st_matrix("e(cptdom)")))) - if ( `: list sizeof dv' != 1 ) | regexm("`dv'", "=") { //multiple dependent variables or no dependent variable where one should be... exit +tokenize `diivs' - display "{err}Invalid equation name specified for " /// - "{cmd:(`dv'`ivlist')} in set {cmd:`one'}." +display _newline "{txt}General dominance statistics: `e(disp_title)'" /// + _newline "{txt}Number of obs{col 27}={res}{col 40}" %12.0f e(N) - local drop_exit = 1 //indicate that constraints will be dropped +display "{txt}Overall Fit Statistic{col 27}={res}{col 36}" /// + %16.4f e(fitstat_o) - continue, break //stop the -while- loop... +if !missing( e(fitstat_a) ) display "{txt}All Sub-models Fit Stat." /// + "{col 27}={res}{col 36}" %16.4f e(fitstat_a) - } +if !missing( e(fitstat_c) ) display "{txt}Constant-only Fit Stat." /// + "{col 27}={res}{col 36}" %16.4f e(fitstat_c) - local ivlist = regexr("`ivlist'", "=", "") //remove the equal sign from the independent variable list +display _newline "{txt}{col 13}{c |}{col 20}Dominance" /// + "{col 35}Standardized{col 53}Ranking" - if ( `: list sizeof ivlist' == 0 ) { //empty independent vatiable list... exit as it breaks process +display "{txt}{col 13}{c |}{col 20}Stat.{col 35}Domin. Stat." - display "{err}Empty set of independent variables specified for " /// - "equation {cmd:`dv'} in set {cmd:`one'}." +display "{txt}{hline 12}{c +}{hline 72}" - local drop_exit = 1 //indicate that constraints will be dropped +local current_eq "" - continue, break //stop the -while- loop... +forvalues iv = 1/`:list sizeof diivs' { - } + if "`current_eq'" != abbrev("`: word `iv' of `eqivs''", 11) /// + display `"{txt}`=abbrev("`: word `iv' of `eqivs''", 11)'{txt}{col 13}{c |}"' - local dv = trim("`dv'") //remove white spaces in dependent variable which can affect putting the string into equations + local current_eq = abbrev("`: word `iv' of `eqivs''", 11) - foreach iv of local ivlist { //loop over the independent variables in the first equation... + local `iv' = abbrev("``iv''", 10) - if substr("`dv'", 1, 1) == "~" { //method for parsing (G)SEM parameters with Stata defaults + display "{txt}{col 2}{lalign 11:``iv''}{c |}{col 14}{res}" /// + %15.4f `gendom'[1,`iv'] "{col 29}" %12.4f /// + `stzd_gendom'[1,`iv'] "{col 53}" %-2.0f `ranks'[1,`iv'] - local dv "/" //(g)sem's DV is a forward slash +} - local iv = subinstr("`iv'", "+", " ", .) //remove all plus signs required for covariance +display "{txt}{hline 12}{c BT}{hline 72}" - if `: list sizeof iv' == 3 /// //formatting if this is a covariance statement (situation where there would be 3 items in the parsed "independent variable") - local iv "`: word 1 of `iv''(`: word 2 of `iv'',`: word 3 of `iv'')" +if `cdltest' { - else local iv "`: word 1 of `iv''(`: word 2 of `iv'')" //formatting if this is a standard parameter + display "{txt}Conditional dominance statistics" _newline "{hline 85}" - } + matrix list e(cdldom), noheader format(%12.4f) - capture constraint free //find an unused constraint number + display "{txt}{hline 85}" - if !_rc { //if an unused constraint can be found... +} - local constr `r(free)' //use the free constraint +if `cpttest' { - constraint `constr' _b[`dv':`iv'] = 0 //establish this constraint as one that -domme- will use + display "{txt}Complete dominance designation" _newline "{hline 85}" - local add "`dv':`iv'" //set up the parameter estimate label + matrix list e(cptdom), noheader - local set`setcount' "`set`setcount'' `add'" //add parameter estimate label to set + display "{txt}{hline 85}" - local cset`setcount' "`cset`setcount'' `constr'" //add constraint to set +} - } +if `=`cpttest'*`cdltest'' { - else { //if no unused constraints remain... + display _newline "{txt}Strongest dominance designations" - display "{err}{cmd:domme} cannot make any more constraints as the " /// - "{help constraint dir} is full (see {help constraint drop})." + tempname bestdom cdl gen decision + if strlen("`e(reverse)'") mata: st_matrix("`bestdom'", /// + st_matrix("e(cptdom)"):*-1) - local drop_exit = 1 //indicate that constraints will be dropped + else matrix `bestdom' = e(cptdom) - continue, break //stop the -while- loop... + forvalues dominator = 1/`=colsof(e(cdldom))-1' { - } + forvalues dominatee = `=`dominator'+1'/`=colsof(e(cdldom))' { - if `drop_exit' continue, break //stop the within equation, -forvalues- loop ... + scalar `cdl' = 0 - } + scalar `gen' = 0 - if `drop_exit' continue, break //stop the -while- loop looking for equations in a set... + mata: st_numscalar("`cdl'", /// + ( sum( st_matrix("e(cdldom)")[`dominator', .] /// + :>st_matrix("e(cdldom)")[`dominatee', .] ) ) /// + :==rows( st_matrix("e(cdldom)") ) ) - } + if !`cdl' mata: /// + st_numscalar("`cdl'", -1*((sum(st_matrix("e(cdldom)")[`dominator', .] /// + :st_matrix("e(b)")[1, `dominatee']) - else { //nothing went wrong, so entire set was processed, thus... + if !`gen' mata: st_numscalar("`gen'", /// + (st_matrix("e(b)")[1, `dominator']" //add constraint set to all other constraints - the chevrons bind sets together; Mata separates them + scalar `decision' = /// + cond(abs(`bestdom'[`dominator', `dominatee']) == 1, /// + `bestdom'[`dominator', `dominatee'], cond(abs(`cdl') == 1, /// + `cdl'*2, cond(abs(`gen') == 1, `gen'*3, 0))) - mata: model_specs.param_list.put("_set:set" + st_local("setcount"), /// - st_local("cset" + st_local("setcount"))) // add the parameters associated with the name to the parameter mapping*/ - - } - - } - - if ( `:list sizeof ivs' < 2 ) & !`drop_exit' { //exit if too few parameter estimates/sets (otherwise prodices cryptic Mata error) - - display "{err}{cmd:domme} requires at least 2 parameter estimates " /// - "or parameter estimate sets." - - local drop_exit = 1 //indicate that constraints will be dropped - - } - - if strlen("`addconstr'") { //if there are custom constraints added to the analysis by number... - - foreach constrnum of numlist `addconstr' { //foreach constraint... - - local ivs "`ivs' _cns:cns`constrnum'" // format the constraint as a domme-able parameter + matrix `bestdom'[`dominator', `dominatee'] = /// + `decision'*`reverse_adj' - local constrs "`constrs' `constrnum'" // add the constraint to the domme-able constraint list + matrix `bestdom'[`dominatee', `dominator'] = /// + -`decision'*`reverse_adj' } } - local nobrkt_constrs = /// - subinstr(subinstr("`constrs'", ">", "", .), "<", "", .) //local macro where chevrons in constraint list is removed to estimate the constant model - - /*finalize setup*/ - //**markouts and sample**// - tempvar touse keep //declare sample marking variables - - tempname obs allfs consfs //declare temporary scalars - - mark `touse' //declare marking variable using Stata norm "touse" - - quietly generate byte `keep' = 1 `if' `in' //generate tempvar that adjusts for "if" and "in" statements + local names `:colfullnames e(b)' - markout `touse' `keep' //mark sample for only really works with the if & in; all other missing-based adjustments must derive from running full model + mata: display(("{txt}", select(vec(tokens(st_local("names"))' /// + :+((st_matrix("`bestdom'"):==1):*" completely dominates ") /// + :+tokens(st_local("names")))', /// + regexm(vec(tokens(st_local("names"))' /// + :+((st_matrix("`bestdom'"):==1):*" completely dominates ") /// + :+tokens(st_local("names")))', /// + "completely dominates")) , /// + select(vec(tokens(st_local("names"))':+((st_matrix("`bestdom'"):==2) /// + :*" conditionally dominates "):+tokens(st_local("names")))', /// + regexm(vec(tokens(st_local("names"))':+((st_matrix("`bestdom'"):==2) /// + :*" conditionally dominates "):+tokens(st_local("names")))', /// + "conditionally dominates")), /// + select(vec(tokens(st_local("names"))':+((st_matrix("`bestdom'"):==3) /// + :*" generally dominates "):+tokens(st_local("names")))', /// + regexm(vec(tokens(st_local("names"))':+((st_matrix("`bestdom'"):==3) /// + :*" {txt}generally dominates{res} "):+tokens(st_local("names")))', /// + "generally dominates")))') - //**inital run to ensure the syntax works and restrict sample based on full model**// - if !`drop_exit' capture `reg' [`weight'`exp'] if `touse', `ropts' constraints(`extraconstr') //run of overall analysis assuming all is well so far; intended to check e(sample) and whether everything works as it should + display "" - if _rc & !`drop_exit' { //exit if model is not estimable or program results in error - return the returned code +} - display "{err}{cmd:`reg'} resulted in an error." +if `=e(setcount)' { - local exit_code = `=_rc' + forvalues set = 1/`=e(setcount)' { - local drop_exit = 1 //indicate that constraints will be dropped + display "{txt}Parameters in set`set': `e(set`set')'" } - //**fitstat processing**// - if regexm("`fitstat'", "^e\(\)") & !`drop_exit' { //if internal fitstat indication is in "fitstat" option; must have "e()" - - gettoken omit fit_opts: fitstat, parse(",") //parse the fitstat to see if "mcf", "est", "aic", or "bic" present - - compute_fitstat `fit_opts' //return fitstat using internal computation program; this is just to determine if the options given work - - if `r(many_few)' { //something's wrong - probably too many or too few fitstat options provided - - local drop_exit = 1 //indicate that constraints will be dropped - - display "{err}Options offered to e() fitstat computation " /// - "resulted in error." - - } - - else { //...otherwise nothing's wrong - indicate to program to use built in approach +} - local fitstat = "r(fitstat)" //placeholder "returned" statistic - changed prior to program completion - - mata: st_numscalar("r(fitstat)", 0) //placeholder "returned" statistic's value +if strlen("`e(all)'") /// + display "{txt}Parameter estimates included in all sub-models: `e(all)'" + +end - local built_in = 1 //indicate that the program is to produce it's own fitstat +**# Mata function adapting Stata input for Mata and initiating the Mata environment +version 15 - local built_in_style "`r(style)'" //indicate the style of built in statistic to ask for +mata: +mata set matastrict on + +void domme_2mata( + string scalar reg, + string scalar fitstat, + string scalar sets, + string scalar all, + string scalar conditional, + string scalar complete, + string scalar ropts, + string scalar reverse, + string scalar weight, + string scalar inif, + string scalar anything + ) { + +/* ~ declarations and initial structure */ + class AssociativeArray scalar model_specs + + real scalar param, index_count, set, set_param, + full_fitstat, all_fitstat, cons_fitstat, + obs, rc + + string scalar dv, marks, bi_type + + string rowvector params, dv_iv_pairs, dv_iv_parsed, ivs, + all_pairs, sets2, set_list, set_params, + dv_iv_cns, set_cns, set_prm_cns, all_cns, + parm_display, eq_display + + transmorphic parser + +/* ~ argument checks ~ */ + if ( strmatch(fitstat, "e()*") ) { + if ( !((fitstat == "e(), mcf") | (fitstat == "e()")) ) { + display(("{err}{opt fitstat()} incorrectly specified." \ + "Must be set up as 'e()' or 'e(), mcf'." \ + "Note that types 'est', 'aic', and 'bic' are disallowed " \ + "from {cmd:domme} version 1.2.")) + exit(198) } - + else bi_type = "mcf" } - - else local built_in = 0 //if not built in, then indicate the estimation model produces it's own fitstat + else bi_type = "" - capture assert `fitstat' != . //is the "fitstat" the user supplied actually returned by the command? - - if _rc & !`drop_exit' { //exit if fitstat can't be found - - display "{err}{cmd:`fitstat'} not returned by {cmd:`reg'} or " /// - "{cmd:`fitstat'} is not scalar valued. See {help return list}." /// - _newline "Alternatively, {err}{cmd:`fitstat'} may have problems " /// - "being computed internally with the information supplied." - - - local drop_exit = 1 //indicate that constraints will be dropped - - } - - capture assert sign(`fitstat') != -1 //what is the sign of the fitstat? -domme- works best with positive ones - warn and proceed - - if _rc & !`drop_exit' { //merely note that -domme- works best with positive ones (that's what's expected) - - display "{err}{cmd:`fitstat'} returned by {cmd:`reg'}" _newline /// - "is negative. {cmd:domme} is programmed to work best" _newline /// - "with positive {opt fitstat()} summary statistics." _newline - - } - - //**count observations**// - if !`drop_exit' { //if nothing's wrong so far... - - quietly replace `touse' = e(sample) //replace sample marking variable with e(sample) from the model run - - if !inlist("`e(wtype)'", "iweight", "fweight") { //if weights don't affect obs (for probability and analytic weights) - - quietly count if `touse' //tally up "touse" - - scalar `obs' = r(N) //pull out the number of observations included - +/* ~ process parm ests and sets ~ */ + /*parse individual parm ests; if empty, return empty vector*/ + if ( strlen(anything) ) { + parser = tokeninit(" ", "", "()") + tokenset(parser, anything) + + /*get all parms; remove parentehses*/ + params = tokengetall(parser) + params = substr(params, 2, strlen(params):-2) + + /*create dummy dv-iv pair vector; create index*/ + dv_iv_pairs = + J(1, sum(ustrwordcount(params)) - 2*length(params), "") + if ( !length(dv_iv_pairs) ) dv_iv_pairs = J(0, 0, "") + index_count = 0 + + /*fill in dv-iv pair vector*/ + for (param=1; param<=length(params); param++) { + + /*parse dv vs ivs on '='*/ + dv_iv_parsed = tokens( params[param], "=") + + /*exit conditions*/ + if ( length(dv_iv_parsed) == 1 ) { + display(("{err}Parameter specification in position " + + strofreal(param) + " missing '='." \ + "DV-IV pairs cannot be parsed.")) + exit(111) + } + if ( (selectindex(dv_iv_parsed:=="=") != 2) | (ustrwordcount(dv_iv_parsed[1]) != 1) ) { + display(("{err}Parameter specification in position " + strofreal(param) + + " has un-parsable DV." \ + "DV-IV pairs cannot be parsed.")) + exit(111) + } + + /*flag ivs and dvs*/ + ivs = tokens(dv_iv_parsed[3..length(dv_iv_parsed)]) + dv = strtrim(dv_iv_parsed[1]) + + /*construct pairs; update index*/ + dv_iv_pairs[index_count+1..index_count+length(ivs)] = + (dv + ":"):+ivs + index_count = index_count + length(ivs) } - - else if inlist("`e(wtype)'", "iweight", "fweight") { //if the weights do affect obs (for frequency and importance weights) - - quietly summarize `=regexr("`e(wexp)'", "=", "")' if `touse' //tally up "touse" by summing weights - - scalar `obs' = r(sum) //pull out the number of observations included - + } + + /*setup all parsing; largely repeats dv-iv processing*/ + if ( strlen(all) ) { + parser = tokeninit(" ", "", "()") + tokenset(parser, all) + + params = tokengetall(parser) + params = substr(params, 2, strlen(params):-2) + + all_pairs = + J(1, sum(ustrwordcount(params)) - 2*length(params), "") + index_count = 0 + for (param=1; param<=length(params); param++) { + dv_iv_parsed = tokens( params[param], "=") + if ( length(dv_iv_parsed) == 1 ) { + display(("{err}All parameter specification in position " + + strofreal(param) + " missing '='." \ + "DV-IV pairs cannot be parsed.")) + exit(111) + } + if ( (selectindex(dv_iv_parsed:=="=") != 2) | (ustrwordcount(dv_iv_parsed[1]) != 1) ) { + display(("{err}All parameter specification in position " + strofreal(param) + + " has un-parsable DV." \ + "DV-IV pairs cannot be parsed.")) + exit(111) + } + ivs = tokens(dv_iv_parsed[3..length(dv_iv_parsed)]) + dv = strtrim(dv_iv_parsed[1]) + all_pairs[index_count+1..index_count+length(ivs)] = + (dv + ":"):+ ivs + index_count = index_count + length(ivs) } - - //**obtain parameters considered to be "part of model" and not adjusted out of fitstat**// - scalar `allfs' = 0 //defining fitstat of "all subsets" parameters as 0 - needed for dominance() function - - if `:list sizeof all' { //if there is something in the "all" option + } + else all_pairs = J(0, 0, "") + + /*record parms implied by -all()- in Stata*/ + st_local("allset", (length(all_pairs) ? invtokens(all_pairs) : "" )) + + /*setup sets parsing; repeats dv-iv processing with extra parsing layer*/ + if ( strlen(sets) ) { + + /*initial parsing layer; binds sets first on brackets*/ + parser = tokeninit(" ", "", "[]") + tokenset(parser, sets) + + /*parse out all sets*/ + sets2 = tokengetall(parser) + sets2 = substr(sets2, 2, strlen(sets2):-2) + + /*list of sets to update with lists of parm ests*/ + set_list = + J(1, length(sets2), "") + + for (set=1; set<=length(set_list); set++) { - local 2 "`all'" //rename the content of "all" to generalize and simplify loop below - - while strlen("`2'") & !`drop_exit' { //process the equation-to-independent variable mapping for "all"... - - gettoken 1 2: 2, bind //parse the equation-to-independent variable mapping to bind all parentheses together and pull out first equation - - if ( substr("`1'", 1, 1) != "(" ) /// - | ( substr("`1'", strlen("`1'"), strlen("`1'")) != ")" ) { //if a left paren does not begin and a right paren does not end the first all equation... - - display "{err}Equation {cmd:`1'} in {cmd:all()} not bound by " /// - "parentheses." - - local drop_exit = 1 //indicate that constraints will be dropped - - continue, break //stop the -while- loop... - - } - - if !regexm("`1'", "=") { //exit if there is no equal sign to make an equation - - display "{err}Equation {cmd:`1'} in {cmd:all()} is missing a " /// - "{cmd:=} to distinguish equation and independent variable names." - - local drop_exit = 1 //indicate that constraints will be dropped - - continue, break //stop the -while- loop... - - } - - local 1 = subinstr("`1'", "(", "", 1) //remove left paren - - local 1 = subinstr("`1'", ")", "", 1) //remove right paren - - gettoken dv ivlist: 1, parse("=") //further parse the focal equation to separate out dependent from independent variables - - if ( `: list sizeof dv' != 1 ) | regexm("`dv'", "=") { //multiple dependent variables or no dependent variable where one should be... exit - - display "{err}Invalid equation name specified for " /// - "{cmd:(`dv'`ivlist')} in {cmd:all()}." - - local drop_exit = 1 //indicate that constraints will be dropped - - continue, break //stop the -while- loop... - + /*second parsing layer; binds parm ests on parentheses*/ + parser = tokeninit(" ", "", "()") + tokenset(parser, sets2[set]) + + set_params = tokengetall(parser) + set_params = substr(set_params, 2, strlen(set_params):-2) + + set_pairs = + J(1, sum(ustrwordcount(set_params)) - 2*length(set_params), "") + index_count = 0 + + /*combined sets of parm ests within a set; note they are combined*/ + for (set_param=1; set_param<=length(set_params); set_param++) { + dv_iv_parsed = tokens( set_params[set_param], "=") + if ( length(dv_iv_parsed) == 1 ) { + display(("{err}Set position " + strofreal(set) + + " parameter specification in position " + strofreal(set_param) + " missing '='." \ + "DV-IV pairs cannot be parsed.")) + exit(111) } - - local ivlist = regexr("`ivlist'", "=", "") //remove the equal sign from the independent variable list - - if ( `: list sizeof ivlist' == 0 ) { //empty independent vatiable list... exit as it breaks process - - display "{err}Empty set of independent variables specified for " /// - "equation {cmd:`dv'} in {cmd:all()}." - - local drop_exit = 1 //indicate that constraints will be dropped - - continue, break //stop the -while- loop... - + if ( (selectindex(dv_iv_parsed:=="=") != 2) | (ustrwordcount(dv_iv_parsed[1]) != 1) ) { + display(("{err}All parameter specification in position " + strofreal(param) + + " has un-parsable DV." \ + "DV-IV pairs cannot be parsed.")) + exit(111) } - - local dv = trim("`dv'") //remove white spaces in dependent variable which can affect putting the string into equations - - foreach iv of local ivlist { //loop over the independent variables in the all equations... - - if substr("`dv'", 1, 1) == "~" { //method for parsing (G)SEM parameters with Stata defaults - - local dv "/" //(g)sem's DV is a forward slash - - local iv = subinstr("`iv'", "+", " ", .) //remove all plus signs required for covariance - - if `: list sizeof iv' == 3 /// //formatting if this is a covariance statement (situation where there would be 3 items in the parsed "independent variable") - local iv "`: word 1 of `iv''(`: word 2 of `iv'',`: word 3 of `iv'')" - - else local iv "`: word 1 of `iv''(`: word 2 of `iv'')" //formatting if this is a standard parameter - - } - - capture constraint free //find an unused constraint number - - if !_rc { //if an unused constraint can be found... - - local constr `r(free)' //use the free constraint - - constraint `constr' _b[`dv':`iv'] = 0 //establish this constraint as one that -domme- will use - - local add "`dv':`iv'" //set up the parameter estimate label - - local allset "`allset' `add'" //add the parameter estimate label to set - - local allcset "`allcset' `constr'" //add constraint to set - - } - - else { //if no unused constraints remain... - - display "{err}{cmd:domme} cannot make any more constraints as the " /// - "{help constraint dir} is full (see {help constraint drop})." - - local drop_exit = 1 //indicate that constraints will be dropped - - continue, break //stop the -while- loop... - - } - - if `drop_exit' continue, break //stop the -forvalues- loop... - - } - + ivs = tokens(dv_iv_parsed[3..length(dv_iv_parsed)]) + dv = strtrim(dv_iv_parsed[1]) + set_pairs[index_count+1..index_count+length(ivs)] = + (dv + ":"):+ ivs + index_count = index_count + length(ivs) } - - } - - //**obtain "constant" model which will adjusted out of fitstat**// - scalar `consfs' = 0 //define constant-only model fitstat as 0 - needed for dominance() function - - quietly `reg' [`weight'`exp'] if `touse', `ropts' /// - constraints(`nobrkt_constrs' `allcset' `extraconstr') //all constraints used - this estimates the "constant" model - - if `built_in' { + set_list[set] = invtokens(set_pairs) + st_local("set" + strofreal(set), invtokens(set_pairs)) - compute_fitstat `fit_opts' iscons //if a built-in fistat desired, estimate it; note constant model + } + + } + else set_list = J(0, 0, "") + + /*record set count in Stata environment*/ + st_local("setcount", strofreal(length(set_list))) + +/* ~ check primary analysis and set estimation sample ~ */ + st_eclear() + + if (inif == " ") inif = "if _n > 0" + + rc = _stata(reg + " [" + weight + "] " + inif + ", " + ropts, 1) + + if ( rc ) { + + display("{err}{cmd:" + reg + "} resulted in an error.") + exit(rc) + + } + + if ( (!strlen(bi_type)) ) { + + if ( !length( st_numscalar( strtrim(fitstat) ) ) ) { - scalar `consfs' = `fitstat' + display("{err}{cmd:" + fitstat + + "} not returned by {cmd:" + reg + "} or {cmd:" + fitstat + + "} is not scalar valued. See {help return list}.") + exit(198) } - - if !`built_in' & !missing(`fitstat') /// - scalar `consfs' = `fitstat' //return constant model's fitstat if user supplied and not missing + } + + if ((bi_type == "mcf") & + (!length(st_numscalar("e(ll)")))) { + + display("{err}{cmd:e(ll} not returned by {cmd:" + reg + "}.") + exit(198) - if !`built_in' & missing(`fitstat') scalar `consfs' = 0 //otherwise assume the constant-only model is 0 - - if strlen("`all'") { //distinguishes "all subsets" from "constant" fitstats/models - - quietly `reg' [`weight'`exp'] if `touse', `ropts' /// - constraints(`nobrkt_constrs' `extraconstr') //all "subsets" and extra constraints used - this estimates the "all" model - - if `built_in' compute_fitstat `fit_opts' consmodel(`=`consfs'') //if a built-in fistat desired, estimate it; supply constant model's value - - scalar `allfs' = `fitstat' //return all subset model's fitstat as distinct from any constant model - - - } - - if `built_in' local cons_no_add = /// - !inlist("`built_in_style'", "mcf", "est") //local macro indicating that the constant model fitstat should be included in the total fitstat - will be 0 for built in "mcf" and "est" and the actual value of "consfs" for non-built in where "consmodel" is not chosen will be 0 - - else local cons_no_add = 1 //...all other cases can be allowed to be 0 - - if strlen("`all'") scalar `allfs' = `allfs' - `consfs'*`cons_no_add' //adjust all subsets fitstat for the constant fitstat in situations where it's not 0 - - /*dominance statistic determination*/ - mata: model_specs.built_in = strtoreal(st_local("built_in")) - mata: model_specs.constraint_list = tokens(st_local("nobrkt_constrs")) - mata: model_specs.built_in_cons = st_numscalar(st_local("consfs")) - - mata: dominance0(model_specs, &domme_call(), /// - st_local("conditional'"), st_local("complete"), /// - st_local("ivs"), /// - st_numscalar(st_local("allfs")), /// - st_numscalar(st_local("consfs"))*strtoreal(st_local("cons_no_add"))) //invoke "dominance()" function in Mata - - /*translate r-class results from me_dominance() into temp results*/ - matrix `gendom' = r(domwgts) //general dominance statistics - - matrix `stzd_gendom' = r(sdomwgts) //standardized general dominance statistics - - matrix `ranks' = r(ranks) //ranks based on general dominance statistics - - if !strlen("`conditional'") matrix `cdldom' = r(cdldom) //conditional dominance statistic matrix - - if !strlen("`complete'") matrix `cptdom' = r(cptdom) //complete dominance designation matrix - - /*processing display results*/ - //**name matrices**// - matrix colnames `gendom' = `ivs' //name the columns of general dominance statistic vector - - if strlen("`reverse'") { //reverse the direction and interpretation of ranked and standardized general dominance statistics - - mata: st_matrix("`stzd_gendom'", /// - ( st_matrix("`gendom'"):*-1 ):/sum( st_matrix("`gendom'"):*-1 ) ) //reverse the signs of the standardized general dominance statistics - - mata: st_matrix("`ranks'", /// - ( ( st_matrix("`ranks'"):-1 ):*-1 ):+cols( st_matrix("`ranks'") ) ) //reverse the sign of the ranked general dominance statistics - - } - - matrix colnames `stzd_gendom' = `ivs' //name the columns of stanadrdized general dominance statistic vector - - matrix colnames `ranks' = `ivs' //name the columns of ranked general dominance statistic vector - - if !strlen("`complete'") { //if the complete dominance matrix was not suppressed... - - local cptivs = subinstr("`ivs'", ":", "_", .) //must remove colons between equation and independent variable; only 1 colon allowed - "dominates?" and "dominated?" are equations for complete dominance matrix - - local cptivs = subinstr("`ivs'", ".", "_", .) //must also remove "."s - causes error in naming - - if strlen("`reverse'") mata: /// - st_matrix("`cptdom'", st_matrix("`cptdom'"):*-1 ) //reverse the sign of the complete dominance designations - - matrix colnames `cptdom' = `cptivs' //name the columns of the complete dominance designations - - mata: st_matrixcolstripe("`cptdom'", /// - ("?":+st_matrixrowstripe("`cptdom'")[,1], st_matrixrowstripe("`cptdom'")[,2])) //add name the equation for the rows ">?" replacement for "dominates?" - - } - - if !strlen("`conditional'") { //if the conditional dominance matrix was not suppressed... - - matrix rownames `cdldom' = `ivs' //name the rows of the conditional dominance matrix - - local colcdl `:colnames `cdldom'' //the columns of the conditional dominance matrix are at defaults "c1 c2 ... cN" - - local colcdl = subinstr("`colcdl'", "c", "", .) //remove the "c"s from all the rownames; keep the values - - matrix colnames `cdldom' = `colcdl' //replace the column names of the conditional dominance matrix with the number of "orders" which matches their counting sequennce - - matrix coleq `cdldom' = #param_ests //equation names for all columns are "#param_ests" - - } - - if strlen("`e(title)'") local title "`e(title)'" //if the estimation command has an "e(title)" returned, save it - - else if !strlen("`e(title)'") & strlen("e(cmd)") local title "`e(cmd)'" //...otherwise save the "e(cmd)" as that's informative too - - else local title "Custom user analysis" //...finally, if none of the options are returned, call it "custom user analysis" - - /*return values*/ - ereturn post `gendom' [`weight'`exp'], obs(`=`obs'') esample(`touse') //primary estimation command returned value command; clears ereturn and returns "gendom" as e(b) - - if strlen("`setcount'") { //if there are sets... - - ereturn hidden scalar setcount = `setcount' //hidden scalar for use in "display" - - forvalues set = 1/`setcount' { //for each set... - - ereturn local set`set' = trim("`set`set''") //make a separate local macro with it's parameter estimate label contents - - } - + /*record the value of the constraint*/ + dv_iv_cns[param] = st_macroexpand("`" + "r(free)" + "'") + + /*map the value of the constraint to a dv-iv pair*/ + stata("constraint " + dv_iv_cns[param] + " _b[" + + dv_iv_pairs[param] + "] = 0") } - - else ereturn hidden scalar setcount = 0 //...otherwise hidden set count is 0 - - ereturn hidden local disp_title "`title'" //hidden title for display (hence, "disp_title") - - ereturn hidden local reverse "`reverse'" //hidden indicator for reverse - for display - - if `:list sizeof all' ereturn local all = strtrim("`allset'") //parameter estimate labels in all subsets - - if strlen("`ropts'") ereturn local ropts `"`ropts'"' //if there were regression command options return them as macro - - ereturn local reg "`reg'" //return command used in -domme- in macro - - if `built_in' local fitstat "e(`built_in_style')" //if a built-in is used, change "r(fitstat)" to the "e()" form with fitstat name - - ereturn local fitstat "`fitstat'" //return the name of the fitstat used - - ereturn local cmd "domme" //this command is -domme-; return that - - ereturn local title `"Dominance analysis for multiple equations"' //The title/description of -domme- - - ereturn local cmdline `"domme `0'"' //full command as read into -domme- the 0 macro is everything after the command's name - - ereturn scalar fitstat_o = r(fs) //overall fitstat value - - if `:list sizeof all' ereturn scalar fitstat_a = `allfs' + /// - `consfs'*`cons_no_add' //all subsets fitstat value - - if sign( `consfs'*`cons_no_add' ) ereturn scalar fitstat_c = /// //constant model fitstat value - `consfs'*`cons_no_add' - - if !strlen("`conditional'") ereturn matrix cdldom `cdldom' //return conditional dominance matrix - - if !strlen("`complete'") ereturn matrix cptdom `cptdom' //return complete dominance designations - - ereturn matrix ranking `ranks' //return ranked general dominance vector - - ereturn matrix std `stzd_gendom' //return standardized general dominance vector - - /*begin display*/ - Display - + } - - /*drop constraints -domme- made; not ones user supplied*/ - if missing(`exit_code') local exit_code = 198 //make the default exit code 198 - - if `drop_exit' { //if there was a problem during the program's execution - - if `: list sizeof constrs' { //if there are constraints used... - - foreach constr of numlist `constrs' { //go through each constraint that was made... - - constraint drop `constr' //drop the constraint - + + /*repeats dv-iv pair process but groups constraints by set*/ + if (length(set_list)) { + + set_cns = J(1, length(set_list), "") + + for (set=1; set<=length(set_list); set++) { + set_params = tokens(set_list[set]) + set_prm_cns = J(1, length(set_params), "") + for (param=1; param<=length(set_prm_cns); param++) { + rc = _stata("constraint free") + if (rc) break + set_prm_cns[param] = st_macroexpand("`" + "r(free)" + "'") + stata("constraint " + set_prm_cns[param] + " _b[" + + set_params[param] + "] = 0") } - + set_cns[set] = invtokens(set_prm_cns) } - - exit `exit_code' //exit using applicable exit code - - } - - foreach constr of numlist `nobrkt_constrs' { //go through each constraint that was made... - - constraint drop `constr' //drop the constraint - - } - -end - - -/*Display program*/ -program define Display - - version 15.1 - - /*set up*/ - tempname gendom stzd_gendom ranks //declare names for temporary data - - matrix `gendom' = e(b) //as in original command "gendom" is general dominance statistics vector; now in e(b) - - matrix `stzd_gendom' = e(std) //as in original command "stzd_gendom" is standardized general dominance statistics vector; now in e(std) - - matrix `ranks' = e(ranking) //as in original command "ranks" is ranked general dominance statistics vector; now in e(ranking) - - local diivs: colnames e(b) //obtain independent variable names - - local eqivs: coleq e(b) //obtain dependent variable/equation names - - mata: st_local("cdltest", strofreal(cols(st_matrix("e(cdldom)")))) //indicator macro for presence of conditional dominance matrix - - mata: st_local("cpttest", strofreal(cols(st_matrix("e(cptdom)")))) //indicator macro for presence of complete dominance matrix - - tokenize "`diivs'" //tokenize list of dependent variables to associate numbers with independent variables - - /*begin displays*/ - display _newline "{txt}General dominance statistics: `e(disp_title)'" /// - _newline "{txt}Number of obs{col 27}={res}{col 40}" %12.0f e(N) - - display "{txt}Overall Fit Statistic{col 27}={res}{col 36}" /// - %16.4f e(fitstat_o) - - if !missing( e(fitstat_a) ) display "{txt}All Subsets Fit Stat." /// - "{col 27}={res}{col 36}" %16.4f e(fitstat_a) - - if !missing( e(fitstat_c) ) display "{txt}Constant-only Fit Stat." /// - "{col 27}={res}{col 36}" %16.4f e(fitstat_c) - - display _newline "{txt}{col 13}{c |}{col 20}Dominance" /// - "{col 35}Standardized{col 53}Ranking" - - display "{txt}{col 13}{c |}{col 20}Stat.{col 35}Domin. Stat." - - display "{txt}{hline 12}{c +}{hline 72}" - - local current_eq "" //for separating equation from independent variable names in -forvalues- loop below - - forvalues iv = 1/`:list sizeof diivs' { //for each entry of the e(b) vector... - - if "`current_eq'" != abbrev("`: word `iv' of `eqivs''", 11) /// - display `"{res}`=abbrev("`: word `iv' of `eqivs''", 11)'{txt}{col 13}{c |}"' //...display equation name only if it changes - - local current_eq = abbrev("`: word `iv' of `eqivs''", 11) //note current equation - truncate to 11 chars - - local `iv' = abbrev("``iv''", 10) //abbreviate independent variable to 10 chars - - display "{txt}{col 2}{lalign 11:``iv''}{c |}{col 14}{res}" /// - %15.4f `gendom'[1,`iv'] "{col 29}" %12.4f /// - `stzd_gendom'[1,`iv'] "{col 53}" %-2.0f `ranks'[1,`iv'] - + } - - display "{txt}{hline 12}{c BT}{hline 72}" - - if `cdltest' { //if conditional dominance matrix exists.... - - display "{txt}Conditional dominance statistics" _newline "{hline 85}" - - matrix list e(cdldom), noheader format(%12.4f) - - display "{txt}{hline 85}" - + + /*repeats dv-iv pair process but groups -all()- params together*/ + if (length(all_pairs)) { + + all_cns = J(1, length(all_pairs), "") + + for (param=1; param<=length(all_pairs); param++) { + rc = _stata("constraint free") + if (rc) break + all_cns[param] = st_macroexpand("`" + "r(free)" + "'") + stata("constraint " + all_cns[param] + " _b[" + + all_pairs[param] + "] = 0") + } + } - - if `cpttest' { //if complete dominance designations exist... - - display "{txt}Complete dominance designation" _newline "{hline 85}" - - matrix list e(cptdom), noheader - - display "{txt}{hline 85}" - + + if (rc) { + display("{err}{cmd:domme} cannot make any more constraints as the " + + "{help constraint dir} is full (see {help constraint drop}).") + + /*drop constraints if generated*/ + exit_constraint((dv_iv_cns, set_cns, all_cns), 198) } + +/* ~ begin collecting effect sizes ~ */ + all_fitstat = 0 - if `=`cpttest'*`cdltest'' { //if _both_ complete and conditional dominance designations exist - determine strongest dominance designation between each pair of parameter estimates - - display _newline "{res}Strongest dominance designations" _newline - - tempname bestdom cdl gen decision //declare temporary names for strongest designation search - - if strlen("`e(reverse)'") mata: st_matrix("`bestdom'", /// //start by determining complete dominance as "best" - if reversed, then reflect values over 0 - st_matrix("e(cptdom)"):*-1) - - else matrix `bestdom' = e(cptdom) //...otherwise take complete dominance values as is - - forvalues dominator = 1/`=colsof(e(cdldom))-1' { //search through all columns save last... - - forvalues dominatee = `=`dominator'+1'/`=colsof(e(cdldom))' { //...as well as all columns save first, dependent on dominator - - scalar `cdl' = 0 //define conditional dominance as 0 or "not" - - scalar `gen' = 0 //define general dominance as 0 or "not" - - mata: st_numscalar("`cdl'", /// - ( sum( st_matrix("e(cdldom)")[`dominator', .] /// //...sum the number of times the values across all columns of the conditional dominance matrix for the row corresponding to the "dominator" - :>st_matrix("e(cdldom)")[`dominatee', .] ) ) /// //...and compred to the values across the same columns of the conditional dominance matrix for the row corresponding to the "dominatee" - :==rows( st_matrix("e(cdldom)") ) ) //...if that's equal to the number of columns - thus every row for the "dominator" is bigger than every row for the "dominatee" - "dominator" conditionally dominates "dominatee" - - if !`cdl' mata: /// - st_numscalar("`cdl'", -1*((sum(st_matrix("e(cdldom)")[`dominator', .] /// - :st_matrix("e(b)")[1, `dominatee']) //compare dominator's general dominance statistic to dominatee's - - if !`gen' mata: st_numscalar("`gen'", /// - (st_matrix("e(b)")[1, `dominator'] 3 local many_few = 1 //too many fitstat options selected - - if !strlen("`ll'") local ll = "e(ll)" - - if strlen("`mcf'`est'") { - - if !strlen("`obs'") & strlen("`est'") local obs = "e(N)" - - if strlen("`iscons'") local fitstat = `ll' //baseline - - else if !strlen("`iscons'") & strlen("`est'") local /// - fitstat = 1 - (`ll'/`consmodel')^(-2*`consmodel'/`obs') //Estrella R2 - - else local fitstat = 1 - `ll'/`consmodel' //McFadden R2 - + + if ( length(all_pairs) & !strlen(bi_type) ) + all_fitstat = all_fitstat - cons_fitstat + + full_fitstat = + full_fitstat - all_fitstat - + (strlen(bi_type)? 0 : cons_fitstat) + +/* ~ invoke dominance ~ */ + model_specs.put("reg", reg) + model_specs.put("fitstat", fitstat) + model_specs.put("weight", weight) + model_specs.put("inif", inif) + model_specs.put("ropts", ropts) + model_specs.put("cns", (dv_iv_cns, set_cns)) + + model_specs.put("reverse", reverse) + + model_specs.put("all_fitstat", all_fitstat) + model_specs.put("cons_fitstat", cons_fitstat) + + model_specs.put("bi_type", bi_type) + + dominance( model_specs, &domme_call(), + (dv_iv_cns, set_cns)', conditional, complete, full_fitstat ) + +/* ~ return values ~ */ + if ( length(set_list) ) dv_iv_pairs = (dv_iv_pairs, "_set:set":+strofreal((1..length(set_list)))) + + parm_display = eq_display = J(1, length(dv_iv_pairs), "") + + for (param=1; param<=length(dv_iv_pairs); param++) { + parm_display[param] = ustrsplit(dv_iv_pairs[param], ":")[2] + eq_display[param] = ustrsplit(dv_iv_pairs[param], ":")[1] } - - else if strlen("`aic'`bic'") { - - if !strlen("`parm'") local parm = "e(rank)" - - if !strlen("`obs'") & strlen("`bic'") local obs = "e(N)" - - if strlen("`aic'") local fitstat = -2*(`ll') + 2*`parm' //AIC - - else local fitstat = -2*(`ll') + `parm'*ln(`obs') //BIC - + + st_matrixcolstripe("r(domwgts)", (eq_display \ parm_display)') + + if ( strlen(reverse) ) + st_matrix("r(cptdom)", + st_matrix("r(cptdom)"):*-1 ) + + if ( !strlen(conditional) ) { + st_matrixrowstripe("r(cdldom)", (eq_display \ parm_display)' ) + st_matrixcolstripe("r(cdldom)", + (J(1, length(dv_iv_pairs), "#param_ests") \ strofreal(1..length(dv_iv_pairs)))' ) } - - else local many_few = 1 - - return local style "`aic'`bic'`mcf'`est'" - - return scalar many_few = `many_few' - - return scalar fitstat = `fitstat' - -end - -**# Container object for domme specs -version 15.1 - -mata: - - mata set matastrict on - - struct domme_specs { - - real scalar built_in, built_in_cons - - string rowvector constraint_list + + if ( !strlen(complete) ) { + st_matrixrowstripe("r(cptdom)", (">?":+eq_display \ parm_display)' ) + st_matrixcolstripe("r(cptdom)", ("=(st_matrix("r(domwgts)")') + ) ) + st_matrixcolstripe("r(ranks)", (eq_display \ parm_display)' ) - class AssociativeArray scalar param_list + st_numscalar("r(N)", obs) + st_local("reg", reg) + st_local("touse", marks[1]) + built_in + st_local("built_in", strofreal(strlen(bi_type)>0)) + st_local("built_in_style", bi_type) + + st_numscalar("r(allfs)", all_fitstat) + st_numscalar("r(consfs)", + (strlen(bi_type) ? 0 : cons_fitstat) ) + st_numscalar("r(fullfs)", full_fitstat) - } +/* mechanism to drop constraints if generated */ + exit_constraint((dv_iv_cns, set_cns, all_cns), 0) + +} end **# Mata function to execute 'domme-flavored' models -version 15.1 +version 15 mata: mata set matastrict on - real scalar domme_call(string scalar Params_in_Model, /// - real scalar all_subsets_fitstat, real scalar constant_model_fitstat, /// - struct domme_specs scalar model_specs) - { - - real scalar fitstat, loop_length, x + real scalar domme_call(string scalar params_to_remove, + class AssociativeArray scalar model_specs) { - string scalar Constraints, remove_list, key + real scalar fitstat, rc - string rowvector Param_List, out_list + string scalar params_in_model - real rowvector which_constraints + string rowvector cns - Param_List = tokens(Params_in_Model) + /*enumerate all constraints generated*/ + cns = tokens(invtokens(model_specs.get("cns"))) - loop_length = length(Param_List) + /*constraints indicate omitted parm est; reverse to get all included*/ + params_in_model = + select(cns, + !colsum((J(length(tokens(params_to_remove)), 1, + strtoreal(cns)):==J(1, length(cns), + strtoreal(tokens(params_to_remove))')))) - remove_list = "" + rc = _stata(model_specs.get("reg") + " [" + model_specs.get("weight") + "] " + + model_specs.get("inif") + ", constraints(" + + invtokens(params_in_model) + ") " + model_specs.get("ropts"), 1) - for (x = 1; x <= loop_length; x++) { + /*drop constraints if generated*/ + if (rc) { - remove_list = remove_list + " " + /// - model_specs.param_list.get(Param_List[x]) + exit_constraint(cns, rc) } + + fitstat = + built_in_fitstat(model_specs.get("fitstat"), + model_specs.get("bi_type"), + model_specs.get("cons_fitstat"), 0) - out_list = tokens(remove_list) + fitstat = + fitstat - + model_specs.get("all_fitstat") - + (strlen(model_specs.get("bi_type")) ? + 0 : model_specs.get("cons_fitstat") ) - which_constraints = /// - !colsum( /// - J(length(out_list), 1, model_specs.constraint_list):==out_list' /// - ) + return(fitstat) + + } + +end + +**# Mata function to drop constraints and exit +version 15 + +mata: - Constraints = invtokens(select(model_specs.constraint_list, which_constraints)) + mata set matastrict on + + void exit_constraint(string rowvector cns_list, numeric scalar exit_num) { + + numeric scalar cns - stata("\`reg' [\`weight'\`exp'] if \`touse', \`ropts'" + /// - " constraints(" + Constraints + " \`extraconstr')", 1) //estimate model/regression in Stata + for (cns=1; cns<=length(cns_list); cns++) { - - if (model_specs.built_in) { + stata("constraint drop " + cns_list[cns]) - stata("compute_fitstat \`fit_opts' consmodel(" /// - + strofreal(model_specs.built_in_cons) + ")", 1) - - fitstat = st_numscalar("r(fitstat)") - } - else fitstat = st_numscalar(st_local("fitstat")) + if (exit_num) exit(exit_num) + + } + +end - fitstat = fitstat - /// - all_subsets_fitstat - /// - constant_model_fitstat +**# Mata program to obtain fit statistics after model estimation - return(fitstat) +version 15 +mata: + + mata set matastrict on + + numeric scalar built_in_fitstat( + string scalar fitstat, string scalar type, numeric scalar constant, + numeric scalar isconstant + ) { + + numeric scalar value + + /*compute McFadden R2*/ + if (type == "mcf") { + + value = (isconstant? + st_numscalar("e(ll)") : + 1 - st_numscalar("e(ll)")/constant) + + } + /*not a built-in, pass the name and assume Stata returns it*/ + else value = st_numscalar( strtrim(fitstat) ) + + return(value) + } end @@ -1070,4 +831,12 @@ end // 1.1.2 - August 14, 2023 - migrated to a sub-method in the -domin- module - call temporary 'dominance0.mata' (v 0.0) for time being until re-designed similar to -domin- to use 'dominance.mata' versions > 0.0 + --- + domme version 1.2.0 - January 12, 2024 + - most of command migrated to Mata + - now dependent on dominance.mata (not dominance0.mata) + - fixed references to 'all subsets' - should be 'all sub-models' + - minimum version sync-ed with -domin- at 15 (not 15.1) + - fixed issue with -all()- option; not consistent with documentation - required full 'all' typed (not 'a()') + - estrella, aic and bic disallowed with built-in; mcfadden remains */ diff --git a/domme.sthlp b/domme.sthlp index 6ae9421..58fea67 100644 --- a/domme.sthlp +++ b/domme.sthlp @@ -1,5 +1,5 @@ {smcl} -{* *! version 1.1.2 August 14, 2023 J. N. Luchman}{...} +{* *! version 1.2.0 January 12, 2024 J. N. Luchman}{...} {cmd:help domme} {title:Title} @@ -45,11 +45,6 @@ sets of parameter estimates (see option {opt sets()} below). Because it is possible to submit only sets of parameter estimates, the initial parameter estimates specification statement is optional. {p_end} -{p 4 6 2} -{cmd:domme} requires installation of the {cmd:moremata} package -(install {stata ssc install moremata:here}). {p_end} - - {title:Table of Contents} {space 4}{help domme##desc: 1. Description} @@ -144,8 +139,6 @@ removed from the {cmd: indepvars} list directly. The way parameter constraints are produced with the {opt all()} and {opt sets()} options is identical to that of the initial statements to {cmd: domme}. -... brief discussion of why _cons not included? ... - {marker opts}{...} {title:3. Options} @@ -181,20 +174,16 @@ command used produces the desired scalar with constraints. {pmore}The second method accommodates Stata commands' tendency to not return pseudo-R-square values with constraints and expands which commands can get a fit statistic using a built-in fit statistic computation. When {opt fitstat()} -is asked for an empty ereturned statistic indicator (i.e., {res:e()}) you must provide a -the three character code as an option to the {opt fitstat()}. Four fit statistic -options are available. These options are the McFadden pseudo-R squared ({res:mcf}), -the Estrella pseudo-R squared ({res:est}), the Akaike information criterion ({res:aic}), -and the Bayesian information criterion ({res:bic}). For example, to ask {cmd:domme} to -compute McFadden's pseudo-R square as a fit statistic, type {res:fitstat(e(), mcf)} -(See Example #1). +is asked for an empty ereturned statistic indicator (i.e., {res:e()}). Using +{opt fitstat(e())} produces the McFadden pseudo-R squared metric. +Prior to {cmd:domme} version 1.2 you had to provide a three character code as +an option to the {opt fitstat()}. {opt fitstat()} still accepts the the +McFadden pseudo-R squared (as {opt fitstat(e(), mcf)}) +(See Example #1). {pmore}Note that {cmd:domme} has no default fit statistic and the user is -required to provide a fit statistic option. In addition, the built-in options -assume the command in {opt reg()} ereturn specific scalars. {res:mcf} only -requires {res: e(ll)}. {res:est} requires both {res: e(ll)} and {res: e(N)}. -{res:aic} requires {res: e(ll)} and {res: e(parm)}. -Finally, {res: bic} requires {res: e(ll)}, {res: e(parm)}, and {res: e(N)}. +required to provide a fit statistic option. The built-in option +assumes the command in {opt reg()} ereturns {res: e(ll)}. {phang}{opt sets([PEset_1] ... [PEset_x])} binds together parameter estimate constraints as a set that are always constrained jointly and act as a @@ -278,30 +267,29 @@ based on overall model fit statistics that decrease with better fit (e.g., AIC, {phang}Example 1: Path analysis/seemingly unrelated regression (SUR) with built in McFadden pseudo-R squared{p_end} {phang} {stata sureg (price = length foreign gear_ratio) (headroom = mpg)} {p_end} -{phang} {stata domme (price = length foreign gear_ratio) (headroom = mpg), reg(sureg (price = length foreign gear_ratio) (headroom = mpg)) fitstat(e(), mcf)} {p_end} +{phang} {stata domme (price = length foreign gear_ratio) (headroom = mpg), reg(sureg (price = length foreign gear_ratio) (headroom = mpg)) fitstat(e())} {p_end} -{phang}Example 2: Zero-inflated Poisson with built in BIC{p_end} +{phang}Example 2: Zero-inflated Poisson{p_end} {phang} {stata generate zi_pr = price*foreign} {p_end} -{phang} {stata zip zi_pr headroom trunk,inflate(gear_ratio turn)} {p_end} -{phang} {stata domme (zi_pr = headroom trunk) (inflate = gear_ratio turn), reg(zip zi_pr headroom trunk) f(e(), bic) ropt(inflate(gear_ratio turn)) reverse} {p_end} +{phang} {stata zip zi_pr headroom trunk, inflate(gear_ratio turn)} {p_end} +{phang} {stata domme (zi_pr = headroom trunk) (inflate = gear_ratio turn), reg(zip zi_pr headroom trunk) f(e()) ropt(inflate(gear_ratio turn))} {p_end} -{phang}Example 3: Path analysis/SUR model with all option {p_end} +{phang}Example 3: Path analysis/SUR model with all option{p_end} {phang} {stata sem (foreign <- headroom) (price <- foreign length weight) (weight <- turn)} {p_end} -{phang} {stata estat ic} {p_end} -{phang} {stata domme (price = length foreign) (foreign = headroom), all((price = weight) (weight = turn)) reg(sem (foreign <- headroom) (price <- foreign length weight) (weight <- turn)) fitstat(e(), aic) reverse} {p_end} +{phang} {stata domme (price = length foreign) (foreign = headroom), all((price = weight) (weight = turn)) reg(sem (foreign <- headroom) (price <- foreign length weight) (weight <- turn)) fitstat(e())} {p_end} -{phang}Example 4: Generalized negative binomial with all and parmeters treated as _cons in the dominance analysis (i.e., _b[price:foreign]) {p_end} +{phang}Example 4: Generalized negative binomial with all and parmeters treated as _cons in the dominance analysis (i.e., _b[price:foreign] using {cmd:e(ll)} as fit statistic {p_end} {phang} {stata gnbreg price foreign weight turn headroom, lnalpha(weight length)} {p_end} -{phang} {stata domme (price = turn headroom) (lnalpha = weight length), reg(gnbreg price foreign weight turn headroom) f(e(), mcf) ropt(lnalpha(weight length)) all( (price = weight) )} {p_end} +{phang} {stata domme (price = turn headroom) (lnalpha = weight length), reg(gnbreg price foreign weight turn headroom) f(e(ll)) ropt(lnalpha(weight length)) all( (price = weight) )} {p_end} {phang}Example 5: Generalized structural equation model with factor variables{p_end} {phang} {stata sysuse nlsw88, clear} {p_end} {phang} {stata gsem (wage <- union hours, regress) (south <- age ib1.race union, logit)} {p_end} -{phang} {stata domme (wage = union hours) (south = age union), reg(gsem (wage <- union hours, regress) (south <- age ib1.race union, logit)) fitstat(e(), mcf) sets([(south = 2.race 3.race)])}{p_end} +{phang} {stata domme (wage = union hours) (south = age union), reg(gsem (wage <- union hours, regress) (south <- age ib1.race union, logit)) fitstat(e()) sets([(south = 2.race 3.race)])}{p_end} {phang}Example 6: Generalized structural equation model with sets to evaluate independent variables{p_end} {phang} {stata gsem (south union <- wage tenure ttl_exp, logit)} {p_end} -{phang} {stata domme, reg(gsem ( south smsa union <- wage tenure ttl_exp, logit)) fitstat(e(), mcf) sets( [(south = wage) (union = wage)] [(south = tenure) (union = tenure)] [(south = ttl_exp) (union = ttl_exp)]) } +{phang} {stata domme, reg(gsem ( south smsa union <- wage tenure ttl_exp, logit)) fitstat(e()) sets( [(south = wage) (union = wage)] [(south = tenure) (union = tenure)] [(south = ttl_exp) (union = ttl_exp)]) } {p_end} {phang}Examples 7: Replicating results from {cmd:domin}{p_end} @@ -314,9 +302,9 @@ based on overall model fit statistics that decrease with better fit (e.g., AIC, {pmore} {stata domin rep78 trunk weight length, reg(ologit) fitstat(e(r2_p)) all(turn)} {p_end} {pmore} {stata domme (rep78 = trunk weight length), reg(ologit rep78 trunk weight length turn) fitstat(e(), mcf) all((rep78 = turn))} {p_end} -{pmore}7c: Poisson regression with BIC fitstat and constant-only comparison using reverse{p_end} -{pmore} {stata domin price mpg rep78 headroom, reg(fitdom, fitstat_fd(r(S)[1,6]) reg_fd(poisson) postestimation(estat ic)) fitstat(e(fitstat)) consmodel reverse} {p_end} -{pmore} {stata domme (price = mpg rep78 headroom), reg(poisson price mpg rep78 headroom) fitstat(e(), bic) reverse} {p_end} +{pmore}7c: Poisson regression with log-likelihood fitstat and constant-only comparison using reverse{p_end} +{pmore} {stata domin price mpg rep78 headroom, reg(poisson) fitstat(e(ll)) consmodel reverse} {p_end} +{pmore} {stata domme (price = mpg rep78 headroom), reg(poisson price mpg rep78 headroom) fitstat(e(ll)) reverse} {p_end} {marker remark}{...} {title:6. Final Remarks} @@ -328,11 +316,8 @@ and considerations apply to {cmd:domme} as well. {pstd}Any parameter estimates in the model's {opt reg()} specification but not in the initial statements, the {opt sets()}, or {opt all()} are considered to be a part of the to a constant-only model (see Examples #4 and #7c). When using {cmd:domme}'s -built-in fit statistics, parameter estimates in the constant-only model -will be used to compute the baseline model for the {res:mcf} and {res:est} -pseudo-R-squares but will be reported as a part of the constant model for the -{res:aic} and {res:bic}. Other fit statistics supplied to {cmd: domme} will, like -{res:aic} and {res:bic}, be reported as a part of the constant model. +built-in fit statistic, the constant-only model will be used to compute the +baseline model. {pstd}Note that {cmd:domme} does not check to ensure that the parameters supplied it are in the model and it is the user's responsibility to ensure that the @@ -372,7 +357,7 @@ Please cite as: {title:Author} {p 4}Joseph N. Luchman{p_end} -{p 4}Principal Scientist{p_end} +{p 4}Research Fellow{p_end} {p 4}Fors Marsh{p_end} {p 4}Arlington, VA{p_end} {p 4}jluchman@forsmarsh.com{p_end} diff --git a/fitdom.sthlp b/fitdom.sthlp index 932b932..ef499de 100644 --- a/fitdom.sthlp +++ b/fitdom.sthlp @@ -134,7 +134,7 @@ matrix name (i.e., {cmd:r(S)}) but also the matrix subscript selecting the AIC i {title:Author} {p 4}Joseph N. Luchman{p_end} -{p 4}Principal Scientist{p_end} +{p 4}Research Fellow{p_end} {p 4}Fors Marsh{p_end} {p 4}Arlington, VA{p_end} {p 4}jluchman@forsmarsh.com{p_end} diff --git a/mi_dom.ado b/mi_dom.ado index 82b5cd0..dfb7035 100644 --- a/mi_dom.ado +++ b/mi_dom.ado @@ -1,4 +1,4 @@ -*! fitdom version 0.0.0 8/14/2023 Joseph N. Luchman +*! mi_dom version 0.0.1 1/12/2024 Joseph N. Luchman program define mi_dom, eclass //history and version information at end of file @@ -45,4 +45,6 @@ end /* programming notes and history - mi_dom version 0.0.0 - August 14, 2023 +// 0.0.1 - January 12, 2024 + - fixed 'which' reference - erroneously noted 'fitdom' as program - thanks to Eric Melse for bug report diff --git a/mi_dom.sthlp b/mi_dom.sthlp index 94606aa..1e06137 100644 --- a/mi_dom.sthlp +++ b/mi_dom.sthlp @@ -1,5 +1,5 @@ {smcl} -{* *! version 0.0.0 August 14, 2023 J. N. Luchman}{...} +{* *! version 0.0.1 January 12, 2024 J. N. Luchman}{...} {cmd:help mi_dom} {hline}{...} @@ -61,7 +61,7 @@ This produces a command structure like {it:mi estimate, miopts: reg_mi} for each {title:Author} {p 4}Joseph N. Luchman{p_end} -{p 4}Principal Scientist{p_end} +{p 4}Research Fellow{p_end} {p 4}Fors Marsh{p_end} {p 4}Arlington, VA{p_end} {p 4}jluchman@forsmarsh.com{p_end} diff --git a/mixdom.sthlp b/mixdom.sthlp index 5ba6af5..565bfe9 100644 --- a/mixdom.sthlp +++ b/mixdom.sthlp @@ -68,7 +68,7 @@ This option was named {opt xtmopt()} in {cmd:mixdom} versions previous to 2.0.0 {title:Author} {p 4}Joseph N. Luchman{p_end} -{p 4}Principal Scientist{p_end} +{p 4}Research Fellow{p_end} {p 4}Fors Marsh{p_end} {p 4}Arlington, VA{p_end} {p 4}jluchman@forsmarsh.com{p_end} diff --git a/mvdom.sthlp b/mvdom.sthlp index f6be8b2..a93f2ad 100644 --- a/mvdom.sthlp +++ b/mvdom.sthlp @@ -57,7 +57,7 @@ It is not recommended for use as an estimation command outside of {cmd:domin}. {title:Author} {p 4}Joseph N. Luchman{p_end} -{p 4}Principal Scientist{p_end} +{p 4}Research Fellow{p_end} {p 4}Fors Marsh{p_end} {p 4}Arlington, VA{p_end} {p 4}jluchman@forsmarsh.com{p_end}