Improved documentation and tests

Author: theogf

docs/create_kernel_plots.jl

@@ -0,0 +1,49 @@
+using Plots; pyplot();
+using Distributions
+using LinearAlgebra
+using KernelFunctions
+# Translational invariants kernels
+
+default(lw=3.0,titlefontsize=28,tickfontsize=18)
+
+x₀ = 0.0; l=0.1
+n_grid = 101
+fill(x₀,n_grid,1)
+xrange = reshape(collect(range(-3,3,length=n_grid)),:,1)
+for k in [SqExponentialKernel,ExponentialKernel]
+    K = kernelmatrix(k(),xrange,obsdim=1)
+    v = rand(MvNormal(K+1e-7I))
+    plot(xrange,v,lab="",title="f(x)",framestyle=:none) |> display
+    savefig(joinpath(@__DIR__,"src","assets","GP_sample_$(k).png"))
+    plot(xrange,kernel.(k(),x₀,xrange),lab="",ylims=(0,1.1),title="k(0,x)") |> display
+    savefig(joinpath(@__DIR__,"src","assets","kappa_function_$(k).png"))
+end
+
+for k in [GammaExponentialKernel(1.0,1.5)]
+    sparse =1
+    while !isposdef(kernelmatrix(k,xrange*sparse,obsdim=1) + 1e-5I); sparse += 1; end
+    v = rand(MvNormal(kernelmatrix(k,xrange*sparse,obsdim=1)+1e-7I))
+    plot(xrange,v,lab="",title="f(x)",framestyle=:none) |> display
+    savefig(joinpath(@__DIR__,"src","assets","GP_sample_GammaExponentialKernel.png"))
+    plot(xrange,kernel.(k,x₀,xrange),lab="",ylims=(0,1.1),title="k(0,x)") |> display
+    savefig(joinpath(@__DIR__,"src","assets","kappa_function_GammaExponentialKernel.png"))
+end
+
+for k in [MaternKernel,Matern32Kernel,Matern52Kernel]
+    K = kernelmatrix(k(),xrange,obsdim=1)
+    v = rand(MvNormal(K+1e-7I))
+    plot(xrange,v,lab="",title="f(x)",framestyle=:none) |> display
+    savefig(joinpath(@__DIR__,"src","assets","GP_sample_$(k).png"))
+    plot(xrange,kernel.(k(),x₀,xrange),lab="",ylims=(0,1.1),title="k(0,x)") |> display
+    savefig(joinpath(@__DIR__,"src","assets","kappa_function_$(k).png"))
+end
+
+
+for k in [RationalQuadraticKernel]
+    K = kernelmatrix(k(),xrange,obsdim=1)
+    v = rand(MvNormal(K+1e-7I))
+    plot(xrange,v,lab="",title="f(x)",framestyle=:none) |> display
+    savefig(joinpath(@__DIR__,"src","assets","GP_sample_$(k).png"))
+    plot(xrange,kernel.(k(),x₀,xrange),lab="",ylims=(0,1.1),title="k(0,x)") |> display
+    savefig(joinpath(@__DIR__,"src","assets","kappa_function_$(k).png"))
+end

docs/make.jl

@@ -9,6 +9,8 @@ makedocs(
             "User Guide" => "userguide.md",
              "Kernel Functions"=>"kernels.md",
              "Transform"=>"transform.md",
+             "Metrics"=>"metrics.md",
+             "Theory"=>"theory.md",
              "API"=>"api.md"]
 )
 

docs/src/api.md

@@ -0,0 +1,67 @@
+# API Library
+
+---
+```@contents
+Pages = ["api.md"]
+```
+
+```@meta
+CurrentModule = KernelFunctions
+```
+
+## Module
+```@docs
+KernelFunctions
+```
+
+## Kernel Functions
+
+```@docs
+SqExponentialKernel
+Exponential
+GammaExponentialKernel
+ExponentiatedKernel
+MaternKernel
+Matern32Kernel
+Matern52Kernel
+LinearKernel
+PolynomialKernel
+RationalQuadraticKernel
+GammaRationalQuadraticKernel
+```
+
+## Kernel Combinations
+
+```@docs
+KernelSum
+KernelProduct
+```
+
+## Transforms
+
+```@docs
+Transform
+IdentityTransform
+ScaleTransform
+LowRankTransform
+FunctionTransform
+```
+
+## Functions
+
+```@docs
+kernelmatrix
+kernelmatrix!
+kerneldiagmatrix
+kerneldiagmatrix!
+transform
+```
+
+
+## Index
+
+```@index
+Pages = ["api.md"]
+Module = ["KernelFunctions"]
+Order = [:type, :function]
+```

docs/src/metrics.md

@@ -0,0 +1,20 @@
+# Metrics
+
+KernelFunctions.jl relies on [Distances.jl]() for computing the pairwise matrix.
+To do so a distance measure is needed for each kernel. Two very common ones can already be used : `SqEuclidean` and `Euclidean`.
+However all kernels do not rely on distances metrics respecting all the definitions. That's why two additional metrics come with the package : `DotProduct` (`<x,y>`) and `Delta` (`δ(x,y)`). If you want to create a new distance just implement the following :
+
+```julia
+struct Delta <: Distances.PreMetric
+end
+
+@inline function Distances._evaluate(::Delta,a::AbstractVector{T},b::AbstractVector{T}) where {T}
+    @boundscheck if length(a) != length(b)
+        throw(DimensionMismatch("first array has length $(length(a)) which does not match the length of the second, $(length(b))."))
+    end
+    return a==b
+end
+
+@inline (dist::Delta)(a::AbstractArray,b::AbstractArray) = Distances._evaluate(dist,a,b)
+@inline (dist::Delta)(a::Number,b::Number) = a==b
+```

docs/src/theory.md

@@ -4,12 +4,14 @@
 There is a more general `Transform`: `FunctionTransform` that uses a function and apply it on each vector via `mapslices`.
 You can also create a pipeline of `Transform` via `TransformChain`. For example `LowRankTransform(rand(10,5))∘ScaleTransform(2.0)`.
 
+One apply a transformation on a matrix or a vector via `transform(t::Transform,v::AbstractVecOrMat)`
+
 ## Transforms :
 
 ```@docs
   IdentityTransform
   ScaleTransform
   LowRankTransform
   FunctionTransform
-  TransformChain
+  ChainTransform
 ```

docs/src/transform.md

@@ -1,5 +1 @@
 # Building kernel and matrices easily!
-
-## Creating a kernel function
-
-## 

docs/src/userguide.md

@@ -1,5 +1,5 @@
 export Transform, IdentityTransform, ScaleTransform, LowRankTransform, FunctionTransform, ChainTransform
-
+export transform
 
 abstract type Transform end
 

src/transform/transform.jl

@@ -40,6 +40,8 @@ k = SqExponentialKernel()
         w1 = 0.4; w2 = 1.2;
         ks2 = KernelSum([k1,k2],weights=[w1,w2])
         @test all(kernelmatrix(ks,A) .== kernelmatrix(k1,A) + kernelmatrix(k2,A))
+        @test all(kernelmatrix(ks+k1,A) .≈ 2*kernelmatrix(k1,A) + kernelmatrix(k2,A))
+        @test all(kernelmatrix(k1+ks,A) .≈ 2*kernelmatrix(k1,A) + kernelmatrix(k2,A))
         @test all(kernelmatrix(ks,A,B) .== kernelmatrix(k1,A,B) + kernelmatrix(k2,A,B))
         @test all(kerneldiagmatrix(ks,A) .== kerneldiagmatrix(k1,A) + kerneldiagmatrix(k2,A))
         @test all(kernelmatrix(ks2,A) .== w1*kernelmatrix(k1,A) + w2*kernelmatrix(k2,A))
@@ -48,9 +50,11 @@ k = SqExponentialKernel()
         k1 = SqExponentialKernel()
         k2 = LinearKernel()
         kp = k1 * k2
-        @test all(kernelmatrix(kp,A) .== kernelmatrix(k1,A) .* kernelmatrix(k2,A))
-        @test all(kernelmatrix(kp,A,B) .== kernelmatrix(k1,A,B) .* kernelmatrix(k2,A,B))
-        @test all(kerneldiagmatrix(kp,A) .== kerneldiagmatrix(k1,A) .* kerneldiagmatrix(k2,A))
-        @test all(kernelmatrix(kp,A) .== kernelmatrix(k1,A) .* kernelmatrix(k2,A))
+        @test all(kernelmatrix(kp,A) .≈ kernelmatrix(k1,A) .* kernelmatrix(k2,A))
+        @test all(kernelmatrix(kp*k1,A) .≈ kernelmatrix(k1,A).^2 .* kernelmatrix(k2,A))
+        @test all(kernelmatrix(k1*kp,A) .≈ kernelmatrix(k1,A).^2 .* kernelmatrix(k2,A))
+        @test all(kernelmatrix(kp,A) .≈ kernelmatrix(k1,A) .* kernelmatrix(k2,A))
+        @test all(kernelmatrix(kp,A,B) .≈ kernelmatrix(k1,A,B) .* kernelmatrix(k2,A,B))
+        @test all(kernelmatrix(kp,A) .≈ kernelmatrix(k1,A) .* kernelmatrix(k2,A))
     end
 end

test/test_kernelmatrix.jl

@@ -188,4 +188,20 @@ x = rand()*2; v1 = rand(3); v2 = rand(3); id = IdentityTransform()
             @test kappa(GammaRationalQuadraticKernel(1.0,a,1.0),x) ≈ kappa(RationalQuadraticKernel(1.0,a),x)
         end
     end
+    @testset "KernelCombinations" begin
+        k1 = LinearKernel()
+        k2 = SqExponentialKernel()
+        X = rand(2,2)
+        @testset "KernelSum" begin
+            k = k1 + k2
+            @test KernelFunctions.metric(k) == [KernelFunctions.DotProduct(),KernelFunctions.SqEuclidean()]
+            @test length(k) == 2
+            @test transform(k) == [transform(k1),transform(k2)]
+            @test transform(k,X) == [transform(k1,X),transform(k2,X)]
+            @test transform(k,X,1) == [transform(k1,X,1),transform(k2,X,1)]
+        end
+        @testset "KernelProduct" begin
+
+        end
+    end
 end