Multidimensional probabilities (#241)

* Increase allowed lower deviation for some DiffEntEsts with negative bias

* Probabilities subtype AbstractArray instead of AbstractVector

* `fasthist` only works on vectors, so adjust docstring

* Update changelog

* Set package version to 2.1

Author: kahaaga

CHANGELOG.md

@@ -6,6 +6,7 @@ Changelog is kept with respect to version 0.11 of Entropies.jl. From version v2.
 
 - Added `Gao` estimator for differential Shannon entropy.
 - Added `Lord` estimator for differential Shannon entropy.
+- `Probabilities` now wraps `AbstractArray{T, N}` instead of `AbstractVector{T}`, so that it can also represent multidimensional probability mass functions. For vectors, it behaves as before.
 
 ## 2.0
 

Project.toml

@@ -2,7 +2,7 @@ name = "ComplexityMeasures"
 uuid = "ab4b797d-85ee-42ba-b621-05d793b346a2"
 authors = "Kristian Agasøster Haaga <kahaaga@gmail.com>, George Datseries <datseris.george@gmail.com>"
 repo = "https://github.com/juliadynamics/ComplexityMeasures.jl.git"
-version = "2.0.0"
+version = "2.1.0"
 
 [deps]
 Combinatorics = "861a8166-3701-5b0c-9a16-15d98fcdc6aa"

src/probabilities.jl

@@ -9,34 +9,36 @@ export outcome_space
 # Types
 ###########################################################################################
 """
-    Probabilities <: AbstractVector
+    Probabilities <: AbstractArray
     Probabilities(x) → p
 
-`Probabilities` is a simple wrapper around `x::AbstractVector{<:Real}` that ensures its
-values sum to 1, so that `p` can be interpreted as probability mass function.
+`Probabilities` is a simple wrapper around `x::AbstractArray{<:Real, N}` that ensures its
+values sum to 1, so that `p` can be interpreted as `N`-dimensional probability mass
+function. In most use cases, `p` will be a vector.
 """
-struct Probabilities{T} <: AbstractVector{T}
-    p::Vector{T}
-    function Probabilities(x::AbstractVector{T}, normed = false) where T <: Real
+struct Probabilities{T, N} <: AbstractArray{T, N}
+    p::Array{T, N}
+    function Probabilities(x::AbstractArray{T, N}, normed = false) where {T <: Real, N}
         if !normed # `normed` is an internal argument that skips checking the sum.
-            s = sum(x)
+            s = sum(x, dims = 1:N)
             if s ≠ 1
                 x = x ./ s
             end
         end
-        return new{T}(x)
+        return new{T, N}(x)
     end
 end
-function Probabilities(x::AbstractVector{<:Integer})
+function Probabilities(x::AbstractArray{<:Integer, N}) where N
     s = sum(x)
     return Probabilities(x ./ s, true)
 end
 
-# extend base Vector interface:
+# extend base Array interface:
 for f in (:length, :size, :eachindex, :eltype, :parent,
     :lastindex, :firstindex, :vec, :getindex, :iterate)
     @eval Base.$(f)(d::Probabilities, args...) = $(f)(d.p, args...)
 end
+
 Base.IteratorSize(::Probabilities) = Base.HasLength()
 # Special extension due to the rules of the API
 @inline Base.sum(::Probabilities{T}) where T = one(T)
@@ -96,7 +98,7 @@ contain `0`s as entries or not depends on the estimator.
 E.g., in [`ValueHistogram`](@ref) `0`s are skipped, while in
 [`SymbolicPermutation`](@ref) `0` are not, because we get them for free.
 
-    probabilities(x::Array_or_Dataset) → p::Probabilities
+    probabilities(x::Vector_or_Dataset) → p::Probabilities
 
 Estimate probabilities by directly counting the elements of `x`, assuming that
 `Ω = sort(unique(x))`, i.e. that the outcome space is the unique elements of `x`.

test/entropies/estimators/alizadeharghami.jl

@@ -13,9 +13,9 @@ ea = entropy(Shannon(), AlizadehArghami(m = 100), rand(npts))
 ea_n = entropy(Shannon(; base = ℯ), AlizadehArghami(m = 100), randn(npts))
 ea_n3 = entropy(Shannon(; base = 3), AlizadehArghami(m = 100), randn(npts))
 
-@test U - max(0.01, U*0.03) ≤ ea ≤ U + max(0.01, U*0.03)
-@test N * 0.98 ≤ ea_n ≤ N * 1.02
-@test N_base3 * 0.98 ≤ ea_n3 ≤ N_base3 * 1.02
+@test U - max(0.02, U*0.03) ≤ ea ≤ U + max(0.01, U*0.03)
+@test N * 0.96 ≤ ea_n ≤ N * 1.02
+@test N_base3 * 0.96 ≤ ea_n3 ≤ N_base3 * 1.02
 
 x = rand(1000)
 @test_throws ArgumentError entropy(Renyi(q = 2), AlizadehArghami(), x)

test/entropies/estimators/correa.jl

@@ -13,9 +13,9 @@ ea = entropy(Shannon(), Correa(m = 100), rand(npts))
 ea_n = entropy(Shannon(; base = ℯ), Correa(m = 100), randn(npts))
 ea_n3 = entropy(Shannon(; base = 3), Correa(m = 100), randn(npts))
 
-@test U - max(0.01, U*0.03) ≤ ea ≤ U + max(0.01, U*0.03)
-@test N * 0.98 ≤ ea_n ≤ N * 1.02
-@test N_base3 * 0.98 ≤ ea_n3 ≤ N_base3 * 1.02
+@test U - max(0.02, U*0.03) ≤ ea ≤ U + max(0.01, U*0.03)
+@test N * 0.96 ≤ ea_n ≤ N * 1.02
+@test N_base3 * 0.96 ≤ ea_n3 ≤ N_base3 * 1.02
 
 x = rand(1000)
 @test_throws ArgumentError entropy(Renyi(q = 2), Correa(), x)

test/entropies/estimators/ebrahimi.jl

@@ -13,9 +13,9 @@ ea = entropy(Shannon(), Ebrahimi(m = 100), rand(npts))
 ea_n = entropy(Shannon(; base = ℯ), Ebrahimi(m = 100), randn(npts))
 ea_n3 = entropy(Shannon(; base = 3), Ebrahimi(m = 100), randn(npts))
 
-@test U - max(0.01, U*0.03) ≤ ea ≤ U + max(0.01, U*0.03)
-@test N * 0.98 ≤ ea_n ≤ N * 1.02
-@test N_base3 * 0.98 ≤ ea_n3 ≤ N_base3 * 1.02
+@test U - max(0.02, U*0.03) ≤ ea ≤ U + max(0.01, U*0.03)
+@test N * 0.96 ≤ ea_n ≤ N * 1.02
+@test N_base3 * 0.96 ≤ ea_n3 ≤ N_base3 * 1.02
 
 x = rand(1000)
 @test_throws ArgumentError entropy(Renyi(q = 2), Ebrahimi(), x)

test/entropies/estimators/lord.jl

@@ -13,9 +13,9 @@ ea = entropy(Shannon(), Lord(k = 20), rand(npts))
 ea_n = entropy(Shannon(; base = ℯ), Lord(k = 20), randn(npts))
 ea_n3 = entropy(Shannon(; base = 3), Lord(k = 20), randn(npts))
 
-@test U - max(0.03, U*0.03) ≤ ea ≤ U + max(0.03, U*0.03)
-@test N * 0.97 ≤ ea_n ≤ N * 1.03
-@test N_base3 * 0.98 ≤ ea_n3 ≤ N_base3 * 1.02
+@test U - max(0.05, U*0.03) ≤ ea ≤ U + max(0.03, U*0.03)
+@test N * 0.96 ≤ ea_n ≤ N * 1.03
+@test N_base3 * 0.96 ≤ ea_n3 ≤ N_base3 * 1.03
 
 x = rand(1000)
 @test_throws ArgumentError entropy(Renyi(q = 2), Lord(), x)

test/entropies/estimators/vasicek.jl

@@ -14,9 +14,9 @@ ea = entropy(Shannon(), Vasicek(m = 100), rand(npts))
 ea_n = entropy(Shannon(; base = ℯ), Vasicek(m = 100), randn(npts))
 ea_n3 = entropy(Shannon(; base = 3), Vasicek(m = 100), randn(npts))
 
-@test U - max(0.01, U*0.03) ≤ ea ≤ U + max(0.01, U*0.03)
-@test N * 0.98 ≤ ea_n ≤ N * 1.02
-@test N_base3 * 0.98 ≤ ea_n3 ≤ N_base3 * 1.02
+@test U - max(0.02, U*0.03) ≤ ea ≤ U + max(0.01, U*0.03)
+@test N * 0.96 ≤ ea_n ≤ N * 1.02
+@test N_base3 * 0.96 ≤ ea_n3 ≤ N_base3 * 1.02
 
 x = rand(1000)
 @test_throws ArgumentError entropy(Renyi(q = 2), Vasicek(), x)