Move to AbstractVector

Project.toml

@@ -20,8 +20,8 @@ Requires = "1.0.1"
 SpecialFunctions = "0.8, 0.9, 0.10"
 StatsBase = "0.32, 0.33"
 StatsFuns = "0.8, 0.9"
-ZygoteRules = "0.2"
 Zygote = "= 0.4.16"
+ZygoteRules = "0.2"
 julia = "1.3"
 
 [extras]

src/KernelFunctions.jl

@@ -36,8 +36,6 @@ using StatsFuns: logtwo
 using InteractiveUtils: subtypes
 using StatsBase
 
-const defaultobs = 2
-
 """
 Abstract type defining a slice-wise transformation on an input matrix
 """

src/basekernels/fbm.jl

@@ -33,45 +33,35 @@ const sqroundoff = 1e-15
 
 _fbm(modX, modY, modXY, h) = (modX^h + modY^h - modXY^h)/2
 
-function kernelmatrix(κ::FBMKernel, X::AbstractMatrix; obsdim::Int = defaultobs)
-    @assert obsdim ∈ [1,2] "obsdim should be 1 or 2 (see docs of kernelmatrix))"
-    modX = sum(abs2, X; dims = feature_dim(obsdim))
-    modXX = pairwise(SqEuclidean(sqroundoff), X, dims = obsdim)
-    return _fbm.(vec(modX), reshape(modX, 1, :), modXX, κ.h)
+_mod(x::AbstractVector{<:Real}) = abs2.(x)
+_mod(x::ColVecs) = vec(sum(abs2, x.X; dims=1))
+_mod(x::RowVecs) = vec(sum(abs2, x.X; dims=2))
+
+function kernelmatrix(κ::FBMKernel, x::AbstractVector)
+    modx = _mod(x)
+    modxx = pairwise(SqEuclidean(sqroundoff), x)
+    return _fbm.(modx, modx', modxx, κ.h)
 end
 
-function kernelmatrix!(K::AbstractMatrix, κ::FBMKernel, X::AbstractMatrix; obsdim::Int = defaultobs)
-    @assert obsdim ∈ [1,2] "obsdim should be 1 or 2 (see docs of kernelmatrix))"
-    modX = sum(abs2, X; dims = feature_dim(obsdim))
-    modXX = pairwise(SqEuclidean(sqroundoff), X, dims = obsdim)
-    K .= _fbm.(vec(modX), reshape(modX, 1, :), modXX, κ.h)
+function kernelmatrix!(K::AbstractMatrix, κ::FBMKernel, x::AbstractVector)
+    modx = _mod(x)
+    modxx = pairwise(SqEuclidean(sqroundoff), x)
+    K .= _fbm.(modx, modx', modxx, κ.h)
     return K
 end
 
-function kernelmatrix(
-    κ::FBMKernel,
-    X::AbstractMatrix,
-    Y::AbstractMatrix;
-    obsdim::Int = defaultobs,
-)
-    @assert obsdim ∈ [1,2] "obsdim should be 1 or 2 (see docs of kernelmatrix))"
-    modX = sum(abs2, X, dims = feature_dim(obsdim))
-    modY = sum(abs2, Y, dims = feature_dim(obsdim))
-    modXY = pairwise(SqEuclidean(sqroundoff), X, Y,dims = obsdim)
-    return _fbm.(vec(modX), reshape(modY, 1, :), modXY, κ.h)
+function kernelmatrix(κ::FBMKernel, x::AbstractVector, y::AbstractVector)
+    modxy = pairwise(SqEuclidean(sqroundoff), x, y)
+    return _fbm.(_mod(x), _mod(y)', modxy, κ.h)
 end
 
 function kernelmatrix!(
     K::AbstractMatrix,
     κ::FBMKernel,
-    X::AbstractMatrix,
-    Y::AbstractMatrix;
-    obsdim::Int = defaultobs,
+    X::AbstractVector,
+    Y::AbstractVector,
 )
-    @assert obsdim ∈ [1,2] "obsdim should be 1 or 2 (see docs of kernelmatrix))"
-    modX = sum(abs2, X, dims = feature_dim(obsdim))
-    modY = sum(abs2, Y, dims = feature_dim(obsdim))
-    modXY = pairwise(SqEuclidean(sqroundoff), X, Y,dims = obsdim)
-    K .= _fbm.(vec(modX), reshape(modY, 1, :), modXY, κ.h)
+    modxy = pairwise(SqEuclidean(sqroundoff), X, Y,dims = obsdim)
+    K .= _fbm.(_mod(x), _mod(y)', modxy, κ.h)
     return K
 end

src/generic.jl

@@ -25,3 +25,16 @@ end
 
 # Fallback implementation of evaluate for `SimpleKernel`s.
 (k::SimpleKernel)(x, y) = kappa(k, evaluate(metric(k), x, y))
+
+# This is type piracy. We should not doing this.
+function Distances.pairwise(d::PreMetric, x::AbstractVector{<:Real})
+    return pairwise(d, reshape(x, :, 1); dims=1)
+end
+
+function Distances.pairwise(
+    d::PreMetric,
+    x::AbstractVector{<:Real},
+    y::AbstractVector{<:Real},
+)
+    return pairwise(d, reshape(x, :, 1), reshape(y, :, 1); dims=1)
+end

src/kernels/kernelproduct.jl

@@ -24,31 +24,18 @@ Base.length(k::KernelProduct) = length(k.kernels)
 
 (κ::KernelProduct)(x, y) = prod(k(x, y) for k in κ.kernels)
 
-hadamard(x,y) = x.*y
+hadamard(x, y) = x .* y
 
-function kernelmatrix(
-    κ::KernelProduct,
-    X::AbstractMatrix;
-    obsdim::Int=defaultobs,
-)
-    reduce(hadamard, kernelmatrix(κ.kernels[i], X, obsdim = obsdim) for i in 1:length(κ))
+function kernelmatrix(κ::KernelProduct, x::AbstractVector)
+    return reduce(hadamard, kernelmatrix(κ.kernels[i], x) for i in 1:length(κ))
 end
 
-function kernelmatrix(
-    κ::KernelProduct,
-    X::AbstractMatrix,
-    Y::AbstractMatrix;
-    obsdim::Int=defaultobs,
-)
-    reduce(hadamard, kernelmatrix(κ.kernels[i], X, Y, obsdim = obsdim) for i in 1:length(κ))
+function kernelmatrix(κ::KernelProduct, x::AbstractVector, y::AbstractVector)
+    return reduce(hadamard, kernelmatrix(κ.kernels[i], x, y) for i in 1:length(κ))
 end
 
-function kerneldiagmatrix(
-    κ::KernelProduct,
-    X::AbstractMatrix;
-    obsdim::Int=defaultobs,
-) #TODO Add test
-    reduce(hadamard, kerneldiagmatrix(κ.kernels[i], X, obsdim = obsdim) for i in 1:length(κ))
+function kerneldiagmatrix(κ::KernelProduct, x::AbstractVector)
+    return reduce(hadamard, kerneldiagmatrix(κ.kernels[i], x) for i in 1:length(κ))
 end
 
 function Base.show(io::IO, κ::KernelProduct)

src/kernels/kernelsum.jl

@@ -48,25 +48,16 @@ Base.length(k::KernelSum) = length(k.kernels)
 
 (κ::KernelSum)(x, y) = sum(κ.weights[i] * κ.kernels[i](x, y) for i in 1:length(κ))
 
-function kernelmatrix(κ::KernelSum, X::AbstractMatrix; obsdim::Int = defaultobs)
-    sum(κ.weights[i] * kernelmatrix(κ.kernels[i], X, obsdim = obsdim) for i in 1:length(κ))
+function kernelmatrix(κ::KernelSum, x::AbstractVector)
+    return sum(κ.weights[i] * kernelmatrix(κ.kernels[i], x) for i in 1:length(κ))
 end
 
-function kernelmatrix(
-    κ::KernelSum,
-    X::AbstractMatrix,
-    Y::AbstractMatrix;
-    obsdim::Int = defaultobs,
-)
-    sum(κ.weights[i] * kernelmatrix(κ.kernels[i], X, Y, obsdim = obsdim) for i in 1:length(κ))
+function kernelmatrix(κ::KernelSum, x::AbstractVector, y::AbstractVector)
+    return sum(κ.weights[i] * kernelmatrix(κ.kernels[i], x, y) for i in 1:length(κ))
 end
 
-function kerneldiagmatrix(
-    κ::KernelSum,
-    X::AbstractMatrix;
-    obsdim::Int = defaultobs,
-)
-    sum(κ.weights[i] * kerneldiagmatrix(κ.kernels[i], X, obsdim = obsdim) for i in 1:length(κ))
+function kerneldiagmatrix(κ::KernelSum, x::AbstractVector)
+    return sum(κ.weights[i] * kerneldiagmatrix(κ.kernels[i], x) for i in 1:length(κ))
 end
 
 function Base.show(io::IO, κ::KernelSum)

src/kernels/scaledkernel.jl

@@ -13,11 +13,42 @@ function ScaledKernel(kernel::Tk, σ²::Tσ²=1.0) where {Tk<:Kernel,Tσ²<:Real
     return ScaledKernel{Tk, Tσ²}(kernel, [σ²])
 end
 
-kappa(k::ScaledKernel, x) = first(k.σ²) * kappa(k.kernel, x)
-
 (k::ScaledKernel)(x, y) = first(k.σ²) * k.kernel(x, y)
 
-metric(k::ScaledKernel) = metric(k.kernel)
+function kernelmatrix(κ::ScaledKernel, x::AbstractVector, y::AbstractVector)
+    return κ.σ² .* kernelmatrix(κ.kernel, x, y)
+end
+
+function kernelmatrix(κ::ScaledKernel, x::AbstractVector)
+    return κ.σ² .* kernelmatrix(κ.kernel, x)
+end
+
+function kerneldiagmatrix(κ::ScaledKernel, x::AbstractVector)
+    return κ.σ² .* kerneldiagmatrix(κ.kernel, x)
+end
+
+function kernelmatrix!(
+    K::AbstractMatrix,
+    κ::ScaledKernel,
+    x::AbstractVector,
+    y::AbstractVector,
+)
+    kernelmatrix!(K, κ, x, y)
+    K .*= κ.σ² 
+    return K
+end
+
+function kernelmatrix!(K::AbstractMatrix, κ::ScaledKernel, x::AbstractVector)
+    kernelmatrix!(K, κ, x)
+    K .*= κ.σ²
+    return K
+end
+
+function kerneldiagmatrix!(K::AbstractVector, κ::ScaledKernel, x::AbstractVector)
+    kerneldiagmatrix!(K, κ, x)
+    K .*= κ.σ²
+    return K
+end
 
 Base.:*(w::Real, k::Kernel) = ScaledKernel(k, w)