From ed49cb53ac7eda308c112e4bab107afed13e670e Mon Sep 17 00:00:00 2001 From: Joseph Luchman Date: Wed, 14 Dec 2022 18:33:23 -0600 Subject: [PATCH 1/2] Major update to conform with compiled Mata function in domin --- domme.ado | 1872 +++++++++++++++++++++++------------------------------ 1 file changed, 823 insertions(+), 1049 deletions(-) diff --git a/domme.ado b/domme.ado index c570989..77c94e0 100644 --- a/domme.ado +++ b/domme.ado @@ -1,716 +1,741 @@ -*! domme version 1.0.1 4/17/2021 Joseph N. Luchman +*! domme version 1.1.0 xx/xx/202x Joseph N. Luchman -program define domme, eclass // ~ history and version information at end of file ~ +program define domme, eclass // ~ history and version information at end of file ~ -version 15 + version 15.1 -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'") { //replay results - error if "by"; domme allows nothing in "anything" - allow it to replay only when there is nothing in options - if ("`e(cmd)'" != "domme") error 301 - - if _by() error 190 - - Display `0' - - exit 0 - -} + if ("`e(cmd)'" != "domme") error 301 -/*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 + if _by() error 190 -/*exit conditions*/ -capture which lmoremata.mlib //is moremata present? + Display `0' -if _rc { //if moremata cannot be found, tell user to install it. + exit 0 - display "{err}Module {cmd:moremata} not found. Install " /// - "{cmd:moremata} here {stata ssc install moremata}." - - exit 198 + } + + /*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 lb_dominance.mlib //is -domin- present? -} -/*general set up*/ -tempname ranks gendom stzd_gendom cdldom cptdom //temporary matrices for results + if _rc { //if -domin- cannot be found, tell user to install it. -local two "`anything'" //rename the input equation-to-independent variable mapping in "anything": done to generalize and simplify loop below + display "{err}Module {cmd:domin} not found. Install " /// + "{cmd:domin} here {stata ssc install domin}." -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 + 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... + 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 + 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... - - } + if !regexm("`one'", "=") { //exit if there is no equal sign to make an equation - local one = regexr("`one'", "[/(]", "") //remove left paren from equation statement - - local one = regexr("`one'", "[/)]", "") //remove right paren from equation statement + display "{err}Equation {cmd:`one'} is missing a {cmd:=} to " /// + "distinguish equation and independent variable names." - gettoken dv ivlist: one, parse("=") //further parse the focal equation to separate out dependent from independent variables + local drop_exit = 1 //indicate that constraints will be dropped + + continue, break //stop the -while- loop... - 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 vatiable 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 - - continue, break //stop the -while- loop... - - } - - 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 one = regexr("`one'", "[/(]", "") //remove left paren from equation statement - local ivs = "`ivs' `add'" //record the parameter estimate for display - - local constrs "`constrs' `constr'" //add the constraint associated with the parameter estimate to the list of constraints - - } - - 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 within equation -forvalues- loop... - } - - if `drop_exit' continue, break //stop the overall -while- loop... - -} + local one = regexr("`one'", "[/)]", "") //remove right paren from equation statement - //**process and obtain sets**// -local two "`sets'" //rename the set-based equation-to-independent variable mapping: done to generalize and simplify loop below + gettoken dv ivlist: one, parse("=") //further parse the focal equation to separate out dependent from independent variables -local setcount = 0 //used if there are sets to start a count and number them + if ( `: list sizeof dv' != 1 ) | regexm("`dv'", "=") { //multiple dependent variables/equations or no dependent variable where one should be... exit -while strlen("`two'") & !`drop_exit' { //if nothing's wrong so far and there are sets... + display "{err}Invalid equation name specified for {cmd:(`dv'`ivlist')}." - gettoken one two: two, bind //parse the sets of equation-to-independent variable mappings in parentheses + local drop_exit = 1 //indicate that constraints will be dropped + + continue, break //stop the -while- loop... - local 2 = trim("`one'") //rename the set to generalize and simplify loop below - - 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... - - display "{err}Set {cmd:`2'} not bound by brackets, {cmd:[ ]}." - - local drop_exit = 1 //indicate that constraints will be dropped - - continue, break //stop the -while- loop... - - } - - 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 - - local 2 = trim(subinstr("`2'", "]", "", 1)) //remove right bracket - - gettoken 1 2: 2, bind //parse the equation-to-independent variable mapping to bind all brackets together and pull out first equation - - 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... - - display "{err}Equation {cmd:`1'} in set {cmd:`one'} not bound by " /// - "parentheses; {cmd:( )}." - - 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 set {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 1 = subinstr("`1'", "(", "", 1) //remove left paren - - local 1 = trim(subinstr("`1'", ")", "", 1)) //remove right paren + 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 - 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 set {cmd:`one'}." - - 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 vatiable list... exit as it breaks process - display "{err}Empty set of independent variables specified for " /// - "equation {cmd:`dv'} in set {cmd:`one'}." - - local drop_exit = 1 //indicate that constraints will be dropped - - continue, break //stop the -while- loop... - + "equation {cmd:`dv'}." + + local drop_exit = 1 //indicate that constraints will be dropped + + continue, break //stop the -while- loop... + } - - 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 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 - + + 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 set`setcount' "`set`setcount'' `add'" //add parameter estimate label to set + 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 + + local constrs "`constrs' `constr'" //add the constraint associated with the parameter estimate to the list of constraints - local cset`setcount' "`cset`setcount'' `constr'" //add constraint to set - + 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... - + + 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... - + + local drop_exit = 1 //indicate that constraints will be dropped + + continue, break //stop the -while- loop... + } - - if `drop_exit' continue, break //stop the within equation, -forvalues- loop ... - + + if `drop_exit' continue, break //stop the within equation -forvalues- loop... } - - if `drop_exit' continue, break //stop the -while- loop looking for equations in a set... - - } - - if `drop_exit' continue, break //stop the -while- overall sets loop... - - else { //nothing went wrong, so entire set was processed, thus... - - local ivs "`ivs' _set:set`setcount'" //add parameter estimate set to all all subsets - the chevrons bind sets together; Mata separates them - - local constrs "`constrs' <`cset`setcount''>" //add constraint set to all other constraints - the chevrons bind sets together; Mata separates them - - } - -} -if ( `:list sizeof ivs' < 2 ) & !`drop_exit' { //exit if too few parameter estimates/sets (otherwise prodices cryptic Mata error) + if `drop_exit' continue, break //stop the overall -while- loop... - display "{err}{cmd:domme} requires at least 2 parameter estimates " /// - "or parameter estimate sets." - - local drop_exit = 1 //indicate that constraints will be dropped + } -} + //**process and obtain sets**// + local two "`sets'" //rename the set-based equation-to-independent variable mapping: done to generalize and simplify loop below -if strlen("`addconstr'") { //if there are custom constraints added to the analysis by number... + local setcount = 0 //used if there are sets to start a count and number them - foreach constrnum of numlist `addconstr' { //foreach constraint... + while strlen("`two'") & !`drop_exit' { //if nothing's wrong so far and there are sets... - local ivs "`ivs' _cns:cns`constrnum'" // format the constraint as a domme-able parameter - - local constrs "`constrs' `constrnum'" // add the constraint to the domme-able constraint list + gettoken one two: two, bind //parse the sets of equation-to-independent variable mappings in parentheses - } - -} + local 2 = trim("`one'") //rename the set to generalize and simplify loop below -local nobrkt_constrs = /// - subinstr(subinstr("`constrs'", ">", "", .), "<", "", .) //local macro where chevrons in constraint list is removed to estimate the constant model + 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... -/*finalize setup*/ - //**markouts and sample**// -tempvar touse keep //declare sample marking variables + display "{err}Set {cmd:`2'} not bound by brackets, {cmd:[ ]}." -tempname obs allfs consfs //declare temporary scalars + local drop_exit = 1 //indicate that constraints will be dropped -mark `touse' //declare marking variable using Stata norm "touse" + continue, break //stop the -while- loop... -quietly generate byte `keep' = 1 `if' `in' //generate tempvar that adjusts for "if" and "in" statements - -markout `touse' `keep' //mark sample for only really works with the if & in; all other missing-based adjustments must derive from running full model + } - //**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 - -if _rc & !`drop_exit' { //exit if model is not estimable or program results in error - return the returned code + local setcount = `setcount' + 1 //increment setcount in the case that the first (and subsequent) set is valid - display "{err}{cmd:`reg'} resulted in an error." + while strlen("`2'") & !`drop_exit' { //mirrors the processing of the equation-to-independent variable mapping for "anything" if they exist within a set... - local exit_code = `=_rc' - - local drop_exit = 1 //indicate that constraints will be dropped + local 2 = subinstr("`2'", "[", "", 1) //remove left bracket -} + local 2 = trim(subinstr("`2'", "]", "", 1)) //remove right bracket - //**fitstat processing**// -if regexm("`fitstat'", "^e\(\)") & !`drop_exit' { //if internal fitstat indication is in "fitstat" option; must have "e()" + gettoken 1 2: 2, bind //parse the equation-to-independent variable mapping to bind all brackets together and pull out first equation - 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 - - local built_in = 1 //indicate that the program is to produce it's own fitstat - - local built_in_style "`r(style)'" //indicate the style of built in statistic to ask for - - } + 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... -} + display "{err}Equation {cmd:`1'} in set {cmd:`one'} not bound by " /// + "parentheses; {cmd:( )}." -else local built_in = 0 //if not built in, then indicate the estimation model produces it's own fitstat + local drop_exit = 1 //indicate that constraints will be dropped -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 + continue, break //stop the -while- loop... - 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 + } -} + if !regexm("`1'", "=") { //exit if there is no equal sign to make an equation -capture assert sign(`fitstat') != -1 //what is the sign of the fitstat? -domme- works best with positive ones - warn and proceed + display "{err}Equation {cmd:`1'} in set {cmd:`one'} is missing a " /// + "{cmd:=} to distinguish equation and independent variable names." -if _rc & !`drop_exit' { //merely note that -domme- works best with positive ones (that's what's expected) + local drop_exit = 1 //indicate that constraints will be dropped - 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 + continue, break //stop the -while- loop... -} + } - //**count observations**// -if !`drop_exit' { //if nothing's wrong so far... + local 1 = subinstr("`1'", "(", "", 1) //remove left paren - quietly replace `touse' = e(sample) //replace sample marking variable with e(sample) from the model run + local 1 = trim(subinstr("`1'", ")", "", 1)) //remove right paren - 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 - - } + gettoken dv ivlist: 1, parse("=") //further parse the focal equation to separate out dependent from independent variables - else if inlist("`e(wtype)'", "iweight", "fweight") { //if the weights do affect obs (for frequency and importance weights) + if ( `: list sizeof dv' != 1 ) | regexm("`dv'", "=") { //multiple dependent variables or no dependent variable where one should be... exit - quietly summarize `=regexr("`e(wexp)'", "=", "")' if `touse' //tally up "touse" by summing weights - - scalar `obs' = r(sum) //pull out the number of observations included - - } + display "{err}Invalid equation name specified for " /// + "{cmd:(`dv'`ivlist')} in set {cmd:`one'}." - //**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 + local drop_exit = 1 //indicate that constraints will be dropped - if `:list sizeof all' { //if there is something in the "all" option - - 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 + continue, break //stop the -while- loop... - 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 + 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 - 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... - - } - - 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... - + "equation {cmd:`dv'} in set {cmd:`one'}." + + local drop_exit = 1 //indicate that constraints will be dropped + + continue, break //stop the -while- loop... + } - - 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 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 - + + 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 - + 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 set`setcount' "`set`setcount'' `add'" //add parameter estimate label to set + + local cset`setcount' "`cset`setcount'' `constr'" //add constraint to set + } - - else { //if no unused constraints remain... - + + 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... - + + 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... - + + if `drop_exit' continue, break //stop the within equation, -forvalues- loop ... + } + + if `drop_exit' continue, break //stop the -while- loop looking for equations in a set... + + } + + if `drop_exit' continue, break //stop the -while- overall sets loop... + + else { //nothing went wrong, so entire set was processed, thus... + + local ivs "`ivs' _set:set`setcount'" //add parameter estimate set to all all subsets - the chevrons bind sets together; Mata separates them + + local constrs "`constrs' <`cset`setcount''>" //add constraint set to all other constraints - the chevrons bind sets together; Mata separates them + 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*/ + } - + } - - //**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' compute_fitstat `fit_opts' iscons //if a built-in fistat desired, estimate it; note constant model - - scalar `consfs' = `fitstat' //return constant model's fitstat - - if strlen("`all'") { //distinguishes "all subsets" from "constant" fitstats/models + 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 + + local constrs "`constrs' `constrnum'" // add the constraint to the domme-able constraint list + + } + + } + + 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 + + markout `touse' `keep' //mark sample for only really works with the if & in; all other missing-based adjustments must derive from running full model + + //**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 + + 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." + + local exit_code = `=_rc' + + local drop_exit = 1 //indicate that constraints will be dropped + + } + + //**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 + + local built_in = 1 //indicate that the program is to produce it's own fitstat + + local built_in_style "`r(style)'" //indicate the style of built in statistic to ask for + + } + + } + + else local built_in = 0 //if not built in, then indicate the estimation model produces it's own fitstat + 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 + + } + + 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 + + } + + //**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 + + 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... + + } + + 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... + + } + + 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... + + } + + } + + } + + //**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' `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 + constraints(`nobrkt_constrs' `allcset' `extraconstr') //all constraints used - this estimates the "constant" model - scalar `allfs' = `fitstat' //return all subset model's fitstat as distinct from any constant model + if `built_in' { - - } - - 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: me_dominance("`constrs'", "`conditional'", "`complete'", `=`allfs'', /// - `=`consfs'', "`mifile'", `built_in', `cons_no_add', "`reg'", "`ropts'", /// - "`extraconstr'", "`weight'`exp'", "`touse'", "`fit_opts'") //invoke "me_dominance()" function in Mata + compute_fitstat `fit_opts' iscons //if a built-in fistat desired, estimate it; note constant model + + scalar `consfs' = `fitstat' + + } -/*translate r-class results from me_dominance() into temp results*/ - matrix `gendom' = r(gendom) //general dominance statistics - - matrix `stzd_gendom' = r(stzd_gendom) //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 + if !`built_in' & !missing(`fitstat') /// + scalar `consfs' = `fitstat' //return constant model's fitstat if user supplied and not missing + + 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: dominance(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 + + matrix coleq `cptdom' = dominated? //name the equation for the columns "dominated?" + + matrix rownames `cptdom' = `cptivs' //name the rows of the complete dominance designations + + matrix roweq `cptdom' = dominates? //name the equation for the columns "dominates?" + + } -/*processing display results*/ - //**name matrices**// - matrix colnames `gendom' = `ivs' //name the columns of general dominance statistic vector + if !strlen("`conditional'") { //if the conditional dominance matrix was not suppressed... - if strlen("`reverse'") { //reverse the direction and interpretation of ranked and standardized general dominance statistics + matrix rownames `cdldom' = `ivs' //name the rows of the conditional dominance matrix - 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 + 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 `stzd_gendom' = `ivs' //name the columns of stanadrdized general dominance statistic vector + matrix colnames `cdldom' = `colcdl' //replace the column names of the conditional dominance matrix with the number of "orders" which matches their counting sequennce - matrix colnames `ranks' = `ivs' //name the columns of ranked general dominance statistic vector + matrix coleq `cdldom' = #param_ests //equation names for all columns are "#param_ests" - 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 + if strlen("`e(title)'") local title "`e(title)'" //if the estimation command has an "e(title)" returned, save it - matrix colnames `cptdom' = `cptivs' //name the columns of the complete dominance designations - - matrix coleq `cptdom' = dominated? //name the equation for the columns "dominated?" - - matrix rownames `cptdom' = `cptivs' //name the rows of the complete dominance designations - - matrix roweq `cptdom' = dominates? //name the equation for the columns "dominates?" - - } + else if !strlen("`e(title)'") & strlen("e(cmd)") local title "`e(cmd)'" //...otherwise save the "e(cmd)" as that's informative too - 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" - - } + else local title "Custom user analysis" //...finally, if none of the options are returned, call it "custom user analysis" - 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 + /*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) - else local title "Custom user analysis" //...finally, if none of the options are returned, call it "custom user analysis" + if strlen("`setcount'") { //if there are sets... -/*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) + ereturn hidden scalar setcount = `setcount' //hidden scalar for use in "display" - if strlen("`setcount'") { //if there are sets... + forvalues set = 1/`setcount' { //for each set... - ereturn hidden scalar setcount = `setcount' //hidden scalar for use in "display" + ereturn local set`set' = trim("`set`set''") //make a separate local macro with it's parameter estimate label contents - 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 - } - - } - else ereturn hidden scalar setcount = 0 //...otherwise hidden set count is 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 disp_title "`title'" //hidden title for display (hence, "disp_title") - ereturn hidden local reverse "`reverse'" //hidden indicator for reverse - for display + 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 `: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 + 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 reg "`reg'" //return command used in -domme- in macro - ereturn local fitstat "`fitstat'" //return the name of the fitstat used + 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 cmd "domme" //this command is -domme-; return that + ereturn local fitstat "`fitstat'" //return the name of the fitstat used - ereturn local title `"Dominance analysis for multiple equations"' //The title/description of -domme- + ereturn local cmd "domme" //this command is -domme-; return that - ereturn local cmdline `"domme `0'"' //full command as read into -domme- the 0 macro is everything after the command's name + ereturn local title `"Dominance analysis for multiple equations"' //The title/description of -domme- - ereturn scalar fitstat_o = r(fs) //overall fitstat value + ereturn local cmdline `"domme `0'"' //full command as read into -domme- the 0 macro is everything after the command's name - if `:list sizeof all' ereturn scalar fitstat_a = `allfs' + /// - `consfs'*`cons_no_add' //all subsets fitstat value + ereturn scalar fitstat_o = r(fs) //overall fitstat value - if sign( `consfs'*`cons_no_add' ) ereturn scalar fitstat_c = /// //constant model fitstat value - `consfs'*`cons_no_add' + if `:list sizeof all' ereturn scalar fitstat_a = `allfs' + /// + `consfs'*`cons_no_add' //all subsets fitstat value - if !strlen("`conditional'") ereturn matrix cdldom `cdldom' //return conditional dominance matrix - - if !strlen("`complete'") ereturn matrix cptdom `cptdom' //return complete dominance designations + 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 - ereturn matrix ranking `ranks' //return ranked general dominance vector + if !strlen("`complete'") ereturn matrix cptdom `cptdom' //return complete dominance designations - ereturn matrix std `stzd_gendom' //return standardized general dominance vector + ereturn matrix ranking `ranks' //return ranked general dominance vector - /*begin display*/ - Display + 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 + /*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 + + } - 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 - } - - } - exit `exit_code' //exit using applicable exit code + exit `exit_code' //exit using applicable exit code -} + } + + foreach constr of numlist `nobrkt_constrs' { //go through each constraint that was made... -foreach constr of numlist `nobrkt_constrs' { //go through each constraint that was made... - constraint drop `constr' //drop the constraint - -} + + } end @@ -718,570 +743,319 @@ end /*Display program*/ program define Display -version 15.1 + version 15.1 -/*set up*/ -tempname gendom stzd_gendom ranks //declare names for temporary data + /*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 `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 `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) + 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 diivs: colnames e(b) //obtain independent variable names -local eqivs: coleq e(b) //obtain dependent variable/equation 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("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 + 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 + 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) + /*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) + 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 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) + 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 _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}{col 13}{c |}{col 20}Stat.{col 35}Domin. Stat." -display "{txt}{hline 12}{c +}{hline 72}" + display "{txt}{hline 12}{c +}{hline 72}" -local current_eq "" //for separating equation from independent variable names in -forvalues- loop below + 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 + forvalues iv = 1/`:list sizeof diivs' { //for each entry of the e(b) vector... - 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'] - -} + 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 -display "{txt}{hline 12}{c BT}{hline 72}" + local current_eq = abbrev("`: word `iv' of `eqivs''", 11) //note current equation - truncate to 11 chars -if `cdltest' { //if conditional dominance matrix exists.... + local `iv' = abbrev("``iv''", 10) //abbreviate independent variable to 10 chars - display "{txt}Conditional dominance statistics" _newline "{hline 85}" - - matrix list e(cdldom), noheader format(%12.4f) - - display "{txt}{hline 85}" - -} + 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'] -if `cpttest' { //if complete dominance designations exist... + } - display "{txt}Complete dominance designation" _newline "{hline 85}" - - matrix list e(cptdom), noheader - - display "{txt}{hline 85}" - -} + display "{txt}{hline 12}{c BT}{hline 72}" -if `=`cpttest'*`cdltest'' { //if _both_ complete and conditional dominance designations exist - determine strongest dominance designation between each pair of parameter estimates + if `cdltest' { //if conditional dominance matrix exists.... - display _newline "{res}Strongest dominance designations" _newline + display "{txt}Conditional dominance statistics" _newline "{hline 85}" + + matrix list e(cdldom), noheader format(%12.4f) + + display "{txt}{hline 85}" - 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']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" - string colvector constr_mat - - real scalar constr_count, ntuples, display, fs, cptsum, comparecount, left_hand_param, /// - right_hand_param, cptdom, constr, param_est, tuple, order, cpt_num, cptorder + if !`cdl' mata: /// + st_numscalar("`cdl'", -1*((sum(st_matrix("e(cdldom)")[`dominator', .] /// + :")) //set up parsing rules - - tokenset(t, constrs) //register the "constrs" matrix as the one to be parsed - - constr_mat = tokengetall(t)' //obtain all constraint sets and constraints - - /*remove characters binding sets together (i.e., "<>")*/ - for ( constr = 1 ; constr <= rows(constr_mat) ; constr++ ) { //loop over all constraint vector entries... - - if ( substr(constr_mat[constr], 1, 1) == "<" ) { //if any entry begins with "<"... - - constr_mat[constr] = substr(constr_mat[constr], 1, /// - strlen(constr_mat[constr]) - 1) //first character removed ("<") - - constr_mat[constr] = substr(constr_mat[constr], 2, /// - strlen(constr_mat[constr])) //last character removed (">") - - } - - } - - /*set-up and compute all n-tuples of predictors and predictor sets*/ - constr_count = rows(constr_mat) //compute total # of constraint sets and constraints - - ntuples = 2^constr_count - 1 //compute total # of models/regressions - - printf("\n{txt}Total of {res}%f {txt}models/regressions\n", ntuples) - - if ( constr_count > 12 ) /// - printf("\n{txt}Computing all parameter estimate combinations\n") + mata: st_numscalar("`gen'", /// + st_matrix("e(b)")[1, `dominator']>st_matrix("e(b)")[1, `dominatee']) //compare dominator's general dominance statistic to dominatee's - indicators = J(constr_count, 2^constr_count, .) //set up matrix to be filled in which will generate all subsets - row "indicates" presence of parameter estimate - - for ( param_est = 1; param_est <= rows(indicators) ; param_est++ ) { //for each row in indicators matrix... - - combin = J( 1 , 2^( param_est-1 ) , 0 ) , /// - J( 1 , 2^( param_est-1 ) , 1 ) //make a binary matrix - start small, a zero and a 1, then 2 0's and 2 1's, etc... - - indicators[param_est, .] = J(1, 2^( constr_count-param_est ), combin) //spread the binary matrix just created across all rows - net effect is staggering all binaries to obtain all subsets in the final matrix - - } - - indicators = indicators[., 1..cols(indicators)-1] //drop the last row where all constraints are included - - indicators = (colsum(indicators) \ indicators)' //create a "counts" column on which to sort - - indicators = sort(indicators, (1..cols(indicators))) //sort, beginning with counts, followed by all other rows - net effect results in same sort order as cvpermute() - - indicators = indicators[|1, 2\ ., .|]' //omit count's column created before - - tuples = indicators:*constr_mat //apply string constraint values to all subsets indicator matrix - - /*all subsets regressions and progress bar syntax if predictors or sets of - predictors is above 5*/ - display = 1 //for the display of dots during estimation - keeps track of where the regressions are - every 5% there is another "." added - - if ( constr_count > 4 ) { - - printf("\n{txt}Progress in running all regression subsets\n{res}0%%" + /// - "{txt}{hline 6}{res}50%%{txt}{hline 6}{res}100%%\n") - - printf(".") - - displayflush() - - } + if !`gen' mata: st_numscalar("`gen'", /// + (st_matrix("e(b)")[1, `dominator'] 4 ) { //display another "." if a sufficient number of models have been estimated - - if ( floor( tuple/ntuples*20 ) > display ) { - - printf(".") - - displayflush() - - display++ - } - + } - preds = tuples[., tuple]' //take the constraint values in column "tuple" and transpose into row - - construse = invtokens(preds) //collpase constraint names into single string separated by spaces - - stata(reg + " [" + weightexp + "] if " + touse + ", " + ropts + /// - " constraints(" + construse + " " + extraconstr + ")", 1) //estimate model/regression in Stata - - if (built_in) stata("compute_fitstat " + fit_opts + " consmodel(" /// - + strofreal(consfs) + ")", 0) //obtain built-in fitstat - - fs = st_numscalar(st_local("fitstat")) - allfs - consfs*cons_no_add //adjust and record fitstat; omit constant and "all" subsets values - - fits = (fits, fs) //add fitstat to vector of fitstats + local names `:colfullnames e(b)' //obtain full names of e(b) vector + + mata: display((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) /// + :*" generally dominates "):+tokens(st_local("names")))', /// + "generally dominates")))') //complex command here - basic idea is that it takes the value in bestdom, displays the correct dominance type based on number, and plugs in the values of the parameter estimate labels + + display "" } - - fits = fits[2..ntuples + 1] //only keep non-empty fitstats (i.e., omit the first empty one) - /*define the incremental prediction matrices and combination rules*/ - include = abs( sign( strlen(tuples) ):-1 ) // matrix indicating whether a parameter estimate was included in a model which can be mapped onto "fits" vector; "-1" is due to constrain removing parameter estimate + if `=e(setcount)' { //if there are sets... - counts = colnonmissing( exp( ln(include) ) ) //number of parameter estimates in each regression - the exp(ln()) combo makes 0's missing + forvalues set = 1/`=e(setcount)' { //foreach set... + + display "{txt}Parameters in set`set': `e(set`set')'" //display the set's parameter estimate labels - noinclude = ( include:-1 ) //matrix indicating whether parameter estimate is not included in the model which can also be mapped into "fits" vector - - combs_at_ord = J(1, ntuples, 1):*comb(constr_count, counts) //matrix indicating the number of combinations at each "order"/# of parameter ests - - combs_less_1 = J(1, ntuples, 1):*comb(constr_count - 1, counts) //matrix indicating the number of combinations at each "order"/# of parameter ests - 1 - - combs_less_1 = (0, combs_less_1[., 2..ntuples]) //add a 0 to # combinations for first value... omit first "." value - - combs_at_ord = combs_at_ord - combs_less_1 //remove # of combinations for the "order" less the value at "order" - 1 - - include = include:*combs_at_ord //put all the adjusted combination counts into matrix when the param est "is included" - - noinclude = noinclude:*combs_less_1 //put all the "order" - 1 combination counts into matrix when the param est "is not included" - - /*compute conditional dominance*/ - if ( strlen(cdlcompu) == 0 ) { //if conditional dominance is to be computed... - - if ( constr_count > 5 ) /// - printf("\n\n{txt}Computing conditional dominance\n") - - cdl = J(constr_count, constr_count, 0) //dummy matrix to hold conditional dominance stats - - /*loop over orders (i.e., # of param ests) to obtain average - incremental prediction within order*/ - for ( order = 1 ; order <= constr_count ; order++ ) { //foreach order/number of params in a model - - cdl_include = include:^-1 //invert the counts for indluded to create denominators for included fitstats: these are a component of the within-order averages - - cdl_noinclude = noinclude:^-1 //invert the counts for non-indluded to create denominators for not included fitstats: these are a component of the within-order averages - - cdl_include = select(cdl_include:*fits, counts:==order) //at the focal order, select the included fitstats and divide by number of combinations - - if ( order > 1 ) { // at all orders (i.e., >1) where the marginal contribution is not the fitstat itself - - cdl_noinclude = select(cdl_noinclude:*fits, counts:==order-1) //at the focal order less one, select the included fitstats and divide by number of combinations - - cdl_sum = rowsum(cdl_include) + rowsum(cdl_noinclude) //now produce the sum of the weighted values, which happens to be the within order marginal average - - } - - else cdl_sum = rowsum(cdl_include) //sum the marginal contributions @ order 1 - - cdl[., order] = cdl_sum //replace the entries in cdl with the current values of cdl_sum, these are the within-order averages - } - - st_matrix("r(cdldom)", cdl) //after all cdl enties are registered, return r-class matrix "cdldom" - + } - - /*define the full design matrix - compute general dominance - (average conditional dominance across orders)*/ - design = (include + noinclude):*constr_count //create matrix that will have positive and negative signs in the correct places to obtain marginal contributions when summed then multiply by number of variables total (between-order average of within-order averages) - - design = design:^-1 //invert design matrix to create within-between order averaging weights - - gendom = colsum((design:*fits)') //general dominance weights created by computing product of weights and fitstats and summing for each IV - fs = rowsum(gendom) + allfs + consfs*cons_no_add //total fitstat is then sum of gen. dom. wgts replacing the constant-only model (when applicable/cons_no_add is "1") and the "all" subsets stat + if strlen("`e(all)'") /// + display "{txt}Parameter estimates included in all subsets: `e(all)'" //if there is an all subsets set, display those parameter estimate labels - st_matrix("r(gendom)", gendom) //return the general dom. wgts as r-class matrix +end - stzd_gendom = gendom:*fs^-1 //generate the standardized gen. dom. wgts - - st_matrix("r(stzd_gendom)", stzd_gendom) //return the stdzd general dom. wgts as r-class matrix - - st_matrix("r(ranks)", mm_ranks(gendom'*-1, 1, 1)') //return the ranks of the general dom. wgts as r-class matrix - st_numscalar("r(fs)", fs) //return overall fit statistic in r-class scalar - - /*compute complete dominance*/ - if (strlen(cptcompu) == 0) { - - if (constr_count > 5) printf("\n{txt}Computing complete dominance\n") +**# Program to compute built-in fitstats +program define compute_fitstat, rclass - cpt = J(constr_count, constr_count, 0) //dummy matrix for complete dominance - - select2cpt = (J(2, 1, 1) \ J(constr_count - 2, 1, 0)) //generate the "base" comparitor matrix for 2 parameter estimates - will be used in cvpermute() - - cvpermcpt = cvpermutesetup(select2cpt) //setup "base" for the permutations (only all pairs of comparisons desired) - - indicator = (1::constr_count) //generate "indicator" matrix for which parameter estimates are being compared - - for (cpt_num = 1; cpt_num <= comb(constr_count, 2); cpt_num++) { - - combincpt = cvpermute(cvpermcpt) //generate a combination of 2 parameter estimates to select and compare... - - rowcol = select(combincpt:*indicator, combincpt:==1) //note the row in which both parameter estimates being comapred are located - - focus = select(!sign(strlen(tuples)), combincpt:==1) //make a selector (1 vs. 0) matrix for pulling out all fitstats, only on focal parameter estimates - - rest = select(!sign(strlen(tuples)), combincpt:==0) //make a selector (1 vs. 0) matrix for pulling out all fitstats, only on non-focal parameter estimates - - cptsum = 0 //used as a index for determining complete dominance for the current comparison of 2 parameter estimatess - - focusfits = focus:*fits //create matrix of fitstats that correspond only to the focal comparisons - - for (cptorder = 1; cptorder <= constr_count - 1; cptorder++) { //for each order (up to # parameter estimates - 1) - - filterfits = select(focusfits, counts:==cptorder) //on the filtered fitstat matrix, pull out comparisons at a specific order - - focusorder = select(focus, counts:==cptorder) //on the indicator matrix, pull out comparisons at a specific order - - focusorder = colsum(focusorder) //on the filtered indicator fit matrix of order "cptorder", enumerate # of parameter estimates in each model - - restorder = select(rest, counts:==cptorder) //on the indicator matrix of non-focal fits, pull out comparisons at a specific order - - comparecount = 1 //counter to keep track of # of comparisons - - within_sel2cpt = (J(cptorder - 1, 1, 1) \ /// //another looped permutation to make all the specific comparisons w/in order - J(constr_count - cptorder - 1, 1, 0)) - - /*make comparisons between fitstat's - matching on predictors*/ - while ((comparecount <= comb(constr_count - 2, cptorder - 1)) & /// //so long as there are > 2 parameter estimates... loop for all comparisons - (constr_count > 2)) { - - if (cptorder == 1) withinfits = select(filterfits, /// //fitstats when only one parameter estimate is in the model (cptorder = 1 per row) - nothing needs to be subtracted to obtain marginal contribution - restorder[comparecount, .]:==0) - - else if (cptorder == 2) { //fitstats when a focal and only 1 other non-focal parameter is in the model - - withinfits = select(filterfits, focusorder:==1) //select the fitstats where only the focal parameters are in the model (i.e., not both parameters) - - withinrest = select(restorder, focusorder:==1) //pull out the associated columns where the non-focal parameter is - - withinfits = select(withinfits, withinrest[comparecount, .]:==1) //then select the fitstats where there are only the focal parameters are in the model with the non-focal paramter (i.e., both focal parameters are not together) + version 15.1 + + syntax, [ll(string) obs(string) parm(string) consmodel(real 0) aic bic mcf est iscons] + + local many_few = 0 + + if !strlen("`aic'`bic'`mcf'`est'") local mcf "mcf" //no fitstat option entry is McFadden by default + + if strlen("`aic'`bic'`mcf'`est'") > 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 - } - - else { //fitstats when >=2 parameter estimates are in the model - - withinfits = select(filterfits, focusorder:==1) //select the fitstats when only the focal parameters are in the model singly (i.e., not both parameters) - - withinrest = select(restorder, focusorder:==1) //pull out the associated columns where the non-focal parameter is - - cvp_withincpt = cvpermutesetup(within_sel2cpt) //set-up cvpermute to capture all possible combinations of a focal variable with multiple non-focal variables - - cvp_withincomb = cvpermute(cvp_withincpt)*10 //activate permutation of specific set of multiple non-focal variables (rescaled by 10 for use in exponentiating which seaparate columns by base-10 place value) - - revind = (constr_count - 2::1) //used for exponentiation below - - within_exp = J(constr_count - 2, 1, 10) //selection matrix to use for selecting fitstats - adjusted below - - cvp_withincomb = cvp_withincomb:^revind*(1/10) //matrix which now indicates location of a variable positionally by # of 0s (re-scaled back down by 10) - - within_exp = within_exp:^revind*(1/10) //obtain a selection matrix which is scaled the same as the combination matrix above - - within_exp = withinrest:*within_exp //rescale the "withinrest" matrix with only the current non-focal parameters are selected - - within_exp = colsum(within_exp) //make within_exp a rowmat so select() can use it - - cvp_withincomb = colsum(cvp_withincomb) //make cvp_withincomb a rowmat so select() can use it - - withinfits = select(withinfits, within_exp:==cvp_withincomb) //obtain only one specific combination of the non-focal parameters for the comparison - - } - - /*here the comparison is actually made and "cptsum" is updated*/ - left_hand_param = rowsum(withinfits[1, .]) //all the fitstats in row 1 are "left_hand_param"'s data - sum them (there should only be 1) - - right_hand_param = rowsum(withinfits[2, .]) //all the fitstats in row 2 are "right_hand_param"'s data - sum them (there should only be 1) - - cptdom = sign(left_hand_param - right_hand_param) //compare left and right parameters - keep sign only - - cptsum = cptsum + cptdom //add sign to current sum - - comparecount++ //increment comparecount and evaluate the while statement above... - - } - - } - - /*determine completely dominate, dominated by or none*/ - if (constr_count == 2) cptsum = sign(rowsum(focusfits[1, .]) - /// //simple situation w/ 2 parameters - rowsum(focusfits[2, .])) - - if (cptsum == 2^(constr_count - 2)) /// //if all the cptdom comarisons were "+" then, there is complete dominance for "left_hand_param" - cpt[rowcol[1, 1], rowcol[2, 1]] = 1 - - else if (cptsum == -2^(constr_count - 2)) /// //if all the cptdom comarisons were "-" then, there is complete dominance for "right_hand_param" - cpt[rowcol[1, 1], rowcol[2, 1]]= -1 - - else cpt[rowcol[1, 1], rowcol[2, 1]] = 0 //otherwise no complete dominance - - } - - cpt = cpt + cpt'*-1 //make cptdom matrix symmetric in what it is telling the user - - st_matrix("r(cptdom)", cpt) //return r-class matrix "cptdom" - } - -} -end + else if strlen("`aic'`bic'") { + if !strlen("`parm'") local parm = "e(rank)" -/*program to compute built-in fitstats for use in domme*/ -program define compute_fitstat, rclass + 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 + + } + + 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 -syntax, [ll(string) obs(string) parm(string) consmodel(real 0) aic bic mcf est iscons] +mata: -local many_few = 0 + mata set matastrict on -if !strlen("`aic'`bic'`mcf'`est'") local mcf "mcf" //no fitstat option entry is McFadden by default + struct domme_specs { -if strlen("`aic'`bic'`mcf'`est'") > 3 local many_few = 1 //too many fitstat options selected + real scalar built_in, built_in_cons -if !strlen("`ll'") local ll = "e(ll)" + string rowvector constraint_list + + class AssociativeArray scalar param_list -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 + } + +end -} +**# Mata function to execute 'domme-flavored' models +version 15.1 -else if strlen("`aic'`bic'") { +mata: - if !strlen("`parm'") local parm = "e(rank)" + mata set matastrict on - 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 + 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 + + string scalar Constraints, remove_list, key + + string rowvector Param_List, out_list + + real rowvector which_constraints + + Param_List = tokens(Params_in_Model) + + loop_length = length(Param_List) + + remove_list = "" + + for (x = 1; x <= loop_length; x++) { + + remove_list = remove_list + " " + /// + model_specs.param_list.get(Param_List[x]) + + } + + out_list = tokens(remove_list) + + which_constraints = /// + !colsum( /// + J(length(out_list), 1, model_specs.constraint_list):==out_list' /// + ) -} + Constraints = invtokens(select(model_specs.constraint_list, which_constraints)) + + stata("\`reg' [\`weight'\`exp'] if \`touse', \`ropts'" + /// + " constraints(" + Constraints + " \`extraconstr')", 1) //estimate model/regression in Stata + -else local many_few = 1 + if (model_specs.built_in) { + + 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")) -return local style "`aic'`bic'`mcf'`est'" + fitstat = fitstat - /// + all_subsets_fitstat - /// + constant_model_fitstat -return scalar many_few = `many_few' + return(fitstat) -return scalar fitstat = `fitstat' + } end /* programming notes and history -- domme version 1.0 - date - July 2, 2019 --base version - // 1.0.1 - April 17, 2021 (initiating new versioning: #major.#minor.#patch) - -update to documentation for SJ article - -bug fix on constraint dropping with all() option and use with xi: + - domme version 1.0 - date - July 2, 2019 + -base version + // 1.0.1 - April 17, 2021 (initiating new versioning: #major.#minor.#patch) + -update to documentation for SJ article + -bug fix on constraint dropping with all() option and use with xi: + // 1.1.0 - xx/xx/xxxx + -leverage dominance() function in -domin-; create own function to function passing; Mata struct to handle input specs + - -domin- now a dependency. + - use an AssociativeArray to map parameters/parameter sets to constraints in a way that is conformable with -domin- */ From 3dd86e8971c459c2b1fba1fc4bbdad370aa2dcae Mon Sep 17 00:00:00 2001 From: Joseph Luchman Date: Tue, 7 Feb 2023 19:10:50 -0600 Subject: [PATCH 2/2] Version 1.1.0 Update --- domme.ado | 82 ++++----- domme.sthlp | 470 ++++++++++++++++++++++++++++++++-------------------- 2 files changed, 330 insertions(+), 222 deletions(-) diff --git a/domme.ado b/domme.ado index 77c94e0..4138da7 100644 --- a/domme.ado +++ b/domme.ado @@ -1,10 +1,10 @@ -*! domme version 1.1.0 xx/xx/202x Joseph N. Luchman +*! domme version 1.1.0 2/7/2023 Joseph N. Luchman -program define domme, eclass // ~ history and version information at end of file ~ +program define domme, eclass // ~ history and version information at end of file ~ version 15.1 - 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'") { //replay results - error if "by"; domme allows nothing in "anything" - allow it to replay only when there is nothing in options if ("`e(cmd)'" != "domme") error 301 @@ -18,16 +18,16 @@ program define domme, eclass // ~ history and version information at end of file /*define syntax*/ syntax [anything(id="equation names" equalok)] [in] [if] /// - [aw pw iw fw], reg(string) Fitstat(string) [Sets(string) /// + [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 + ROPts(string) ADDConstr(numlist)] //addconstr() and extraconstr() undocumented - for use in possible extensions /*exit conditions*/ - capture which lb_dominance.mlib //is -domin- present? + capture which lb_dominance.mlib //is -domin- present? - if _rc { //if -domin- cannot be found, tell user to install it. + if _rc { //if -domin- cannot be found, tell user to install it. display "{err}Module {cmd:domin} not found. Install " /// "{cmd:domin} here {stata ssc install domin}." @@ -37,56 +37,56 @@ program define domme, eclass // ~ history and version information at end of file } /*general set up*/ - mata: model_specs = domme_specs() //initiate instance of domme_specs() structure + mata: model_specs = domme_specs() //initiate instance of domme_specs() structure - tempname ranks gendom stzd_gendom cdldom cptdom //temporary matrices for results + 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 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 + 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... + 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 + 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 + 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 + local drop_exit = 1 //indicate that constraints will be dropped - continue, break //stop the -while- loop... + continue, break //stop the -while- loop... } - local one = regexr("`one'", "[/(]", "") //remove left paren from equation statement + local one = regexr("`one'", "[/(]", "") //remove left paren from equation statement - local one = regexr("`one'", "[/)]", "") //remove right 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 + 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 + 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 + local drop_exit = 1 //indicate that constraints will be dropped - continue, break //stop the -while- loop... + continue, break //stop the -while- loop... } - local ivlist = regexr("`ivlist'", "=", "") //remove the equal sign from the independent variable list + 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 + 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 + local drop_exit = 1 //indicate that constraints will be dropped - continue, break //stop the -while- loop... + continue, break //stop the -while- loop... } @@ -425,7 +425,7 @@ program define domme, eclass // ~ history and version information at end of file 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 - + 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"... @@ -545,16 +545,16 @@ program define domme, eclass // ~ history and version information at end of file if `built_in' { - compute_fitstat `fit_opts' iscons //if a built-in fistat desired, estimate it; note constant model + compute_fitstat `fit_opts' iscons //if a built-in fistat desired, estimate it; note constant model scalar `consfs' = `fitstat' } if !`built_in' & !missing(`fitstat') /// - scalar `consfs' = `fitstat' //return constant model's fitstat if user supplied and not missing + scalar `consfs' = `fitstat' //return constant model's fitstat if user supplied and not missing - if !`built_in' & missing(`fitstat') scalar `consfs' = 0 //otherwise assume the constant-only model is 0 + 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 @@ -584,12 +584,12 @@ program define domme, eclass // ~ history and version information at end of file 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 + 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 `gendom' = r(domwgts) //general dominance statistics - matrix `stzd_gendom' = r(sdomwgts) //standardized general dominance statistics + matrix `stzd_gendom' = r(sdomwgts) //standardized general dominance statistics matrix `ranks' = r(ranks) //ranks based on general dominance statistics @@ -625,12 +625,14 @@ program define domme, eclass // ~ history and version information at end of file 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 - - matrix coleq `cptdom' = dominated? //name the equation for the columns "dominated?" - + + mata: st_matrixcolstripe("`cptdom'", /// + ("?":+st_matrixrowstripe("`cptdom'")[,1], st_matrixrowstripe("`cptdom'")[,2])) //add name the equation for the rows ">?" replacement for "dominates?" } @@ -733,7 +735,7 @@ program define domme, eclass // ~ history and version information at end of file foreach constr of numlist `nobrkt_constrs' { //go through each constraint that was made... - constraint drop `constr' //drop the constraint + constraint drop `constr' //drop the constraint } @@ -1054,8 +1056,10 @@ end // 1.0.1 - April 17, 2021 (initiating new versioning: #major.#minor.#patch) -update to documentation for SJ article -bug fix on constraint dropping with all() option and use with xi: - // 1.1.0 - xx/xx/xxxx + // 1.1.0 - February 7, 2023 -leverage dominance() function in -domin-; create own function to function passing; Mata struct to handle input specs - -domin- now a dependency. - use an AssociativeArray to map parameters/parameter sets to constraints in a way that is conformable with -domin- + - extensive documentation update + - fixed complete dominance table display; retains equation name and adds informative prefix */ diff --git a/domme.sthlp b/domme.sthlp index b9a6b63..449d95c 100644 --- a/domme.sthlp +++ b/domme.sthlp @@ -1,207 +1,275 @@ {smcl} -{* *! version 1.0.1 April 17, 2021 J. N. Luchman}{...} +{* *! version 1.1.0 February 7, 2023 J. N. Luchman}{...} {cmd:help domme} -{hline}{...} {title:Title} {pstd} -{ul on}Dom{ul off}inance analysis for {ul on}m{ul off}ulitple {ul on}e{ul off}quation models{p_end} +{bf:domme} {hline 2} {ul on}Dom{ul off}inance analysis for {ul on}m{ul off}ulitple {ul on}e{ul off}quation models{p_end} + {title:Syntax} -{phang} +{p 8 16 2} {cmd:domme} [{cmd:(}{it:eqname1 = parmnamelist1}{cmd:)} {cmd:(}{it:eqname2 = parmnamelist2}{cmd:)} ... {cmd:(}{it:eqnameN = parmnamelistN}{cmd:)}] -{ifin} {weight}{cmd:,} {opt r:eg(full_estcmd)} -{opt f:itstat(returned_scalar | built_in_options)} -[{it:options}] +{ifin} {weight}{cmd:,} [{it:options}] +{synoptset 35 tabbed}{...} +{synopthdr} {synoptline} -{phang}{cmd:pweight}s, {cmd:aweight}s, {cmd:iweight}s, and {cmd:fweight}s are allowed but must -be able to be used by the command in {opt reg()}, see help {help weights:weights}. -{help Time series operators} are also allowed for commands in {opt reg()} that accept them. -Finally, {help Factor variables} are also allowed, but like weights and time series operators, -must be accepted by the command in {opt reg()}. - -{phang}{cmd:domme} requires installation of Ben Jann's {cmd:moremata} package -(install here: {stata ssc install moremata}). Users are strongly encouraged to install -{stata ssc install domin:domin} as well and read over its help file for basic information on -dominance analysis. - -{title:Development Webpage} +{syntab:Model} +{synopt :{opt r:eg(full_estimation_command)}}preditive model command to call{p_end} +{synopt :{opt rop:ts(command_options)}}options to command in {opt reg()}{p_end} +{synopt :{opt f:itstat(fitstat_opts)}}fit statistic returned by {opt reg()} or computed using built-in method{p_end} +{synopt :{opt s:ets([PEset_1] ... [PEset_x])}}sets of indepdendent variables{p_end} +{synopt :{opt a:ll(PEall)}}indepdendent variables included in all subets{p_end} + +{syntab:Reporting} +{synopt :{opt nocon:ditional}}suppresses computation of conditional dominance statistics{p_end} +{synopt :{opt nocom:plete}}suppresses computation of complete dominance designations{p_end} +{synopt :{opt rev:erse}}reverses interpretation for statistics that decrease with better fit{p_end} +{synoptline} +{p 4 6 2} +Command in {opt reg()} must accept {help constraint}s as a command option. {p_end} -{phang} Additional discussion of results, options, and conceptual issues on: +{p 4 6 2} +{cmd:aweight}s, {cmd:fweight}s, {cmd:iweight}s, and {cmd:pweight}s are +allowed; see {help weight}. Weight use is restricted to commands in {opt reg()} +that accept them.{p_end} -{phang}{browse "http://github.com/jluchman/domme/blob/master/README.md"} +{p 4 6 2} +Note that {cmd:domme} requires at least two parameter estimates or +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} -{phang} Please report bugs, requests for features, and contribute to as well as follow on-going development of {cmd:domin} on: +{p 4 6 2} +{cmd:domme} requires installation of the {cmd:domin} package +(install {stata ssc install domin:here}). {p_end} -{phang}{browse "http://github.com/jluchman/domme"} {title:Table of Contents} -{help domme##desc: 1. Description} -{help domme##setup: 2. Set-up} -{help domme##opts: 3. Options} -{help domme##remark: 4. Final Remarks} -{help domme##examp: 5. Examples} -{help domme##sav: 6. Saved Results} -{help domme##refs: 7. References} +{space 4}{help domme##desc: 1. Description} +{space 4}{help domme##setup: 2. Set-up} +{space 4}{help domme##opts: 3. Options} +{space 4}{help domme##sav: 4. Saved Results} +{space 4}{help domme##examp: 5. Examples} +{space 4}{help domme##remark: 6. Final Remarks} +{space 4}{help domme##refs: 7. References} + {marker desc}{...} {title:1. Description} {pstd} -Dominance analysis for mulitple equation models is an extention of standard dominance analysis -(see {help domin}) which focuses on finding the relative importance -of parameter estimates in an estimation model based on contribution of each parameter -estimate to an overall model fit statistic (see Luchman, 2019 for a discussion). -As an extension of standard dominance analysis, it is recommended that the user -familiarize themselves with standard dominance analysis before attempting to use the -multiple equation version of the methodology. +Dominance analysis (DA) is a methodology for determining the relative +importance of independent variables or parameter estimates in a predictive model. +The {cmd: domme} approach extends on {cmd: domin} by being able to accommodate more +predictive models such as those with multiple equations/dependent variables +(see Luchman, Lei, and Kaplan, 2020 for a discussion). As an extension of single +equation dominance analysis, it is recommended that the user familiarize +themselves with single equation dominance analysis (i.e., see {help domin:domin}) +before attempting to use the multiple equation version of the methodology. {pstd} -Dominance analysis for mulitple equation models differs from standard dominance analysis -primarily in how the ensemble of fit metrics are collected. Standard dominance analysis -obtains the ensemble of fit metrics to compute dominance statistics by including or -excluding independent variables from a statistical model. Dominance analysis for mulitple -equation models obtains the ensemble of fit metrics to compute dominance statistics by -using {help constraint}s which permit each parameter estimate to be estimated from the -data or constrained to zero in a given statistical model. Constraining a parameter estimate to -zero effectively omits the parameter from estimation and it cannot contribute to model fit. +This multiple equation DA implementation differs from the implementation of +the single equation version of DA in how it includes or excludes parameter +estimates from the model. Multiple equation DA uses {help constraint}s to exclude +parameter estimates from the model by constraining their values to be 0. When constrained +to be 0, the parameter estimate cannot affect prediction and, thus, is effectively +excluded. This is why {cmd: domme} can only be used with commands that accept +constraints. + {marker setup}{...} {title:2. Set-up} {pstd} -{cmd:domme} requires that all parameters to be dominance analyzed are written out in the -initial {res:(eqname = parmnamelist)} statements. {cmd:domme} will use the {res:(eqname = parmnamelist)} -statements (similar to those of commands like {help sureg}) to create parameter statements that -it will produce {help constraint:constraints} from. Each entry in {res:parmnamelist} is given a -separate constraint with the associated {res:eqname}. For example, the statement: +This implementation of multiple equatuion DA must be provided with the components +of the parameter estimates from which it will create "parameter equals 0" +constraints. How these constraints are constructed follows from the way in which +Stata names each parameter estimate. Consider, for example, the following logistic +regression model (estimated on the {cmd: sysuse auto} data). {pstd} -{res:(price = mpg turn trunk foreign)} +{cmd:logit foreign price mpg turn trunk} {pstd} -will create the series of four parameters: +Following the estimation of this model, the user can ask for the names of all +the parameters in the {cmd: e(b)} or coefficient matrix by using: + +{pstd} +{cmd: display "`: colfullnames e(b)'"} + +{pstd} +which produces + +{pstd} +{res:foreign:price foreign:mpg foreign:turn foreign:trunk foreign:_cons} + +{pstd} +This series of names are the parameter names for all the coefficients in the +logit model. In order to use DA on this model, the user needs +to supply {cmd: domme} with the names of all the parameters that will be used. +One way to supply these parameter names is to use the initial +{res:(eqname = parmnamelist)} statements. For instance: + +{pstd} +{cmd: domme (foreign = price mpg turn trunk), ...} + +{pstd} +implies four constraints: + +{phang} +{res: _b[foreign:price] = 0} + +{phang} +{res: _b[foreign:mpg] = 0} + +{phang} +{res: _b[foreign:turn] = 0} + +{phang} +{res: _b[foreign:trunk] = 0} {pstd} -{res:_b[price:mpg] _b[price:turn] _b[price:trunk] _b[price:foreign]} +{cmd: domme} uses these constraints to "remove" parameters by constraining their +value to 0. This mimics {cmd: domin}'s method where the parameter's name is +removed from the {cmd: indepvars} list directly. {pstd} -such parameters would be produced by a model like {cmd:glm price mpg turn trunk foreign} +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}. -{pstd}Note that the current version of {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 lists -supplied are valid parameters in the estimated model. +... brief discussion of why _cons not included? ... {marker opts}{...} {title:3. Options} -{phang}{opt reg(full_estimation_command)} refers {cmd:domme} to a command that accepts {help constraint}s, -uses {help ml} to estimate parameters, and that can produce the scalar in the {opt fitstat()} option. -{cmd:domme} is quite flexible and can be applied to any built-in or user-written -{help program}. - -{pmore}The {it:full_estcmd} is the full estimation command, not including options following the -comma, as would be submitted to Stata. The {opt reg()} option has no default and the user is -required to provide a valid statistical model. - -{phang}{opt f:itstat(returned_scalar | built_in_options)} refers {cmd:domme} to ascalar valued -model fit summary statistic used to compute all dominance statistics. The scalar in -{opt fitstat()} can be any {help return:returned}, {help ereturn:ereturned}, or other -{help scalar:scalar} produced by the estimation command in {opt reg()}. - -{pmore}In addition to fit statistics produced by the estimation command in {opt reg()}, {cmd:domme} -also allows several built-in model fit statistics to be computed using the model log-likelihood and -degrees of freedom. Four fit statistics are available using the built-in options for {cmd:domme}. -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}). - -{pmore}To instruct {cmd:domme} to compute a built-in fit statistic, supply the {opt fitstat()} option -with an empty ereturned statistic indicator (i.e., {res:e()}) and provide the three character code -for the desired fit statistic. For example, to ask {cmd:domme} to compute McFadden's pseuod-R -square as a fit statistic, type {res:fitstat(e(), mcf)}. Note that {cmd:domme} has no default -and the user is required to provide a valid fit statistic. - -{phang}{opt sets([(eqname1_set1 = parmnamelist1_set1) ... (eqnameR_set1 = parmnamelistR_set1)] ... [ ... (eqnameR_setN = parmnamelistR_setN)])} -binds together parameter estimates as a set in the all possible combinations -ensemble. Hence, all parameter estimates in a set will always appear together and are considered a -single parameter estimate in the all possible combinations ensemble. - -{pmore}{opt sets()} are generated in a way similar to that of the initial statements to {cmd:domme} -in that a series of {res:(eqname = parmnamelist)} statements must be provided and then bound together -to produce a set. Any set of {res:(eqname = paramlist)} in a single set must be bound by brackets -"{res:[]}". For example, consider again the model {cmd:glm price mpg turn trunk foreign}. To -produce two sets of parameters, one that includes {it:mpg} and {it:turn} as well as a second that -includes {it:trunk} and {it:foreign}, the {opt sets()} type {res:sets( [(price = mpg turn)]} -{res:[(price = trunk foreign)] )}. - -{pmore}This above {opt sets()} statement is rather simple and refers to single equations within -a model. A single set can include parameters from multiple equations - in fact, doing so is -how independent variable dominance statistics can be computed in {cmd:domme} ... note re: independent variables ... - -{phang}{opt all((eqname1_all = parmnamelist1_all) ... (eqnameR_all = parmnamelistR_all))} -defines a set of parameter estimates to be included in all the combinations in the -ensemble. Thus, all parameter estimates included in the {opt all()} option are effectively used -as of covariates which are to be included in the model fit metric, but for which dominance -statistics will not be computed. Thus, the magnitude of the overall fit statistic associated -with the set of parameter in the {opt all()} option are subtracted from overall fit metric prior -to the computation of dominance statistics for all the parameter estimiates to be dominance -analyzed. ...note this is how to take them out of the constant... - -{pmore}The {opt all()} statements are set up in a way identical to that of the initial statments -in a {res:(eqname = parmnamelist)} format. - -{phang}{opt ropts(command_options)} supplies the command in {opt reg()} with any relevant estimation options. -Any options formally following the comma in standard Stata syntax, besides {opt constraints()}, -can be supplied to the statisical model this way. +{dlgtab:Model} + +{phang}{opt reg(full_estimation_command)} refers {cmd:domme} to a command that accepts +{help constraint}s, uses {help ml} to estimate parameters, and that can produce the +scalar-valued statistic referenced in the {opt fitstat()} option. {cmd:domme} can be +applied to any built-in or user-written {help program} that meets these criteria. + +{pmore}{it:full_estimation_command} is the full estimation command, not including +a comma or options following the comma, as would be submitted to Stata. The +{opt reg()} option has no default and the user is required to provide a +{cmd:domme}-compatible statistical model. + +{phang}{opt ropts(command_options)} supplies the command in {opt reg()} with any +relevant estimation options. Any options normally following the comma in standard +Stata syntax can be supplied to the statisical model this way. The only exception +to is the use of {opt constraints()}; {cmd:domme} cannot, at current, accept +constraints other than those it creates. + +{phang}{opt f:itstat(fitstat_opts)} the scalar-valued model fit summary statistic +used in the dominance analysis. There are two ways {cmd:domme} points to fit +statistics. + +{pmore}The first method is identical to {cmd:domin}'s approach. {cmd:domme} +accepts any {help return:returned}, {help ereturn:ereturned}, or other +{help scalar:scalar} produced by the estimation command in {opt reg()}. Note +that some Stata commands change their list of ereturn-ed results when +constraints are applied (e.g., {cmd:logit}, {cmd:poisson}). Ensure that the +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). + +{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)}. + +{phang}{opt sets([PEset_1] ... [PEset_x])} binds together parameter estimate +constraints as a set that are always constrained jointly and act as a +single parameter estimate. + +{pmore}Each {it:PEset} is put together in the same way as the initial statements in +that they are constructed from a series of {res:(eqname = parmnamelist)} statements. +All {it:PEset}s must be bound by brackets "{res:[]}". For example, consider again +the model {cmd:logit foreign price mpg turn trunk}. To produce two sets of parameters, +one that includes {it:price} and {it:mpg} as well as a second that includes {it:turn} +and {it:trunk}, the {opt sets()} type {res:sets( [(foreign = price mpg)]} +{res:[(foreign = turn trunk)] )}. + +{pmore}Note that a single set can include parameters from multiple equations +(see Example #6). + +{phang}{opt all(PEall)} defines a set of parameter estimate constraints that +are allowed to explain the fit metric with a higher priority than the parameter +estimates in the initial statements or the {opt sets()} option (see Example #3). +In effect, the parameter estimates defined in the {opt all()} option are used +like covariates. + +{pmore}The {it:PEall} statement is set up in a way similar to the {it:PE_set}s +in a {res:(eqname = parmnamelist)} format and can accept parameters from multiple +equations. + +{dlgtab:Reporting} {phang}{opt noconditional} suppresses the computation and display of of the conditional dominance -statistics. Suppressing the computation of the conditional dominance statistics can save -computation time when conditional dominance statistics are not desired. Suppressing the -computation of conditional dominance statistics also suppresses the -"strongest dominance designations" list. +statistics which can save computation time when conditional dominance statistics +are not desired. Suppressing the computation of conditional dominance statistics +also suppresses the "strongest dominance designations" list. -{phang}{opt nocomplete} suppresses the computation of the complete dominance designations. -Suppressing the computation of the complete dominancedesignations can save computation time -when complete dominance designations are not desired. Suppressing the computation of -complete dominance designations also suppresses the "strongest dominance designations" list. +{phang}{opt nocomplete} suppresses the computation of the complete dominance designations +which can save computation time when complete dominance designations are not desired. +Suppressing the computation of complete dominance designations also suppresses +the "strongest dominance designations" list. {phang}{opt reverse} reverses the interpretation of all dominance statistics in the -{cmd:e(ranking)} vector, {cmd:e(cptdom)} matrix, fixes the computation of the -{cmd:e(std)} vector, and the "strongest dominance designations" list. -{cmd:domin} assumes by default that higher values on overall fit statistics constitute -better fit, as dominance analysis has historically been based on the explained-variance R2 metric. -However, dominance analysis can be applied to any model fit statistic -(see Azen, Budescu, & Reiser, 2001 for other examples). {opt reverse} is then useful -for the interpetation of dominance statistics based on overall model fit statistics -that decrease with better fit (e.g., the built in AIC, BIC statistics). +{cmd:e(ranking)} vector, {cmd:e(cptdom)} matrix, and corrects the computation of the +{cmd:e(std)} vector as well as the "strongest dominance designations" list. +{cmd:domme} assumes by default that higher values on overall fit statistics constitute +better fit. {opt reverse} is useful for the interpetation of dominance statistics +based on overall model fit statistics that decrease with better fit (e.g., AIC, BIC). + +{marker sav}{...} +{title:4. Saved Results} + +{phang}{cmd:domme} saves the following results to {cmd: e()}: + +{synoptset 16 tabbed}{...} +{p2col 5 15 19 2: scalars}{p_end} +{synopt:{cmd:e(N)}}number of observations{p_end} +{synopt:{cmd:e(fitstat_o)}}overall fit statistic value{p_end} +{synopt:{cmd:e(fitstat_a)}}fit statistic value associated with variables in {opt all()}{p_end} +{synopt:{cmd:e(fitstat_c)}}fit statistic value computed by default when the constant model is non-zero{p_end} +{p2col 5 15 19 2: macros}{p_end} +{synopt:{cmd:e(cmdline)}}command as typed{p_end} +{synopt:{cmd:e(title)}}{cmd:Dominance analysis for multiple equations}{p_end} +{synopt:{cmd:e(cmd)}}{cmd:domme}{p_end} +{synopt:{cmd:e(fitstat)}}contents of the {opt fitstat()} option{p_end} +{synopt:{cmd:e(reg)}}contents of the {opt reg()} option{p_end} +{synopt:{cmd:e(ropts)}}contents of the {opt ropts()} option{p_end} +{synopt:{cmd:e(properties)}}{cmd:b}{p_end} +{synopt:{cmd:e(set{it:#})}}parameters included in {opt set(#)}{p_end} +{synopt:{cmd:e(all)}}parameters included in {opt all()}{p_end} +{p2col 5 15 19 2: matrices}{p_end} +{synopt:{cmd:e(b)}}general dominance statistics vector{p_end} +{synopt:{cmd:e(std)}}general dominance standardized statistics vector{p_end} +{synopt:{cmd:e(ranking)}}rank ordering based on general dominance statistics vector{p_end} +{synopt:{cmd:e(cdldom)}}conditional dominance statistics matrix{p_end} +{synopt:{cmd:e(cptdom)}}complete dominance designation matrix{p_end} +{p2col 5 15 19 2: functions}{p_end} +{synopt:{cmd:e(sample)}}marks estimation sample{p_end} -{marker remark}{...} -{title:4. Final Remarks} - -{pstd}Any parameter estimates in the model but not in either the initial syntax for -to dominance analyze, the {opt sets()} option, or the {opt all()} option are assumed to act as -a part of the model constant and are in some cases ignored in computing the model fit statistic. -Thus, any parameter estimate not included in some modeling statement will be treated like a -regression constant in most regression models; that is, as a baseline against which the -full model is compared in terms of the log likelihood. "Constant" parameters are omitted entirely -from fit statistic computations for the built in {res:mcf} and {res:est} options but are -reported as a part of the constant model fit statistic for the {res:aic} and {res:bic} -options in {opt fitstat()}. - -{pstd}When not using the built-in options, it is the responsibility of the user to supply -{cmd:domme} with an overall fit statistic that can be validly dominance analyzed. Non-R2 -overall fit statistics can be used however {cmd:domme} assumes that the fit statistic supplied -{it:acts} like an R2 statistic. Thus, {cmd:domin} assumes that better model fit is associated -with increases to the fit statistic and all marginal contributions can be obtained by subtraction. -For model fit statistics that decrease with better fit (i.e., AIC, BIC, deviance), the -interpretation of the dominance relationships need to be reversed (see Example #2). {marker examp}{...} {title:5. Examples} @@ -229,42 +297,54 @@ interpretation of the dominance relationships need to be reversed (see Example # {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 2.race 3.race), reg(gsem (wage <- union hours, regress) (south <- age ib1.race union, logit)) fitstat(e(), mcf)}{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}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)]) } {p_end} -{marker sav}{...} -{title:6. Saved Results} +{phang}Examples 7: Replicating results from {cmd:domin}{p_end} +{pmore}7a: Logit model with factor varaible{p_end} +{pmore} {stata sysuse auto, clear} {p_end} +{pmore} {stata domin foreign price mpg turn trunk, reg(logit) fitstat(e(r2_p)) sets((i.rep78))} {p_end} +{pmore} {stata domme (foreign = price mpg turn trunk), reg(logit foreign price mpg turn trunk ib1.rep78) fitstat(e(), mcf) sets([(foreign = 3.rep78 4.rep78)])} {p_end} -{phang}{cmd:domme} saves the following results to {cmd: e()}: +{pmore}7b: Ordered logit model with covariate{p_end} +{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} + +{marker remark}{...} +{title:6. Final Remarks} + +{pstd}See {stata help domin:domin's help file} for an extensive discussion of the role of dominance +analysis as a postestimation method and caveats about its use. All these notes +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. + +{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 +parameter estimates created in the initial statements, as well as those created by +{opt sets()} and {opt all()}, are valid parameters in the estimated model. +{cmd:domme} also attempts to clean-up parameter constraints it creates but +under certain circumstances, when {cmd:domme} fails to execute in a way that it does +not capture, parameter constraints will remain in memory. {cmd:domme} also will never +overwrite existing parameter estimate constraints and, if there are insufficient parameter +constraints in memory, {cmd:domme} will fail with an error noting insufficient free +constraints. Use {help constraint dir} to list all defined constraints in memory. -{synoptset 16 tabbed}{...} -{p2col 5 15 19 2: scalars}{p_end} -{synopt:{cmd:e(N)}}number of observations{p_end} -{synopt:{cmd:e(fitstat_o)}}overall fit statistic value{p_end} -{synopt:{cmd:e(fitstat_a)}}fit statistic value associated with variables in {opt all()}{p_end} -{synopt:{cmd:e(fitstat_c)}}fit statistic value computed by default when the constant model is non-zero{p_end} -{p2col 5 15 19 2: macros}{p_end} -{synopt:{cmd:e(cmdline)}}command as typed{p_end} -{synopt:{cmd:e(title)}}{cmd:Dominance analysis for multiple equations}{p_end} -{synopt:{cmd:e(cmd)}}{cmd:domme}{p_end} -{synopt:{cmd:e(fitstat)}}contents of the {opt fitstat()} option{p_end} -{synopt:{cmd:e(reg)}}contents of the {opt reg()} option{p_end} -{synopt:{cmd:e(ropts)}}contents of the {opt ropts()} option{p_end} -{synopt:{cmd:e(properties)}}{cmd:b}{p_end} -{synopt:{cmd:e(set{it:#})}}parameters included in {opt set(#)}{p_end} -{synopt:{cmd:e(all)}}parameters included in {opt all()}{p_end} -{p2col 5 15 19 2: matrices}{p_end} -{synopt:{cmd:e(b)}}general dominance statistics vector{p_end} -{synopt:{cmd:e(std)}}general dominance standardized statistics vector{p_end} -{synopt:{cmd:e(ranking)}}rank ordering based on general dominance statistics vector{p_end} -{synopt:{cmd:e(cdldom)}}conditional dominance statistics matrix{p_end} -{synopt:{cmd:e(cptdom)}}complete dominance designation matrix{p_end} -{p2col 5 15 19 2: functions}{p_end} -{synopt:{cmd:e(sample)}}marks estimation sample{p_end} {marker refs}{...} {title:7. References} @@ -272,10 +352,34 @@ interpretation of the dominance relationships need to be reversed (see Example # {p 4 8 2}Luchman, J. N., Lei, X., and Kaplan, S. A. (2020). Relative importance analysis with multivariate models: Shifting the focus from independent variables to parameter estimates. {it:Journal of Applied Structural Equation Modeling, 4(2)}, 40–59.{p_end} +{title:Development Webpage} + +{phang} Additional discussion of results, options, and conceptual issues on: + +{phang}{browse "http://github.com/jluchman/domme/blob/master/README.md"} + +{phang} Please report bugs, requests for features, and contribute to as well as follow on-going development of {cmd:domme} on: + +{phang}{browse "http://github.com/jluchman/domme"} + +{title:Article} + +Please cite as: + +{p 4 8 2}Luchman, J. N. (2021). Determining relative importance in Stata using dominance analysis: domin and domme. {it:The Stata Journal, 21(2)}, 510–538. https://doi.org/10.1177/1536867X211025837{p_end} + + {title:Author} {p 4}Joseph N. Luchman{p_end} -{p 4}Senior Scientist{p_end} -{p 4}Fors Marsh Group LLC{p_end} +{p 4}Principal Scientist{p_end} +{p 4}Fors Marsh{p_end} {p 4}Arlington, VA{p_end} -{p 4}jluchman@forsmarshgroup.com{p_end} +{p 4}jluchman@forsmarsh.com{p_end} + + +{title:See Also} + +{browse "https://CRAN.R-project.org/package=domir":R package domir}, {browse "https://cran.r-project.org/web/packages/domir/vignettes/domir_basics.html":Detailed description of Dominance Analysis} + +