st: desparate attempt at penalized likelihood with d2

["Mark Anders" <markandothers@gmail.com>]

I am trying to code a penalized likelihood estimator. See Gould,
Pitblado, and Sribney (2006)'s d2 evaluator example for ordinary
probit below (my actual likelihood function is a lot messier, so I am
posting this example instead). The penalty, as Firth (1993) suggests,
is a function of the information matrix added to the likelihood
function, or exactly .5*ln(det(negH)). My problem is how to factor the
penalty back into the likelihood function using the Hessian in the d2
evaluator. I tried adding

qui replace `llf' = 'llf'+.5*ln(det(`negH'))

to the end of the evaluator and use it as a d0 evaluator, but that
didn't work. I think lf evaluator is not an option either because of
the matrix and vector calculations needed. I don't know how to proceed
at this point. Thanks for your help in advance, and happy New Year to
all!

Mark

program myprobit_d2
	args todo b lnf g negH g1
	tempvar xb lj
	mleval `xb' = `b'

	qui {
    	gen double `lj' = normal( `xb')   if $ML_y1 == 1
    	replace `lj' = normal(-`xb')  if $ML_y1 == 0
    	mlsum `lnf' = ln(`lj')‏
	if (`todo'==0 | `lnf'==.) exit

	replace `g1' = normalden(`xb')/`lj' if $ML_y1 == 1
	replace `g1' = -normalden(`xb')/`lj'   if $ML_y1 == 0
	mlvecsum `lnf' `g' = `g1', eq(1)‏
      	if (`todo'==1 | `lnf'==.) exit

	mlmatsum `lnf' `negH' = `g1'*(`g1'+`xb'), eq(1,1)‏
	}

end

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