[Octave-bug-tracker] [bug #31671] betapdf fails for large parameter values

[Christos Dimitrakakis]

URL:
  <http://savannah.gnu.org/bugs/?31671>

                 Summary: betapdf fails for large parameter values
                 Project: GNU Octave
            Submitted by: olethros
            Submitted on: Tue Nov 16 12:06:23 2010
                Category: Libraries
                Severity: 3 - Normal
                Priority: 5 - Normal
              Item Group: Incorrect Result
                  Status: None
             Assigned to: None
         Originator Name: Christos Dimitrakakis
        Originator Email: 
             Open/Closed: Open
         Discussion Lock: Any
                 Release: 3.2.3
        Operating System: GNU/Linux

    _______________________________________________________

Details:

Example failure:

octave-3.2.3:6> betapdf(0.5, 1000,1000)
warning: division by zero
ans = NaN

this is due to a naive implementation of this function and prohibits the use
of betapdf for Bayesian statistics. Here is a version which works for all
values:

function p=betapdf(x, alpha, beta)
  z = lgamma(alpha + beta) - lgamma(alpha) - lgamma(beta);
  ya = log(x) .* (alpha - 1);
  ya(x==0) = 0;
  yb = log(1 - x) .* (beta - 1);
  p = exp(z + ya + yb);
  p(x==0) = 0;
  p(x==1) = 0;
endfunction


Example success:
octave-3.2.3:10> betapdf(0.5, 1000,1000)
ans =  35.678






    _______________________________________________________

Reply to this item at:

  <http://savannah.gnu.org/bugs/?31671>

_______________________________________________
  Message sent via/by Savannah
  http://savannah.gnu.org/