error in ml stat

From: <dph_at_email.domain.hidden>
Date: Wed, 12 Dec 2012 11:30:05 -0500
We have found an error in the maximum likelihood statistic in ISIS,
"ml".

Formula 7.51 (in the version of the manual dated 6/2/11) is incorrect.
Instead of ln Gamma[Mi+1], it should be ln Gamma[Ci+1]. It should also
have factor of 2 to make it distributed as chisq.  (Thanks to Herman
Marshall for providing details.  And for being the first user,
apparently, of "ml".)

I have tested this by defining user statistics as printed in the manual
(accidentally, but which did indicate that that is what is implemented
in ISIS) and with the correct formula (which gives good results).  I fit
a well-exposed emission line with a Gaussian plus polynomial, and I made
confidence contours in line center and area using "chisqr", "cash", and
"myml" stats, with similar results.  "myml" and "cash" gave identical
contours.

E.g.

> cat myml_stat.sl

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

require( "gsl" );

define myml_stat( y, fx, w )
{
    variable v = 2 * ( fx + lngamma( y + 1 ) - y * log( fx ) );
    return( v, sum( v ) );
}

define myml_report ( stat, npts, nvpars )
{
    variable s = sprintf ("  My ML = %0.4g / %d = %0.4g\n", stat, npts-nvpars, stat/(npts-nvpars) );
    return s;
}

add_slang_statistic( "myml", &myml_stat, &myml_report; delta_is_chisqr ) ; 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%



-- Dave

David Huenemoerder  617-253-4283 (o); -253-8084 (f); http://space.mit.edu/home/dph
MIT Kavli Institute for Astrophysics and Space Research
One Hampshire Street, Room NE80-6065
Cambridge, MA  02139
[Admin. Asst.: Elaine Tirrell, 617-253-7480, egt_at_email.domain.hidden
----
You received this message because you are
subscribed to the isis-users list.
To unsubscribe, send a message to
isis-users-request_at_email.domain.hiddenwith the first line of the message as:
unsubscribe
Received on Wed Dec 12 2012 - 11:30:29 EST

This archive was generated by hypermail 2.3.0 : Fri May 02 2014 - 08:35:47 EDT