Module: Distribution::Beta::Ruby_

Extended by:
Math
Defined in:
lib/distribution/beta/ruby.rb

Constant Summary

Constants included from MathExtension

MathExtension::EULER, MathExtension::LN2, MathExtension::LNPI, MathExtension::LOG_FLOAT_MIN, MathExtension::ROOT3_FLOAT_EPSILON, MathExtension::ROOT3_FLOAT_MIN, MathExtension::ROOT4_FLOAT_EPSILON, MathExtension::ROOT4_FLOAT_MIN, MathExtension::ROOT5_FLOAT_EPSILON, MathExtension::ROOT5_FLOAT_MIN, MathExtension::ROOT6_FLOAT_EPSILON, MathExtension::ROOT6_FLOAT_MIN, MathExtension::SQRT2, MathExtension::SQRTPI

Class Method Summary collapse

Methods included from Math

beta, binomial_coefficient, binomial_coefficient_gamma, combinations, erfc_e, exact_regularized_beta, factorial, fast_factorial, gammp, gammq, incomplete_beta, incomplete_gamma, lbeta, logbeta, loggamma, permutations, regularized_beta, rising_factorial, unnormalized_incomplete_gamma

Methods included from MathExtension

#beta, #binomial_coefficient, #binomial_coefficient_gamma, #binomial_coefficient_multiplicative, #erfc_e, #exact_regularized_beta, #exp_err, #factorial, #fast_factorial, #gammq, #incomplete_beta, #incomplete_gamma, #lbeta, #logbeta, #loggamma, #permutations, #regularized_beta, #rising_factorial, #unnormalized_incomplete_gamma

Class Method Details

.cdf(x, a, b) ⇒ Object

Gamma cumulative distribution function Translated from GSL-1.9: cdf/beta.c gsl_cdf_beta_P



32
33
34
35
36
 
# File 'lib/distribution/beta/ruby.rb', line 32

def cdf(x,a,b)
  return 0.0 if x <= 0.0
  return 1.0 if x >= 1.0
  Math::IncompleteBeta.axpy(1.0, 0.0, a,b,x)
end

.p_value(p, a, b, rmin = 0, rmax = 1) ⇒ Object

Inverse of the beta distribution function



39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
 
# File 'lib/distribution/beta/ruby.rb', line 39

def p_value(p,a,b, rmin=0, rmax=1)
  raise "a <= 0" if a <= 0
  raise "b <= 0" if b <= 0
  raise "rmin == rmax" if rmin == rmax
  raise "p <= 0" if p <= 0
  raise "p > 1" if p > 1

  precision=8.88e-016
  max_iterations=256

  ga = 0
  gb = 2

  i = 1
  while ((gb - ga) > precision) && (i < max_iterations)
    guess = (ga + gb) / 2.0
    result = cdf(guess, a, b)

    if (result == p) || (result == 0)
      gb = ga
    elsif result > p
      gb = guess
    else
      ga = guess
    end

    raise 'No value' if i == max_iterations

    i+=1
  end

  rmin + guess * (rmax - rmin)
end

.pdf(x, a, b) ⇒ Object

Beta distribution probability density function

Adapted from GSL-1.9 (apparently by Knuth originally), found in randist/beta.c

Form: p(x) dx = (Gamma(a + b)/(Gamma(a) Gamma(b))) x^(a-1) (1-x)^(b-1) dx

References



16
17
18
19
20
21
22
23
24
25
26
27
28
 
# File 'lib/distribution/beta/ruby.rb', line 16

def pdf(x,a,b)
  return 0 if x < 0 || x > 1

  gab = Math.lgamma(a+b).first
  ga  = Math.lgamma(a).first
  gb  = Math.lgamma(b).first

  if x == 0.0 || x == 1.0
    Math.exp(gab - ga - gb) * x**(a-1) * (1-x)**(b-1)
  else
    Math.exp(gab - ga - gb + Math.log(x)*(a-1) + Math::Log.log1p(-x)*(b-1))
  end
end