class SimpleRandom
Constants
- C_32_BIT
- F_32_BIT
- GAMMA_NAUGHT
- GAMMA_VALUES
Public Class Methods
new()
click to toggle source
# File lib/simple-random/simple_random.rb, line 7 def initialize @m_w = 521288629 @m_z = 362436069 end
Public Instance Methods
beta(a, b)
click to toggle source
# 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)
click to toggle source
# 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)
click to toggle source
# 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)
click to toggle source
# 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)
click to toggle source
Get exponential random sample with specified mean
# 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)
click to toggle source
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.
# 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)
click to toggle source
# 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)
click to toggle source
# 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)
click to toggle source
# 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)
click to toggle source
Sample normal distribution with given mean and standard deviation
# 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)
click to toggle source
# 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)
click to toggle source
# 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)
click to toggle source
Get triangular random sample with specified lower limit, mode, upper limit
# 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)
click to toggle source
Produce a uniform random sample from the open interval (lower, upper).
# 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)
click to toggle source
# 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
Private Instance Methods
gamma_function(x)
click to toggle source
# File lib/simple-random/simple_random.rb, line 176 def gamma_function(x) return 1e308 if x > 171.0 if x.to_f == x.to_i return unless x > 0 return 1 if x.to_i == 1 (1...x).inject(&:*) else z = if x.abs > 1.0 x.abs - x.abs.to_i else x end gr = GAMMA_VALUES.inject(GAMMA_NAUGHT) do |sum, g| sum * z + g end r = if x.abs > 1 (1..(x.abs.to_i)).inject(1.0) { |prod, i| prod * (x.abs - i) } else 1.0 end if x < 0 && x.abs > 1 -Math::PI * gr * z / (x * r * Math.sin(Math::PI * x)) else r / (gr * z) end end end
generate_temporal_seed(timestamp = Time.now)
click to toggle source
# File lib/simple-random/simple_random.rb, line 170 def generate_temporal_seed(timestamp = Time.now) x = (timestamp.to_f * 1000000).to_i [x >> 16, x % 4294967296] end
get_unsigned_int()
click to toggle source
This is the heart of the generator. It uses George Marsaglia’s MWC algorithm to produce an unsigned integer. See www.bobwheeler.com/statistics/Password/MarsagliaPost.txt
# File lib/simple-random/simple_random.rb, line 150 def get_unsigned_int @m_z = 36969 * (@m_z & 65535) + (@m_z >> 16); @m_w = 18000 * (@m_w & 65535) + (@m_w >> 16); ((@m_z << 16) + (@m_w & 65535)) % 4294967296 end
validate_seeds!(*args)
click to toggle source
# File lib/simple-random/simple_random.rb, line 156 def validate_seeds!(*args) return true if args.compact.empty? unless args[0].to_f.abs > 0 fail InvalidSeedArgument, 'Seeds must be strictly positive' end unless args[1].nil? || args[1].to_f.abs > 0 fail InvalidSeedArgument, 'Seeds must be strictly positive' end true end