[Top][All Lists]

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

[Help-gsl] gamma distribution

From: Sam Mason
Subject: [Help-gsl] gamma distribution
Date: Thu, 2 May 2013 18:02:02 +0100

I've been using the GSL for a while and have just realised that I'm
working around bugs in the Gamma distribution's RNG and PDF.  Some
parameters (i.e. small alpha and beta) cause the RNG to return zero.
As far as I understand this is incorrect mathematically and is infact
due to rounding errors.  I'm OK with this "feature", but when you feed
this zero back to gsl_ran_gamma_pdf() (with the same parameters) you
get zero probability–which feels wrong.

The implementation in R returns +Inf when alpha < 1, which seems more useful.

I'm currently working around this by adding DBL_MIN onto the value I
get back from gsl_ran_gamma(), but this doesn't seem very nice.

Would a patch that adjusts the behaviour of gsl_ran_gamma_pdf(), i.e.
changing the second "if" statement to:

  else if (x == 0)
      if (a == 1)
        return 1/b ;
      else if (a < 1)
        return GSL_POSINF ;
        return 0 ;

do the right thing in everybody's eyes?



reply via email to

[Prev in Thread] Current Thread [Next in Thread]