Merge branch 'interface-fix' of https://github.com/jaksle/KernelDensi…

…ty.jl into interface-fix

Project.toml

@@ -4,15 +4,10 @@ authors = ["Simon Byrne and various contributors"]
 version = "0.6.8"
 
 [deps]
-BenchmarkTools = "6e4b80f9-dd63-53aa-95a3-0cdb28fa8baf"
 Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f"
 DocStringExtensions = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae"
 FFTW = "7a1cc6ca-52ef-59f5-83cd-3a7055c09341"
 Interpolations = "a98d9a8b-a2ab-59e6-89dd-64a1c18fca59"
-Loess = "4345ca2d-374a-55d4-8d30-97f9976e7612"
-LsqFit = "2fda8390-95c7-5789-9bda-21331edee243"
-Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
-Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"
 StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91"
 
 [compat]

src/interp.jl

@@ -1,6 +1,6 @@
 import Interpolations: interpolate, scale
 
-struct InterpKDE{K,I} <: AbstractKDE
+mutable struct InterpKDE{K,I} <: AbstractKDE
     kde::K
     itp::I
     InterpKDE{K,I}(kde::K, itp::I) where {K,I} = new{K,I}(kde, itp)
@@ -23,10 +23,20 @@ function InterpKDE(kde::BivariateKDE, opts...)
 end
 InterpKDE(kde::BivariateKDE) = InterpKDE(kde::BivariateKDE, BSpline(Quadratic(Line(OnGrid()))))
 
-pdf(k::UnivariateKDE,x) = pdf(InterpKDE(k),x)
+
+# pdf interface implementation
+# it should be consistent with Distributions.pdf
+
+# any dimension
 pdf(ik::InterpKDE,x::Real...) = ik.itp(x...)
-pdf(ik::InterpKDE,xs::AbstractVector) = [ik.itp(x) for x in xs]
-Base.broadcast(::typeof(pdf),k::UnivariateKDE,xs) = InterpKDE(k).itp.(xs)
+pdf(ik::InterpKDE, V::AbstractVector) = ik.itp(V...)
+pdf(ik::InterpKDE, M::AbstractArray{<:Real, N}) where N = pdf.(ik,eachslice(M, dims=ntuple(i->i+1, N-1)) )
+
+# 1 dimension
+pdf(k::UnivariateKDE,x) = pdf(InterpKDE(k),x)
+Base.Broadcast.broadcasted(::typeof(pdf),k::UnivariateKDE,xs) = Base.Broadcast.broadcasted(InterpKDE(k).itp, xs)
 
+# 2 dimensions
 pdf(k::BivariateKDE,x,y) = pdf(InterpKDE(k),x,y)
-pdf(ik::InterpKDE,xs::AbstractVector,ys::AbstractVector) = [ik.itp(x,y) for x in xs, y in ys]
+pdf(ik::InterpKDE,xs::AbstractVector,ys::AbstractVector) = ik.itp.(xs,ys')
+pdf(k::BivariateKDE, M) = pdf(InterpKDE(k),M)

src/univariate.jl

@@ -13,7 +13,7 @@ sum(density) * step(x) ≈ 1
 
 $(FIELDS)
 """
-struct UnivariateKDE{R<:AbstractRange} <: AbstractKDE
+mutable struct UnivariateKDE{R<:AbstractRange} <: AbstractKDE
     "Gridpoints for evaluating the density."
     x::R
     "Kernel density at corresponding gridpoints `x`."

test/interp.jl

@@ -1,25 +1,43 @@
 using Test
 using KernelDensity
 
-X = randn(100)
-Y = randn(100)
-
-k = kde(X)
-@test pdf(k, k.x) ≈ k.density
-
-k = kde((X,Y))
-@test pdf(k, k.x, k.y) ≈ k.density
-
-# Try to evaluate the KDE outside the interpolation domain
-# The KDE is allowed to be zero, but it should not be greater than the exact solution
-k = kde([0.0], bandwidth=1.0)
-@test pdf(k, k.x) ≈ k.density
-@test pdf(k, -10.0) ≤ pdf(Normal(), -10.0)
-@test pdf(k, +10.0) ≤ pdf(Normal(), +10.0)
-
-k = kde(([0.0],[0.0]), bandwidth=(1.0, 1.0))
-@test pdf(k, k.x, k.y) ≈ k.density
-@test pdf(k, -10.0, 0.0) ≤ pdf(MvNormal(2, 1.0), [-10.0, 0.0])
-@test pdf(k, +10.0, 0.0) ≤ pdf(MvNormal(2, 1.0), [+10.0, 0.0])
-@test pdf(k, 0.0, -10.0) ≤ pdf(MvNormal(2, 1.0), [0.0, -10.0])
-@test pdf(k, 0.0, +10.0) ≤ pdf(MvNormal(2, 1.0), [0.0, +10.0])
+@testset "interpolation computation" begin
+    X = randn(100)
+    Y = randn(100)
+
+    k = kde(X)
+    @test pdf.(k, k.x) ≈ k.density
+
+    k = kde((X,Y))
+    @test pdf(k, k.x, k.y) ≈ k.density
+
+    # Try to evaluate the KDE outside the interpolation domain
+    # The KDE is allowed to be zero, but it should not be greater than the exact solution
+    k = kde([0.0], bandwidth=1.0)
+    @test pdf.(k, k.x) ≈ k.density
+    @test pdf(k, -10.0) ≤ pdf(Normal(), -10.0)
+    @test pdf(k, +10.0) ≤ pdf(Normal(), +10.0)
+
+    k = kde(([0.0],[0.0]), bandwidth=(1.0, 1.0))
+    @test pdf(k, k.x, k.y) ≈ k.density
+    @test pdf(k, -10.0, 0.0) ≤ pdf(MvNormal(2, 1.0), [-10.0, 0.0])
+    @test pdf(k, +10.0, 0.0) ≤ pdf(MvNormal(2, 1.0), [+10.0, 0.0])
+    @test pdf(k, 0.0, -10.0) ≤ pdf(MvNormal(2, 1.0), [0.0, -10.0])
+    @test pdf(k, 0.0, +10.0) ≤ pdf(MvNormal(2, 1.0), [0.0, +10.0])
+end
+
+@testset "pdf method interface" begin
+    k = kde([-1., 1.])
+    ik = InterpKDE(k)
+
+    @test pdf.(k, [0., 1.]) ≈ [pdf(k, 0.), pdf(k, 1.)] ≈ pdf.(k,[0., 1.]) ≈ [pdf(ik, 0.), pdf(ik, 1.)]
+    @test all( pdf.(k, (0., 1.)) .≈ (pdf(k, 0.), pdf(k, 1.)) )
+    @test pdf.(k, [0. 1.; 2. -1.]) ≈ [pdf(k, 0.) pdf(k, 1.); pdf(k, 2.) pdf(k, -1.)]
+
+    k2d = kde(([-1., 1.], [0., 1.]))
+    ik2d = InterpKDE(k2d)
+
+    @test pdf(k2d, [0.5, 0.1]) ≈ pdf(k2d, [0.5; 0.1]) ≈ pdf(k2d, 0.5, 0.1) ≈ pdf(ik2d, 0.5, 0.1)
+    @test pdf(k2d, [0.5 1.; 0.1 1.]) ≈ [pdf(ik2d, 0.5, 0.1), pdf(ik2d, 1., 1.)]
+    @test pdf(k2d, [0.5; 1. ;;; 0.1; 1.]) ≈ [pdf(ik2d, 0.5, 1.) pdf(ik2d, 0.1, 1.)]
+end