Re: st: MLE of a function with complex coefficients

[jpitblado@stata.com (Jeff Pitblado, StataCorp LP)]

Pete Aling <pete@unblinkinginsight.com> is having trouble getting -ml- to
iterate with his likelihood evaluator:

> I'm using STATA to do MLE and I'm using the lf method.
> 
> I have only 1 variable, interest rate and it is explained only as a function
> of its previous value. The coefficients however are complex and non-linear
> thus I need to hard code the lag into the likelihood function. I'm not sure
> whether I've done this correctly since I get the following result when I run
> my do file.

> initial:       log likelihood =     -<inf>  (could not be evaluated)
> 
> could not find feasible values
> r(491);
> 
> The process that the interest rate follows is this
> 
> tbrate = exp(beta)*L.tbrate + (alpha/beta)*(exp(beta)-1) + e
> 
> where e is an error term. The likelihood function is in the paper and is
> based on the form of tbrate and the properties of the error term.

> Here is a copy of my code. If anyone has more experience with this and might
> be able to help I'd appreciate it.  tbrate is the interest rate variable and
> ltbrate is its lagged values.

> capture program drop myest
> program define myest
> 	version 6
> 	
> 	args lnf theta1 theta2 theta3 theta4
> 	tempvar a b c d e
> 	quietly gen double `a'= `theta3'/(2*`theta2')
> 	quietly gen double `b'= exp(2*`theta2') - 1
> 	quietly gen double `c'=	2*`theta4'
> 	quietly gen double `d'= exp(`theta2')
> 	quietly gen double `e'= `theta1'/`theta2'
> 	quietly replace `lnf' = (ln(`a'*`b'*(ltbrate)^`c')) + (($ml_y1-`d'*ltbrate-`e'*(`d'-1))^2)/(`a'*`b'*(ltbrate)^`c')
> end
> ml model lf myest (tbrate) /alpha /beta /sigma2 /gamma
> ml maximize

Pete does not say which Stata release he is using, and although he has a
-version 6- in his likelihood evaluator I will assume Stata 8.2.

I found the following issues with Pete's code:

1.  The -ml model- statements specifies 5 equations, but -myest- only
identifies 4.  I'll assume that '/alpha' was meant to be the first equation,
thus the -ml model- statement should be:

	ml model lf myest (alpha: tbrate) /beta /sigma2 /gamma

2.  The dependent variable is not specified correctly; as coded above -ml-
thinks -tbrate- is an independent variable, thus the -ml model- statement
should be further changed to:

	ml model lf myest (alpha: tbrate =) /beta /sigma2 /gamma

3.  -myest- refers to '$ml_y1' when it should be referring to '$ML_y1'
(capital M and capital L).

4.  The log-likelihood formula is up-side-down, and appears to be normal
density function (without the square root of 2*pi and missing a sqrt or 0.5 in
the first -ln()- term).  When possible I prefer to code with Stata's built-in
functions to simplify my expressions, this also keeps me from making mistakes
in my formulas.  In this case, I would use the -normden()- function
(also known as -normalden()- in Stata 9).

5.  -ml- does not get to see -ltbrate-, thus it cannot exclude the 1 expected
missing value from the estimation sample.  This is easily fixed with an -if-
statement in the call to -ml-.

Since I do not have Pete's dataset, I used the auto data to verify that my
changes would result in a converged model fit.  Here is my do-file:

***** BEGIN:
clear

capture program drop myest
program define myest
        version 8.2

        args lnf alpha beta lns2 gamma
        tempvar a b c d e m s
        quietly gen double `a' = exp(`lns2')/(2*`beta')
        quietly gen double `b' = exp(2*`beta') - 1
        quietly gen double `c' = 2*`gamma'
        quietly gen double `d' = exp(`beta')
        quietly gen double `e' = `alpha'/`beta'

        quietly gen double `m' = `d'*ltbrate-`e'*(`d'-1)
        quietly gen double `s' = `a'*`b'*(ltbrate)^`c'

        quietly replace `lnf' = ln(normalden($ML_y1,`m',`s'))
end

        * surrogate for Pete's data
        sysuse auto
        sort mpg
        gen t = _n
        tsset t
        rename mpg tbrate
        gen ltbrate = l.tbrate

ml model lf myest (alpha: tbrate=) /beta /lnsigma2 /gamma if ltbrate < .
ml init /beta = 1 /gamma = 1
ml maximize
***** END:

Note that I added my own initial values for '/beta' and '/gamma', this was to
prevent -ml- from having to use random initial values caused by starting out
with a missing log-likelihood (-ml-'s default initial values are all 0's
resulting in missing values in tempvars `a' and `e').

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