The problem here is that you’ve got a Product of two different types of distributions. If you look at the type of the vector being passed into the constructor for Product, you’ll see the following:
julia> typeof([Exponential(1.0), Normal(0.0, 1.0)])
Vector{Distribution{Univariate, Continuous}} (alias for Array{Distribution{ArrayLikeVariate{0}, Continuous}, 1})
Replacing [Exponential(1.0), Normal(0.0, 1.0)] with [Normal(2.0, 1.0), Normal(0.0, 1.0)] in your example above yields the following results on my machine when I run your code:
BenchmarkTools.Trial: 616 samples with 1 evaluation per sample.
Range (min … max): 7.258 ms … 10.340 ms ┊ GC (min … max): 0.00% … 10.45%
Time (median): 8.138 ms ┊ GC (median): 5.06%
Time (mean ± σ): 8.123 ms ± 395.103 μs ┊ GC (mean ± σ): 4.04% ± 3.37%
▆▇▄ ▃█▂ ▄▂▃
▂▁▁▁▁▁▁▁▁▄███▆▅▃▄▂▃▃▄████████▇▄▃▃▁▃▁▂▂▂▃▁▁▂▂▂▁▂▁▃▂▁▂▃▂▂▂▁▂▃ ▃
7.26 ms Histogram: frequency by time 9.54 ms <
Memory estimate: 15.26 MiB, allocs estimate: 3.
BenchmarkTools.Trial: 622 samples with 1 evaluation per sample.
Range (min … max): 6.856 ms … 10.599 ms ┊ GC (min … max): 0.00% … 32.92%
Time (median): 7.948 ms ┊ GC (median): 11.17%
Time (mean ± σ): 8.035 ms ± 529.012 μs ┊ GC (mean ± σ): 12.23% ± 4.55%
▁▇▃▃▁▄▄ ▂▁▅▅█▆▄▅▃▄▁▂ ▁ ▁ ▁ ▁ ▁
▂▂▃▂▂▃▂▁▁▆███████▇████████████▇▆▇█████▇▅█▇█▄█▆▃▃▃▂▄▃▁▁▁▁▂▂▂ ▅
6.86 ms Histogram: frequency by time 9.51 ms <
Memory estimate: 45.78 MiB, allocs estimate: 12.
i.e. it continues to take very slightly longer to run, but is almost as fast.
Now, clearly, this isn’t the use-case you have in mind (I’m assuming that you actually really do want to sample from this pair of distributions). I don’t know exactly what your use-case is, but you could sample the same set of random numbers doing something like
julia> N = 1_000_000;
julia> @benchmark vcat(rand(Exponential(1.0), $N)', rand(Normal(0.0, 1.0), $N)')
BenchmarkTools.Trial: 729 samples with 1 evaluation per sample.
Range (min … max): 5.989 ms … 9.828 ms ┊ GC (min … max): 0.00% … 34.27%
Time (median): 6.642 ms ┊ GC (median): 10.00%
Time (mean ± σ): 6.855 ms ± 450.783 μs ┊ GC (mean ± σ): 12.89% ± 5.28%
▇▇█▅▃
▂▃▃▂▃▁▂▁▂▁▂▁▁▂▆█████▇▅▄▄▃▃▃▃▄▃▂▃▃▂▂▃▅▆▅▅▄▃▂▂▃▄▃▅▄▄▂▃▃▂▂▂▁▁▂ ▃
5.99 ms Histogram: frequency by time 8 ms <
Memory estimate: 30.52 MiB, allocs estimate: 9.
which seems to give reasonable performance.
I do basically agree with what others have said though: if you really do need a distribution over a pair of distributions whose types are different, product_distribution probably isn’t the tool for you.
edit: having said all of this, your implementation of rand for Product distributions might be a better choice in general than the one that is currently in Distributions.jl. I’m sure the maintainers would appreciate you filing an issue highlighting this performance problem, and pointing out that your implementation of rand_prod seems to yield quite substantial improvements.