SquaredExponentialKernel becomes SqExponentialKernel

Author: theogf

benchmark/kernelmatrix.jl

@@ -2,7 +2,7 @@ using KernelFunctions
 
 SUITE["KernelFunctions"] = BenchmarkGroup()
 
-kernelnames = ["SquaredExponentialKernel"]
+kernelnames = ["SqExponentialKernel"]
 kerneltypes = ["ARD","ISO"]
 kernels=Dict{String,Dict{String,KernelFunctions.Kernel}}()
 for k in kernelnames

dev/debugAD.jl

@@ -12,39 +12,39 @@ testfunction(k,A,B) = det(kernelmatrix(k,A,B))
 testfunction(k,A) = sum(kernelmatrix(k,A))
 k = MaternKernel(vl)
 KernelFunctions.kappa(k,3)
-testfunction(SquaredExponentialKernel(vl),A)
+testfunction(SqExponentialKernel(vl),A)
 testfunction(MaternKernel(vl),A)
 @which kernelmatrix(MaternKernel(vl),A,B)
 #For debugging
 @info "Running Zygote gradients"
 Zygote.refresh()
 ## Zygote
-Zygote.gradient(x->testfunction(SquaredExponentialKernel(x),A),vl)
+Zygote.gradient(x->testfunction(SqExponentialKernel(x),A),vl)
 Zygote.gradient(x->testfunction(MaternKernel(x),A),vl)
-Zygote.gradient(x->testfunction(SquaredExponentialKernel(x),A,B),vl)[1]
+Zygote.gradient(x->testfunction(SqExponentialKernel(x),A,B),vl)[1]
 Zygote.gradient(x->testfunction(MaternKernel(x),A,B),vl)[1]
-Zygote.gradient(x->testfunction(SquaredExponentialKernel(x),A,B),l)
+Zygote.gradient(x->testfunction(SqExponentialKernel(x),A,B),l)
 Zygote.gradient(x->testfunction(MaternKernel(x),A,B),l)
-Zygote.gradient(x->testfunction(SquaredExponentialKernel(x),A),l)
+Zygote.gradient(x->testfunction(SqExponentialKernel(x),A),l)
 Zygote.gradient(x->testfunction(MaternKernel(x),A),l)
 Zygote.gradient(x->testfunction(MaternKernel(x),A),l)
 Zygote.gradient(x->kernelmatrix(MaternKernel(x,1.0),A)[1],l)
 @info "Running Tracker gradients"
 ## Tracker
-# Tracker.gradient(x->testfunction(SquaredExponentialKernel(vl),x,B),A)
-# Tracker.gradient(x->testfunction(SquaredExponentialKernel(l),x[:,:]),A)
-# # Tracker.gradient(x->testfunction(SquaredExponentialKernel(x),A,B),vl)
-# Tracker.gradient(x->testfunction(SquaredExponentialKernel(x),A),vl)
-# Tracker.gradient(x->testfunction(SquaredExponentialKernel(x),A,B),l)
-# Tracker.gradient(x->testfunction(SquaredExponentialKernel(x),A),l)
+# Tracker.gradient(x->testfunction(SqExponentialKernel(vl),x,B),A)
+# Tracker.gradient(x->testfunction(SqExponentialKernel(l),x[:,:]),A)
+# # Tracker.gradient(x->testfunction(SqExponentialKernel(x),A,B),vl)
+# Tracker.gradient(x->testfunction(SqExponentialKernel(x),A),vl)
+# Tracker.gradient(x->testfunction(SqExponentialKernel(x),A,B),l)
+# Tracker.gradient(x->testfunction(SqExponentialKernel(x),A),l)
 
 @info "Running ForwardDiff gradients"
 ## ForwardDiff
-ForwardDiff.gradient(x->testfunction(SquaredExponentialKernel(x),A,B),vl) #✓
+ForwardDiff.gradient(x->testfunction(SqExponentialKernel(x),A,B),vl) #✓
 ForwardDiff.gradient(x->testfunction(MaternKernel(x),A,B),vl) #✓
-ForwardDiff.gradient(x->testfunction(SquaredExponentialKernel(x),A),vl) #✓
+ForwardDiff.gradient(x->testfunction(SqExponentialKernel(x),A),vl) #✓
 ForwardDiff.gradient(x->testfunction(MaternKernel(x),A),vl) #✓
-ForwardDiff.gradient(x->testfunction(SquaredExponentialKernel(x[1]),A,B),[l])
+ForwardDiff.gradient(x->testfunction(SqExponentialKernel(x[1]),A,B),[l])
 ForwardDiff.gradient(x->testfunction(MaternKernel(x[1]),A,B),[l])
-ForwardDiff.gradient(x->testfunction(SquaredExponentialKernel(x[1]),A),[l])
+ForwardDiff.gradient(x->testfunction(SqExponentialKernel(x[1]),A),[l])
 ForwardDiff.gradient(x->testfunction(MaternKernel(x[1]),A),[l])

dev/matrixvsvectors.jl

@@ -14,7 +14,7 @@ timekf = similar(timestheno); memkf = similar(timestheno)
     B = randn(D,1001)
 
     # Standardised eq kernel with length-scale 0.1.
-    medkf = median(@benchmark KernelFunctions.kernelmatrix(SquaredExponentialKernel(0.01),$A,$B,obsdim=2))
+    medkf = median(@benchmark KernelFunctions.kernelmatrix(SqExponentialKernel(0.01),$A,$B,obsdim=2))
     timekf[i] = medkf.time/1e6; memkf[i] = medkf.memory/2^20
     medstheno = median(@benchmark pw(eq(; l=0.1), ColsAreObs($A), ColsAreObs($B)))
     timestheno[i] = medstheno.time/1e6; memstheno[i] = medstheno.memory/2^20

src/KernelFunctions.jl

@@ -1,7 +1,7 @@
 module KernelFunctions
 
 export kernelmatrix, kernelmatrix!, kerneldiagmatrix, kerneldiagmatrix!, kappa
-export Kernel, SquaredExponentialKernel, MaternKernel, Matern32Kernel, Matern52Kernel
+export Kernel, SqExponentialKernel, MaternKernel, Matern32Kernel, Matern52Kernel
 
 export Transform, ScaleTransform
 
@@ -18,7 +18,7 @@ include("utils.jl")
 include("transform/transform.jl")
 include("kernelmatrix.jl")
 
-kernels = ["squaredexponential","matern"]
+kernels = ["sqexponential","matern"]
 for k in kernels
     include(joinpath("kernels",k*".jl"))
 end

src/generic.jl

@@ -1,11 +1,12 @@
 
 @inline metric(κ::Kernel) = κ.metric
 kernels =
-for k in [:SquaredExponentialKernel,:MaternKernel,:Matern32Kernel,:Matern52Kernel]
+for k in [:SqExponentialKernel,:MaternKernel,:Matern32Kernel,:Matern52Kernel]
     eval(quote
         @inline (κ::$k)(d::Real) = kappa(κ,d)
         @inline (κ::$k)(x::AbstractVector{T},y::AbstractVector{T}) where {T} = kernel(κ,evaluate(κ.(metric),x,y))
         @inline (κ::$k)(x::AbstractMatrix{T},y::AbstractMatrix{T},obsdim::Integer=defaultobs) where {T} = kernelmatrix(κ,x,y,obsdim=obsdim)
+            @inline (κ::$k)(x::AbstractMatrix{T},obsdim::Integer=defaultobs) where {T} = kernelmatrix(κ,x,obsdim=obsdim)
     end)
 end
 ### Transform generics

src/kernels/sqexponential.jl

@@ -0,0 +1,53 @@
+"""
+    SqExponentialKernel([α=1])
+
+The squared exponential kernel is an isotropic Mercer kernel given by the formula:
+
+```
+    κ(x,y) = exp(-‖x-y‖²)
+```
+
+See also [`ExponentialKernel`](@ref) for a
+related form of the kernel or [`GammaExponentialKernel`](@ref) for a generalization.
+
+# Examples
+
+```jldoctest; setup = :(using KernelFunctions)
+julia> SqExponentialKernel()
+SqExponentialKernel{Float64,Float64}(1.0)
+
+julia> SqExponentialKernel(2.0f0)
+SqExponentialKernel{Float32,Float32}(2.0)
+
+julia> SqExponentialKernel([2.0,3.0])
+SqExponentialKernel{Float64,Array{Float64}}(1.0)
+```
+"""
+struct SqExponentialKernel{T,Tr<:Transform} <: Kernel{T,Tr}
+    transform::Tr
+    metric::SemiMetric
+    function SqExponentialKernel{T,Tr}(transform::Tr) where {T,Tr<:Transform}
+        return new{T,Tr}(transform,SqEuclidean())
+    end
+end
+
+function SqExponentialKernel(α::T=1.0) where {T<:Real}
+    SqExponentialKernel{T,ScaleTransform{T}}(ScaleTransform(α))
+end
+
+function SqExponentialKernel(α::A) where {A<:AbstractVector{<:Real}}
+    SqExponentialKernel{eltype(A),ScaleTransform{A}}(ScaleTransform(α))
+end
+
+function SqExponentialKernel(t::T) where {T<:Transform}
+    SqExponentialKernel{eltype(t),T}(t)
+end
+
+@inline kappa(κ::SqExponentialKernel, d²::Real) where {T} = exp(-d²)
+
+# function convert(
+#         ::Type{K},
+#         κ::SqExponentialKernel
+#     ) where {K>:SqExponentialKernel{T,A} where {T,A}}
+#     return SqExponentialKernel{T}(T.(κ.α))
+# end

src/kernels/squaredexponential.jl

@@ -3,23 +3,30 @@ using KernelFunctions, Test, Distances
 # test type conversion
 l = 2.0
 vl = [l,l]
+s = ScaleTransform(3.0)
 
-## SquaredExponentialKernel
-@testset "SquaredExponentialKernel" begin
-    @test KernelFunctions.metric(SquaredExponentialKernel(l)) == SqEuclidean()
-    @test KernelFunctions.transform(SquaredExponentialKernel(l)) == ScaleTransform(l)
-    @test KernelFunctions.transform(SquaredExponentialKernel(vl)) == ScaleTransform(vl)
+## SqExponentialKernel
+@testset "SqExponentialKernel" begin
+    @test KernelFunctions.metric(SqExponentialKernel(l)) == SqEuclidean()
+    @test KernelFunctions.transform(SqExponentialKernel(l)) == ScaleTransform(l)
+    @test KernelFunctions.transform(SqExponentialKernel(vl)) == ScaleTransform(vl)
+    @test KernelFunctions.transform(SqExponentialKernel(s)) == s
 end
 
+## MaternKernel
+
 @testset "MaternKernel" begin
     @test KernelFunctions.metric(MaternKernel(l)) == Euclidean()
-    @test KernelFunctions.metric(MaternKernel(l,1.5)) == Euclidean()
-    @test KernelFunctions.metric(MaternKernel(l,2.5)) == Euclidean()
+    @test KernelFunctions.metric(Matern32Kernel(l)) == Euclidean()
+    @test KernelFunctions.metric(Matern52Kernel(l)) == Euclidean()
     @test KernelFunctions.transform(MaternKernel(l)) == ScaleTransform(l)
+    @test KernelFunctions.transform(Matern32Kernel(l)) == ScaleTransform(l)
+    @test KernelFunctions.transform(Matern52Kernel(l)) == ScaleTransform(l)
     @test KernelFunctions.transform(MaternKernel(vl)) == ScaleTransform(vl)
-    @test isa(MaternKernel(),Matern32Kernel)
-    @test isa(MaternKernel(1.0,1.0),MaternKernel)
-    @test isa(MaternKernel(1.0,1.5),Matern32Kernel)
-    @test isa(MaternKernel(1.0,2.5),Matern52Kernel)
-    @test isa(MaternKernel(1.0,Inf),SquaredExponentialKernel)
+    @test KernelFunctions.transform(Matern32Kernel(vl)) == ScaleTransform(vl)
+    @test KernelFunctions.transform(Matern52Kernel(vl)) == ScaleTransform(vl)
+    @test KernelFunctions.transform(MaternKernel(s)) == s
+    @test KernelFunctions.transform(Matern32Kernel(s)) == s
+    @test KernelFunctions.transform(Matern52Kernel(s)) == s
+
 end

test/constructors.jl

@@ -7,7 +7,7 @@ A = rand(dims...)
 B = rand(dims...)
 K = [zeros(dims[1],dims[1]),zeros(dims[2],dims[2])]
 Kdiag = [zeros(dims[1]),zeros(dims[2])]
-kernels = [SquaredExponentialKernel(),MaternKernel(1.0,1.0),Matern32Kernel(),Matern52Kernel()]
+kernels = [SqExponentialKernel(),MaternKernel(),Matern32Kernel(),Matern52Kernel()]
 @testset "Inplace Kernel Matrix" begin
     for k in kernels
         @testset "$k" begin

test/kernelmatrix.jl

@@ -7,7 +7,7 @@ dims = [10,5]
 A = rand(dims...)
 B = rand(dims...)
 K = [zeros(dims[1],dims[1]),zeros(dims[2],dims[2])]
-kernels = [SquaredExponentialKernel,MaternKernel]
+kernels = [SqExponentialKernel,MaternKernel,Matern32Kernel,Matern52Kernel]
 l = 2.0
 vl = l*ones(dims[1])
 testfunction(k,A,B) = det(kernelmatrix(k,A,B))