Statalist The Stata Listserver


[Date Prev][Date Next][Thread Prev][Thread Next][Date index][Thread index]

Re: st: ml init


From   jpitblado@stata.com (Jeff Pitblado, StataCorp LP)
To   statalist@hsphsun2.harvard.edu
Subject   Re: st: ml init
Date   Fri, 16 Mar 2007 13:54:01 -0500

Martin Heissel <martin.heissel@post.au.dk> asks for advice on how to use
-ml init-:

> I have a problem using the ml init command
> 
> I want to estimate a model with a restriction one a coefficient not being
> negative, using a procedure I found on
> http://www.stata.com/support/faqs/stat/intconst.html.
> 
> But i run into trouble during the 2nd estimation (setting separate equations
> for the coefficients) with a log-likelihood = -inf
> 
> I thought of using the residuals of an intial regression as starting values
> for sigma, but I cannot figure out how to use the ml init command. Can
> anyone help?

Martin could use the natural log of the standard deviation of his depvar as an
intial value for -/lnsigma-.

	sum mpg
	local lns = ln(r(sd))
	...
	ml init /lnsigma = `lns'

In the following, I organized the Stata code Martin's sent into a single
do-file and used the auto data so that you can copy it into the do-file editor
and run it.

***** BEGIN:
* xmpl.do

capture program drop mynormal_d0
capture program drop mynormal_b

program mynormal_d0
	version 9.2
	args todo b lnf
	tempname lnsigma sigma
	tempvar mu 
	mleval `mu' = `b', eq(1)
	mleval `lnsigma' = `b', eq(2) scalar
	quietly {
		scalar `sigma' = exp(`lnsigma')
		mlsum `lnf' = ln( normalden($ML_y1,`mu',`sigma') )
	}
end

sysuse auto

quietly sum mpg
local lns = ln(r(sd))
ml model d0 mynormal_d0 ///
	(xb: mpg = turn trunk) ///
	(lnsigma:), ///
	diparm(lnsigma, exp label(sigma))  
ml init /lnsigma = `lns'			// <-- NEW CODE
ml maximize

program mynormal_b
	version 9.2
	args todo b lnf
	tempname a lnsigma sigma
	tempvar xb mu 
	mleval `a' = `b', eq(1) scalar
	mleval `xb' = `b', eq(2)
	mleval `lnsigma' = `b', eq(3) scalar
	quietly {
		generate double `mu' = `xb' +`a'*$x2
		scalar `sigma' = exp(`lnsigma')
		mlsum `lnf' = ln( normalden($ML_y1,`mu',`sigma') )
	}
end
  
global x2 gear
  
ml model d0 mynormal_b ///
	(a: ) ///
	(xb: mpg = turn trunk) ///
	(lnsigma:), ///
	diparm(lnsigma, exp label(sigma))
ml init /lnsigma = `lns'			// <-- NEW CODE
ml maximize

* end: xmpl.do
***** END:

--Jeff
jpitblado@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/



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