Re: st: marginal effects for user defined programs

[May Boggess <mboggess@stata.com>]

On Wednesday morning Steve asked about defining your own
predict for a user-written command:

> What would the mypredict program look like if we want to generalise 
> and have be able to get the marginal effects for the
> predicted probability of every outcome?

To make a predict command that can predict more than one thing,
it will need options. In my example today I am going
to do my own multinomial logit model using -ml-. Then I'll compare
my answers to that given by -mfx- following -mlogit-.
First I'll define my estimation command:

cap program drop maymlogit
program maymlogit
        args lnf theta3 theta4
        tempvar g
	gen double `g'=1+exp(`theta3')+exp(`theta4')
        qui replace `lnf'=        -ln(`g') if $ML_y1==2
        qui replace `lnf'=`theta3'-ln(`g') if $ML_y1==3
        qui replace `lnf'=`theta4'-ln(`g') if $ML_y1==4
end

 cap program drop mymlogit
 program mymlogit, eclass
 ml model lf maymlogit /*
        */ (theta3: rep = mpg for) /*
        */ (theta4: rep = mpg for)
 ml max
 eret local cmd="mymlogit"
 eret local predict="mypredict"
 end

Now I need to define mypredict program:

 cap program drop mypredict
 program mypredict
 syntax newvarlist(min=1 max=1) [if] [in], [ p2 p3 p4]
  tempvar xb3 xb4 touse
  mark `touse'
  _predict double `xb3' if `touse', xb equation(#1)
  _predict double `xb4' if `touse', xb equation(#2)
  if "`p2'"!=""{
    gen double `varlist'=1/(1+exp(`xb3')+exp(`xb4'))
  }
  else if "`p3'"!=""{
    gen double `varlist'=exp(`xb3')/(1+exp(`xb3')+exp(`xb4'))
  }
  else if "`p4'"!=""{
    gen double `varlist'=exp(`xb4')/(1+exp(`xb3')+exp(`xb4'))
  }
  label var `varlist'
 end

As Steve said, I am using the fact that I know the categories of
my response variable and what the base category is going to be.
It is possible of course to write a more general
estimator and predict, but if I know how many cateories that are, then
it is certainly much easier to code something for that specific case.

There are a couple of things to note about my predict that are important
when using -mfx-. First is that I am using type double. This is 
crucially important for numerical accuracy. If you take the double out 
and let it use float, you'll see a noticable reduction in accuracy
of the derivatives. Second is the [if] and [in]. Inside the
.ado code for -mfx- is uses if and in, so if you don't have them there
it will definitely error out.

Now I am ready to run my example.

clear
sysuse auto
replace rep=2 if rep==1 | rep==5| rep==.
mymlogit
mfx, predict(p2) nonlinear
mfx, predict(p3) nonlinear
mfx, predict(p4) nonlinear

mlogit rep mpg for, base(2)
mfx, predict(outcome(2))
mfx, predict(outcome(3))
mfx, predict(outcome(4))


You'll notice I used the option -nonlinear- on -mfx-.
This is because the default method for obtaining derivatives
in -mfx- is the linear method. After an offical Stata estimation
command, this is fine as the default, since we know that the information
that -mfx- is relying on the detect whether or not the linear
method is still appropriate is all there. After a user-written command
though, it is wise to play safe and use the nonlinear method.
If the linear method is appropriate, there will be no difference between
the answers given by the two methods. The only difference would
be that the linear method will be faster on models with
many variables.


-- May
mmb@stata.com


*
*   For searches and help try:
*   http://www.stata.com/support/faqs/res/findit.html
*   http://www.stata.com/support/statalist/faq
*   http://www.ats.ucla.edu/stat/stata/