Bookmark and Share

Notice: On March 31, it was announced that Statalist is moving from an email list to a forum. The old list will shut down at the end of May, and its replacement, statalist.org is already up and running.


[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

RE: st: AW: Logit maximum likelihood in mata code


From   nick bungy <nickbungystata@hotmail.co.uk>
To   "statalist@hsphsun2.harvard.edu" <statalist@hsphsun2.harvard.edu>
Subject   RE: st: AW: Logit maximum likelihood in mata code
Date   Fri, 19 Apr 2013 15:52:13 +0100

Hi Klaus,
Thanks for the reply. I intend to have a look through the support materials you linked. My aim here was to derive some basic MLEs in the optimizer to help get to grips with what it can do (as it seems a
powerful tool)
I cleaned up my code and the optimizer runs with just one indep variable (see first set of code) but the moment I try to add a second indep variable I just get invalid subscript. I am presumably making an error somewhere but struggling to find it.
One variable - does work:
sysuse auto, clear
mata
st_view(X=0,.,("weight"))
st_view(y=0,.,("foreign"))
void logistic(todo, p, y, X, lf, g, H)
{
b = p[1,1::1]
lf = y'*ln(invlogit(X*b)):+(1:-y)'*ln(1:-invlogit(X*b))
}
S = optimize_init()
optimize_init_evaluator(S, &logistic())
optimize_init_evaluatortype(S, "gf0")
optimize_init_params(S, 0)
optimize_init_argument(S, 1, y)
optimize_init_argument(S, 2, X)
optimize_init_which(S,"max")
p = optimize(S)
--------------
Two variables - subscript invalid:
sysuse auto, clear
gen one = 1
mata
st_view(X=0,.,("weight", "one"))
st_view(y=0,.,("foreign"))
void logistic(todo, p, y, X, lf, g, H)
{
b = p[1, (1::cols(X))]
lf = y'*ln(invlogit(X*b)):+(1:-y)'*ln(1:-invlogit(X*b))
}
S = optimize_init()
optimize_init_evaluator(S, &logistic())
optimize_init_evaluatortype(S, "gf0")
optimize_init_params(S, (0))
optimize_init_argument(S, 1, y)
optimize_init_argument(S, 2, X)
optimize_init_which(S,"max")
p = optimize(S)
> From: kpforr@googlemail.com
> To: statalist@hsphsun2.harvard.edu
> Subject: st: AW: Logit maximum likelihood in mata code
> Date: Fri, 19 Apr 2013 12:45:07 +0200
> 
> <>
> Dear Nicholas,
> 
> I do not understand your code, but your problem is extensively discussed in
> this Stata-Press-book:
> http://www.statapress.com/books/maximum-likelihood-estimation-stata/
> You find a working mata-code for your problem here.
> http://www.statapress.com/data/ml4.html
> 
> best
> 
> Klaus
> 
> 
> __________________________________
> 
> Dr. Klaus Pforr
> GESIS -- Leibniz Institut für Sozialwissenschaft
> B2,1
> Postfach 122155
> D - 68072 Mannheim
> Tel: +49 621 1246 298
> Fax: +49 621 1246 100 
> E-Mail: klaus.pforr@gesis.org
> __________________________________
> 
> -----Ursprüngliche Nachricht-----
> Von: owner-statalist@hsphsun2.harvard.edu
> [mailto:owner-statalist@hsphsun2.harvard.edu] Im Auftrag von nick bungy
> Gesendet: Freitag, 19. April 2013 09:34
> An: statalist@hsphsun2.harvard.edu
> Betreff: st: Logit maximum likelihood in mata code
> 
> Dear Statalist,
> I'm having trouble setting up the optimizer to find the maximum likelihood
> of a logit model. I'm still learning the syntax and set up of mata, though
> I've tried to follow 'help mf_optimize' as best I can.
> Could someone assist in where I am going wrong with the code?
> ----------
> sysuse auto, clear
> mata
> st_view(x=0,.,("weight"))
> st_view(y=0,.,("foreign"))
> x=x,J(rows(x),1,1)
> void logistic(todo, p, y, x, llf, g, H)
> {
> b = p[1,1::2]
> One =J(rows(x),1,1)
> llf = y'*log(One:/One+exp(-(x*b)))+(One-y)'*log(One:/One+exp(-(x*b)))
> }
> S = optimize_init()
> optimize_init_evaluator(S, &logistic()) optimize_init_params(S, 0)
> optimize_init_argument(S, 1, y) optimize_init_argument(S, 2, x)
> optimize_init_which(S,"max")
> p = optimize(S)
> -----------
> Many thanks,
> Nicholas 
> *
> * For searches and help try:
> * http://www.stata.com/help.cgi?search
> * http://www.stata.com/support/faqs/resources/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/faqs/resources/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/faqs/resources/statalist-faq/
*   http://www.ats.ucla.edu/stat/stata/


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