RE: st: RE: simulated maximum likelihood estimation

["Jun Xu" <mystata@hotmail.com>]

Nick,

Thanks a lot as always! What I am trying to do is to estimate a simulated 
maximum likelihood  logit models. In logit models,

eta = x'b , where eta = ln(p/1-p), and based on soc theories, I suspect
eta = x'b + error, where error is a positively valued disturbance term (say 
expotentially distributed). Because there is probably no closed form 
solution, I need to turn to simulated maximum likelihood estimation. That 
is, in the codes I provided, I try to add an error (expotentially 
distributed with mean of 5) to x'b, repeat this process for 100 times, then 
averge them out and we get simulated likelihood (use expectation to 
approximate the integral, or to integrate out that disturbance term). Then I 
use ml to maximize the function.

I tried to dummy up mvprobit_ll.ado file, in which they kind of created # of 
replication times of variables like:

		forval i = 1/$S_MLE_M {
			gen double `d`i'' = 0
			gen double `sp`i'' = 0
			gen double `arg`i'' = 0
			gen double `k`i'' = 2*``i''-1
		}
...
...
...

replace `sp`i'' = normprob((`arg`i'')/`c`i'`i'')*`sp`j''
...


Jun Xu
Ph.D. Candidate
Department of Sociology
Indiana University at Bloomington
http://mypage.iu.edu/~junxu/home




> From: "Nick Cox" <n.j.cox@durham.ac.uk>
> Reply-To: statalist@hsphsun2.harvard.edu
> To: <statalist@hsphsun2.harvard.edu>
> Subject: RE: st: RE: simulated maximum likelihood estimation
> Date: Mon, 9 May 2005 21:18:25 +0100
>
> As I understand it, you don't need 200
> temporary variables to do what you
> want, as 2 will suffice. Thus your
> code, at a hasty hack, seems to boil down
> to this. I would still recommend building
> up the log likelihood as a sum of logs, not
> as a product to be logged later. I don't
> understand what you are trying to do,
> so can't help further.
>
> set trace off
> set more off
> cap program drop my_ll
> program my_ll
>      version 8.2
>      args lnf xb
>      tempname b sp0
>      sca `b' = 5
>
>      tempvar ed sp
>      gen double `sp' = 1
>      gen double `ed' = 1
>
>      qui forval i = 1/100 {
>          replace `ed' = -ln(uniform()) * `b'
>          replace `sp'= `sp' * (invlogit(`xb' + `ed')) if $ML_y1 == 1
>          replace `sp'= `sp' * (invlogit(-(`xb' + `ed'))) if $ML_y1 == 0
>      }
>
>      qui replace `lnf' = ln(`sp')
> end
>
> Nick
> n.j.cox@durham.ac.uk
>
> Jun Xu
>
> > Sorry for that mistake. The following is a revised version.
> > Still encounter
> > difficulty calculating numerical derivatives. I don't think it's data
> > problem. I tried to look into mvprob_ll.ado by Dr. Cappellari and Dr.
> > Jenkins, and I am stuck there. No one is responsible for
> > solving my problem,
> > except myself; however, I do need some even slight hint to
> > get me through.
> >
> > ********************************
> > set trace off
> > set more off
> > cap program drop my_ll
> > program my_ll
> >     version 8.2
> >     args lnf xb
> >     tempname b sp0
> >     sca `b' = 5
> >
> >     forval i = 1/100 {
> >         tempname ed`i' like`j' sp`i'
> >         gen double `sp`i'' = 0
> >     }
> >
> >     gen double `sp0' = 1
> >
> >     qui replace `lnf' = 0
> >     forval i = 1/100 {
> >         loc j = `i' - 1
> >         qui gen double `ed`i'' = -ln(uniform())*`b'
> >         qui replace `sp`i''= `sp`j''*(invlogit(
> > `xb'+`ed`i''))      if
> > $ML_y1 == 1
> >         qui replace `sp`i''=
> > `sp`j''*(invlogit(-(`xb'+`ed`i'')))    if
> > $ML_y1 == 0
> >     }
> >
> >     qui replace `lnf' = ln(`sp100')
> > end
> >
> > use binlfp2.dta, clear
> > ml model lf my_ll (lfp = k5 k618 age), technique(nr bhhh dfp bfgs)
> > ml search
> > ml maximize
> >
>
> *
> *   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/
_________________________________________________________________
Don�t just search. Find. Check out the new MSN Search! 
http://search.msn.click-url.com/go/onm00200636ave/direct/01/

*
*   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/