add tests, fix a few small errors

Author: simonbyrne

.travis.yml

@@ -15,3 +15,4 @@ before_install:
 script:
   - julia -e 'Pkg.init(); run(`ln -s $(pwd()) $(Pkg.dir("KDE"))`); Pkg.pin("KDE"); Pkg.resolve()'
   - julia -e 'using KDE; @assert isdefined(:KDE); @assert typeof(KDE) === Module'
+  - julia ./runtests.jl

runtests.jl

@@ -0,0 +1,12 @@
+tests = [
+         "univariate",
+         "bivariate",
+         ]
+
+println("Running tests:")
+
+for t in tests
+    test_fn = joinpath("test", "$t.jl")
+    println(" * $test_fn")
+    include(test_fn)
+end

src/bivariate.jl

@@ -58,14 +58,19 @@ function conv(k::BivariateKDE, dist::(UnivariateDistribution,UnivariateDistribut
     # Convolve fft with characteristic function of kernel
     cx = -twoπ/(step(k.x)*Kx)
     cy = -twoπ/(step(k.y)*Ky)
-    for j = 1:size(ft,2)
-        for i = 1:size(ft,1)
-            ft[i,j] *= cf(distx,(i-1)*cx)*cf(disty,min(j-1,Ky-j+1)*cy)
+    for j = 0:size(ft,2)-1
+        for i = 0:size(ft,1)-1
+            ft[i+1,j+1] *= cf(distx,i*cx)*cf(disty,min(j,Ky-j)*cy)
         end
     end
+    dens = irfft(ft, Kx)
+
+    for i = 1:length(dens)
+        dens[i] = max(0.0,dens[i])
+    end
 
     # Invert the Fourier transform to get the KDE
-    BivariateKDE(k.x, k.y, irfft(ft, Kx))
+    BivariateKDE(k.x, k.y, dens)
 end
 
 typealias BivariateDistribution Union(MultivariateDistribution,(UnivariateDistribution,UnivariateDistribution))

src/univariate.jl

@@ -16,6 +16,7 @@ kernel_dist{D}(::Type{D},w::Real) = (s = w/std(D(0.0,1.0)); D(0.0,s))
 function default_bandwidth(data::Vector{Float64}, alpha::Float64 = 0.9)
     # Determine length of data
     ndata = length(data)
+    ndata <= 1 && return alpha
 
     # Calculate width using variance and IQR
     var_width = std(data)
@@ -103,12 +104,18 @@ function conv(k::UnivariateKDE, dist::UnivariateDistribution)
     #  = e^{i*t*a} \sum_{k=0}^K e^{-2pi*i*k*(-t*s*K/2pi)/K} N_k / N
     #  = A * fft(N_k/N)[-t*s*K/2pi + 1]
     c = -twoπ/(step(k.x)*K)
-    for j = 1:length(ft)
-        ft[j] *= cf(dist,(j-1)*c)
+    for j = 0:length(ft)-1
+        ft[j+1] *= cf(dist,j*c)
+    end
+
+    dens = irfft(ft, K)
+    # fix rounding error.
+    for i = 1:K
+        dens[i] = max(0.0,dens[i])
     end
 
     # Invert the Fourier transform to get the KDE
-    UnivariateKDE(k.x, irfft(ft, K))
+    UnivariateKDE(k.x, dens)
 end
 
 # main kde interface methods

test/bivariate.jl

@@ -0,0 +1,60 @@
+using Base.Test
+using Distributions
+using KDE
+
+import KDE: kernel_dist, default_bandwidth, kde_boundary, kde_range, tabulate
+
+for D in [(Normal,Normal),(Uniform,Uniform),(Logistic,Logistic)]
+    d = kernel_dist(D,(0.5,0.5))
+    @test isa(d,D)
+    dx,dy = d
+    @test mean(dx) == 0.0
+    @test mean(dy) == 0.0
+    @test_approx_eq std(dx) 0.5
+    @test_approx_eq std(dy) 0.5
+end
+
+r = kde_range((-2.0,2.0), 128)
+@test step(r) > 0
+
+for X in ([0.0], [0.0,0.0], [0.0,0.5], [-0.5:0.1:0.5])
+    w = default_bandwidth(X)
+    @test w > 0
+    lo, hi = kde_boundary(X,w)
+    @test lo < hi    
+    kr = kde_range((lo,hi), 10)
+    @test step(kr) > 0    
+    
+    for D in (Normal, )
+        k1 = tabulate((X,X),(r,r))
+        @test isa(k1,BivariateKDE)
+        @test size(k1.density) == (length(k1.x), length(k1.y))
+        @test all(k1.density .>= 0.0)
+        @test_approx_eq sum(k1.density)*step(k1.x)*step(k1.y) 1.0
+
+        k2 = conv(k1,kernel_dist((D,D),(0.1,0.1)))
+        @test isa(k2,BivariateKDE)
+        @test size(k2.density) == (length(k2.x), length(k2.y))
+        @test all(k2.density .>= 0.0)
+        @test_approx_eq sum(k2.density)*step(k2.x)*step(k2.y) 1.0
+
+        k3 = kde((X,X);kernel=D)
+        @test isa(k3,BivariateKDE)
+        @test size(k3.density) == (length(k3.x), length(k3.y))
+        @test all(k3.density .>= 0.0)
+        @test_approx_eq sum(k3.density)*step(k3.x)*step(k3.y) 1.0
+        
+        k4 = kde((X,X),(r,r);kernel=D)
+        @test isa(k4,BivariateKDE)
+        @test size(k4.density) == (length(k4.x), length(k4.y))
+        @test all(k4.density .>= 0.0)
+        @test_approx_eq sum(k4.density)*step(k4.x)*step(k4.y) 1.0
+
+        k5 = kde([X X];kernel=D)
+        @test isa(k5,BivariateKDE)
+        @test size(k5.density) == (length(k5.x), length(k5.y))
+        @test all(k5.density .>= 0.0)
+        @test_approx_eq sum(k5.density)*step(k5.x)*step(k5.y) 1.0
+
+    end
+end

test/univariate.jl

@@ -0,0 +1,51 @@
+using Base.Test
+using Distributions
+using KDE
+
+import KDE: kernel_dist, default_bandwidth, kde_boundary, kde_range, tabulate
+
+for D in [Normal,Uniform,Logistic]
+    d = kernel_dist(D,0.5)
+    @test isa(d,D)
+    @test mean(d) == 0.0
+    @test_approx_eq std(d) 0.5
+end
+
+r = kde_range((-2.0,2.0), 128)
+@test step(r) > 0
+
+for X in ([0.0], [0.0,0.0], [0.0,0.5], [-0.5:0.1:0.5])
+    w = default_bandwidth(X)
+    @test w > 0
+    lo, hi = kde_boundary(X,w)
+    @test lo < hi    
+    kr = kde_range((lo,hi), 10)
+    @test step(kr) > 0    
+    
+    for D in (Normal, )
+        k1 = tabulate(X,r)
+        @test isa(k1,UnivariateKDE)
+        @test length(k1.density) == length(k1.x)
+        @test all(k1.density .>= 0.0)
+        @test_approx_eq sum(k1.density)*step(k1.x) 1.0
+
+        k2 = conv(k1,kernel_dist(D,0.1))
+        @test isa(k2,UnivariateKDE)
+        @test length(k2.density) == length(k2.x)
+        @test all(k2.density .>= 0.0)
+        @test_approx_eq sum(k2.density)*step(k2.x) 1.0
+
+        k3 = kde(X;kernel=D)
+        @test isa(k3,UnivariateKDE)
+        @test length(k3.density) == length(k3.x)
+        @test all(k3.density .>= 0.0)
+        @test_approx_eq sum(k3.density)*step(k3.x) 1.0
+        
+        k4 = kde(X,r;kernel=D)
+        @test isa(k4,UnivariateKDE)
+        @test length(k4.density) == length(k4.x)
+        @test all(k4.density .>= 0.0)
+        @test_approx_eq sum(k4.density)*step(k4.x) 1.0
+
+    end
+end