I work with a lot of custom distribution functions and it is often convenient and more performant to pass an array to the logpdf function rather than dot vectorization. The reason is that I can reuse intermediate calculations for each data point rather than recomputing them. I encountered a deprecation warning for certain array types, but not others. I want to determine whether this will affect my work, but its not clear to me why a the deprecation warning appears in some cases but not others. Here is an example:
using Distributions
import Distributions: logpdf
mutable struct mydist{T1} <: ContinuousUnivariateDistribution
p1::T1
end
#The following fails
function logpdf(dist::mydist,data::Array{Real,1})
#Do stuff here...
return 0.0
end
logpdf(mydist(0.0),[1.0,2.0,3.0])
┌ Warning: `logpdf(d::UnivariateDistribution, X::AbstractArray)` is deprecated, use `logpdf.(d, X)` instead.
│ caller = top-level scope at none:0
â”” @ Core none:0
#The following works
function logpdf(dist::mydist,data::Array{Float64,1})
#Do stuff here...
return 0.0
end
logpdf(mydist(0.0),[1.0,2.0,3.0])
#The following works
function logpdf(dist::mydist,data::AbstractArray{T}) where T
#Do stuff here...
return 0.0
end
logpdf(mydist(0.0),[1.0,2.0,3.0])
In your first two examples, you’re not actually calling the methods you’ve defined.
If you do @which logpdf(mydist(0.0),[1.0,2.0,3.0]) after defining only the first method, you will see that it dispatches to a method inside Distributions rather than your implementation. The reason is that [1.0,2.0,3.0] is an Array{Float64,1} which is not a subtype of Array{Real,1}. If you want to accept Arrays of subtypes of Real, you should declare your method as logpdf(dist::mydist,data::Array{T,1}) where T<:Real
or with the shorthand logpdf(dist::mydist,data::Array{<:Real,1}).
By giving up dot-syntax, you’re also giving up the possibility of fusion later down the road. You might be able to get the best of both worlds by intercepting the broadcast and doing that intermediate work up-front. I also highly recommend using a non-mutable struct if at all possible if you’re trying to eek out the best possible performance.
struct MyDist ...
end
function Broadcast.broadcasted(::typeof(logpdf), dist::MyDist, data)
return Broadcast.Broadcasted(logpdf, (preplogpdf(dist, data), data))
end
Then you can have preplogpdf return another type that holds those intermediate calculations. This is getting down into the internals of extending broadcast — which is a rather complicated system as a whole —but you don’t really need to know much of it to hook in at this level. You’re just intercepting the system that builds up the broadcast syntax tree as its being built. The broadcasted function is that entry point and the Broadcasted type is the building block of this tree.
This isn’t something I’ve considered doing before, but it should be totally supported and may be a generally useful technique for others, too.
A simpler alternative could be to pre-compute values you might need for logpdf in your distribution’s constructor and store them in the struct. Some of the distributions defined in Distributions do just that. For example, MvNormal uses a specialized matrix type for storing the covariance in its Cholesky form, which is computed on construction.
Thank you yha and mbauman for your helpful replies.
I did not realize I was calling a method in Distributions. The methods in my actual code do not seem to be going through Distributions AbstractArray method. I think I just encountered it on accident because I passed an argument of a different type than I intended. I’m glad to hear that Arrays in general aren’t being deprecated. I don’t want to rewrite all that code!
mbauman, great catch. I have been in the habit of using mutable structs because that is usually what I need. But you are right. It is not necessary and I might be able to get some increase in performance without sacrificing any flexibility.
Thanks for the ideas about using dot vectorization or precomputing values in the constructor. I was wondering what advantages fusion might provide above and beyond my current approach where I pass an array to logpdf?