Dev (#37)

* some refactoring

* pretty printing

* some updates for resettablerng

* deepcopy on iterate

* Move ResettableRNG to MeasureTheory

* Make Pretty a const

* drop extra spaces

* formatting

* formatting

* formatting

* add DensityInterface dependency

* faster `rootmeasure`

* updates

* primitive measure docs

* get tests passing

* cleaning up

* formatting

* drop some unneeded methods

* update for DensityInterface

* updating densities to DensityInterface approach

* update domains

* moving things around

* updates to IntegerBounds

* more domain mucking

* fix some exports

* update deps

* update `using`

* working on tests

* Working toward tests passing

* some refactoring

* working on tests

* cleaning up

* get tests to pass

* tests passing

* update logdensity_def for pointwise product

* PrettyPrinting + tests

* speed up `rootmeasure`

* oops didn't mean to include that

* tile(::FactoredBase)

* Update Half and FactoredBase

* drop exp.jl

* simplify basemeasure

* drop old `include`

* drop redundant method

* update compat Returns

* add AbstractDensity

* Move Affine to MeasureTheory

* drop some old `For` code

* update counting measure

* add testvalue(::Type{T})

* bugfix

* some dispatch adjustments

* simplify show

* formatting

* transformed measures

* ZeroSets

* using LinearAlgebra, Statistics

* working on MeasureTheory tests

* updates

* updates

* typo

* update default to mimic Base

* fix tile(::Lebesgue)

* Add test_interface function

* small doc update

* DensityKind(::Likelihood)

* fixing show(::Likelihood)

* adding some docs

* update `rand` method

* drop old integration code

* add `rebase`

* export rebase

* law for ⊗

* typo

* Make Likelihood more flexible

* update kernel methods

* add Likelihood method (avoid stack overflow)

* refactoring

* compat

* Maybe Comat just works?

* refactoring

* some new stuff

* working on tests

* update powermeasure combinator

* bugfix

* comment out debugging lines

* more refactoring

* test @inferred basemeasure_depth(μ)

* drop `constructor` (just use ConstructionBase.constructorof`)

* debugging

* update help

* update interface

* make tests harder

* fixes

* Dirac bugfix

* formatting

* improve type inference

* working on type inference

* update interface

* udpates

* get test passing

* @test !isabstracttype(typejoin(...))

* work on show methods

* update CI

* remove old code

* update productmeasure

* prettyprinting stuff

* Drop te @constprop :aggressive stuff (maybe don't need it?)

* nerline

* dropping some old code

* update tbasemeasure_type(::PowerMeasure)

* moar tests

* update SpikeMixture

* update superpose type parameter name

* drop old tests

* func_string

* more updates

* getting closer

* almost there!

* generated function for type stability

* tests passing!

* newline

* more fixes

* exports and bugfix

* insupport(μ::Counting{T}, x) where {T<:Type}

* working on MeasureTheory tests

* MeasureTheory tests passing

* drop some old code

* inlining

* improve inference

* update `tile(::For)`

* tighten down infrerence

* update basemeasure(::For) for generators

* loosen type bound on instance_type

* drop debugging code

* small update for Likelihood, and a test

* fixing up likelihoods

* improve `basemeasure_depth` dispatch

* still some trouble with inferred basemeasure_depth

* clean up `For` dispatch

* simplify _logdensityof

* optimize for Returns{True} case

* rework basemeasure_depth

* aggressive tests passing!!

* drop type-level stuff

* drop help

* license

* affero

* copyright notice

* merge

* Drop Create Commons

* cleanup after merge

* update support computations

* insupport(d::SuperpositionMeasure, x)

* dorp ParamWeighted

* insupport(d::FactoredBase, x)

* export unsafe_logdensityof

* call promote_type instead of promote_rule

* logdensity_def for named tuple product measures

* type annotation for now

* debugging

* drop shows

* speed up mapped arrays

* throw an error for `Union{}` types

* MT tests passing

* updates

* get tests passing

* MIT license for MeasureBase

* bump version

Author: cscherrer

LICENSE

@@ -1,13 +1,14 @@
 name = "MeasureBase"
 uuid = "fa1605e6-acd5-459c-a1e6-7e635759db14"
 authors = ["Chad Scherrer <[email protected]> and contributors"]
-version = "0.5.1"
+version = "0.6.0"
 
 [deps]
 Compat = "34da2185-b29b-5c13-b0c7-acf172513d20"
 ConstructionBase = "187b0558-2788-49d3-abe0-74a17ed4e7c9"
 DensityInterface = "b429d917-457f-4dbc-8f4c-0cc954292b1d"
 FillArrays = "1a297f60-69ca-5386-bcde-b61e274b549b"
+IfElse = "615f187c-cbe4-4ef1-ba3b-2fcf58d6d173"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 LogExpFunctions = "2ab3a3ac-af41-5b50-aa03-7779005ae688"
 MappedArrays = "dbb5928d-eab1-5f90-85c2-b9b0edb7c900"
@@ -24,11 +25,12 @@ Compat = "3.35"
 ConstructionBase = "1.3"
 DensityInterface = "0.4"
 FillArrays = "0.12"
+IfElse = "0.1"
 LogExpFunctions = "0.3"
 MappedArrays = "0.4"
 PrettyPrinting = "0.3"
 Reexport = "1"
-Static = "0.4"
+Static = "0.5"
 Tricks = "0.1"
 julia = "1.3"
 

Project.toml

@@ -45,6 +45,7 @@ gentype(μ::AbstractMeasure) = typeof(testvalue(μ))
 
 # gentype(μ::AbstractMeasure) = gentype(basemeasure(μ))
 
+import IfElse: ifelse
 export logdensity_def
 export basemeasure
 export basekleisli
@@ -64,6 +65,7 @@ function logdensity_def end
 
 using Compat
 
+include("splat.jl")
 include("partial-static.jl")
 include("proxies.jl")
 include("kleisli.jl")
@@ -85,14 +87,12 @@ include("combinators/factoredbase.jl")
 include("combinators/weighted.jl")
 include("combinators/superpose.jl")
 include("combinators/product.jl")
-include("combinators/for.jl")
 include("combinators/power.jl")
 include("combinators/spikemixture.jl")
 include("combinators/likelihood.jl")
 include("combinators/pointwise.jl")
 include("combinators/restricted.jl")
 include("combinators/smart-constructors.jl")
-include("combinators/conditional.jl")
 
 include("rand.jl")
 
@@ -102,11 +102,13 @@ module Interface
 using Reexport
 using MeasureBase
 using MeasureBase:basemeasure_depth, proxy
+using MeasureBase: insupport
 @reexport using Test
 
 export test_interface
 export basemeasure_depth
 export proxy
+export insupport
 
 include("interface.jl")
 end # module Interface

src/MeasureBase.jl

@@ -1,6 +0,0 @@
-struct ConditionalMeasure{M,C} <: AbstractMeasure
-    parent::M 
-    constraint::C
-end
-
-Base.:|(μ::AbstractMeasure, constraint) = ConditionalMeasure(μ, constraint)

src/combinators/conditional.jl

@@ -1,25 +1,21 @@
 export FactoredBase
 
-struct FactoredBase{R,C,V,B} <: AbstractMeasure
-    inbounds::R
+struct FactoredBase{C,V,B} <: AbstractMeasure
     constℓ::C
     varℓ::V
     base::B
 end
 
-@inline function logdensity_def(d::FactoredBase{Returns{True}}, x)
-    d.constℓ + d.varℓ()
-end
-
 @inline function logdensity_def(d::FactoredBase, x)
-    d.inbounds(x) || return -Inf
     d.constℓ + d.varℓ()
 end
 
 function Pretty.tile(fb::FactoredBase)
     result = Pretty.literal("FactoredBase")
-    result *= Pretty.list_layout(Pretty.tile.([fb.inbounds, fb.constℓ, fb.varℓ, fb.base]))
+    result *= Pretty.list_layout(Pretty.tile.([fb.constℓ, fb.varℓ, fb.base]))
     result
 end
 
 basemeasure(d::FactoredBase) = d.base
+
+insupport(d::FactoredBase, x) = insupport(d.base, x)

src/combinators/factoredbase.jl

@@ -1,193 +0,0 @@
-
-export For
-using Random
-import Base
-
-struct For{T, F, I} <: AbstractProductMeasure
-    f::F
-    inds::I
-
-    @inline function For(f::F, inds::I) where {F,I<:Tuple}
-        T = Core.Compiler.return_type(f, Tuple{_eltype.(inds)...})
-        new{T,instance_type(f),I}(f, inds)
-    end
-
-    For{T,F,I}(f::F, inds::I) where {T,F,I} = new{T,F,I}(f,inds)
-end
-
-# For(f, gen::Base.Generator) = ProductMeasure(Base.Generator(f ∘ gen.f, gen.iter))
-
-@inline function logdensity_def(d::For{T,F,I}, x::AbstractVector{X}) where {X,T,F,I<:Tuple{<:AbstractVector}}
-    ℓ = zero(float(Core.Compiler.return_type(logdensity_def, Tuple{T,X})))
-    @inbounds for j in eachindex(x)
-        ℓ += logdensity_def(d.f(j), x[j])
-    end
-    ℓ
-end
-
-function logdensity_def(d::For, x::AbstractVector) 
-    sum(eachindex(x)) do i
-        @inbounds logdensity_def(d.f(getindex.(d.inds,i)...), x[i])
-    end
-end
-
-function logdensity_def(d::For{T,F,I}, x::AbstractArray{X}) where {T,F,I,X}
-    ℓ = zero(float(Core.Compiler.return_type(logdensity_def, Tuple{T,X})))
-
-    @inbounds for j in CartesianIndices(x)
-        i = (getindex(ind, j) for ind in d.inds)
-        ℓ += logdensity_def(d.f(i...), x[j])
-    end
-    ℓ
-end
-
-function logdensity_def(d::For{T,F,I}, x) where {N,T,F,I<:NTuple{N,<:Base.Generator}}
-    sum(zip(x, d.inds...)) do (xⱼ, dⱼ...)
-        logdensity_def(d.f(dⱼ...), xⱼ)
-    end
-end
-
-function logdensity_def(d::For{T,F,I}, x::AbstractVector) where {N,T,F,I<:NTuple{N,<:Base.Generator}}
-    sum(zip(x, d.inds...)) do (xⱼ, dⱼ...)
-        logdensity_def(d.f(dⱼ...), xⱼ)
-    end
-end
-
-function marginals(d::For{T,F,I}) where {T,F,I}
-    f(x...) = d.f(x...)::T
-    mappedarray(f, d.inds...)
-end
-
-function marginals(d::For{T,F,I}) where {N,T,F,I<:NTuple{N,<:Base.Generator}}
-    f(x...) = d.f(x...)::T
-    Iterators.map(f, d.inds...)
-end
-
-@inline function basemeasure(d::For{T,F,I}) where {T,F,I}
-    B = Core.Compiler.return_type(basemeasure, Tuple{T})
-    _basemeasure(d, B, static(Base.issingletontype(B)))
-end
-
-@inline function _basemeasure(d::For{T,F,I}, ::Type{B}, ::True) where {T,F,I,B}
-    instance(B) ^ axes(d.inds)
-end
-
-@inline function _basemeasure(d::For{T,F,I}, ::Type{B}, ::False) where {T,F,I,B<:AbstractMeasure}
-    f = basekleisli(d.f)
-    _For(B, f, d.inds)
-end
-
-@inline function _basemeasure(d::For{T,F,I}, ::Type{B}, ::False) where {T,F,I,B}
-    productmeasure(basemeasure.(marginals(d)))
-end
-
-function _basemeasure(d::For{T,F,I}, ::Type{B}, ::True) where {N,T<:AbstractMeasure,F,I<:NTuple{N,<:Base.Generator},B}
-    return instance(B) ^ minimum(length, d.inds)
-end
-
-function _basemeasure(d::For{T,F,I}, ::Type{B}, ::False) where {N,T<:AbstractMeasure,F,I<:NTuple{N,<:Base.Generator},B}
-    f = basekleisli(d.f)
-    _For(B, f, d.inds)
-end
-
-function Pretty.tile(d::For{T}) where {T}
-    result = Pretty.literal("For{")
-    result *= Pretty.tile(T)
-    result *= Pretty.literal("}")
-    result *= Pretty.list_layout(
-        [
-            Pretty.literal(func_string(d.f, Tuple{_eltype.(d.inds)...})),
-            Pretty.tile.(d.inds)...
-        ]
-    )
-end
-
-function _For(::Type{T}, f::F, inds::I) where {T,F,I}
-    For{T,F,I}(f,inds)
-end
-
-
-"""
-    For(f, base...)
-
-`For` provides a convenient way to construct a `ProductMeasure`. There are
-several options for the `base`. With Julia's `do` notation, this can look very
-similar to a standard `for` loop, while maintaining semantics structure that's
-easier to work with.
-
-------------
-
-# `For(f, base::Int...)`
-
-When one or several `Int` values are passed for `base`, the result is treated as
-depending on `CartesianIndices(base)`. 
-
-```
-julia> For(3) do λ Exponential(λ) end |> marginals
-3-element mappedarray(MeasureBase.var"#17#18"{var"#15#16"}(var"#15#16"()), ::CartesianIndices{1, Tuple{Base.OneTo{Int64}}}) with eltype Exponential{(:λ,), Tuple{Int64}}:
- Exponential(λ = 1,)
- Exponential(λ = 2,)
- Exponential(λ = 3,)
-```
-
-```
-julia> For(4,3) do μ,σ Normal(μ,σ) end |> marginals
-4×3 mappedarray(MeasureBase.var"#17#18"{var"#11#12"}(var"#11#12"()), ::CartesianIndices{2, Tuple{Base.OneTo{Int64}, Base.OneTo{Int64}}}) with eltype Normal{(:μ, :σ), Tuple{Int64, Int64}}:
- Normal(μ = 1, σ = 1)  Normal(μ = 1, σ = 2)  Normal(μ = 1, σ = 3)
- Normal(μ = 2, σ = 1)  Normal(μ = 2, σ = 2)  Normal(μ = 2, σ = 3)
- Normal(μ = 3, σ = 1)  Normal(μ = 3, σ = 2)  Normal(μ = 3, σ = 3)
- Normal(μ = 4, σ = 1)  Normal(μ = 4, σ = 2)  Normal(μ = 4, σ = 3)
-```
-
--------
-
-# `For(f, base::AbstractArray...)``
-
-In this case, `base` behaves as if the arrays are `zip`ped together before
-applying the map.
-
-```
-julia> For(randn(3)) do x Exponential(x) end |> marginals
-3-element mappedarray(x->Main.Exponential(x), ::Vector{Float64}) with eltype Exponential{(:λ,), Tuple{Float64}}:
- Exponential(λ = -0.268256,)
- Exponential(λ = 1.53044,)
- Exponential(λ = -1.08839,)
-```
-
-```
-julia> For(1:3, 1:3) do μ,σ Normal(μ,σ) end |> marginals
-3-element mappedarray((:μ, :σ)->Main.Normal(μ, σ), ::UnitRange{Int64}, ::UnitRange{Int64}) with eltype Normal{(:μ, :σ), Tuple{Int64, Int64}}:
- Normal(μ = 1, σ = 1)
- Normal(μ = 2, σ = 2)
- Normal(μ = 3, σ = 3)
-```
-
-----
-
-# `For(f, base::Base.Generator)`
-
-For `Generator`s, the function maps over the values of the generator:
-
-```
-julia> For(eachrow(rand(4,2))) do x Normal(x[1], x[2]) end |> marginals |> collect
-4-element Vector{Normal{(:μ, :σ), Tuple{Float64, Float64}}}:
- Normal(μ = 0.255024, σ = 0.570142)
- Normal(μ = 0.970706, σ = 0.0776745)
- Normal(μ = 0.731491, σ = 0.505837)
- Normal(μ = 0.563112, σ = 0.98307)
-```
-
-"""
-For(f, inds...) = For(f, inds)
-For(f, n::Integer) = For(f, Base.OneTo(n))
-For(f, inds::Integer...) = For(i -> f(Tuple(i)...), Base.CartesianIndices(inds))
-# For(f, inds::Base.Generator) = productmeasure(mymap(f, inds))
-
-function Random.rand!(rng::AbstractRNG, d::For, x) 
-    mar = marginals(d)
-    @inbounds for (dⱼ, j) in zip(marginals(d), eachindex(x))
-        x[j] = rand(rng,dⱼ)
-    end
-    return x
-end
-

src/combinators/for.jl

@@ -14,15 +14,20 @@ unhalf(μ::Half) = μ.parent
 isnonnegative(x) = x ≥ 0.0
 
 @inline function basemeasure(μ::Half)
-    inbounds = isnonnegative
     constℓ = static(logtwo)
     varℓ = Returns(static(0.0))
     base = basemeasure(unhalf(μ))
-    FactoredBase(inbounds, constℓ, varℓ, base)
+    FactoredBase(constℓ, varℓ, base)
 end
 
-function Base.rand(rng::AbstractRNG, T::Type, μ::Half)
+function Base.rand(rng::AbstractRNG, ::Type{T}, μ::Half) where {T}
     return abs(rand(rng, T, unhalf(μ)))
 end
 
 logdensity_def(μ::Half, x) = logdensity_def(unhalf(μ), x)
+
+@inline function insupport(d::Half, x)
+    ifelse(isnonnegative(x), insupport(unhalf(d), x), false)
+end
+
+testvalue(::Half) = 1.0

src/combinators/half.jl

@@ -10,8 +10,6 @@ struct PointwiseProductMeasure{M,L} <: AbstractMeasure
     end
 end
 
-Base.size(μ::PointwiseProductMeasure) = size(μ.data)
-
 function Base.show(io::IO, μ::PointwiseProductMeasure)
     io = IOContext(io, :compact => true)
     print(io, μ.measure, " ⊙ ", μ.likelihood)
@@ -29,9 +27,7 @@ function Base.show_unquoted(io::IO, μ::PointwiseProductMeasure, indent::Int, pr
     return nothing
 end
 
-Base.length(m::PointwiseProductMeasure{T}) where {T} = length(m.data)
-
-⊙(args...) = pointwiseproduct(args...)
+⊙(μ, ℓ) = pointwiseproduct(μ, ℓ)
 
 @inline function logdensity_def(d::PointwiseProductMeasure, x)
     logdensity_def(d.measure, x) + logdensity_def(d.likelihood, x)

src/combinators/pointwise.jl

@@ -65,11 +65,41 @@ params(d::PowerMeasure) = params(first(marginals(d)))
 end
 
 @inline function logdensity_def(d::PowerMeasure{M}, x) where {M}
-    T = eltype(x)
-    ℓ = zero(float(Core.Compiler.return_type(logdensity_def, Tuple{M,T})))
+    ℓ = 0.0
+    # ℓ = zero(typeintersect(AbstractFloat,Core.Compiler.return_type(logdensity_def, Tuple{M,T})))
     parent = d.parent
-    @inbounds for xj in x
-        ℓ += logdensity_def(parent, xj)
+    @simd for xj in x
+        Δℓ = logdensity_def(parent, xj)
+        ℓ += Δℓ
     end
     ℓ
 end
+
+@generated function logdensity_def(d::PowerMeasure{M, Tuple{Base.OneTo{StaticInt{N}}}}, x) where {M,N}
+    quote
+        $(Expr(:meta, :inline))
+        ℓ = 0.0
+        parent = d.parent
+        @inbounds @simd for j in 1:$N
+            Δℓ = logdensity_def(parent, x[j])
+            ℓ += Δℓ
+        end
+        ℓ
+    end
+end
+
+@inline function insupport(μ::PowerMeasure, x)
+    p = μ.parent
+    all(x) do xj
+        # https://github.com/SciML/Static.jl/issues/36
+        dynamic(insupport(p, xj))
+    end
+end
+
+@inline function insupport(μ::PowerMeasure, x::AbstractArray)
+    p = μ.parent
+    all(x) do xj
+        # https://github.com/SciML/Static.jl/issues/36
+        dynamic(insupport(p, xj))
+    end
+end