Class: SimpleRandom

Inherits:
Object
  • Object
show all
Defined in:
lib/simple-random/simple_random.rb

Direct Known Subclasses

MultiThreadedSimpleRandom

Defined Under Namespace

Classes: InvalidSeedArgument

Constant Summary collapse

C_32_BIT = 
4294967296
F_32_BIT = 
4294967296.0

Instance Method Summary collapse

Constructor Details

#initialize ⇒ SimpleRandom

Returns a new instance of SimpleRandom.



7
8
9
10
 
# File 'lib/simple-random/simple_random.rb', line 7

def initialize
  @m_w = 521288629
  @m_z = 362436069
end

Instance Method Details

#beta(a, b) ⇒ Object 



102
103
104
105
106
107
 
# File 'lib/simple-random/simple_random.rb', line 102

def beta(a, b)
  fail ArgumentError, "Parameters must be strictly positive" unless a > 0 && b > 0
  u = gamma(a, 1)
  v = gamma(b, 1)
  u / (u + v)
end

#cauchy(median, scale) ⇒ Object 



115
116
117
118
119
 
# File 'lib/simple-random/simple_random.rb', line 115

def cauchy(median, scale)
  fail ArgumentError, 'Scale must be positive' unless scale > 0

  median + scale * Math.tan(Math::PI * (uniform - 0.5))
end

#chi_square(degrees_of_freedom) ⇒ Object 



94
95
96
 
# File 'lib/simple-random/simple_random.rb', line 94

def chi_square(degrees_of_freedom)
  gamma(0.5 * degrees_of_freedom, 2.0)
end

#dirichlet(*parameters) ⇒ Object 



139
140
141
142
143
 
# File 'lib/simple-random/simple_random.rb', line 139

def dirichlet(*parameters)
  sample = parameters.map { |a| gamma(a, 1) }
  sum = sample.inject(0.0) { |sum, g| sum + g }
  sample.map { |g| g / sum }
end

#exponential(mean = 1) ⇒ Object

Get exponential random sample with specified mean



42
43
44
45
46
 
# File 'lib/simple-random/simple_random.rb', line 42

def exponential(mean = 1)
  fail ArgumentError, "Mean must be strictly positive" unless mean > 0

  -1.0 * mean * Math.log(uniform)
end

#gamma(shape, scale) ⇒ Object

Implementation based on “A Simple Method for Generating Gamma Variables” by George Marsaglia and Wai Wan Tsang. ACM Transactions on Mathematical Software Vol 26, No 3, September 2000, pages 363-372.



66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
 
# File 'lib/simple-random/simple_random.rb', line 66

def gamma(shape, scale)
  fail ArgumentError, 'Shape must be strictly positive' unless shape > 0

  base = if shape < 1
    gamma(shape + 1.0, 1.0) * uniform ** -shape
  else
    d = shape - 1 / 3.0
    c = (9 * d) ** -0.5

    begin
      z = normal

      condition1 = z > (-1.0 / c)
      condition2 = false

      if condition1
        u = uniform
        v = (1 + c * z) ** 3
        condition2 = Math.log(u) < (0.5 * (z ** 2) + d * (1.0 - v + Math.log(v)))
      end
    end while !condition2

    d * v
  end

  scale * base
end

#inverse_gamma(shape, scale) ⇒ Object 



98
99
100
 
# File 'lib/simple-random/simple_random.rb', line 98

def inverse_gamma(shape, scale)
  1.0 / gamma(shape, 1.0 / scale)
end

#laplace(mean, scale) ⇒ Object 



127
128
129
130
131
132
133
 
# File 'lib/simple-random/simple_random.rb', line 127

def laplace(mean, scale)
  u_1 = uniform(-0.5, 0.5)
  u_2 = uniform

  sign = u_1 / u_1.abs
  mean + sign * scale * Math.log(1 - u_2)
end

#log_normal(mu, sigma) ⇒ Object 



135
136
137
 
# File 'lib/simple-random/simple_random.rb', line 135

def log_normal(mu, sigma)
  Math.exp(normal(mu, sigma))
end

#normal(mean = 0.0, standard_deviation = 1.0) ⇒ Object

Sample normal distribution with given mean and standard deviation



35
36
37
38
39
 
# File 'lib/simple-random/simple_random.rb', line 35

def normal(mean = 0.0, standard_deviation = 1.0)
  fail ArgumentError, 'Standard deviation must be strictly positive' unless standard_deviation > 0

  mean + standard_deviation * ((-2.0 * Math.log(uniform)) ** 0.5) * Math.sin(2.0 * Math::PI * uniform)
end

#set_seed(*args) ⇒ Object 



12
13
14
15
16
17
18
19
20
21
22
23
24
25
 
# File 'lib/simple-random/simple_random.rb', line 12

def set_seed(*args)
  validate_seeds!(*args)

  @m_w, @m_z = if args.size > 1
    args[0..1].map(&:to_i)
  elsif args.first.is_a?(Numeric)
    [@m_w, args.first.to_i]
  else
    generate_temporal_seed(args.first || Time.now)
  end

  @m_w %= C_32_BIT
  @m_z %= C_32_BIT
end

#student_t(degrees_of_freedom) ⇒ Object 



121
122
123
124
125
 
# File 'lib/simple-random/simple_random.rb', line 121

def student_t(degrees_of_freedom)
  fail ArgumentError, 'Degrees of freedom must be strictly positive' unless degrees_of_freedom > 0

  normal / ((chi_square(degrees_of_freedom) / degrees_of_freedom) ** 0.5)
end

#triangular(lower, mode, upper) ⇒ Object

Get triangular random sample with specified lower limit, mode, upper limit



49
50
51
52
53
54
55
56
57
58
59
60
61
 
# File 'lib/simple-random/simple_random.rb', line 49

def triangular(lower, mode, upper)
  fail ArgumentError, 'Upper bound must be greater than lower bound.' unless lower < upper
  fail ArgumentError, 'Mode must lie between the upper and lower limits' if mode > upper || mode < lower

  f_c = (mode - lower) / (upper - lower)
  uniform_rand_num = uniform

  if uniform_rand_num < f_c
    lower + Math.sqrt(uniform_rand_num * (upper - lower) * (mode - lower))
  else
    upper - Math.sqrt((1 - uniform_rand_num) * (upper - lower) * (upper - mode))
  end
end

#uniform(lower = 0, upper = 1) ⇒ Object

Produce a uniform random sample from the open interval (lower, upper).



28
29
30
31
32
 
# File 'lib/simple-random/simple_random.rb', line 28

def uniform(lower = 0, upper = 1)
  fail ArgumentError, 'Upper bound must be greater than lower bound.' unless lower < upper

  ((get_unsigned_int + 1) * (upper - lower) / F_32_BIT) + lower
end

#weibull(shape, scale) ⇒ Object 



109
110
111
112
113
 
# File 'lib/simple-random/simple_random.rb', line 109

def weibull(shape, scale)
  fail ArgumentError, 'Shape and scale must be positive' unless shape > 0 && scale > 0

  scale * ((-Math.log(uniform)) ** (1.0 / shape))
end