Output distribution of rand(Float32) and rand(Float64), thread 2

Internals & Design
,
45

No, that’s not the case now and won’t be the case in an improved implementation:

julia> function test_rand(b::Bernoulli, N)
           hits = 0
           misses = 0
           for i in 1:N
              sample(b) ? (hits += 1) : (misses += 1)
           end
           @show N, hits, hits / N, b.p
           @show N, misses, misses / N, 1 - b.p
           return nothing
       end
test_rand (generic function with 1 method)

julia> test_rand(Bernoulli(prevfloat(1f0)), typemax(Int32))
(N, hits, hits / N, b.p) = (2147483647, 2147483529, 0.9999999450519681, 0.99999994f0)
(N, misses, misses / N, 1 - b.p) = (2147483647, 118, 5.4948031927900404e-8, 5.9604645f-8)

The current issue on the low end is that the sample space of rand is much too sparse to reflect the precision with which one can specify a small frequency p::Float32. On the upper end, the available precision for p is lower but perfectly matched by the sample space, so the Bernoulli random variable is weighted according to specification.

1 Like