eltype and rand with Univariate/Continuous distributions

[halleysfifthinc]

I'm not sure how much more endemic this issue is (I've only tested with a Normal distribution), but here are a couple examples of counter-intuitive (and presumably wrong) behavior.

If I have a distribution (eg. normal), I would expect the output from sampling to match the types parameterizing the distribution.

julia> n = Normal(1.0f0, 0.1f0)
Normal{Float32}(μ=1.0f0, σ=0.1f0)

julia> typeof(rand(n)) # I would expect this to be Float32
Float64

In the case of a Normal dist., that could be fixed with this patch:

diff --git a/src/univariate/continuous/normal.jl b/src/univariate/continuous/normal.jl
index 8d6b3d3c..563696e2 100644
--- a/src/univariate/continuous/normal.jl
+++ b/src/univariate/continuous/normal.jl
@@ -87,7 +87,7 @@ cf(d::Normal, t::Real) = exp(im * t * d.μ - d.σ^2/2 * t^2)

 #### Sampling

-rand(rng::AbstractRNG, d::Normal) = d.μ + d.σ * randn(rng)
+rand(rng::AbstractRNG, d::Normal{T}) where T = d.μ + d.σ * randn(rng, T)


 #### Fitting

But the following example would still produce Float64:

julia> rand(n, 4) # I would expect an Array{Float32,1}
4-element Array{Float64,1}:
 0.9081488251686096
 1.1430985927581787
 1.0427565574645996
 0.9090480804443359

because

julia> eltype(n) # this is hard coded in common.jl#51
Float64

and

rand(rng::AbstractRNG, s::Sampleable{Univariate}, dims::Dims) =
    rand!(rng, sampler(s), Array{eltype(s)}(undef, dims))
# at src/univariates.jl#160

Presumably the second issue would affect most/all Univariate distributions, however, I am not sure how typical of an example the code for the Normal dist. is as far as the type signature and returned type. (Unfortunately I have neither the familiarity with the Distributions.jl source nor the time to do a more thorough investigation and/or make a PR.)

Unrelated to the main issue, I also wonder (again, not being familiar with the code) whether in the second example the sampler(s) is necessary (and/or a no-op) given s::Sampleable.

[matbesancon]

yes it's indeed a legacy issue we are trying to erase one step at a time. See #882 for instance.

The Univariate case should be easy to fix, there is this generic fallback in commons.jl, but you can specialize it for all distributions, so for example:

eltype(::Normal{T}) where {T} = T

[matbesancon]

PR very welcome on that :)