Finding Mata code

[Richard Williams]

I rarely use Mata, so this is my first ever post in this forum.

heckprobit.ado includes the lines

`vv' mata: heckprob_init("inits","`seldep'", ///
"`depvar'","`touse'")

Assuming it can be done, I would like to see the code for heckprob_init. How do I do that? I assume it would be part of a file with a .mata extension, but there seem to be zillions of those.

[Richard Williams]

OK, I think I found the command:

. mata: mata which heckprob_init
  heckprob_init():  lmataado

But, the only file I can find is lmataado.mlib. I take it that means I am out of luck with regards to seeing the source code?

[Niels Henrik Bruun]

Hi Richard
You're not the first to get stuck on that: http://www.stata.com/statalist/archi.../msg01054.html

[Niels Henrik Bruun]

The closest I could get was

mata describe using lmataado

. But that I actually rather useless, sorry

[Jeff Pitblado (StataCorp)]

The Mata source code that goes into lmataado.mlib is not public. However if Richard is curious about the purpose of
the Mata function heckprob_init(), it is a subroutine that puts information into a Mata structure that the heckprobit
likelihood evaluator needs for computing the log-likelihood value and its derivatives. The information is basically identifying
the selected and not-selected observations in the estimation sample.

logit, ologit, and mlogit all do similar things with Mata. For example, in logit.ado you will see the line

mata: mopt__pl_init("perfect")

and the source for mopt__pl_init() can be viewed via

viewsource moptimize_evaluators.mata

There you will see that mopt__pl_init() creates an external Mata object and returns its name in the local macro
named perfect. This object is just a 4 column row vector that the logit likelihood evaluator will fill with
information about completely determined outcomes. The completely determined outcomes information is posted
by Mata function mopt__pl_post(), which gives you a clue what the elements of the 4 column row vector
mean.

[Richard Williams]

Thanks Jeff. Based on this presentation by Jeffrey Wooldridge,

http://www.stata.com/meeting/chicago...wooldridge.pdf

I have written a couple of beta programs that work with fractional variables, i.e. vars that range between 0 and 1. As Wooldridge notes, this basically involved relaxing Stata's insistence that dependent vars be coded 0/1. Somebody asked if this could be done with heckprob too. I suspect it can but the hidden MAta code was causing me grief. I'll see if I can figure it out now.

[Sophie McKee]

Hello Dr. Williams, I would like to estimate a model with a fractional dependent variable with sample selection. I have tried tweaking heckprobit by replacing calls to "probit" by calls to "fracreg probit" when applicable, but the edited program crashes after the call to`vv' mata: heckprob_init("inits","`seldep'", "`depvar'","`touse'") or after the call to capture noi ml model lf2 heckprob_lf2() [...].

From what I read in your post, it seems you have managed to solve the issue of using heckprobit with a fractional dependent variable, and I would be very grateful if you could share the code.
Sincerely

[FernandoRios]

Hi Sophie
This is not exactly via Mata, but you could use -ml- to program your own heckprobit and heckfracprobit.
Below an example:

** This first two programs define the MLE for heckman probit, and heckman fracprobit 
capture program drop mle_heckprobit
program mle_heckprobit, 
    args lnf xb zg arho
    qui {
        tempvar rho
        gen double `rho'=tanh(`arho') 
        replace `lnf'=log(binormal(`xb',`zg',`rho')) if $ML_y1!=0 & $ML_y2!=0 
        replace `lnf'=log(binormal(-`xb',`zg',-`rho')) if $ML_y1==0 & $ML_y2!=0
        replace `lnf'=log(normal(-`zg')) if $ML_y2==0
    }
end

capture program drop mle_heckfprobit
program mle_heckfprobit, 
    args lnf xb zg arho
    *arho
    qui {
        tempvar rho
        gen double `rho'=tanh(`arho') 
        replace `lnf'=$ML_y1*log(binormal(`xb',`zg',`rho'))+(1-$ML_y1)*log(binormal(-`xb',`zg',-`rho')) if $ML_y2!=0 
        replace `lnf'=log(normal(-`zg')) if $ML_y2==0
    }
end

* here an example that compares both procedures

 webuse school, clear
 ml model lf mle_heckprobit (private:private=years logptax) (vote:vote=years loginc logptax) /arho , missing maximize technique(nr bhhh)  nolog
 ml display
 heckprobit private years logptax, sel(vote=years loginc logptax)
 predict priv_xb, xb
 predict vote_xb, xbsel
* and an example using simulated data
 expand 10
 drawnorm u1 u2, corr(1 -.6994978 1) cstorage(lower)
 gen vote_hat = (vote_xb+u2)>0
 gen private_hat = normal(priv_xb+u1)
 
 ml model lf mle_heckfprobit (private:private_hat=years logptax) (vote:vote_hat=years loginc logptax) /arho , missing maximize technique(nr bhhh) 
 ml display

HTH

[David Roodman]

cmp may be able to do what you need. It offers the fractional probit model and lets you mix in selection equations, instrumenting equations, etc.

[Alfonso Sánchez-Peñalver]

Following what Jeff Pitblado (StataCorp) said in #5 above, I would like to suggest that all classes and methods compiled into a Mata library should be documented in help files. A class should have all its elements, public private and protected, and methods well documented, and the methods clearly which arguments are inputs, which are outputs, its purpose and limitations, and what it returns. Yes, it is a lot of work, but it will be extremely useful for cases where Stata doesn't provide exactly what we need. Particularly with classes that can be extended, for example if you had a fixed parameters Logit class, and you wanted to do fractional response Logit, as Richard Williams is trying to do here. You could simply override the method(s) that limit the dependent variable to be 0/1 and check that all values are between 0 and 1 (inclusive), and set the default variance matrix to be the robust one because it is quasi-maximum likelihood. If the class is well designed you would have the different methods that calculate the different variance matrices programmed, and don't need to do it again.

[Richard Williams]

Yes, it is a lot of work, but it will be extremely useful for cases where Stata doesn't provide exactly what we need. Particularly with classes that can be extended, for example if you had a fixed parameters Logit class, and you wanted to do fractional response Logit, as Richard Williams is trying to do here.

Yes, for a couple of my programs, I copied Stata's code, found where it was doing a 0/1 check, and changed it to doing a check for being in the range 0 to 1. It was pretty easy. Stata now has fracreg but I'm surprised it doesn't have even more programs that work with proportions rather than 0/1.

My most recent handout on this is at

https://www3.nd.edu/~rwilliam/stats3...onseModels.pdf