RE: st: Changes in Stata's ml routine d0? Stata 8.2 vs. Stata 11.2

[<Jean-Francois.Bertrand@fin.gc.ca>]

Hello,

My two cents worth, when I use ml init for at starting value at 0.01, I typically avoid the Cholesky factorization problem.  It may occur later, but it typically works.

I also had difficulty replicating some papers result with Stata 11.

Cheers 


Jean-François Bertrand
Senior Economist | Économiste principal

-----Original Message-----
From: owner-statalist@hsphsun2.harvard.edu [mailto:owner-statalist@hsphsun2.harvard.edu] On Behalf Of Phil1899@gmx.de
Sent: March 5, 2012 8:53 AM
To: statalist@hsphsun2.harvard.edu
Subject: Re: st: Changes in Stata's ml routine d0? Stata 8.2 vs. Stata 11.2

Thanks for the replies!

I should have mentioned that the program “mdraws” (findit mdraws) from Cappellari and Jenkins has to be installed in order to run the code below. 

I was not aware of the possibility to run the program in Stata 8.2 using version control. Unfortunately, the problem is still there when doing that. The warning “not concave” still appears at each iteration and at iteration 2, the Cholesky factorization is not possible. 

I am surprised by these problems as I just copied the code from the Stata Journal article into a do-file using the same data as the authors of the article. Any suggestions for an explanation?

Best,
Phil


-------- Original-Nachricht --------
> Datum: Mon, 5 Mar 2012 20:52:31 +1100
> Von: Phil Clayton <philclayton@internode.on.net>
> An: statalist@hsphsun2.harvard.edu
> Betreff: Re: st: Changes in Stata\'s ml routine d0? Stata 8.2 vs. 
> Stata 11.2

> Have you tried running it under version control? eg by putting:
> version 8.2
> at the top of your do-file.
> 
> Phil
> 
> On 05/03/2012, at 8:43 PM, Phil1899@gmx.de wrote:
> 
> > Thanks for the hint Nick!
> > 
> > I looked through the file, however, as far as I can understand none 
> > of
> the fixes / changes to ml can explain why a code working for Stata 8.2 
> does not work for Stata 11.2.  Does anyone have an idea why the code 
> below no longer works?
> > 
> > Best,
> > Phil
> > 
> > -------- Original-Nachricht --------
> >> Datum: Mon, 5 Mar 2012 09:08:43 +0000
> >> Von: Nick Cox <njcoxstata@gmail.com>
> >> An: statalist@hsphsun2.harvard.edu
> >> Betreff: Re: st: Changes in Stata\'s ml routine d0? Stata 8.2 vs. 
> >> Stata
> 11.2
> > 
> >> Changes in -ml- (and much else) are documented in the whatsnew* 
> >> help files. Start with -help whatsnew- and scroll to the bottom of 
> >> the viewer for an index.
> >> 
> >> Nick
> >> 
> >> On Mon, Mar 5, 2012 at 8:20 AM,  <Phil1899@gmx.de> wrote:
> >> 
> >>> I would like to ask whether there have been any major changes in
> >> Stata’s ml routine during the last six years (in particular in case
> of d0)? I am
> >> wondering about that as I am trying to replicate results from a 
> >> Stata Journal article from the year 2006 describing a code for 
> >> Stata’s ml
> routine.
> >> I am using the same code and the same data as in the article, 
> >> however,
> the
> >> program no longer converges (not concave). I am working with Stata
> 11.2,
> >> the authors of the Stata Journal article used Stata 8.2.
> >>> 
> >>> The article I am referring to is the following: Peter Haan & Arne
> >> Uhlendorff, 2006. "Estimation of multinomial logit models with
> unobserved
> >> heterogeneity using maximum simulated likelihood," Stata Journal, vol.
> 6(2), pages
> >> 229-245.
> >>> 
> >>> The data used in the article is available at:
> >> http://www.gllamm.org/jspmix.dat
> >>> 
> >>> I present the code below. I would appreciate any advice in this
> matter.
> >>> 
> >>> Best,
> >>> Phil
> >>> 
> >>> 
> >>> clear
> >>> clear matrix
> >>> clear mata
> >>> set mem 400m
> >>> set more off
> >>> cd C:\temp
> >>> set matsize 5000
> >>> 
> >>> cd ""
> >>> 
> >>> matrix p = (7, 11)
> >>> global draws "50"
> >>> infile scy3 id sex stag ravi fry3 tby using jspmix.dat, clear keep 
> >>> scy3 sort scy3
> >>> 
> >>> by scy3: keep if _n==1
> >>> mdraws, neq(2) dr($draws) prefix(c) burn(10) prime (p)
> >>> 
> >>> local repl=${draws}
> >>> local r=1
> >>> while `r' <= `repl' {
> >>> by scy3: gen random_1`r'=invnorm(c1_`r') by scy3: gen 
> >>> random_2`r'=invnorm(c2_`r') local r=`r'+1 } sort scy3 save 
> >>> mdraws_${draws}, replace
> >>> 
> >>> clear
> >>> infile scy3 id sex stag ravi fry3 tby using jspmix.dat, clear sort 
> >>> scy3 merge m:1 scy3 using mdraws_${draws}.dta, update drop _merge 
> >>> sort scy3
> >>> 
> >>> * starting values
> >>> mlogit tby sex, base(1)
> >>> matrix init= e(b)
> >>> matrix list init
> >>> local c1 = init[1,3]
> >>> local c2 = init[1,4]
> >>> local c3 = init[1,5]
> >>> local c4 = init[1,6]
> >>> matrix C = (`c1', `c2', `c3', `c4') matrix list C
> >>> 
> >>> gen a1=0
> >>> gen a2=0
> >>> gen a3=0
> >>> replace a1=1 if tby==1
> >>> replace a2=1 if tby==2
> >>> replace a3=1 if tby==3
> >>> sort scy3
> >>> 
> >>> program drop _all
> >>> program define mlogit_sim_d0
> >>> args todo b lnf
> >>> tempvar etha2 etha3 random1 random2 lj pi1 pi2 pi3 sum lnpi L1 L2 
> >>> last tempname lnsig1 lnsig2 atrho12 sigma1 sigma2 cov12 mleval 
> >>> `etha2' = `b', eq(1) mleval `etha3' = `b', eq(2) mleval `lnsig1' = 
> >>> `b', eq(3) scalar mleval `lnsig2' = `b', eq(4) scalar mleval 
> >>> `atrho12' = `b', eq(5) scalar
> >>> 
> >>> qui {
> >>> scalar `sigma1'=(exp(`lnsig1'))^2
> >>> scalar `sigma2'=(exp(`lnsig2'))^2
> >>> scalar
> >>
> `cov12'=[exp(2*`atrho12')-1]/[exp(2*`atrho12')+1]*(exp(`lnsig2'))*(exp
> (`lnsig1'))
> >>> gen double `random1' = 0
> >>> gen double `random2' = 0
> >>> gen double `lnpi'=0
> >>> gen double `sum'=0
> >>> gen double `L1'=0
> >>> gen double `L2'=0
> >>> by scy3: gen byte `last'=(_n==_N)
> >>> gen double `pi1'= 0
> >>> gen double `pi2'= 0
> >>> gen double `pi3'= 0
> >>> }
> >>> matrix W = ( `sigma1' , `cov12' \ `cov12' , `sigma2') capture 
> >>> matrix L=cholesky(W) if _rc != 0 { di "Warning: cannot do Cholesky 
> >>> factorization of rho matrix"
> >>> }
> >>> local l11=L[1,1]
> >>> local l21=L[2,1]
> >>> local l22=L[2,2]
> >>> 
> >>> local repl=${draws}
> >>> local r=1
> >>> while `r' <= `repl' {
> >>> qui {
> >>> replace `random1' = random_1`r'*`l11'
> >>> replace `random2' = random_2`r'*`l22' + random_1`r'*`l21'
> >>> 
> >>> replace `pi1'= 1/(1 + 
> >>> exp(`etha2'+`random1')+exp(`etha3'+`random2'))
> >>> replace `pi2'= exp(`etha2'+`random1')*`pi1'
> >>> replace `pi3'= exp(`etha3'+`random2')*`pi1'
> >>> 
> >>> replace `lnpi'=ln(`pi1'*a1+`pi2'*a2+`pi3'*a3)
> >>> 
> >>> by scy3: replace `sum'=sum(`lnpi') by scy3: replace `L1' 
> >>> =exp(`sum'[_N]) if _n==_N by scy3: replace `L2'=`L2'+`L1' if 
> >>> _n==_N
> >>> 
> >>> }
> >>> local r=`r'+1
> >>> }
> >>> qui gen `lj'=cond(!`last',0, ln(`L2'/`repl')) qui mlsum `lnf'=`lj'
> >>> if (`todo'==0|`lnf'>=.) exit
> >>> end
> >>> 
> >>> 
> >>> ml model d0 mlogit_sim_d0 ( tby = sex) ( tby = sex) /lnsig1 
> >>> /lnsig2
> >> /atsig12
> >>> ml init C 0.5 0.5 0.5, copy
> >>> ml maximize
> >> 
> >> *
> >> *   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/
> > 
> > -- 
> > NEU: FreePhone 3-fach-Flat mit kostenlosem Smartphone!                  
>                
> > Jetzt informieren: 
> > http://mobile.1und1.de/?ac=OM.PW.PW003K20328T7073a
> > *
> > *   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/
> 
> 
> *
> *   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/

-- 
NEU: FreePhone 3-fach-Flat mit kostenlosem Smartphone!                                  
Jetzt informieren: http://mobile.1und1.de/?ac=OM.PW.PW003K20328T7073a
*
*   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/

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