From e63d9fa39e4ba22e3b595d1d110a114dcf93eb63 Mon Sep 17 00:00:00 2001
From: Kristian Haaga <kahaaga@gmail.com>
Date: Mon, 29 Aug 2022 12:35:46 +0200
Subject: [PATCH 01/20] Approximate entropy

---
 docs/make.jl                                  |  3 +-
 docs/src/ApproximateEntropy.md                | 83 +++++++++++++++++++
 src/Entropies.jl                              |  1 +
 .../approximate_entropy.jl                    | 80 ++++++++++++++++++
 test/runtests.jl                              | 19 +++++
 5 files changed, 185 insertions(+), 1 deletion(-)
 create mode 100644 docs/src/ApproximateEntropy.md
 create mode 100644 src/approximate_entropy/approximate_entropy.jl

diff --git a/docs/make.jl b/docs/make.jl
index 98e2dc343..19a3b68d8 100644
--- a/docs/make.jl
+++ b/docs/make.jl
@@ -43,7 +43,8 @@ PAGES = [
         "DispersionEntropy.md",
         "NearestNeighbors.md",
         "NaiveKernel.md",
-        "TimeScaleMODWT.md"
+        "TimeScaleMODWT.md",
+        "ApproximateEntropy.md"
     ],
     "Non-exported" => "nonexported.md"
 ]
diff --git a/docs/src/ApproximateEntropy.md b/docs/src/ApproximateEntropy.md
new file mode 100644
index 000000000..c60c94cfe
--- /dev/null
+++ b/docs/src/ApproximateEntropy.md
@@ -0,0 +1,83 @@
+# Approximate entropy
+
+```@docs
+approx_entropy
+```
+
+## Example
+
+Here, we reproduce the Henon map example with ``R=0.8`` from Pincus (1991),
+comparing our values with relevant values from table 1 in Pincus (1991).
+
+We use `DiscreteDynamicalSystem` from `DynamicalSystems.jl` to represent the map,
+and use the `trajectory` function from the same package to iterate the map
+for different initial conditions, for multiple time series lengths.
+
+Finally, we summarize our results in box plots and compare the values to those
+obtained by Pincus (1991).
+
+```@example
+using DynamicalSystems, StatsBase, PyPlot
+
+# Equation 13 in Pincus (1991)
+function eom_henon(x, p, n)
+    R = p[1]
+    x, y = (x...,)
+    dx = R*y + 1 - 1.4*x^2
+    dy = 0.3*R*x
+
+    return SVector{2}(dx, dy)
+end
+
+function henon(; u₀ = rand(2), R = 0.8)
+    DiscreteDynamicalSystem(eom_henon, u₀, [R])
+end
+
+ts_lengths = [300, 1000, 2000, 3000]
+nreps = 100
+apens_08 = [zeros(nreps) for i = 1:length(ts_lengths)]
+
+# For some initial conditions, the Henon map as specified here blows up,
+# so we need to check for infinite values.
+containsinf(x) = any(isinf.(x))
+
+for (i, L) in enumerate(ts_lengths)
+    k = 1
+    while k <= nreps
+        sys = henon(u₀ = rand(2), R = 0.8)
+        t = trajectory(sys, L, Ttr = 5000)
+
+        if !any([containsinf(tᵢ) for tᵢ in t])
+            x, y = columns(t)
+            apen = approx_entropy(x, r = 0.05, m = 2)
+            apens_08[i][k] = apen
+            k += 1
+        end
+    end
+end
+
+f = figure(figsize = (6, 5))
+
+# First subplot is an example time series
+sys = henon(u₀ = [0.5, 0.1], R = 0.8)
+x, y = columns(trajectory(sys, 100, Ttr = 500))
+
+subplot(211)
+plot(x, label = "x")
+plot(y, label = "y")
+xlabel("Time (t)")
+ylabel("Value")
+
+# Second subplot contains the approximate entropy values
+subplot(212)
+boxplot(apens_08, positions = ts_lengths, widths = [150, 150, 150, 150], notch = true)
+scatter(ts_lengths, [0.337, 0.385, NaN, 0.394], label = "Pincus (1991)")
+xlabel("Time series length (L)")
+ylabel("ApEn(m = 2, r = 0.05)")
+
+legend()
+tight_layout()
+savefig("approx_entropy_pincus.png") # hide
+```
+
+![Approximate entropy](approx_entropy_pincus.png)
diff --git a/src/Entropies.jl b/src/Entropies.jl
index 854ae31d0..fe1be197b 100644
--- a/src/Entropies.jl
+++ b/src/Entropies.jl
@@ -8,4 +8,5 @@ module Entropies
     include("wavelet/wavelet.jl")
     include("nearest_neighbors/nearest_neighbors.jl")
     include("dispersion/dispersion_entropy.jl")
+    include("approximate_entropy/approximate_entropy.jl")
 end
diff --git a/src/approximate_entropy/approximate_entropy.jl b/src/approximate_entropy/approximate_entropy.jl
new file mode 100644
index 000000000..759972df5
--- /dev/null
+++ b/src/approximate_entropy/approximate_entropy.jl
@@ -0,0 +1,80 @@
+using DelayEmbeddings
+using Neighborhood
+using StatsBase
+using Distances
+
+export approx_entropy
+
+"""
+    approx_entropy(x, m = 3, r = 0.5 * StatsBase.std(x), base = MathConstants.e) → ApEn
+
+Compute the approximate entropy (ApEn; Pincus, 1991)[^Pincus1991] of the scalar-valued
+time series `x`, using embedding dimension (pattern length) `m` and tolerance radius `r`,
+using logarithms to the given `base`.
+
+[^Pincus1991]: Pincus, S. M. (1991). Approximate entropy as a measure of system complexity. Proceedings of the National Academy of Sciences, 88(6), 2297-2301.
+"""
+function approx_entropy(x; m = 2, r = 0.2 * StatsBase.std(x), base = MathConstants.e)
+    m >= 1 || throw(ArgumentError("m must be >= 1. Got m=$m."))
+    N = length(x)
+    ϕᵐ = ϕmr_tree(x, m, r, N, base = base)
+    ϕᵐ⁺¹ = ϕmr_tree(x, m + 1, r, N, base = base)
+    return ϕᵐ - ϕᵐ⁺¹
+end
+
+"""
+    ϕmr_tree(x::AbstractVector{T}, k, r, N; base = MathConstants.e)
+
+Construct the embedding
+
+```math
+u = \\{{\\bf u}_n \\}_{n = 1}^{N - k + 1} =
+\\{[x(i), x(i + 1), \\ldots, x(i + k - 1)]\\}_{n = 1}^{N - k + 1}
+```
+
+and use a tree-and-nearest-neighbor search approach to compute
+
+
+```math
+\\phi^k(r) = \\dfrac{1}{N - k + 1} \\sum_{i}^{N - k + 1} \\log_{b}{(C_i^k(r))},
+```
+
+taking logarithms to `base` ``b``, and where
+
+```math
+C_i^k(r) = \\textrm{number of } j \\textrm{ such that } d({\\bf u}_i, {\\bf u}_j) < r,
+```
+
+where ``d`` is the maximum (Chebyshev) distance, `r` is the tolerance, and `N` is the
+length of the original scalar-valued time series `x`.
+
+## Examples
+
+```jldoctest; setup = :(using Entropies)
+julia> x = repeat([2, 3, 4, 8, 8, 9, 3, 6, 4, 5, 8, 7], 100);
+
+julia> approx_entropy(x, m = 2, r = 3)
+0.42020216510849995
+```
+"""
+function ϕmr_tree(x::AbstractVector{T}, k, r, N::Int;
+        base = MathConstants.e) where T <: Real
+    pts = genembed(x, 0:(k - 1))
+
+    # For each xᵢ ∈ pts, find the indices of neighbors within radius `r` according to
+    # the Chebyshev (maximum) metric. Self-inclusion is wanted, so no Theiler window
+    # is specified. This tree-approach is more than an order of magnitude
+    # faster than checking point-by-point using nested loops, and also allocates
+    # more than an order of magnitude less.
+    tree = KDTree(pts, Chebyshev())
+    idxs_of_neighbors = bulkisearch(tree, pts, WithinRange(r))
+
+    q = N - k + 1
+    ϕ = 0.0
+    for nn_idxsᵢ in idxs_of_neighbors
+        # Self-inclusion during the neighbor search means that there are always neighbors,
+        # so we never encounter log(0) here.
+        ϕ += log(base, length(nn_idxsᵢ) / q)
+    end
+    return ϕ / q
+end
diff --git a/test/runtests.jl b/test/runtests.jl
index a9392b602..60cb66568 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -3,6 +3,7 @@ using Entropies
 using DelayEmbeddings
 using Wavelets
 using StaticArrays
+using StatsBase
 using Neighborhood: KDTree, BruteForce
 
 @testset "Histogram estimation" begin
@@ -376,4 +377,22 @@ end
         @test typeof(de) <: Real
         @test de >= 0.0
     end
+
+    @testset "Approximate entropy" begin
+        x = rand(100)
+        N = length(x)
+        r = StatsBase.std(x)
+        m = 2 # embedding dimension
+        base = MathConstants.e
+
+        # Counting function
+        ϕᵐ = Entropies.ϕmr_tree(x, m, r, N, base = base)
+        ϕᵐ⁺¹ =  Entropies.ϕmr_tree(x, m + 1, r, N, base = base)
+        @test ϕᵐ isa T where T <: Real
+        @test ϕᵐ⁺¹ isa T where T <: Real
+
+        # Approximate entropy
+        ap_en = approx_entropy(x, m = m, r = r, base = base)
+        @test ap_en isa T where T <: Real
+    end
 end

From 14058b2b2c315d541855a0aa7e6c1172d172a1a5 Mon Sep 17 00:00:00 2001
From: Kristian Haaga <kahaaga@gmail.com>
Date: Mon, 5 Sep 2022 08:21:33 +0200
Subject: [PATCH 02/20] Add estimator and use api

---
 docs/src/ApproximateEntropy.md                |  1 +
 .../approximate_entropy.jl                    | 49 +++++++++++++++----
 test/runtests.jl                              |  4 ++
 3 files changed, 45 insertions(+), 9 deletions(-)

diff --git a/docs/src/ApproximateEntropy.md b/docs/src/ApproximateEntropy.md
index c60c94cfe..c37f42df4 100644
--- a/docs/src/ApproximateEntropy.md
+++ b/docs/src/ApproximateEntropy.md
@@ -1,6 +1,7 @@
 # Approximate entropy
 
 ```@docs
+ApproximateEntropy
 approx_entropy
 ```
 
diff --git a/src/approximate_entropy/approximate_entropy.jl b/src/approximate_entropy/approximate_entropy.jl
index 759972df5..00bdfa315 100644
--- a/src/approximate_entropy/approximate_entropy.jl
+++ b/src/approximate_entropy/approximate_entropy.jl
@@ -3,25 +3,57 @@ using Neighborhood
 using StatsBase
 using Distances
 
+export ApproximateEntropy
 export approx_entropy
 
 """
-    approx_entropy(x, m = 3, r = 0.5 * StatsBase.std(x), base = MathConstants.e) → ApEn
+    ApproximateEntropy(; r = 0.1, m::Int = 2)
 
-Compute the approximate entropy (ApEn; Pincus, 1991)[^Pincus1991] of the scalar-valued
-time series `x`, using embedding dimension (pattern length) `m` and tolerance radius `r`,
-using logarithms to the given `base`.
+Compute the approximate entropy (ApEn; Pincus, 1991)[^Pincus1991] of a scalar-valued
+time series, using embedding dimension (pattern length) `m` and tolerance radius `r`.
+
+The tolerance radius `r` should be determined from the input data, for example as
+some fraction of the standard deviation of the input.
+
+!!! info
+    This estimator is only available for entropy estimation.
+    Probabilities cannot be obtained directly.
 
 [^Pincus1991]: Pincus, S. M. (1991). Approximate entropy as a measure of system complexity. Proceedings of the National Academy of Sciences, 88(6), 2297-2301.
 """
-function approx_entropy(x; m = 2, r = 0.2 * StatsBase.std(x), base = MathConstants.e)
-    m >= 1 || throw(ArgumentError("m must be >= 1. Got m=$m."))
+Base.@kwdef struct ApproximateEntropy <: EntropyEstimator
+    r::Real = 0.1
+    m::Int = 2
+end
+
+function genentropy(x::AbstractVector, est::ApproximateEntropy; base = MathConstants.e)
+    est.m >= 1 || throw(ArgumentError("m must be >= 1. Got m=$(est.m)."))
     N = length(x)
-    ϕᵐ = ϕmr_tree(x, m, r, N, base = base)
-    ϕᵐ⁺¹ = ϕmr_tree(x, m + 1, r, N, base = base)
+    ϕᵐ = ϕmr_tree(x, est.m, est.r, N, base = base)
+    ϕᵐ⁺¹ = ϕmr_tree(x, est.m + 1, est.r, N, base = base)
     return ϕᵐ - ϕᵐ⁺¹
 end
 
+function genentropy(x::AbstractDataset, est::ApproximateEntropy)
+    throw(
+        ArgumentError("Approximate entropy is currently not defined for multivariate data.")
+    )
+end
+
+
+"""
+    approx_entropy(x::AbstractVector ; m = 2, r = 0.2*StatsBase.std(x), base = MathConstants.e)
+Shorthand for `genentropy(x, ApproximateEntropy(m = m, r = r); base = base)` which
+calculates the approximate entropy of for dimension `m` with tolerance radius `r`.
+
+If `r` is not specified, `r` is estimated to a recommended value based on the standard
+deviation of `x`.
+"""
+function approx_entropy(x; m = 2, r = 0.2*StatsBase.std(x), base = MathConstants.e)
+    est = ApproximateEntropy(r, m)
+    genentropy(x, est; base = base)
+end
+
 """
     ϕmr_tree(x::AbstractVector{T}, k, r, N; base = MathConstants.e)
 
@@ -34,7 +66,6 @@ u = \\{{\\bf u}_n \\}_{n = 1}^{N - k + 1} =
 
 and use a tree-and-nearest-neighbor search approach to compute
 
-
 ```math
 \\phi^k(r) = \\dfrac{1}{N - k + 1} \\sum_{i}^{N - k + 1} \\log_{b}{(C_i^k(r))},
 ```
diff --git a/test/runtests.jl b/test/runtests.jl
index 60cb66568..eab0611ed 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -392,7 +392,11 @@ end
         @test ϕᵐ⁺¹ isa T where T <: Real
 
         # Approximate entropy
+        est = ApproximateEntropy(r = 0.1, m = 2)
+        @test genentropy(x, est; base = base) isa T where T <: Real
         ap_en = approx_entropy(x, m = m, r = r, base = base)
         @test ap_en isa T where T <: Real
+
+        @test_throws ArgumentError genentropy(Dataset(rand(100, 3)), est)
     end
 end

From 568c26c4a0c5b9298f6c4b1b260f3ad093cd04ee Mon Sep 17 00:00:00 2001
From: Kristian Haaga <kahaaga@gmail.com>
Date: Fri, 30 Sep 2022 17:38:15 +0200
Subject: [PATCH 03/20] New API

---
 docs/src/complexity_measures.md               |   6 +
 .../approximate_entropy.jl                    | 111 ------------------
 .../approximate_entropy.jl                    |  91 ++++++++++++++
 .../complexity_measures.jl                    |   1 +
 test/complexity_measures/approx_entropy.jl    |  16 +++
 .../complexity_measures.jl                    |   1 +
 6 files changed, 115 insertions(+), 111 deletions(-)
 delete mode 100644 src/approximate_entropy/approximate_entropy.jl
 create mode 100644 src/complexity_measures/approximate_entropy.jl
 create mode 100644 test/complexity_measures/approx_entropy.jl

diff --git a/docs/src/complexity_measures.md b/docs/src/complexity_measures.md
index c063d28c1..e8c35f110 100644
--- a/docs/src/complexity_measures.md
+++ b/docs/src/complexity_measures.md
@@ -6,3 +6,9 @@
 reverse_dispersion
 distance_to_whitenoise
 ```
+
+## Approximate entropy
+
+```@docs
+approx_entropy
+```
diff --git a/src/approximate_entropy/approximate_entropy.jl b/src/approximate_entropy/approximate_entropy.jl
deleted file mode 100644
index 00bdfa315..000000000
--- a/src/approximate_entropy/approximate_entropy.jl
+++ /dev/null
@@ -1,111 +0,0 @@
-using DelayEmbeddings
-using Neighborhood
-using StatsBase
-using Distances
-
-export ApproximateEntropy
-export approx_entropy
-
-"""
-    ApproximateEntropy(; r = 0.1, m::Int = 2)
-
-Compute the approximate entropy (ApEn; Pincus, 1991)[^Pincus1991] of a scalar-valued
-time series, using embedding dimension (pattern length) `m` and tolerance radius `r`.
-
-The tolerance radius `r` should be determined from the input data, for example as
-some fraction of the standard deviation of the input.
-
-!!! info
-    This estimator is only available for entropy estimation.
-    Probabilities cannot be obtained directly.
-
-[^Pincus1991]: Pincus, S. M. (1991). Approximate entropy as a measure of system complexity. Proceedings of the National Academy of Sciences, 88(6), 2297-2301.
-"""
-Base.@kwdef struct ApproximateEntropy <: EntropyEstimator
-    r::Real = 0.1
-    m::Int = 2
-end
-
-function genentropy(x::AbstractVector, est::ApproximateEntropy; base = MathConstants.e)
-    est.m >= 1 || throw(ArgumentError("m must be >= 1. Got m=$(est.m)."))
-    N = length(x)
-    ϕᵐ = ϕmr_tree(x, est.m, est.r, N, base = base)
-    ϕᵐ⁺¹ = ϕmr_tree(x, est.m + 1, est.r, N, base = base)
-    return ϕᵐ - ϕᵐ⁺¹
-end
-
-function genentropy(x::AbstractDataset, est::ApproximateEntropy)
-    throw(
-        ArgumentError("Approximate entropy is currently not defined for multivariate data.")
-    )
-end
-
-
-"""
-    approx_entropy(x::AbstractVector ; m = 2, r = 0.2*StatsBase.std(x), base = MathConstants.e)
-Shorthand for `genentropy(x, ApproximateEntropy(m = m, r = r); base = base)` which
-calculates the approximate entropy of for dimension `m` with tolerance radius `r`.
-
-If `r` is not specified, `r` is estimated to a recommended value based on the standard
-deviation of `x`.
-"""
-function approx_entropy(x; m = 2, r = 0.2*StatsBase.std(x), base = MathConstants.e)
-    est = ApproximateEntropy(r, m)
-    genentropy(x, est; base = base)
-end
-
-"""
-    ϕmr_tree(x::AbstractVector{T}, k, r, N; base = MathConstants.e)
-
-Construct the embedding
-
-```math
-u = \\{{\\bf u}_n \\}_{n = 1}^{N - k + 1} =
-\\{[x(i), x(i + 1), \\ldots, x(i + k - 1)]\\}_{n = 1}^{N - k + 1}
-```
-
-and use a tree-and-nearest-neighbor search approach to compute
-
-```math
-\\phi^k(r) = \\dfrac{1}{N - k + 1} \\sum_{i}^{N - k + 1} \\log_{b}{(C_i^k(r))},
-```
-
-taking logarithms to `base` ``b``, and where
-
-```math
-C_i^k(r) = \\textrm{number of } j \\textrm{ such that } d({\\bf u}_i, {\\bf u}_j) < r,
-```
-
-where ``d`` is the maximum (Chebyshev) distance, `r` is the tolerance, and `N` is the
-length of the original scalar-valued time series `x`.
-
-## Examples
-
-```jldoctest; setup = :(using Entropies)
-julia> x = repeat([2, 3, 4, 8, 8, 9, 3, 6, 4, 5, 8, 7], 100);
-
-julia> approx_entropy(x, m = 2, r = 3)
-0.42020216510849995
-```
-"""
-function ϕmr_tree(x::AbstractVector{T}, k, r, N::Int;
-        base = MathConstants.e) where T <: Real
-    pts = genembed(x, 0:(k - 1))
-
-    # For each xᵢ ∈ pts, find the indices of neighbors within radius `r` according to
-    # the Chebyshev (maximum) metric. Self-inclusion is wanted, so no Theiler window
-    # is specified. This tree-approach is more than an order of magnitude
-    # faster than checking point-by-point using nested loops, and also allocates
-    # more than an order of magnitude less.
-    tree = KDTree(pts, Chebyshev())
-    idxs_of_neighbors = bulkisearch(tree, pts, WithinRange(r))
-
-    q = N - k + 1
-    ϕ = 0.0
-    for nn_idxsᵢ in idxs_of_neighbors
-        # Self-inclusion during the neighbor search means that there are always neighbors,
-        # so we never encounter log(0) here.
-        ϕ += log(base, length(nn_idxsᵢ) / q)
-    end
-    return ϕ / q
-end
diff --git a/src/complexity_measures/approximate_entropy.jl b/src/complexity_measures/approximate_entropy.jl
new file mode 100644
index 000000000..25123ef51
--- /dev/null
+++ b/src/complexity_measures/approximate_entropy.jl
@@ -0,0 +1,91 @@
+using DelayEmbeddings
+using Neighborhood
+using Distances
+
+export approx_entropy
+
+# So we don't need to import StatsBase
+_std(x) = sqrt(sum((x .- mean(x)) .^ 2) / (length(x) - 1))
+
+"""
+    approx_entropy(x; m = 2, τ = 1, r = 0.2 * StatsBase.std(x), base = MathConstants.e)
+
+Compute the approximate entropy (ApEn; Pincus, 1991)[^Pincus1991] of a univariate
+timeseries `x`, using embedding dimension (pattern length) `m`, embedding lag `τ`,
+and tolerance radius `r`.
+
+The tolerance radius `r` should be determined from the input data, for example as
+some fraction of the standard deviation of the input.
+
+[^Pincus1991]: Pincus, S. M. (1991). Approximate entropy as a measure of system complexity.
+    Proceedings of the National Academy of Sciences, 88(6), 2297-2301.
+"""
+function approx_entropy(x::AbstractVector{T}; m::Int = 2, τ::Int = 1,
+        r = 0.2 * _std(x), base = MathConstants.e) where T
+    m >= 1 || throw(ArgumentError("m must be >= 1. Got m=$(m)."))
+
+    # Definition in https://www.ncbi.nlm.nih.gov/pmc/articles/PMC7515030/
+    if m == 1
+       return compute_ϕ(x, r; m, τ, base)
+    else
+        ϕᵐ = compute_ϕ(x, r; m, τ, base)
+        ϕᵐ⁺¹ = compute_ϕ(x, r; m = m + 1, τ, base)
+
+        return ϕᵐ - ϕᵐ⁺¹
+    end
+end
+
+approx_entropy(x; kwargs...) =
+    throw(ArgumentError(
+        "`approx_entropy` not implemented for input data of type $(typeof(x))"
+    ))
+
+"""
+    compute_ϕ(x, r; m::Int = 2, τ::Int = 1, base = MathConstants.e)
+
+Construct the embedding
+
+```math
+u = \\{{\\bf u}_n \\}_{n = 1}^{N - k + 1} =
+\\{[x(i), x(i + 1), \\ldots, x(i + k - 1)]\\}_{n = 1}^{N - k + 1}
+```
+
+and use a tree-and-nearest-neighbor search approach to compute
+
+```math
+\\phi^k(r) = \\dfrac{1}{N - k + 1} \\sum_{i}^{N - k + 1} \\log_{b}{(C_i^k(r))},
+```
+
+taking logarithms to `base` ``b``, and where
+
+```math
+C_i^k(r) = \\textrm{number of } j \\textrm{ such that } d({\\bf u}_i, {\\bf u}_j) < r,
+```
+
+where ``d`` is the maximum (Chebyshev) distance, `r` is the tolerance, and `N` is the
+length of the original scalar-valued time series `x`.
+"""
+function compute_ϕ(x::AbstractVector{T}, r; m::Int = 2, τ::Int = 1,
+        base = MathConstants.e) where T <: Real
+
+    τs = 0:-τ:-(m - 1)*τ
+    pts = genembed(x, τs)
+    f = length(pts)
+
+    # For each xᵢ ∈ pts, find the indices of neighbors within radius `r` according to
+    # the Chebyshev (maximum) metric. Self-inclusion is wanted, so no Theiler window
+    # is specified. This tree-approach is more than an order of magnitude
+    # faster than checking point-by-point using nested loops, and also allocates
+    # more than an order of magnitude less.
+    tree = KDTree(pts, Chebyshev())
+    idxs_of_neighbors = bulkisearch(tree, pts, WithinRange(r))
+
+    ϕ = 0.0
+    for nn_idxsᵢ in idxs_of_neighbors
+        # Self-inclusion during the neighbor search means that there are always neighbors,
+        # so we never encounter log(0) here.
+        Cᵢᵐ = length(nn_idxsᵢ)
+        ϕ += log(base, Cᵢᵐ / f)
+    end
+    return ϕ / f
+end
diff --git a/src/complexity_measures/complexity_measures.jl b/src/complexity_measures/complexity_measures.jl
index a3719a54e..8820d4f91 100644
--- a/src/complexity_measures/complexity_measures.jl
+++ b/src/complexity_measures/complexity_measures.jl
@@ -1 +1,2 @@
 include("reverse_dispersion_entropy.jl")
+include("approximate_entropy.jl")
diff --git a/test/complexity_measures/approx_entropy.jl b/test/complexity_measures/approx_entropy.jl
new file mode 100644
index 000000000..9018e1d9d
--- /dev/null
+++ b/test/complexity_measures/approx_entropy.jl
@@ -0,0 +1,16 @@
+
+@test_throws ArgumentError approx_entropy(Dataset(rand(100, 3)))
+
+# "Analytical tests" - compare with results from the `approximateEntropy` function
+# in MATLAB (which uses natural logarithms)
+x = [0, 1, 0, 1, 0, 1, 0] # A regular signal that should have approx entropy close to zero
+m, τ, r  = 2, 2, 0.5
+res_x_matlab = -0.036497498714443
+res_x_ent = approx_entropy(x; r, m, τ, base = MathConstants.e)
+@test round(res_x_ent, digits = 5) == round(res_x_matlab, digits = 5)
+
+y = repeat([collect(-0.5:0.1:0.5); collect(0.4:-0.1:-0.4)], 5)
+m, τ, r = 2, 1, 0.3
+res_y_matlab = 0.195753687224351
+res_y_ent = approx_entropy(y; r, m, τ, base = MathConstants.e)
+@test round(res_y_ent, digits = 5) == round(res_y_matlab, digits = 5)
diff --git a/test/complexity_measures/complexity_measures.jl b/test/complexity_measures/complexity_measures.jl
index 21d193bd0..1be298fe6 100644
--- a/test/complexity_measures/complexity_measures.jl
+++ b/test/complexity_measures/complexity_measures.jl
@@ -2,4 +2,5 @@ using Entropies, Test
 
 @testset "Complexity measures" begin
     testfile("./reverse_dispersion.jl")
+    testfile("./approx_entropy.jl")
 end

From 566f787fb7d14d8c88d97a20b78b1cf31295d887 Mon Sep 17 00:00:00 2001
From: Kristian Haaga <kahaaga@gmail.com>
Date: Fri, 30 Sep 2022 18:59:08 +0200
Subject: [PATCH 04/20] Update example

---
 docs/src/ApproximateEntropy.md | 84 ----------------------------------
 docs/src/examples.md           | 80 +++++++++++++++++++++++++++++++-
 2 files changed, 79 insertions(+), 85 deletions(-)
 delete mode 100644 docs/src/ApproximateEntropy.md

diff --git a/docs/src/ApproximateEntropy.md b/docs/src/ApproximateEntropy.md
deleted file mode 100644
index c37f42df4..000000000
--- a/docs/src/ApproximateEntropy.md
+++ /dev/null
@@ -1,84 +0,0 @@
-# Approximate entropy
-
-```@docs
-ApproximateEntropy
-approx_entropy
-```
-
-## Example
-
-Here, we reproduce the Henon map example with ``R=0.8`` from Pincus (1991),
-comparing our values with relevant values from table 1 in Pincus (1991).
-
-We use `DiscreteDynamicalSystem` from `DynamicalSystems.jl` to represent the map,
-and use the `trajectory` function from the same package to iterate the map
-for different initial conditions, for multiple time series lengths.
-
-Finally, we summarize our results in box plots and compare the values to those
-obtained by Pincus (1991).
-
-```@example
-using DynamicalSystems, StatsBase, PyPlot
-
-# Equation 13 in Pincus (1991)
-function eom_henon(x, p, n)
-    R = p[1]
-    x, y = (x...,)
-    dx = R*y + 1 - 1.4*x^2
-    dy = 0.3*R*x
-
-    return SVector{2}(dx, dy)
-end
-
-function henon(; u₀ = rand(2), R = 0.8)
-    DiscreteDynamicalSystem(eom_henon, u₀, [R])
-end
-
-ts_lengths = [300, 1000, 2000, 3000]
-nreps = 100
-apens_08 = [zeros(nreps) for i = 1:length(ts_lengths)]
-
-# For some initial conditions, the Henon map as specified here blows up,
-# so we need to check for infinite values.
-containsinf(x) = any(isinf.(x))
-
-for (i, L) in enumerate(ts_lengths)
-    k = 1
-    while k <= nreps
-        sys = henon(u₀ = rand(2), R = 0.8)
-        t = trajectory(sys, L, Ttr = 5000)
-
-        if !any([containsinf(tᵢ) for tᵢ in t])
-            x, y = columns(t)
-            apen = approx_entropy(x, r = 0.05, m = 2)
-            apens_08[i][k] = apen
-            k += 1
-        end
-    end
-end
-
-f = figure(figsize = (6, 5))
-
-# First subplot is an example time series
-sys = henon(u₀ = [0.5, 0.1], R = 0.8)
-x, y = columns(trajectory(sys, 100, Ttr = 500))
-
-subplot(211)
-plot(x, label = "x")
-plot(y, label = "y")
-xlabel("Time (t)")
-ylabel("Value")
-
-# Second subplot contains the approximate entropy values
-subplot(212)
-boxplot(apens_08, positions = ts_lengths, widths = [150, 150, 150, 150], notch = true)
-scatter(ts_lengths, [0.337, 0.385, NaN, 0.394], label = "Pincus (1991)")
-xlabel("Time series length (L)")
-ylabel("ApEn(m = 2, r = 0.05)")
-
-legend()
-tight_layout()
-savefig("approx_entropy_pincus.png") # hide
-```
-
-![Approximate entropy](approx_entropy_pincus.png)
diff --git a/docs/src/examples.md b/docs/src/examples.md
index a8376b9cf..446adbc0a 100644
--- a/docs/src/examples.md
+++ b/docs/src/examples.md
@@ -225,4 +225,82 @@ end
 ```
 
 You see that while the direct entropy values of the chaotic and noisy signals change massively with `N` but they are almost the same for the normalized version.
-For the regular signals, the entropy decreases nevertheless because the noise contribution of the Fourier computation becomes less significant.
\ No newline at end of file
+For the regular signals, the entropy decreases nevertheless because the noise contribution of the Fourier computation becomes less significant.
+
+## Approximate entropy
+
+Here, we reproduce the Henon map example with ``R=0.8`` from Pincus (1991),
+comparing our values with relevant values from table 1 in Pincus (1991).
+
+We use `DiscreteDynamicalSystem` from `DynamicalSystems.jl` to represent the map,
+and use the `trajectory` function from the same package to iterate the map
+for different initial conditions, for multiple time series lengths.
+
+Finally, we summarize our results in box plots and compare the values to those
+obtained by Pincus (1991).
+
+```@example
+using Entropies, DynamicalSystems, CairoMakie
+
+# Equation 13 in Pincus (1991)
+function eom_henon(x, p, n)
+    R = p[1]
+    x, y = (x...,)
+    dx = R*y + 1 - 1.4*x^2
+    dy = 0.3*R*x
+
+    return SVector{2}(dx, dy)
+end
+
+function henon(; u₀ = rand(2), R = 0.8)
+    DiscreteDynamicalSystem(eom_henon, u₀, [R])
+end
+
+ts_lengths = [300, 1000, 2000, 3000]
+nreps = 100
+apens_08 = [zeros(nreps) for i = 1:length(ts_lengths)]
+
+# For some initial conditions, the Henon map as specified here blows up,
+# so we need to check for infinite values.
+containsinf(x) = any(isinf.(x))
+
+for (i, L) in enumerate(ts_lengths)
+    k = 1
+    while k <= nreps
+        sys = henon(u₀ = rand(2), R = 0.8)
+        t = trajectory(sys, L, Ttr = 5000)
+
+        if !any([containsinf(tᵢ) for tᵢ in t])
+            x, y = columns(t)
+            apen = approx_entropy(x, r = 0.05, m = 2)
+            apens_08[i][k] = apen
+            k += 1
+        end
+    end
+end
+
+fig = Figure()
+
+# Example time series
+a1 = Axis(fig[1,1]; xlabel = "Time (t)", ylabel = "Value")
+sys = henon(u₀ = [0.5, 0.1], R = 0.8)
+x, y = columns(trajectory(sys, 100, Ttr = 500))
+lines!(a1, 1:length(x), x, label = "x")
+lines!(a1, 1:length(y), y, label = "y")
+
+# Approximate entropy values, compared to those of the original paper (black dots).
+a2 = Axis(fig[2, 1]; 
+    xlabel = "Time series length (L)", 
+    ylabel = "ApEn(m = 2, r = 0.05)")
+
+# hacky boxplot, but this seems to be how it's done in Makie at the moment
+n = length(ts_lengths)
+for i = 1:n
+    boxplot!(a2, fill(ts_lengths[i], n), apens_08[i];
+        width = 200)
+end
+
+scatter!(a2, ts_lengths, [0.337, 0.385, NaN, 0.394]; 
+    label = "Pincus (1991)", color = :black)
+fig
+```

From 582f127a4583385a8442a10734dd484cbc30a0c6 Mon Sep 17 00:00:00 2001
From: Kristian Haaga <kahaaga@gmail.com>
Date: Fri, 30 Sep 2022 19:02:11 +0200
Subject: [PATCH 05/20] Just use Statistics

---
 src/complexity_measures/approximate_entropy.jl | 8 +++-----
 1 file changed, 3 insertions(+), 5 deletions(-)

diff --git a/src/complexity_measures/approximate_entropy.jl b/src/complexity_measures/approximate_entropy.jl
index 25123ef51..be6c8c54f 100644
--- a/src/complexity_measures/approximate_entropy.jl
+++ b/src/complexity_measures/approximate_entropy.jl
@@ -1,14 +1,12 @@
 using DelayEmbeddings
 using Neighborhood
 using Distances
+using Statistics
 
 export approx_entropy
 
-# So we don't need to import StatsBase
-_std(x) = sqrt(sum((x .- mean(x)) .^ 2) / (length(x) - 1))
-
 """
-    approx_entropy(x; m = 2, τ = 1, r = 0.2 * StatsBase.std(x), base = MathConstants.e)
+    approx_entropy(x; m = 2, τ = 1, r = 0.2 * Statistics.std(x), base = MathConstants.e)
 
 Compute the approximate entropy (ApEn; Pincus, 1991)[^Pincus1991] of a univariate
 timeseries `x`, using embedding dimension (pattern length) `m`, embedding lag `τ`,
@@ -21,7 +19,7 @@ some fraction of the standard deviation of the input.
     Proceedings of the National Academy of Sciences, 88(6), 2297-2301.
 """
 function approx_entropy(x::AbstractVector{T}; m::Int = 2, τ::Int = 1,
-        r = 0.2 * _std(x), base = MathConstants.e) where T
+        r = 0.2 * std(x), base = MathConstants.e) where T
     m >= 1 || throw(ArgumentError("m must be >= 1. Got m=$(m)."))
 
     # Definition in https://www.ncbi.nlm.nih.gov/pmc/articles/PMC7515030/

From f48e1b623f8670435bd0d37c9ff8716e3ef62642 Mon Sep 17 00:00:00 2001
From: Kristian Haaga <kahaaga@gmail.com>
Date: Mon, 3 Oct 2022 09:43:32 +0200
Subject: [PATCH 06/20] Test approx_entropy using Pincus' example

---
 test/complexity_measures/approx_entropy.jl | 50 ++++++++++++++++------
 1 file changed, 38 insertions(+), 12 deletions(-)

diff --git a/test/complexity_measures/approx_entropy.jl b/test/complexity_measures/approx_entropy.jl
index 9018e1d9d..1be8e8454 100644
--- a/test/complexity_measures/approx_entropy.jl
+++ b/test/complexity_measures/approx_entropy.jl
@@ -1,16 +1,42 @@
+using DynamicalSystemsBase
+using Statistics
 
 @test_throws ArgumentError approx_entropy(Dataset(rand(100, 3)))
 
-# "Analytical tests" - compare with results from the `approximateEntropy` function
-# in MATLAB (which uses natural logarithms)
-x = [0, 1, 0, 1, 0, 1, 0] # A regular signal that should have approx entropy close to zero
-m, τ, r  = 2, 2, 0.5
-res_x_matlab = -0.036497498714443
-res_x_ent = approx_entropy(x; r, m, τ, base = MathConstants.e)
-@test round(res_x_ent, digits = 5) == round(res_x_matlab, digits = 5)
+# Here, we try to reproduce Pincus' results within reasonable accuracy
+# for the Henon map. He doesn't give initial conditions, so we just check that our
+#  results +- 1σ approaches what he gets for this system for time series length 1000).
+# ---------------------------------------------------------------------------
+# Equation 13 in Pincus (1991)
+function eom_henon(u, p, n)
+    R = p[1]
+    x, y = (u...,)
+    dx = R*y + 1 - 1.4*x^2
+    dy = 0.3*R*x
 
-y = repeat([collect(-0.5:0.1:0.5); collect(0.4:-0.1:-0.4)], 5)
-m, τ, r = 2, 1, 0.3
-res_y_matlab = 0.195753687224351
-res_y_ent = approx_entropy(y; r, m, τ, base = MathConstants.e)
-@test round(res_y_ent, digits = 5) == round(res_y_matlab, digits = 5)
+    return SVector{2}(dx, dy)
+end
+henon(; u₀ = rand(2), R = 0.8) = DiscreteDynamicalSystem(eom_henon, u₀, [R])
+
+# For some initial conditions, the Henon map as specified here blows up,
+# so we need to check for infinite values.
+containsinf(x) = any(isinf.(x))
+
+function calculate_hs(; nreps = 50, L = 1000)
+    # Calculate approx entropy for 50 different initial conditions
+    hs = zeros(nreps)
+    k = 1
+    while k <= nreps
+        sys = henon(u₀ = rand(2), R = 0.8)
+        t = trajectory(sys, L, Ttr = 5000)
+
+        if !any([containsinf(tᵢ) for tᵢ in t])
+            x = t[:, 1]
+            hs[k] = approx_entropy(x, r = 0.05, m = 2)
+            k += 1
+        end
+    end
+    return hs
+end
+hs = calculate_hs()
+@test mean(hs) - std(hs) <= 0.385 <= mean(hs) + std(hs)

From 80b735e051d25550be07bf50bb8cb70d791a5183 Mon Sep 17 00:00:00 2001
From: Kristian Haaga <kahaaga@gmail.com>
Date: Mon, 3 Oct 2022 09:44:23 +0200
Subject: [PATCH 07/20] Use Project.toml for tests

This is the standard for Julia 1.3 an onwards.
---
 test/Project.toml | 5 +++++
 1 file changed, 5 insertions(+)
 create mode 100644 test/Project.toml

diff --git a/test/Project.toml b/test/Project.toml
new file mode 100644
index 000000000..4e6f0385b
--- /dev/null
+++ b/test/Project.toml
@@ -0,0 +1,5 @@
+[deps]
+DynamicalSystemsBase = "6e36e845-645a-534a-86f2-f5d4aa5a06b4"
+Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
+Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"
+Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"

From 8331e34a0b29ab59c74783e5391b788f6e964c61 Mon Sep 17 00:00:00 2001
From: Kristian Haaga <kahaaga@gmail.com>
Date: Mon, 3 Oct 2022 09:44:34 +0200
Subject: [PATCH 08/20] Dynamical system notation

---
 docs/src/examples.md | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/docs/src/examples.md b/docs/src/examples.md
index 446adbc0a..23627e1e0 100644
--- a/docs/src/examples.md
+++ b/docs/src/examples.md
@@ -243,9 +243,9 @@ obtained by Pincus (1991).
 using Entropies, DynamicalSystems, CairoMakie
 
 # Equation 13 in Pincus (1991)
-function eom_henon(x, p, n)
+function eom_henon(u, p, n)
     R = p[1]
-    x, y = (x...,)
+    x, y = u
     dx = R*y + 1 - 1.4*x^2
     dy = 0.3*R*x
 

From 432e246437820c078dc22886ea9cfc48c4b52bd5 Mon Sep 17 00:00:00 2001
From: Kristian Haaga <kahaaga@gmail.com>
Date: Mon, 3 Oct 2022 09:53:15 +0200
Subject: [PATCH 09/20] Update docstring with description

---
 .../approximate_entropy.jl                    | 81 ++++++++++++++-----
 1 file changed, 60 insertions(+), 21 deletions(-)

diff --git a/src/complexity_measures/approximate_entropy.jl b/src/complexity_measures/approximate_entropy.jl
index be6c8c54f..da612c4bb 100644
--- a/src/complexity_measures/approximate_entropy.jl
+++ b/src/complexity_measures/approximate_entropy.jl
@@ -9,11 +9,49 @@ export approx_entropy
     approx_entropy(x; m = 2, τ = 1, r = 0.2 * Statistics.std(x), base = MathConstants.e)
 
 Compute the approximate entropy (ApEn; Pincus, 1991)[^Pincus1991] of a univariate
-timeseries `x`, using embedding dimension (pattern length) `m`, embedding lag `τ`,
-and tolerance radius `r`.
+timeseries `x`, by first embedding `x` using embedding dimension `m` and embedding lag `τ`,
+then searching for similar vectors within tolerance radius `r`, using the estimator
+described below, with logarithms to the given `base` (natural logarithm is used in
+Pincus, 1991).
 
-The tolerance radius `r` should be determined from the input data, for example as
-some fraction of the standard deviation of the input.
+`r` should be determined from the input data, for example as some fraction of the
+standard deviation of the input.
+
+## Description
+
+Approximate entropy is defined as
+
+```math
+ApEn(m ,r) = \\lim_{N \\to \\infty} \\left[ \\phi(x, m, r) - \\phi(x, m + 1, r) \\right].
+```
+
+For a finite-length timeseries `x`, an estimator for ``ApEn(m ,r)`` is
+
+```math
+ApEn(m, r, N) = \\phi(x, m, r, N) -  \\phi(x, m + 1, r, N),
+```
+
+where `N = length(x)` and
+
+```math
+\\phi(x, k, r, N) =
+\\dfrac{1}{N-(k-1)\\tau} \\sum_{i=1}^{N - (k-1)\\tau}
+\\log{\\left(
+    \\sum_{j = 1}^{N-(k-1)\\tau} \\dfrac{\\theta(d({\\bf x}_i^m, {\\bf x}_j^m) \\leq r)}{N-(k-1)\\tau}
+    \\right)}.
+```
+
+Here, ``\\theta(\\cdot)`` returns 1 if the argument is true and 0 otherwise,
+ ``d({\\bf x}_i, {\\bf x}_j)`` returns the Chebyshev distance between vectors
+ ``{\\bf x}_i`` and ``{\\bf x}_j``, and the `k`-dimensional embedding vectors are
+constructed from the input timeseries ``x(t)`` as
+
+```math
+{\\bf x}_i^k = (x(i), x(i+τ), x(i+2τ), \\ldots, x(i+(k-1)\\tau)).
+```
+
+Note: in the original paper, they fix `τ = 1`; in our implementation, the normalization
+constant is modified to account for embeddings with `τ != 1`.
 
 [^Pincus1991]: Pincus, S. M. (1991). Approximate entropy as a measure of system complexity.
     Proceedings of the National Academy of Sciences, 88(6), 2297-2301.
@@ -24,10 +62,10 @@ function approx_entropy(x::AbstractVector{T}; m::Int = 2, τ::Int = 1,
 
     # Definition in https://www.ncbi.nlm.nih.gov/pmc/articles/PMC7515030/
     if m == 1
-       return compute_ϕ(x, r; m, τ, base)
+       return compute_ϕ(x; k = m, r, τ, base)
     else
-        ϕᵐ = compute_ϕ(x, r; m, τ, base)
-        ϕᵐ⁺¹ = compute_ϕ(x, r; m = m + 1, τ, base)
+        ϕᵐ = compute_ϕ(x; k = m, r, τ, base)
+        ϕᵐ⁺¹ = compute_ϕ(x; k = m + 1, r, τ, base)
 
         return ϕᵐ - ϕᵐ⁺¹
     end
@@ -39,7 +77,8 @@ approx_entropy(x; kwargs...) =
     ))
 
 """
-    compute_ϕ(x, r; m::Int = 2, τ::Int = 1, base = MathConstants.e)
+    compute_ϕ(x::AbstractVector{T}; r = 0.2 * Statistics.std(x), k::Int = 2,
+        τ::Int = 1, base = MathConstants.e) where T <: Real
 
 Construct the embedding
 
@@ -51,7 +90,7 @@ u = \\{{\\bf u}_n \\}_{n = 1}^{N - k + 1} =
 and use a tree-and-nearest-neighbor search approach to compute
 
 ```math
-\\phi^k(r) = \\dfrac{1}{N - k + 1} \\sum_{i}^{N - k + 1} \\log_{b}{(C_i^k(r))},
+\\phi^k(r) = \\dfrac{1}{N - kτ + 1} \\sum_{i}^{N - kτ + 1} \\log_{b}{(C_i^k(r))},
 ```
 
 taking logarithms to `base` ``b``, and where
@@ -63,12 +102,15 @@ C_i^k(r) = \\textrm{number of } j \\textrm{ such that } d({\\bf u}_i, {\\bf u}_j
 where ``d`` is the maximum (Chebyshev) distance, `r` is the tolerance, and `N` is the
 length of the original scalar-valued time series `x`.
 """
-function compute_ϕ(x::AbstractVector{T}, r; m::Int = 2, τ::Int = 1,
-        base = MathConstants.e) where T <: Real
-
-    τs = 0:-τ:-(m - 1)*τ
+function compute_ϕ(x::AbstractVector{T}; r = 0.2 * Statistics.std(x), k::Int = 2,
+        τ::Int = 1, base = MathConstants.e) where T <: Real
+    τs = 0:τ:(k - 1)*τ
     pts = genembed(x, τs)
-    f = length(pts)
+
+    N = length(x)
+    # The original normalization uses `τ = 1`. This needs to be adjusted to account for
+    # `τ != 1` in the embedding.
+    f = N - (k-1)*τ
 
     # For each xᵢ ∈ pts, find the indices of neighbors within radius `r` according to
     # the Chebyshev (maximum) metric. Self-inclusion is wanted, so no Theiler window
@@ -78,12 +120,9 @@ function compute_ϕ(x::AbstractVector{T}, r; m::Int = 2, τ::Int = 1,
     tree = KDTree(pts, Chebyshev())
     idxs_of_neighbors = bulkisearch(tree, pts, WithinRange(r))
 
-    ϕ = 0.0
-    for nn_idxsᵢ in idxs_of_neighbors
-        # Self-inclusion during the neighbor search means that there are always neighbors,
-        # so we never encounter log(0) here.
-        Cᵢᵐ = length(nn_idxsᵢ)
-        ϕ += log(base, Cᵢᵐ / f)
-    end
+    # Self-inclusion during the neighbor search means that there are always neighbors,
+    # so we never encounter log(0) here.
+    ϕ = sum(log(base, length(neighborsᵢ) / f) for neighborsᵢ in idxs_of_neighbors)
+
     return ϕ / f
 end

From e27372565ef13f73aa7d51366fca0acb510c4ba9 Mon Sep 17 00:00:00 2001
From: Kristian Haaga <kahaaga@gmail.com>
Date: Sat, 22 Oct 2022 15:42:37 +0200
Subject: [PATCH 10/20] New api

Shorten
---
 docs/src/complexity_measures.md               |  10 +-
 docs/src/examples.md                          |   5 +-
 src/complexity.jl                             |  13 ++-
 .../approximate_entropy.jl                    | 106 ++++++++++++------
 test/complexity_measures/approx_entropy.jl    |  12 +-
 5 files changed, 101 insertions(+), 45 deletions(-)

diff --git a/docs/src/complexity_measures.md b/docs/src/complexity_measures.md
index 510cdd9fd..2eb8d6a90 100644
--- a/docs/src/complexity_measures.md
+++ b/docs/src/complexity_measures.md
@@ -1,4 +1,4 @@
-# [Complexity measures API](@id complexity_measures)
+## [Complexity measures API](@id complexity_measures)
 
 ```@docs
 complexity
@@ -14,6 +14,14 @@ distance_to_whitenoise
 
 ## Approximate entropy
 
+```@docs
+ApproxEntropy
+```
+
+## Convenience functions
+
+We provide a few convenience functions for widely used "entropy-like" complexity measures, such as "approximate entropy". Other arbitrary specialized convenience functions can easily be defined in a couple lines of code.
+
 ```@docs
 approx_entropy
 ```
diff --git a/docs/src/examples.md b/docs/src/examples.md
index a82c00c56..6ae6fd990 100644
--- a/docs/src/examples.md
+++ b/docs/src/examples.md
@@ -320,6 +320,8 @@ apens_08 = [zeros(nreps) for i = 1:length(ts_lengths)]
 # so we need to check for infinite values.
 containsinf(x) = any(isinf.(x))
 
+c = ApproxEntropy(r = 0.05, m = 2)
+
 for (i, L) in enumerate(ts_lengths)
     k = 1
     while k <= nreps
@@ -328,8 +330,7 @@ for (i, L) in enumerate(ts_lengths)
 
         if !any([containsinf(tᵢ) for tᵢ in t])
             x, y = columns(t)
-            apen = approx_entropy(x, r = 0.05, m = 2)
-            apens_08[i][k] = apen
+            apens_08[i][k] = complexity(c, x)
             k += 1
         end
     end
diff --git a/src/complexity.jl b/src/complexity.jl
index 55ddaf193..a3fa113cf 100644
--- a/src/complexity.jl
+++ b/src/complexity.jl
@@ -16,8 +16,13 @@ Estimate the complexity measure `c` for [input data](@ref input_data) `x`, where
 be any of the following measures:
 
 - [`ReverseDispersion`](@ref).
+- [`ApproxEntropy`](@ref).
 """
-function complexity end
+function complexity(c::C, x) where C <: ComplexityMeasure
+    T = typeof(x)
+    msg = "`complexity` not implemented for $C and input data of type $T."
+    throw(ArgumentError(msg))
+end
 
 """
     complexity_normalized(c::ComplexityMeasure, x) → m ∈ [a, b]
@@ -25,8 +30,10 @@ function complexity end
 Estimate the [`complexity`](@ref) measure `c` for [input data](@ref input_data) `x`,
 normalized to the interval `[a, b]`, where `[a, b]` depends on `c`.
 """
-function complexity_normalized(c::C, args...; kwargs...) where {C <: ComplexityMeasure}
-    throw(ArgumentError("complexity_normalized not implemented for $C."))
+function complexity_normalized(c::C, x) where {C <: ComplexityMeasure}
+    T = typeof(x)
+    msg = "`complexity_normalized` not implemented for $C and input data of type $T."
+    throw(ArgumentError(msg))
 end
 
 include("complexity_measures/complexity_measures.jl")
diff --git a/src/complexity_measures/approximate_entropy.jl b/src/complexity_measures/approximate_entropy.jl
index da612c4bb..9022ee2ed 100644
--- a/src/complexity_measures/approximate_entropy.jl
+++ b/src/complexity_measures/approximate_entropy.jl
@@ -2,20 +2,28 @@ using DelayEmbeddings
 using Neighborhood
 using Distances
 using Statistics
+using Neighborhood.NearestNeighbors: inrangecount
 
+export ApproxEntropy
 export approx_entropy
 
 """
-    approx_entropy(x; m = 2, τ = 1, r = 0.2 * Statistics.std(x), base = MathConstants.e)
+    ApproxEntropy([x]; r = 0.2std(x), kwargs...)
+
+An estimator for the approximate entropy (ApEn; Pincus, 1991)[^Pincus1991] complexity
+measure, used with [`complexity`](@ref).
 
-Compute the approximate entropy (ApEn; Pincus, 1991)[^Pincus1991] of a univariate
-timeseries `x`, by first embedding `x` using embedding dimension `m` and embedding lag `τ`,
-then searching for similar vectors within tolerance radius `r`, using the estimator
-described below, with logarithms to the given `base` (natural logarithm is used in
-Pincus, 1991).
+The keyword argument `r` is mandatory if an input timeseries `x` is not provided.
 
-`r` should be determined from the input data, for example as some fraction of the
-standard deviation of the input.
+## Keyword arguments
+- `r::Real`: The radius used when querying for nearest neighbors around points. Its value
+    should be determined from the input data, for example as some proportion of the
+    standard deviation of the data.
+- `m::Int = 2`: The embedding dimension.
+- `τ::Int = 1`: The embedding lag.
+- `base::Real = MathConstants.e`: The base to use for the logarithm. Pincus (1991) uses the
+    natural logarithm.
+- `metric`: The metric used to compute distances.
 
 ## Description
 
@@ -25,7 +33,12 @@ Approximate entropy is defined as
 ApEn(m ,r) = \\lim_{N \\to \\infty} \\left[ \\phi(x, m, r) - \\phi(x, m + 1, r) \\right].
 ```
 
-For a finite-length timeseries `x`, an estimator for ``ApEn(m ,r)`` is
+Approximate entropy is estimated for a timeseries `x`, by first embedding `x` using
+embedding dimension `m` and embedding lag `τ`, then searching for similar vectors within
+tolerance radius `r`, using the estimator described below, with logarithms to the given
+`base` (natural logarithm is used in Pincus, 1991).
+
+Specifically, for a finite-length timeseries `x`, an estimator for ``ApEn(m ,r)`` is
 
 ```math
 ApEn(m, r, N) = \\phi(x, m, r, N) -  \\phi(x, m + 1, r, N),
@@ -50,15 +63,39 @@ constructed from the input timeseries ``x(t)`` as
 {\\bf x}_i^k = (x(i), x(i+τ), x(i+2τ), \\ldots, x(i+(k-1)\\tau)).
 ```
 
-Note: in the original paper, they fix `τ = 1`; in our implementation, the normalization
-constant is modified to account for embeddings with `τ != 1`.
+!!! note "Flexible embedding lag"
+    In the original paper, they fix `τ = 1`. In our implementation, the normalization
+    constant is modified to account for embeddings with `τ != 1`.
 
 [^Pincus1991]: Pincus, S. M. (1991). Approximate entropy as a measure of system complexity.
     Proceedings of the National Academy of Sciences, 88(6), 2297-2301.
 """
-function approx_entropy(x::AbstractVector{T}; m::Int = 2, τ::Int = 1,
-        r = 0.2 * std(x), base = MathConstants.e) where T
-    m >= 1 || throw(ArgumentError("m must be >= 1. Got m=$(m)."))
+Base.@kwdef struct ApproxEntropy{I, M, B, R} <: ComplexityMeasure
+    m::I = 2
+    τ::I = 1
+    metric::M = Chebyshev()
+    base::B = MathConstants.e
+    r::R
+
+    function ApproxEntropy(m::I, τ::I, r::R, base::B, metric::M) where {I, R, M, B}
+        m >= 1 || throw(ArgumentError("m must be >= 1. Got m=$(m)."))
+        r > 0 || throw(ArgumentError("r must be > 0. Got r=$(r)."))
+        new{I, M, B, R}(m, τ, metric, base, r)
+    end
+end
+
+function ApproxEntropy(x::AbstractVector{T}; m::Int = 2, τ::Int = 1, metric = Chebyshev(),
+        base = MathConstants.e) where T
+    r = 0.2 * Statistics.std(x)
+    ApproxEntropy(m, τ, r, base, metric)
+end
+
+function ApproxEntropy(; r, m::Int = 2, τ::Int = 1, metric = Chebyshev())
+    ApproxEntropy(m, τ, r, base, metric)
+end
+
+function complexity(c::ApproxEntropy, x::AbstractVector{T}) where T
+    (; m, τ, r, base) = c
 
     # Definition in https://www.ncbi.nlm.nih.gov/pmc/articles/PMC7515030/
     if m == 1
@@ -66,16 +103,10 @@ function approx_entropy(x::AbstractVector{T}; m::Int = 2, τ::Int = 1,
     else
         ϕᵐ = compute_ϕ(x; k = m, r, τ, base)
         ϕᵐ⁺¹ = compute_ϕ(x; k = m + 1, r, τ, base)
-
         return ϕᵐ - ϕᵐ⁺¹
     end
 end
 
-approx_entropy(x; kwargs...) =
-    throw(ArgumentError(
-        "`approx_entropy` not implemented for input data of type $(typeof(x))"
-    ))
-
 """
     compute_ϕ(x::AbstractVector{T}; r = 0.2 * Statistics.std(x), k::Int = 2,
         τ::Int = 1, base = MathConstants.e) where T <: Real
@@ -106,23 +137,28 @@ function compute_ϕ(x::AbstractVector{T}; r = 0.2 * Statistics.std(x), k::Int =
         τ::Int = 1, base = MathConstants.e) where T <: Real
     τs = 0:τ:(k - 1)*τ
     pts = genembed(x, τs)
-
-    N = length(x)
-    # The original normalization uses `τ = 1`. This needs to be adjusted to account for
-    # `τ != 1` in the embedding.
-    f = N - (k-1)*τ
-
-    # For each xᵢ ∈ pts, find the indices of neighbors within radius `r` according to
-    # the Chebyshev (maximum) metric. Self-inclusion is wanted, so no Theiler window
-    # is specified. This tree-approach is more than an order of magnitude
-    # faster than checking point-by-point using nested loops, and also allocates
-    # more than an order of magnitude less.
     tree = KDTree(pts, Chebyshev())
-    idxs_of_neighbors = bulkisearch(tree, pts, WithinRange(r))
 
-    # Self-inclusion during the neighbor search means that there are always neighbors,
-    # so we never encounter log(0) here.
-    ϕ = sum(log(base, length(neighborsᵢ) / f) for neighborsᵢ in idxs_of_neighbors)
+    # Account for `τ != 1` in the normalization constant.
+    f = length(x) - (k - 1)*τ
+
+    # `inrangecount` counts the query point itself, which is wanted for approximate entropy,
+    # because there are always neighbors and thus log(0) is never encountered.
+    ϕ = sum(log(base, inrangecount(tree, pᵢ, r) / f) for pᵢ in pts)
 
     return ϕ / f
 end
+
+"""
+    approx_entropy(x; m = 2, τ = 1, r = 0.2 * Statistics.std(x), base = MathConstants.e)
+
+Convenience syntax for computing the approximate entropy (Pincus, 1991) for timeseries `x`.
+
+This is just a wrapper for `complexity(ApproxEntropy(; m, τ, r, base), x)` (see
+also [`ApproxEntropy`](@ref)).
+"""
+function approx_entropy(x; m = 2, τ = 1, r = 0.2 * Statistics.std(x),
+         base = MathConstants.e)
+    c = ApproxEntropy(; m, τ, r, base)
+    return complexity(c, x)
+end
diff --git a/test/complexity_measures/approx_entropy.jl b/test/complexity_measures/approx_entropy.jl
index 1be8e8454..bd6305dfe 100644
--- a/test/complexity_measures/approx_entropy.jl
+++ b/test/complexity_measures/approx_entropy.jl
@@ -1,7 +1,7 @@
 using DynamicalSystemsBase
 using Statistics
 
-@test_throws ArgumentError approx_entropy(Dataset(rand(100, 3)))
+@test_throws ArgumentError complexity(ApproxEntropy(), Dataset(rand(100, 3)))
 
 # Here, we try to reproduce Pincus' results within reasonable accuracy
 # for the Henon map. He doesn't give initial conditions, so we just check that our
@@ -25,6 +25,7 @@ containsinf(x) = any(isinf.(x))
 function calculate_hs(; nreps = 50, L = 1000)
     # Calculate approx entropy for 50 different initial conditions
     hs = zeros(nreps)
+    hs_conv = zeros(nreps)
     k = 1
     while k <= nreps
         sys = henon(u₀ = rand(2), R = 0.8)
@@ -32,11 +33,14 @@ function calculate_hs(; nreps = 50, L = 1000)
 
         if !any([containsinf(tᵢ) for tᵢ in t])
             x = t[:, 1]
-            hs[k] = approx_entropy(x, r = 0.05, m = 2)
+            hs[k] = complexity(ApproxEntropy(r = 0.05, m = 2), x)
+            hs_conv[k] = approx_entropy(x, r = 0.05, m = 2)
             k += 1
         end
     end
-    return hs
+    return hs, hs_conv
 end
-hs = calculate_hs()
+hs, hs_conv = calculate_hs()
+
 @test mean(hs) - std(hs) <= 0.385 <= mean(hs) + std(hs)
+@test mean(hs_conv) - std(hs_conv) <= 0.385 <= mean(hs_conv) + std(hs_conv)

From df8729a520e8fc2c18b1190cfe91c6c784e40c5d Mon Sep 17 00:00:00 2001
From: Kristian Haaga <kahaaga@gmail.com>
Date: Sat, 22 Oct 2022 16:02:37 +0200
Subject: [PATCH 11/20] Update tests

---
 test/complexity_measures/approx_entropy.jl | 3 ++-
 1 file changed, 2 insertions(+), 1 deletion(-)

diff --git a/test/complexity_measures/approx_entropy.jl b/test/complexity_measures/approx_entropy.jl
index bd6305dfe..9982089d3 100644
--- a/test/complexity_measures/approx_entropy.jl
+++ b/test/complexity_measures/approx_entropy.jl
@@ -1,7 +1,8 @@
 using DynamicalSystemsBase
 using Statistics
 
-@test_throws ArgumentError complexity(ApproxEntropy(), Dataset(rand(100, 3)))
+@test_throws UndefKeywordError ApproxEntropy()
+@test_throws ArgumentError complexity(ApproxEntropy(r = 0.2), Dataset(rand(100, 3)))
 
 # Here, we try to reproduce Pincus' results within reasonable accuracy
 # for the Henon map. He doesn't give initial conditions, so we just check that our

From b239afc687e770675d174e0d751a842f4d355115 Mon Sep 17 00:00:00 2001
From: Kristian Haaga <kahaaga@gmail.com>
Date: Sat, 22 Oct 2022 16:32:12 +0200
Subject: [PATCH 12/20] Correct error

---
 test/complexity_measures/approx_entropy.jl | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/test/complexity_measures/approx_entropy.jl b/test/complexity_measures/approx_entropy.jl
index 9982089d3..e8dac57f3 100644
--- a/test/complexity_measures/approx_entropy.jl
+++ b/test/complexity_measures/approx_entropy.jl
@@ -2,7 +2,7 @@ using DynamicalSystemsBase
 using Statistics
 
 @test_throws UndefKeywordError ApproxEntropy()
-@test_throws ArgumentError complexity(ApproxEntropy(r = 0.2), Dataset(rand(100, 3)))
+@test_throws MethodError complexity(ApproxEntropy(r = 0.2), Dataset(rand(100, 3)))
 
 # Here, we try to reproduce Pincus' results within reasonable accuracy
 # for the Henon map. He doesn't give initial conditions, so we just check that our

From fc86eac207fa40843c16513fe507b743775687b9 Mon Sep 17 00:00:00 2001
From: Kristian Haaga <kahaaga@gmail.com>
Date: Sat, 22 Oct 2022 17:46:43 +0200
Subject: [PATCH 13/20] Use new version of Neighborhood.jl

---
 Project.toml                                   | 2 +-
 src/complexity_measures/approximate_entropy.jl | 3 +--
 2 files changed, 2 insertions(+), 3 deletions(-)

diff --git a/Project.toml b/Project.toml
index 448ddb1ea..02d403029 100644
--- a/Project.toml
+++ b/Project.toml
@@ -23,7 +23,7 @@ Wavelets = "29a6e085-ba6d-5f35-a997-948ac2efa89a"
 DelayEmbeddings = "2.3.4"
 Distances = "0.9, 0.10"
 FFTW = "1"
-Neighborhood = "0.1, 0.2"
+Neighborhood = "0.2.4"
 QuadGK = "2"
 Scratch = "1"
 SpecialFunctions = "0.10, 1.0, 2"
diff --git a/src/complexity_measures/approximate_entropy.jl b/src/complexity_measures/approximate_entropy.jl
index 9022ee2ed..df3055228 100644
--- a/src/complexity_measures/approximate_entropy.jl
+++ b/src/complexity_measures/approximate_entropy.jl
@@ -1,8 +1,7 @@
 using DelayEmbeddings
-using Neighborhood
+using Neighborhood: inrangecount
 using Distances
 using Statistics
-using Neighborhood.NearestNeighbors: inrangecount
 
 export ApproxEntropy
 export approx_entropy

From 74bde4d8ad7456f3548e56c5d4ede5a353026418 Mon Sep 17 00:00:00 2001
From: Kristian Haaga <kahaaga@gmail.com>
Date: Sat, 22 Oct 2022 22:40:21 +0200
Subject: [PATCH 14/20] Correct argument order

---
 src/complexity_measures/approximate_entropy.jl | 6 +++---
 1 file changed, 3 insertions(+), 3 deletions(-)

diff --git a/src/complexity_measures/approximate_entropy.jl b/src/complexity_measures/approximate_entropy.jl
index df3055228..2410a4750 100644
--- a/src/complexity_measures/approximate_entropy.jl
+++ b/src/complexity_measures/approximate_entropy.jl
@@ -76,7 +76,7 @@ Base.@kwdef struct ApproxEntropy{I, M, B, R} <: ComplexityMeasure
     base::B = MathConstants.e
     r::R
 
-    function ApproxEntropy(m::I, τ::I, r::R, base::B, metric::M) where {I, R, M, B}
+    function ApproxEntropy(m::I, τ::I, metric::M, base::B, r::R) where {I, R, M, B}
         m >= 1 || throw(ArgumentError("m must be >= 1. Got m=$(m)."))
         r > 0 || throw(ArgumentError("r must be > 0. Got r=$(r)."))
         new{I, M, B, R}(m, τ, metric, base, r)
@@ -86,11 +86,11 @@ end
 function ApproxEntropy(x::AbstractVector{T}; m::Int = 2, τ::Int = 1, metric = Chebyshev(),
         base = MathConstants.e) where T
     r = 0.2 * Statistics.std(x)
-    ApproxEntropy(m, τ, r, base, metric)
+    ApproxEntropy(m, τ, metric, base, r)
 end
 
 function ApproxEntropy(; r, m::Int = 2, τ::Int = 1, metric = Chebyshev())
-    ApproxEntropy(m, τ, r, base, metric)
+    ApproxEntropy(m, τ, metric, base, r)
 end
 
 function complexity(c::ApproxEntropy, x::AbstractVector{T}) where T

From 02cc5cc7573b1327ce3909e30171dc48d05de555 Mon Sep 17 00:00:00 2001
From: Kristian Haaga <kahaaga@gmail.com>
Date: Sat, 22 Oct 2022 22:41:53 +0200
Subject: [PATCH 15/20] Fix tests

---
 test/complexity_measures/approx_entropy.jl | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/test/complexity_measures/approx_entropy.jl b/test/complexity_measures/approx_entropy.jl
index e8dac57f3..9982089d3 100644
--- a/test/complexity_measures/approx_entropy.jl
+++ b/test/complexity_measures/approx_entropy.jl
@@ -2,7 +2,7 @@ using DynamicalSystemsBase
 using Statistics
 
 @test_throws UndefKeywordError ApproxEntropy()
-@test_throws MethodError complexity(ApproxEntropy(r = 0.2), Dataset(rand(100, 3)))
+@test_throws ArgumentError complexity(ApproxEntropy(r = 0.2), Dataset(rand(100, 3)))
 
 # Here, we try to reproduce Pincus' results within reasonable accuracy
 # for the Henon map. He doesn't give initial conditions, so we just check that our

From 35d381918e191f58aea4c7eb51cb6eaa52fb91d0 Mon Sep 17 00:00:00 2001
From: Kristian Haaga <kahaaga@gmail.com>
Date: Mon, 24 Oct 2022 00:03:45 +0200
Subject: [PATCH 16/20] Rename ApproxEntropy to ApproximateEntropy

---
 docs/src/complexity_measures.md               |  2 +-
 docs/src/examples.md                          |  2 +-
 src/complexity.jl                             |  2 +-
 .../approximate_entropy.jl                    | 24 +++++++++----------
 test/complexity_measures/approx_entropy.jl    |  6 ++---
 5 files changed, 18 insertions(+), 18 deletions(-)

diff --git a/docs/src/complexity_measures.md b/docs/src/complexity_measures.md
index 2eb8d6a90..68a57004c 100644
--- a/docs/src/complexity_measures.md
+++ b/docs/src/complexity_measures.md
@@ -15,7 +15,7 @@ distance_to_whitenoise
 ## Approximate entropy
 
 ```@docs
-ApproxEntropy
+ApproximateEntropy
 ```
 
 ## Convenience functions
diff --git a/docs/src/examples.md b/docs/src/examples.md
index 6ae6fd990..d1815e5bf 100644
--- a/docs/src/examples.md
+++ b/docs/src/examples.md
@@ -320,7 +320,7 @@ apens_08 = [zeros(nreps) for i = 1:length(ts_lengths)]
 # so we need to check for infinite values.
 containsinf(x) = any(isinf.(x))
 
-c = ApproxEntropy(r = 0.05, m = 2)
+c = ApproximateEntropy(r = 0.05, m = 2)
 
 for (i, L) in enumerate(ts_lengths)
     k = 1
diff --git a/src/complexity.jl b/src/complexity.jl
index a3fa113cf..08314bfdd 100644
--- a/src/complexity.jl
+++ b/src/complexity.jl
@@ -16,7 +16,7 @@ Estimate the complexity measure `c` for [input data](@ref input_data) `x`, where
 be any of the following measures:
 
 - [`ReverseDispersion`](@ref).
-- [`ApproxEntropy`](@ref).
+- [`ApproximateEntropy`](@ref).
 """
 function complexity(c::C, x) where C <: ComplexityMeasure
     T = typeof(x)
diff --git a/src/complexity_measures/approximate_entropy.jl b/src/complexity_measures/approximate_entropy.jl
index 2410a4750..c0c6ad9a9 100644
--- a/src/complexity_measures/approximate_entropy.jl
+++ b/src/complexity_measures/approximate_entropy.jl
@@ -3,11 +3,11 @@ using Neighborhood: inrangecount
 using Distances
 using Statistics
 
-export ApproxEntropy
+export ApproximateEntropy
 export approx_entropy
 
 """
-    ApproxEntropy([x]; r = 0.2std(x), kwargs...)
+    ApproximateEntropy([x]; r = 0.2std(x), kwargs...)
 
 An estimator for the approximate entropy (ApEn; Pincus, 1991)[^Pincus1991] complexity
 measure, used with [`complexity`](@ref).
@@ -69,31 +69,31 @@ constructed from the input timeseries ``x(t)`` as
 [^Pincus1991]: Pincus, S. M. (1991). Approximate entropy as a measure of system complexity.
     Proceedings of the National Academy of Sciences, 88(6), 2297-2301.
 """
-Base.@kwdef struct ApproxEntropy{I, M, B, R} <: ComplexityMeasure
+Base.@kwdef struct ApproximateEntropy{I, M, B, R} <: ComplexityMeasure
     m::I = 2
     τ::I = 1
     metric::M = Chebyshev()
     base::B = MathConstants.e
     r::R
 
-    function ApproxEntropy(m::I, τ::I, metric::M, base::B, r::R) where {I, R, M, B}
+    function ApproximateEntropy(m::I, τ::I, metric::M, base::B, r::R) where {I, R, M, B}
         m >= 1 || throw(ArgumentError("m must be >= 1. Got m=$(m)."))
         r > 0 || throw(ArgumentError("r must be > 0. Got r=$(r)."))
         new{I, M, B, R}(m, τ, metric, base, r)
     end
 end
 
-function ApproxEntropy(x::AbstractVector{T}; m::Int = 2, τ::Int = 1, metric = Chebyshev(),
+function ApproximateEntropy(x::AbstractVector{T}; m::Int = 2, τ::Int = 1, metric = Chebyshev(),
         base = MathConstants.e) where T
     r = 0.2 * Statistics.std(x)
-    ApproxEntropy(m, τ, metric, base, r)
+    ApproximateEntropy(m, τ, metric, base, r)
 end
 
-function ApproxEntropy(; r, m::Int = 2, τ::Int = 1, metric = Chebyshev())
-    ApproxEntropy(m, τ, metric, base, r)
+function ApproximateEntropy(; r, m::Int = 2, τ::Int = 1, metric = Chebyshev())
+    ApproximateEntropy(m, τ, metric, base, r)
 end
 
-function complexity(c::ApproxEntropy, x::AbstractVector{T}) where T
+function complexity(c::ApproximateEntropy, x::AbstractVector{T}) where T
     (; m, τ, r, base) = c
 
     # Definition in https://www.ncbi.nlm.nih.gov/pmc/articles/PMC7515030/
@@ -153,11 +153,11 @@ end
 
 Convenience syntax for computing the approximate entropy (Pincus, 1991) for timeseries `x`.
 
-This is just a wrapper for `complexity(ApproxEntropy(; m, τ, r, base), x)` (see
-also [`ApproxEntropy`](@ref)).
+This is just a wrapper for `complexity(ApproximateEntropy(; m, τ, r, base), x)` (see
+also [`ApproximateEntropy`](@ref)).
 """
 function approx_entropy(x; m = 2, τ = 1, r = 0.2 * Statistics.std(x),
          base = MathConstants.e)
-    c = ApproxEntropy(; m, τ, r, base)
+    c = ApproximateEntropy(; m, τ, r, base)
     return complexity(c, x)
 end
diff --git a/test/complexity_measures/approx_entropy.jl b/test/complexity_measures/approx_entropy.jl
index 9982089d3..9de2b5028 100644
--- a/test/complexity_measures/approx_entropy.jl
+++ b/test/complexity_measures/approx_entropy.jl
@@ -1,8 +1,8 @@
 using DynamicalSystemsBase
 using Statistics
 
-@test_throws UndefKeywordError ApproxEntropy()
-@test_throws ArgumentError complexity(ApproxEntropy(r = 0.2), Dataset(rand(100, 3)))
+@test_throws UndefKeywordError ApproximateEntropy()
+@test_throws ArgumentError complexity(ApproximateEntropy(r = 0.2), Dataset(rand(100, 3)))
 
 # Here, we try to reproduce Pincus' results within reasonable accuracy
 # for the Henon map. He doesn't give initial conditions, so we just check that our
@@ -34,7 +34,7 @@ function calculate_hs(; nreps = 50, L = 1000)
 
         if !any([containsinf(tᵢ) for tᵢ in t])
             x = t[:, 1]
-            hs[k] = complexity(ApproxEntropy(r = 0.05, m = 2), x)
+            hs[k] = complexity(ApproximateEntropy(r = 0.05, m = 2), x)
             hs_conv[k] = approx_entropy(x, r = 0.05, m = 2)
             k += 1
         end

From f7b9669f13f4aa6f6cfb6378b6f5deb1e00d4e0d Mon Sep 17 00:00:00 2001
From: Kristian Haaga <kahaaga@gmail.com>
Date: Mon, 24 Oct 2022 00:39:24 +0200
Subject: [PATCH 17/20] Fix constructor, which fixes doc example

---
 docs/src/examples.md                           |  3 ++-
 src/complexity_measures/approximate_entropy.jl | 13 +++++--------
 2 files changed, 7 insertions(+), 9 deletions(-)

diff --git a/docs/src/examples.md b/docs/src/examples.md
index d1815e5bf..cc816e88f 100644
--- a/docs/src/examples.md
+++ b/docs/src/examples.md
@@ -225,10 +225,11 @@ des = zeros(length(windows))
 pes = zeros(length(windows))
 
 m, c = 2, 6
+est_rd = ReverseDispersion(symbolization = GaussianSymbolization(c), m = m, τ = 1)
 est_de = Dispersion(symbolization = GaussianSymbolization(c), m = m, τ = 1)
 
 for (i, window) in enumerate(windows)
-    rdes[i] = reverse_dispersion(y[window], est_de; normalize = true)
+    rdes[i] = complexity_normalized(est_rd, y[window])
     des[i] = entropy_normalized(Renyi(), y[window], est_de)
 end
 
diff --git a/src/complexity_measures/approximate_entropy.jl b/src/complexity_measures/approximate_entropy.jl
index c0c6ad9a9..5fd491a46 100644
--- a/src/complexity_measures/approximate_entropy.jl
+++ b/src/complexity_measures/approximate_entropy.jl
@@ -81,17 +81,14 @@ Base.@kwdef struct ApproximateEntropy{I, M, B, R} <: ComplexityMeasure
         r > 0 || throw(ArgumentError("r must be > 0. Got r=$(r)."))
         new{I, M, B, R}(m, τ, metric, base, r)
     end
+    function ApproximateEntropy(x::AbstractVector{T}; m::Int = 2, τ::Int = 1,
+            metric = Chebyshev(), base = MathConstants.e) where T
+        r = 0.2 * Statistics.std(x)
+        ApproximateEntropy(m, τ, metric, base, r)
+    end
 end
 
-function ApproximateEntropy(x::AbstractVector{T}; m::Int = 2, τ::Int = 1, metric = Chebyshev(),
-        base = MathConstants.e) where T
-    r = 0.2 * Statistics.std(x)
-    ApproximateEntropy(m, τ, metric, base, r)
-end
 
-function ApproximateEntropy(; r, m::Int = 2, τ::Int = 1, metric = Chebyshev())
-    ApproximateEntropy(m, τ, metric, base, r)
-end
 
 function complexity(c::ApproximateEntropy, x::AbstractVector{T}) where T
     (; m, τ, r, base) = c

From 92226812589ab8cc186f813691a2c64806c0dbe8 Mon Sep 17 00:00:00 2001
From: Kristian Haaga <kahaaga@gmail.com>
Date: Mon, 24 Oct 2022 01:16:08 +0200
Subject: [PATCH 18/20] Missing parenthesis

---
 docs/src/examples.md | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/docs/src/examples.md b/docs/src/examples.md
index d19c36f8a..cffc79b3b 100644
--- a/docs/src/examples.md
+++ b/docs/src/examples.md
@@ -359,7 +359,7 @@ for i = 1:n
 end
 
 scatter!(a2, ts_lengths, [0.337, 0.385, NaN, 0.394]; 
-    label = "Pincus (1991)", color = :black
+    label = "Pincus (1991)", color = :black)
 fig
 ```
 

From 69f90914da7710228765bfcfce6ba6561a337f51 Mon Sep 17 00:00:00 2001
From: Kristian Haaga <kahaaga@gmail.com>
Date: Mon, 24 Oct 2022 01:22:19 +0200
Subject: [PATCH 19/20] Fix failing doc example due to erroneous SampleEntropy
 constructor

---
 src/complexity_measures/sample_entropy.jl | 18 +++++++++++-------
 1 file changed, 11 insertions(+), 7 deletions(-)

diff --git a/src/complexity_measures/sample_entropy.jl b/src/complexity_measures/sample_entropy.jl
index b8a4f186f..c2c286a05 100644
--- a/src/complexity_measures/sample_entropy.jl
+++ b/src/complexity_measures/sample_entropy.jl
@@ -72,14 +72,18 @@ Base.@kwdef struct SampleEntropy{I, R, M} <: ComplexityMeasure
     τ::I = 1
     metric::M = Chebyshev()
     r::R
-end
-function SampleEntropy(x::AbstractVector; m::Int = 2, τ::Int = 1, metric = Chebyshev())
-    r = 0.2 * Statistics.std(x)
-    SampleEntropy(; m, τ, metric, r)
-end
 
-function SampleEntropy(; r, m::Int = 2, τ::Int = 1, metric = Chebyshev())
-    SampleEntropy(; m, τ, metric, r)
+    function SampleEntropy(m::I, τ::I, metric::M, r::R) where {I, R, M, B}
+        m >= 1 || throw(ArgumentError("m must be >= 1. Got m=$(m)."))
+        r > 0 || throw(ArgumentError("r must be > 0. Got r=$(r)."))
+        new{I, R, M}(m, τ, metric, r)
+    end
+
+    function SampleEntropy(x::AbstractVector{T}; m::Int = 2, τ::Int = 1,
+            metric = Chebyshev()) where T
+        r = 0.2 * Statistics.std(x)
+        SampleEntropy(m, τ, metric, r)
+    end
 end
 
 # See comment in https://github.com/JuliaDynamics/Entropies.jl/pull/71 for why

From e461492dd87bd31b985800ef17de9d1f8903f245 Mon Sep 17 00:00:00 2001
From: Kristian Haaga <kahaaga@gmail.com>
Date: Mon, 24 Oct 2022 01:42:15 +0200
Subject: [PATCH 20/20] Remove rogue variable from example

---
 docs/src/examples.md | 1 -
 1 file changed, 1 deletion(-)

diff --git a/docs/src/examples.md b/docs/src/examples.md
index cffc79b3b..604b07437 100644
--- a/docs/src/examples.md
+++ b/docs/src/examples.md
@@ -408,7 +408,6 @@ x_periodic .= (x_periodic .- mean(x_periodic)) ./ std(x_periodic)
 rs = 10 .^ range(-1, 0, length = 30)
 base = 2
 m = 2
-c = 
 hs_U = [complexity_normalized(SampleEntropy(m = m, r = r), x_U) for r in rs]
 hs_N = [complexity_normalized(SampleEntropy(m = m, r = r), x_N) for r in rs]
 hs_periodic = [complexity_normalized(SampleEntropy(m = m, r = r), x_periodic) for r in rs]