*! version 1.0.0 PR 14jan2009 * Based on Stata 10 stepwise.ado version 6.1.1 16jul2007 program define mimstepwise, eclass byable(onecall) version 9 local version : di "version " string(_caller()) ", missing:" if _by() { local BY `"by `_byvars'`_byrc0':"' } if _caller() < 9 { di as err "not supported for versions earlier than 9" exit 198 * `version' `BY' sw_8 `0' * exit } capture _on_colon_parse `0' if c(rc) | `"`s(after)'"' == "" { gettoken old : 0, parse(" ,") if "`old'" != "," { capture which `old' if !c(rc) { `version' `BY' sw_8 `0' exit } } if (_by()) error 190 // !! PR: line replaced by mim equivalent if ("`e(cmd)'" == "") error 301 if ("`e(MIM_prefix2)'" == "") error 301 Display `0' exit } `version' `BY' StepWise `0' ereturn local cmdline `"mimstepwise `0'"' end program Display * local cmd = cond("`e(cmd2)'" == "", "`e(cmd)'", "`e(cmd2)'") mim `0' end program StepWise, eclass byable(recall) version 9 local version : di "version " string(_caller()) ":" // : _on_colon_parse `0' local command `"`s(after)'"' local 0 `"`s(before)'"' // !! PR: lr option disabled syntax [fw iw pw aw] [if] [in] [, /// pr(numlist max=1 >0 <1) /// pe(numlist max=1 >0 <1) /// FORWard /// HIERarchical /// LOCkterm1 /// Level(passthru) /// ] local hier `hierarchical' if "`weight'" != "" { local wgt [`weight'`exp'] } if `"`pr'`pe'"' == "" { di as err "at least one of pr(#) and pe(#) must be specified" exit 198 } if "`pr'" != "" & "`pe'" != "" { if "`hier'" != "" { di as err "hierarchical pe(#) pr(#) invalid" exit 198 } if `pr' <= `pe' { di as err "pr(`pr') <= pe(`pe') invalid" exit 198 } } if "`forward'" != "" & "`pe'" == "" { di as err "option forward invalid or option pe(#) is missing" exit 198 } // parse the command and check for conflicts `version' _prefix_command stepwise /// `wgt' `if' `in' , `level' : `command' local version `"`s(version)'"' local cmdname `"`s(cmdname)'"' local termlist `"`s(anything0)'"' local wgt `"`s(wgt)'"' local weight `"`s(wtype)'"' local exp `"`s(wexp)'"' local if `"`s(if)'"' local in `"`s(in)'"' local cmdopts `"`s(options)'"' local rest `"`s(rest)'"' CheckProps `cmdname' `lr' Check4Robust "`lr'" `wgt', `cmdopts' local constant `"`s(constant)'"' // Display options local diopts `"`s(efopt)' level(`s(level)')"' // build up part of the command that follows the varlist local cmdrest `"`wgt'"' if `"`cmdopts'"' != "" { local cmdrest `"`cmdrest', `cmdopts'"' } if `"`rest'"' != "" { local cmdrest `"`cmdrest' `rest'"' } // note: marksample looks at `varlist' `if' `in' and [`weight'`exp'] // and takes into account if I am being called with the -by:- prefix marksample touse, novarlist // !! PR: remove _mj==0 from model quietly replace `touse' = 0 if _mj == 0 // parse termlist, and build the list of variable names ParseTerms depvar termlist /// : `touse' `cmdname' "`wgt'" "`constant'" /// : `termlist' TermMacros term k curr : `lockterm1' : `termlist' // !! PR: fit full model using mim (was `version' `cmdname' `depvar' `curr' if `touse' `cmdrest') quietly `version' mim, storebv:`cmdname' `depvar' `curr' if `touse' `cmdrest' local colna : colname e(b) local varlist : subinstr local colna "_cons" "", all local dropped : list curr - varlist // check for missing standard errors tempname v matrix `v' = e(V) if diag0cnt(`v') { local dim = colsof(`v') forval i = 1/`dim' { local vname : word `i' of `colna' if `v'[`i',`i'] == 0 & "`vname'" != "_cons" { local dropped `dropped' `vname' } } } // check for dropped variables if "`dropped'" != "" { gettoken term termlist : termlist, bind match(par) while `"`term'"' != "" { local term : list term - dropped if `"`term'"' != "" { local tlist `"`tlist' (`term')"' } gettoken term termlist : termlist, bind match(par) } local termlist `"`tlist'"' foreach dvar of local dropped { di as txt "note: `dvar' dropped because of estimability" } } // check for dropped observations quietly count if `touse' & !e(sample) if r(N) { di as txt "note: `r(N)' obs. dropped because of estimability" quietly replace `touse' = 0 if !e(sample) } // !! PR: `cmdname' is the Stata command name. Needs to be prefixed with mim:. local cmdname mim:`cmdname' if "`pr'" != "" { if "`pe'" == "" { BackSel `pr' /// "`hier'" /// "`lr'" /// "`lockterm1'" /// "`version'" /// `cmdname' /// "`depvar'" /// "`termlist'" /// `touse' /// "`cmdrest'" } else if "`forward'" == "" { BackStep `pr' /// `pe' /// "`lr'" /// "`lockterm1'" /// "`version'" /// `cmdname' /// "`depvar'" /// "`termlist'" /// `touse' /// "`cmdrest'" } else { ForStep `pr' /// `pe' /// "`lr'" /// "`lockterm1'" /// "`version'" /// `cmdname' /// "`depvar'" /// "`termlist'" /// `touse' /// "`cmdrest'" } } else { ForSel `pe' /// "`hier'" /// "`lr'" /// "`lockterm1'" /// "`version'" /// `cmdname' /// "`depvar'" /// "`termlist'" /// `touse' /// "`cmdrest'" } ereturn local stepwise stepwise Display, `diopts' end program BackSel args pr hier lr lock version cmdname depvar termlist touse cmdrest TermMacros term k curr : `lock' : `termlist' if "`lr'" != "" { local LRtest "LR test" } di in smcl as txt "{p2colset 0 23 32 2}{...}" di as txt "{p2col :`LRtest'}begin with full model{p_end}" // note: the "full" model has already been fit, and e() should already // contain it's results if "`lr'" != "" { tempname ll1 df1 ll0 df0 GetLR `ll1' `df1' } local prfmt : display %6.4f `pr' local done 0 local k0 `k' while !`done' { local drop_i 0 local p_i 0 if "`hier'" != "" & `k' > 0 { local start `k' } else local start 1 forval i = `start'/`k' { if "`lr'" != "" { local curr1 : subinstr local curr "`term`i''" "", word quietly `version' `cmdname' `depvar' `curr1' if `touse' `cmdrest' quietly LRtest `ll1' `df1' } else quietly mim: testparm `term`i'' if r(p) > `p_i' { local p_i = r(p) local drop_i `i' if "`lr'" != "" { GetLR `ll0' `df0' } } } if `p_i' >= `pr' { local p : display %6.4f `p_i' di as txt /// "{p2col :p = {res:`p'} >= `prfmt'}removing {res:`term`drop_i''}{p_end}" local newcurr : list curr - term`k' local curr : /// subinstr local newcurr "`term`drop_i''" "`term`k''", word local term`drop_i' `term`k'' local term`k' local --k } else { if "`hier'" != "" & `k' > 0 { local p : display %6.4f `p_i' di as txt /// "{p2col :p = {res:`p'} <{space 2}`prfmt'}" /// "keeping{space 2}{res:`term`drop_i''}{p_end}" } local done 1 } if "`lr'" != "" { scalar `ll1' = `ll0' scalar `df1' = `df0' } else { quietly `version' `cmdname' `depvar' `curr' if `touse' `cmdrest' } } if "`lr'" != "" { quietly `version' `cmdname' `depvar' `curr' if `touse' `cmdrest' } if `k' == `k0' { di as txt /// "{p2col :p < `prfmt'}for all terms in model{p_end}" } di in smcl as txt "{p2colreset}{...}" end program BackStep args pr pe lr lock version cmdname depvar termlist touse cmdrest TermMacros term k curr : `lock' : `termlist' if "`lr'" != "" { local LRtest "LR test" } di in smcl as txt "{p2colset 0 23 32 2}{...}" di as txt "{p2col :`LRtest'}begin with full model{p_end}" // note: the "full" model has already been fit, and e() should already // contain it's results if "`lr'" != "" { tempname ll1 df1 ll0 df0 GetLR `ll1' `df1' } local prfmt : display %6.4f `pr' local pefmt : display %6.4f `pe' local done 0 local k0 `k' local kp1 = `k' + 1 local first 1 while !`done' { // backward step local drop_i 0 local p_i 0 forval i = 1/`k' { if "`lr'" != "" { local curr1 : subinstr local curr "`term`i''" "", word quietly `version' `cmdname' `depvar' `curr1' if `touse' `cmdrest' quietly LRtest `ll1' `df1' } else quietly mim: testparm `term`i'' if r(p) > `p_i' { local p_i = r(p) local drop_i `i' if "`lr'" != "" { GetLR `ll0' `df0' } } } if `p_i' >= `pr' { local p : display %6.4f `p_i' di as txt /// "{p2col :p = {res:`p'} >= `prfmt'}removing {res:`term`drop_i''}{p_end}" local newcurr : list curr - term`k' local curr : /// subinstr local newcurr "`term`drop_i''" "`term`k''", word local tmp `term`k'' local term`k' `term`drop_i'' local term`drop_i' `tmp' local --k local --kp1 if `first' & "`lr'" != "" { scalar `ll1' = `ll0' scalar `df1' = `df0' } } else { local done 1 if "`lr'" != "" { scalar `ll0' = `ll1' scalar `df0' = `df1' } } if `first' { if "`lr'" == "" { quietly `version' `cmdname' `depvar' `curr' if `touse' `cmdrest' } local first 0 continue } // forward step local add_i 0 local p_i 2 forval i = `kp1'/`k0' { quietly `version' `cmdname' `depvar' `curr' `term`i'' if `touse' `cmdrest' if "`lr'" != "" { quietly LRtest `ll0' `df0' } else quietly mim: testparm `term`i'' if r(p) < `p_i' { local p_i = r(p) local add_i `i' if "`lr'" != "" { GetLR `ll1' `df1' } } } if `p_i' < `pe' { local p : display %6.4f `p_i' di as txt /// "{p2col :p = {res:`p'} <{space 2}`pefmt'}" /// "adding{space 3}{res:`term`add_i''}{p_end}" local curr `curr' `term`add_i'' if "`hier'" != "" { local ++start } else { local tmp `term`kp1'' local term`kp1' `term`add_i'' local term`add_i' `tmp' local ++k local ++kp1 } local done 0 } else { if "`lr'" != "" { scalar `ll1' = `ll0' scalar `df1' = `df0' } } if "`lr'" == "" { quietly `version' `cmdname' `depvar' `curr' if `touse' `cmdrest' } } quietly `version' `cmdname' `depvar' `curr' if `touse' `cmdrest' if `k' == `k0' { di as txt /// "{p2col :p < `prfmt'}for all terms in model{p_end}" } di in smcl as txt "{p2colreset}{...}" end program ForSel args pe hier lr lock version cmdname depvar termlist touse cmdrest TermMacros term k curr : `lock' : `termlist' local curr `term0' // is empty if `"`lock'"' == "" if "`lr'" != "" { local LRtest "LR test" } di in smcl as txt "{p2colset 0 23 32 2}{...}" if "`lock'" == "" { di as txt "{p2col :`LRtest'}begin with empty model{p_end}" } else di as txt "{p2col :`LRtest'}begin with term 1 model{p_end}" if "`lr'" != "" { // fit the empty model quietly `version' `cmdname' `depvar' if `touse' `cmdrest' tempname ll0 df0 ll1 df1 GetLR `ll0' `df0' } local start 1 local pefmt : display %6.4f `pe' local done 0 while !`done' { local add_i 0 local p_i 2 if "`hier'" != "" { local last `start' } else local last `k' forval i = `start'/`last' { quietly `version' `cmdname' `depvar' `curr' `term`i'' if `touse' `cmdrest' if "`lr'" != "" { quietly LRtest `ll0' `df0' } else quietly mim: testparm `term`i'' if r(p) < `p_i' { local p_i = r(p) local add_i `i' if "`lr'" != "" { GetLR `ll1' `df1' } } } if `p_i' < `pe' { local p : display %6.4f `p_i' di as txt /// "{p2col :p = {res:`p'} <{space 2}`pefmt'}" /// "adding{space 2}{res:`term`add_i''}{p_end}" local curr `curr' `term`add_i'' if "`hier'" != "" { local ++start if `start' > `k' { local done 1 } } else { local term`add_i' `term`k'' local term`k' local --k } } else { if "`hier'" != "" { local p : display %6.4f `p_i' di as txt /// "{p2col :p = {res:`p'} >= `pefmt'}testing {res:`term`add_i''}{p_end}" } local done 1 } if "`lr'" != "" { scalar `ll0' = `ll1' scalar `df0' = `df1' } } if "`curr'" == "`term0'" { di as txt /// "{p2col :p >= `pefmt'}for all terms in model{p_end}" } quietly `version' `cmdname' `depvar' `curr' if `touse' `cmdrest' di in smcl as txt "{p2colreset}{...}" end program ForStep args pr pe lr lock version cmdname depvar termlist touse cmdrest TermMacros term k curr : `lock' : `termlist' local curr `term0' if "`lr'" != "" { local LRtest "LR test" } di in smcl as txt "{p2colset 0 23 32 2}{...}" if "`lock'" == "" { di as txt "{p2col :`LRtest'}begin with empty model{p_end}" } else di as txt "{p2col :`LRtest'}begin with term 1 model{p_end}" if "`lr'" != "" { // fit the empty model quietly `version' `cmdname' `depvar' `curr' if `touse' `cmdrest' tempname ll0 df0 ll1 df1 GetLR `ll0' `df0' } local prfmt : display %6.4f `pr' local pefmt : display %6.4f `pe' local done 0 local k0 `k' local k 0 local kp1 1 local first 1 while !`done' { // forward step local add_i 0 local p_i 2 forval i = `kp1'/`k0' { quietly `version' `cmdname' `depvar' `curr' `term`i'' if `touse' `cmdrest' if "`lr'" != "" { quietly LRtest `ll0' `df0' } else quietly mim: testparm `term`i'' if r(p) < `p_i' { local p_i = r(p) local add_i `i' if "`lr'" != "" { GetLR `ll1' `df1' } } } if `p_i' < `pe' { local p : display %6.4f `p_i' di as txt /// "{p2col :p = {res:`p'} <{space 2}`pefmt'}" /// "adding{space 3}{res:`term`add_i''}{p_end}" local curr `curr' `term`add_i'' if "`hier'" != "" { local ++start } else { local tmp `term`kp1'' local term`kp1' `term`add_i'' local term`add_i' `tmp' local ++k local ++kp1 } if `first' & "`lr'" != "" { scalar `ll0' = `ll1' scalar `df0' = `df1' } } else { if "`lr'" != "" { scalar `ll1' = `ll0' scalar `df1' = `df0' } local done 1 } if `first' { local first 0 continue } if "`lr'" == "" { quietly `version' `cmdname' `depvar' `curr' if `touse' `cmdrest' } // backward step local drop_i 0 local p_i 0 forval i = 1/`k' { if "`lr'" != "" { local curr1 : subinstr local curr "`term`i''" "", word quietly `version' `cmdname' `depvar' `curr1' if `touse' `cmdrest' quietly LRtest `ll1' `df1' } else quietly mim: testparm `term`i'' if r(p) > `p_i' { local p_i = r(p) local drop_i `i' if "`lr'" != "" { GetLR `ll0' `df0' } } } if `p_i' >= `pr' { local p : display %6.4f `p_i' di as txt /// "{p2col :p = {res:`p'} >= `prfmt'}removing {res:`term`drop_i''}{p_end}" local newcurr : list curr - term`k' local curr : /// subinstr local newcurr "`term`drop_i''" "`term`k''", word local tmp `term`k'' local term`k' `term`drop_i'' local term`drop_i' `tmp' local --k local --kp1 local done 0 } else { if "`lr'" != "" { scalar `ll0' = `ll1' scalar `df0' = `df1' } } } if "`curr'" == "`term0'" { di as txt /// "{p2col :p >= `pefmt'}for all terms in model{p_end}" } quietly `version' `cmdname' `depvar' `curr' if `touse' `cmdrest' di in smcl as txt "{p2colreset}{...}" end // Syntax: : "[weight]" : // // Parse the , removing observations with missing values, and // collinear variables. program ParseTerms _on_colon_parse `0' local c_depvar : word 1 of `s(before)' local c_terms : word 2 of `s(before)' _on_colon_parse `s(after)' local touse : word 1 of `s(before)' local cmdname : word 2 of `s(before)' local wgt : word 3 of `s(before)' local constant : word 4 of `s(before)' local termlist `"`s(after)'"' local stcmds streg stcox local is_st : list cmdname in stcmds // there is no for -stcox- and -streg-, // otherwise it is required if !`is_st' { gettoken depvar termlist : termlist, bind match(par) if "`par'" != "" { di as err "depvar cannot be part of a group" exit 198 } unab depvar : `depvar' if "`cmdname'" != "intreg" { markout `touse' `depvar' } if "`cmdname'" == "intreg" { gettoken depvar2 termlist : termlist, bind match(par) if "`par'" != "" { di as err "depvar cannot be part of a group" exit 198 } unab depvar2 : `depvar2' if "`cmdname'" != "intreg" { markout `touse' `depvar2' } local depvar `depvar' `depvar2' } } local k 0 while `"`:list retok termlist'"' != "" { local ++k gettoken term termlist : termlist, bind match(par) if "`par'" == "" { unab term : `term' if `:word count `term'' != 1 { gettoken term rest : term local termlist `"`rest' `termlist'"' } } markout `touse' `term' quietly _rmcoll `term' if `touse' `wgt', `constant' local term`k' `r(varlist)' local xvars0 `xvars0' `term`k'' } if `k' == 0 { di as err "no independent variables" exit 198 } quietly _rmcoll `xvars0' if `touse' `wgt', `constant' local dups0 : list dups xvars0 local xvars `r(varlist)' local colvar : list xvars0 - xvars if `"`colvar'"' != "" { local colvar : word 1 of `colvar' di as err "between-term collinearity, variable `colvar'" exit 459 } forval i = 1/`k' { local dups : list term`i' & dups0 if "`dups'" != "" { local term`i' : list term`i' - dups local dups0 : list dups0 - dups } if "`term`i''" != "" { local terms `"`terms' (`term`i'')"' } } // saved results c_local `c_depvar' `depvar' c_local `c_terms' `"`:list retok terms'"' end // Syntax: : // // Build local macros for the caller that contain each individual term, how // many terms there are, and a list of all the current term variables. program TermMacros _on_colon_parse `0' local c_stub : word 1 of `s(before)' local c_k : word 2 of `s(before)' local c_curr : word 3 of `s(before)' _on_colon_parse `s(after)' local lock `s(before)' local termlist `"`s(after)'"' local k 0 while `"`termlist'"' != "" { if "`lock'" != "" { local lock } else local ++k gettoken term termlist : termlist, bind match(par) c_local `c_stub'`k' `term' local curr `"`curr' `term'"' } c_local `c_k' `k' c_local `c_curr' `"`:list retok curr'"' end program CheckProps args cmdname lr local props : properties `cmdname' local swlist sw swml local swprop : list swlist & props if "`swprop'" == "" { di as err "`cmdname' is not supported by stepwise" exit 199 } if "`lr'" != "" { local swlist swml if `"`:list swlist & props'"' == "" { di as err /// "`cmdname' does not allow the lr option of stepwise" exit 199 } } end program Check4Robust, sclass gettoken lr 0 : 0 if ("`lr'" == "") { syntax [fw pw aw iw] [, noConstant * ] sreturn local constant `constant' exit } syntax [fw pw aw iw] [, Robust CLuster(passthru) vce(string) /// noConstant * ] local lvce : length local vce if "`weight'" == "pweight" | "`robust'`cluster'`svy'" != "" /// | `"`vce'"' == substr("robust",1,max(1,`lvce')) { if "`svy'" != "" { di as err "option `lr' is not allowed with svy" } else di as err /// "option `lr' is not allowed with pweights, robust, or cluster()" exit 198 } sreturn local constant `constant' end program LRtest, rclass args ll1 df1 tempname ll0 df0 GetLR `ll0' `df0' return scalar p = chiprob(abs(`df1' - `df0'), abs(2*(`ll1'-`ll0'))) end program GetLR args ll df capture confirm number `e(ll)' if c(rc) { di as err /// "option lr requires that e(ll) contains the " /// "log likelihood value for the model fit" exit 322 } capture { confirm integer number `e(df_m)' assert e(df_m) >= 0 } if c(rc) { di as err /// "option lr requires that e(df_m) contains the " /// "model degress of freedom" exit 322 } scalar `ll' = e(ll) scalar `df' = e(df_m) end exit