Statalist


[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

st: AW: swboot.ado question


From   "Martin Weiss" <martin.weiss1@gmx.de>
To   <statalist@hsphsun2.harvard.edu>
Subject   st: AW: swboot.ado question
Date   Thu, 26 Feb 2009 19:58:21 +0100

<> 


Well, you have to pass the "options" for -stepwise- as "options" for
swboot.ado (findit swboot if anyone wants to weigh in. Just pasting the code
does not replace proper source...). So do not try to impose the syntax
diagram for the -sw- prefix on a user-written -program-. See the example:


************
sysuse auto, clear
swboot  mpg price rep78 headroom trunk weight length turn displacement
gear_ratio foreign, reps(50) pe(0.05) pr(0.1)
// this leads to an error
swboot, reps(50) pe(0.05) pr(0.1)  mpg price rep78 headroom trunk weight
length turn displacement gear_ratio foreign
************



HTH
Martin

-----Ursprüngliche Nachricht-----
Von: owner-statalist@hsphsun2.harvard.edu
[mailto:owner-statalist@hsphsun2.harvard.edu] Im Auftrag von ALEXANDER J
SHACKMAN
Gesendet: Donnerstag, 26. Februar 2009 18:26
An: statalist@hsphsun2.harvard.edu
Cc: ajshackman@gmail.com
Betreff: st: swboot.ado question

Hi,

I'm new to Stata, so apologies for asking what is probably a naive question
about swboot.ado (designed to permit bootstraps of forward/backward
regressions; appended below).

Basically, I want to run: 

stepwise, pe(.05) pr(.20) forward lockterm1: regress anxtemp (age
female0male) peak_1 peak_2 peak_3 ... peak_157

I have verified that this command works properly, and now I want to adapt
swboot.ado to use it.

The key part of the .ado seems to be this bit:

       if "`forward'"=="forward" {
         quietly sw regress `yvar' `xvars', pe(`pe') pr(`pr') forward
         }
 
Unfortunately, you cannot just replace that with: 

         quietly sw, pe(`pe') pr(`pr') forward lockterm1: regress `yvar'
`xvars'

That doesn't work. Any suggestions for a newcomer?

Thanks,
Alex Shackman

----------------- swboot.ado ----------------------------

*! version 2.2.0 JMGarrett 13Aug00
/* Bootstrap sw linear or logistic regression; count variables selected  */
/* Form:  swboot yvar xvarlist, reps(#) options                          */
/* Options allowed:  reps(#), pe(#), pr(#), forward, N(#) model roc gof  */

program define swboot
  version 6.0
  #delimit ;
    syntax varlist (min=2) [if] [in], [Reps(int 1) N(int 0) PE(real .05)
    PR(real .1) FORward MODel ROC GOF] ;
  #delimit cr
  marksample touse
  markout `touse'
  preserve
  quietly keep if `touse'
  keep `varlist'
  local yvar="`1'"
  local varlbly : variable label `yvar'
  capture assert `yvar'==1 | `yvar'==0
  if _rc==0 {local regtype "log"}
  if _rc~=0 {local regtype "lin"}

* Read in x variables
  tokenize "`varlist'"
  local i=1
  while "`2'"~=""  {
    local x`i'="`2'"
    * quietly drop if `x`i''==.  /* do not need this? */
    local xvars `xvars' `x`i''
    local count`i'=0
    local i=`i'+1
    macro shift
    }
  local numx=`i'-1
  if `n'~=0 {local sample=`n'}
  if `n'==0 {local sample=_N}

* select bootstrap sample and run stepwise regressions
  if "`forward'"=="forward" {local forprnt="(Forward Selection)"}
  local numreps=1
  tempfile tempds
  disp "  "
  if "`regtype'"=="lin" {
    #delimit ;
       disp _n(1) in gre "Model:" in gre "    Stepwise " in blue "Linear "
         in gre "Regression `forprnt'" ;
       disp in gre "Outcome:" in gre "  `varlbly' -- " in yel "`yvar'";
       disp _n(1) in gre  "Options:  pr=" in yel "`pr'" in green "  pe="
          in yel "`pe'" in green "  n=" in yel "`sample'" ;
    #delimit cr
    disp "  "
    }
  if "`regtype'"=="log" {
    #delimit ;
       disp _n(1) in gre "Model:" in gre "    Stepwise " in blue "Logistic "
         in gre "Regression `forprnt'" ;
       disp in gre "Outcome:" in gre "  `varlbly' -- " in yel "`yvar'";
       disp _n(1) in gre  "Options:  pr=" in yel "`pr'" in green "  pe="
          in yel "`pe'" in green "  n=" in yel "`sample'" ;
    #delimit cr
    disp "  "
    }
    if `n'>_N  {
      disp in red "There are only " in yel _N in red " available " /*
        */ "observations for this model;"
      }
  while `numreps'<=`reps'  {
    quietly save `tempds', replace
    bsample `sample'
    if "`regtype'"=="lin" {
       if "`forward'"=="" {
         quietly sw regress `yvar' `xvars', pe(`pe') pr(`pr')
         }
       if "`forward'"=="forward" {
         quietly sw regress `yvar' `xvars', pe(`pe') pr(`pr') forward
         }
       }
    if "`regtype'"=="log" {
       if "`forward'"=="" {
         quietly sw logistic `yvar' `xvars', pe(`pe') pr(`pr')
         }
       if "`forward'"=="forward" {
         quietly sw logistic `yvar' `xvars', pe(`pe') pr(`pr') forward
         }
       if "`roc'"=="roc" {
          quietly lroc, nograph
          local arearoc=round(r(area),.0001)
          }
       if "`gof'"=="gof" {
          quietly lfit, group(10)
          local pval=round(chiprob(r(df),r(chi2)),.0001)
          }
       }
    mat b=e(b)
    local keepX : colnames(b)
    tokenize `keepX'
    local i 1
    local xlist ""
    while "`1'"~="_cons"  {
      local keepx`i'="`1'"
      local xlist `xlist' `keepx`i''
      local i=`i'+1
      macro shift
      }
    local numkeep=`i'-1

    * If model is requested, create dataset to hold and print estimates
    if "`model'"=="model"  {
      drop _all
      mat b=b'
      mat v=e(V)
      mat v=vecdiag(v)
      mat v=v'
      mat bv=b,v
      set more off
      quietly svmat bv
      set more on
      quietly gen str8 varname=" "
      quietly drop if _n==_N
      quietly gen se=sqrt(bv2)
      quietly rename bv1 coef
      local i 1
        while `i'<=`numkeep' {
          quietly replace varname="`keepx`i''" if _n==`i'
          local i=`i'+1
          }
      if "`regtype'"=="lin"  {
         quietly gen t=coef/se
         quietly gen p=tprob(`sample'-2,t)
         format coef se p t %10.4f
         disp " "
         disp in blue "Rep `numreps':"
         list varname coef se p, noob 
        }
      if "`regtype'"=="log"  {
         quietly gen OR=exp(coef)
         quietly gen lower=exp(coef-1.96*se)
         quietly gen upper=exp(coef+1.96*se)
         quietly gen z=abs(coef/se)
         quietly gen p=2*(1-normprob(z))
         format OR lower upper %6.2f
         format se p %7.4f
         disp _n(1) in blue "Rep `numreps':"
         list varname OR se p lower upper, noob
         disp " "
         if "`roc'"=="roc" {
           disp in gr "      Area under ROC curve     = " in yel `arearoc'
           }
         if "`gof'"=="gof" {
           disp in gr "      Goodness-of-fit test:  p = " in yel `pval'
           }
         }
      }

    * If model is not selected, print a summary of variables for each rep
    if "`model'"==""  {
       if "`roc'"=="roc" | "`gof'"=="gof"  {disp "  "}
       disp in gr "Rep `numreps':" in yel "  `xlist'"
       if "`roc'"=="roc" & "`regtype'"=="log" {
          disp in gr "          Area under ROC curve   =   " in yel
`arearoc'
          }
       if "`gof'"=="gof" {
          disp _s(1) in gr "         Goodness-of-fit test:  p = " in yel
`pval'
          }
       }

    * Count the number of times each variable is selected
    local i 1
    while `i'<=`numx'  {
      local j 1
      while `j'<=`numkeep'  {
        if "`x`i''"=="`keepx`j''"  {local count`i'=`count`i''+1}
        local j=`j'+1
        }
      local i=`i'+1
      }
    quietly use `tempds', clear
    local numreps=`numreps'+1
    }
  
* if last rep
  local i=1
  disp _n(2) in green "Summary: (Number of times each variable is selected)"
  disp in green "----------------------------------------------------"
  while `i'<=`numx' {
    disp "   `x`i'':" _col(15) "`count`i''"
    local i=`i'+1
    }
end


-- 
Alexander J. Shackman, Ph.D.
Laboratory for Affective Neuroscience
Waisman Laboratory for Brain Imaging & Behavior
University of Wisconsin-Madison
1202 West Johnson Street
Madison, Wisconsin 53706

Telephone: +1 (608) 358-5025
Fax: +1 (608) 265-2875
Email: shackman@wisc.edu
http://psyphz.psych.wisc.edu/~shackman
*
*   For searches and help try:
*   http://www.stata.com/help.cgi?search
*   http://www.stata.com/support/statalist/faq
*   http://www.ats.ucla.edu/stat/stata/


*
*   For searches and help try:
*   http://www.stata.com/help.cgi?search
*   http://www.stata.com/support/statalist/faq
*   http://www.ats.ucla.edu/stat/stata/



© Copyright 1996–2014 StataCorp LP   |   Terms of use   |   Privacy   |   Contact us   |   What's new   |   Site index