LoadError: ArgumentError: Sampler for this object is not defined

New to Julia
, ,
1

Hi All,

I’m learning Julia and coming from R and a bit of Python. To help me figure Julia out (in a way that’s directly related to my current work), I’ve decided to write a module based on the Vose/Walker Alias method for sampling randomly from arbitrary discrete distributions. I realize there are some relevant functions in the StatsBase and Distributions packages, but I find learning through errors is helpful. I have encountered an error, though, that’s proven to be hard for me to resolve on my own and I’m hoping someone can point me in the right direction. I’m sure this problem boils down to confusion about Types and how to work with them properly.

Here’s the code:

#= implement the Vose algorithm in Julia http://www.keithschwarz.com/darts-dice-coins/
=#

using Random

import Base.rand
import Random: Sampler, SamplerSimple, SamplerTrivial, Repetition, eltype

rng = MersenneTwister()

n = 5
a = collect(1:n)
b = collect(1:n)
c = rand(Float32,n)

# discrete distribution struct
struct DiscreteDistribution
    prob::Vector
    support::Vector
end

struct AliasTable
    prob::Vector
    support::Vector
    accept::Vector
    alias::Vector
end

dd = DiscreteDistribution(c,a)

# build alias table
function vose_method(dd::DiscreteDistribution)
    n_prob = length(dd.prob)
    n_support = length(dd.support)
    sum_prob = sum(dd.prob)
    prob_scaled = dd.prob * n_prob
    accept = Vector(undef, n_prob)
    alias = Vector{typeof(dd.support[1])}(undef, n_prob)
    small = collect(1:n_prob)[prob_scaled .< 1]
    large = collect(1:n_prob)[prob_scaled .>= 1]
    while length(large) > 0 && length(small) > 0
        s = popfirst!(small)
        l = popfirst!(large)
        accept[s] = prob_scaled[s]
        alias[s] = dd.support[l]
        prob_scaled[l] = (prob_scaled[l] + prob_scaled[s]) - 1
        prob_scaled[l] < 1 ? push!(small, l) : push!(large, l)
    end
    while length(large) > 0
        l = popfirst!(large)
        accept[l] = 1
    end
    while length(small) > 0
        s = popfirst!(small)
        accept[s] = 1
    end
    return AliasTable(dd.prob, dd.support, accept, alias)
end

# sampling method
function sample_vose(rng, aliastable::AliasTable)
    j = rand(rng, 1:length(aliastable.support))
    outcome = rand(rng) <= aliastable.prob[j] ? aliastable.support[j] : aliastable.alias[j]
end

# new sampler
function Random.Sampler(::Type{<:AbstractRNG}, d::DiscreteDistribution, ::Repetition)
    SamplerSimple(d, vose_method(d))
end

Random.eltype(::Type{<:DiscreteDistribution}) = Int

sp = Random.Sampler(rng, dd)

isa(sp, SamplerSimple{<:DiscreteDistribution})

# extending rand
function Random.rand(rng::Type{<:AbstractRNG}, sp::SamplerSimple{<:DiscreteDistribution})
    sample_vose(rng, sp.data)
end

rand(rng, sp)

The last function call returns the error,

LoadError: ArgumentError: Sampler for this object is not defined

But, isa(sp, SamplerSimple{<:DiscreteDistribution}) returns true, so I’m not sure where the dispatch is running into problems.

Any information would be greatly appreciated. I did read the docs in the Standard Library for the Random module and was able to produce a similar error when the example function make_alias_table is defined to return a custom type (without that, if it returns a vector for example, the dispatch seems to work as expected). But, I can’t see why that’s not working either. In case it helps, here’s the code for the example straight out of the docs with the modification to make_alias_table just described,

using Random

import Base.rand
import Random: Sampler, SamplerSimple, Repetition, eltype

rng = MersenneTwister()

struct DiscreteDistribution
    probabilities::Vector
end

struct OtherType
    b::Vector
end

Random.eltype(::Type{<:DiscreteDistribution}) = Int

dd = DiscreteDistribution([0.2,0.2,0.2,0.2,0.2])

function make_alias_table(p)
    OtherType(p)
end

function Random.Sampler(::Type{<:AbstractRNG}, distribution::DiscreteDistribution, ::Repetition)
    SamplerSimple(distribution, make_alias_table(distribution.probabilities))
end

sp = Sampler(rng, dd)

function draw_number(rng, sp)
    rand(rng, sp.data)
end

function rand(rng::AbstractRNG, sp::SamplerSimple{<:DiscreteDistribution})
    draw_number(rng, sp)
end

Many thanks for any help!
Chris

2

I looked only at the simplified second example, the problem is that in your draw_number function, you call rand on sp.data which is of type OtherType, but you didn’t define what this should do. When OtherType is replaced with Vector, the Random module already knows how to sample from that. So you can for example add

Random.Sampler(::Type{RNG}, ot::OtherType, n::Repetition) where {RNG<:AbstractRNG} =
    SamplerSimple(ot, Sampler(RNG, ot.b, n))

rand(rng::AbstractRNG, sp::SamplerSimple{<:OtherType}) = rand(rng, sp.data)
3

Ah of course. I think that’s probably it. I’ll double check and then mark this as the solution. Thanks!

4

Okay, so you were right about the simplified example. Once corrected, I was then able to track down the type mismatch in the signature for the call to rand in the first example. Apparently, julia did not like the variable rng having the type constraint of ‘Type{<:AbstractRNG}’. After changing it to simply ‘rng::AbstractRNG’ (as it appeared in the simplified example), the first example worked. I don’t understand why yet. But, for those interested and maybe working through the documentation as I was, here’s the revised code for the first example that appears to be working:

#= implement the Vose algorithm in Julia http://www.keithschwarz.com/darts-dice-coins/
=#

using Random

import Base.rand
import Random: Sampler, SamplerSimple, SamplerTrivial, Repetition, eltype

rng = MersenneTwister()

n = 5
a = collect(1:n)
b = collect(1:n)
c = rand(Float32,n)

# discrete distribution struct
struct DiscreteDistribution
    prob::Vector
    support::Vector
end

struct AliasTable
    prob::Vector
    support::Vector
    accept::Vector
    alias::Vector
end

dd = DiscreteDistribution(c,a)

# build alias table
function vose_method(dd::DiscreteDistribution)
    n_prob = length(dd.prob)
    n_support = length(dd.support)
    sum_prob = sum(dd.prob)
    prob_scaled = dd.prob * n_prob
    accept = Vector(undef, n_prob)
    alias = Vector{typeof(dd.support[1])}(undef, n_prob)
    small = collect(1:n_prob)[prob_scaled .< 1]
    large = collect(1:n_prob)[prob_scaled .>= 1]
    while length(large) > 0 && length(small) > 0
        s = popfirst!(small)
        l = popfirst!(large)
        accept[s] = prob_scaled[s]
        alias[s] = dd.support[l]
        prob_scaled[l] = (prob_scaled[l] + prob_scaled[s]) - 1
        prob_scaled[l] < 1 ? push!(small, l) : push!(large, l)
    end
    while length(large) > 0
        l = popfirst!(large)
        accept[l] = 1
    end
    while length(small) > 0
        s = popfirst!(small)
        accept[s] = 1
    end
    return AliasTable(dd.prob, dd.support, accept, alias)
end

# sampling method
function sample_vose(rng, aliastable::AliasTable)
    j = rand(rng, 1:length(aliastable.support))
    outcome = rand(rng) <= aliastable.prob[j] ? aliastable.support[j] : aliastable.alias[j]
end

# new sampler
function Random.Sampler(::Type{<:AbstractRNG}, dd::DiscreteDistribution, ::Repetition)
    SamplerSimple(dd, vose_method(dd))
end

sp = Random.Sampler(rng, dd)

# extending rand
function rand(rng::AbstractRNG, sp::SamplerSimple{<:DiscreteDistribution})
    sample_vose(rng, sp.data)
end

rand(rng, sp)

If anyone can explain the difference, then, between rng::Type{<:AbstractRNG} and rng::AbstractRNG to me, I’d appreciate it. And, why does it seem to work in the call to Random.Sampler, but not in the call to rand?

5

rand expects an AbstractRNG object, not a type. It’s like the difference between 1 and Integer; 1 is an instance of an Int <: Integer, whereas Integer is a type. Hopefully this example helps illustrate:

julia> f(a::Integer) = a
f (generic function with 1 method)

julia> f(::Type{<:Integer}) = "f(Integer) was called"
f (generic function with 2 methods)

julia> f(1)
1

julia> f(Int)
"f(Integer) was called"

rand wants an AbstractRNG object, whereas Random.Sampler just wants the type. Not sure why Random.Sampler just needs the type though.

6

Of course, I see it now. Thanks!

7

Generating a random value has two stages, calling Sampler and then calling rand. Sampler is not meant to consume entropy bits (only rand is), so the rng type is enough.