From 864bfab507cd1ed0b4777e0fa7dcb0540219e262 Mon Sep 17 00:00:00 2001 From: David Widmann Date: Sun, 17 Jan 2021 17:42:30 +0100 Subject: [PATCH 01/12] Rename `TensorProduct` and implement `TensorCore.tensor` --- Project.toml | 4 +- docs/src/kernels.md | 2 +- src/KernelFunctions.jl | 8 +++- src/deprecated.jl | 5 +++ src/kernels/kernelproduct.jl | 25 ----------- src/kernels/kernelsum.jl | 22 --------- src/kernels/overloads.jl | 22 +++++++++ src/kernels/tensorproduct.jl | 76 ++++++++++++++++++++++---------- src/matrix/kernelkroneckermat.jl | 6 +-- src/utils.jl | 2 - test/kernels/tensorproduct.jl | 29 +++++++++--- 11 files changed, 115 insertions(+), 86 deletions(-) create mode 100644 src/kernels/overloads.jl diff --git a/Project.toml b/Project.toml index ee7323fc2..b5a93b983 100644 --- a/Project.toml +++ b/Project.toml @@ -1,6 +1,6 @@ name = "KernelFunctions" uuid = "ec8451be-7e33-11e9-00cf-bbf324bd1392" -version = "0.8.16" +version = "0.8.17" [deps] Compat = "34da2185-b29b-5c13-b0c7-acf172513d20" @@ -13,6 +13,7 @@ Requires = "ae029012-a4dd-5104-9daa-d747884805df" SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b" StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91" StatsFuns = "4c63d2b9-4356-54db-8cca-17b64c39e42c" +TensorCore = "62fd8b95-f654-4bbd-a8a5-9c27f68ccd50" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" ZygoteRules = "700de1a5-db45-46bc-99cf-38207098b444" @@ -24,5 +25,6 @@ Requires = "1.0.1" SpecialFunctions = "0.8, 0.9, 0.10, 1" StatsBase = "0.32, 0.33" StatsFuns = "0.8, 0.9" +TensorCore = "0.1" ZygoteRules = "0.2" julia = "1.3" diff --git a/docs/src/kernels.md b/docs/src/kernels.md index 9ed837c20..7c6844013 100644 --- a/docs/src/kernels.md +++ b/docs/src/kernels.md @@ -268,7 +268,7 @@ The [`KernelProduct`](@ref) is defined as a product of kernels: ### Tensor Product -The [`TensorProduct`](@ref) is defined as: +The [`KernelTensorProduct`](@ref) is defined as a tensor product of kernels: ```math k(x,x';\{k_i\}) = \prod_i k_i(x_i,x'_i) diff --git a/src/KernelFunctions.jl b/src/KernelFunctions.jl index 7885b0917..07c89f2dd 100644 --- a/src/KernelFunctions.jl +++ b/src/KernelFunctions.jl @@ -33,9 +33,8 @@ export LinearKernel, PolynomialKernel export RationalQuadraticKernel, GammaRationalQuadraticKernel export GaborKernel, PiecewisePolynomialKernel export PeriodicKernel, NeuralNetworkKernel -export KernelSum, KernelProduct +export KernelSum, KernelProduct, KernelTensorProduct export TransformedKernel, ScaledKernel -export TensorProduct export Transform, SelectTransform, @@ -56,6 +55,9 @@ export ColVecs, RowVecs export MOInput export IndependentMOKernel, LatentFactorMOKernel +# Reexports +export tensor, ⊗ + using Compat using Requires using Distances, LinearAlgebra @@ -65,6 +67,7 @@ using ZygoteRules: @adjoint, pullback using StatsFuns: logtwo using InteractiveUtils: subtypes using StatsBase +using TensorCore abstract type Kernel end abstract type SimpleKernel <: Kernel end @@ -105,6 +108,7 @@ include(joinpath("matrix", "kernelmatrix.jl")) include(joinpath("kernels", "kernelsum.jl")) include(joinpath("kernels", "kernelproduct.jl")) include(joinpath("kernels", "tensorproduct.jl")) +include(joinpath("kernels", "overloads.jl")) include(joinpath("approximations", "nystrom.jl")) include("generic.jl") diff --git a/src/deprecated.jl b/src/deprecated.jl index 5c3ca3131..eeae1eb4c 100644 --- a/src/deprecated.jl +++ b/src/deprecated.jl @@ -7,3 +7,8 @@ @deprecate PiecewisePolynomialKernel{V}(A::AbstractMatrix{<:Real}) where {V} transform( PiecewisePolynomialKernel{V}(size(A, 1)), LinearTransform(cholesky(A).U) ) + +@deprecate TensorProduct(kernels) KernelTensorProduct(kernels) +@deprecate TensorProduct(kernel::Kernel, kernels::Kernel...) KernelTensorProduct( + kernel, kernels... +) diff --git a/src/kernels/kernelproduct.jl b/src/kernels/kernelproduct.jl index 8201815e0..7f2800077 100644 --- a/src/kernels/kernelproduct.jl +++ b/src/kernels/kernelproduct.jl @@ -41,31 +41,6 @@ end @functor KernelProduct -Base.:*(k1::Kernel, k2::Kernel) = KernelProduct(k1, k2) - -function Base.:*( - k1::KernelProduct{<:AbstractVector{<:Kernel}}, - k2::KernelProduct{<:AbstractVector{<:Kernel}}, -) - return KernelProduct(vcat(k1.kernels, k2.kernels)) -end - -function Base.:*(k1::KernelProduct, k2::KernelProduct) - return KernelProduct(k1.kernels..., k2.kernels...) -end - -function Base.:*(k::Kernel, ks::KernelProduct{<:AbstractVector{<:Kernel}}) - return KernelProduct(vcat(k, ks.kernels)) -end - -Base.:*(k::Kernel, kp::KernelProduct) = KernelProduct(k, kp.kernels...) - -function Base.:*(ks::KernelProduct{<:AbstractVector{<:Kernel}}, k::Kernel) - return KernelProduct(vcat(ks.kernels, k)) -end - -Base.:*(kp::KernelProduct, k::Kernel) = KernelProduct(kp.kernels..., k) - Base.length(k::KernelProduct) = length(k.kernels) (κ::KernelProduct)(x, y) = prod(k(x, y) for k in κ.kernels) diff --git a/src/kernels/kernelsum.jl b/src/kernels/kernelsum.jl index 176e870ea..017642d95 100644 --- a/src/kernels/kernelsum.jl +++ b/src/kernels/kernelsum.jl @@ -41,28 +41,6 @@ end @functor KernelSum -Base.:+(k1::Kernel, k2::Kernel) = KernelSum(k1, k2) - -function Base.:+( - k1::KernelSum{<:AbstractVector{<:Kernel}}, k2::KernelSum{<:AbstractVector{<:Kernel}} -) - return KernelSum(vcat(k1.kernels, k2.kernels)) -end - -Base.:+(k1::KernelSum, k2::KernelSum) = KernelSum(k1.kernels..., k2.kernels...) - -function Base.:+(k::Kernel, ks::KernelSum{<:AbstractVector{<:Kernel}}) - return KernelSum(vcat(k, ks.kernels)) -end - -Base.:+(k::Kernel, ks::KernelSum) = KernelSum(k, ks.kernels...) - -function Base.:+(ks::KernelSum{<:AbstractVector{<:Kernel}}, k::Kernel) - return KernelSum(vcat(ks.kernels, k)) -end - -Base.:+(ks::KernelSum, k::Kernel) = KernelSum(ks.kernels..., k) - Base.length(k::KernelSum) = length(k.kernels) (κ::KernelSum)(x, y) = sum(k(x, y) for k in κ.kernels) diff --git a/src/kernels/overloads.jl b/src/kernels/overloads.jl new file mode 100644 index 000000000..b3ba76c7f --- /dev/null +++ b/src/kernels/overloads.jl @@ -0,0 +1,22 @@ +for (M, op, T) in ( + (:Base, :+, :KernelSum), + (:Base, :*, :KernelProduct), + (:TensorCore, :tensor, :KernelTensorProduct), +) + @eval begin + $M.$op(k1::Kernel, k2::Kernel) = $T(k1, k2) + + $M.$op(k1::$T, k2::$T) = $T(k1.kernels..., k2.kernels...) + function $M.$op( + k1::$T{<:AbstractVector{<:Kernel}}, k2::$T{<:AbstractVector{<:Kernel}} + ) + return $T(vcat(k1.kernels, k2.kernels)) + end + + $M.$op(k::Kernel, ks::$T) = $T(k, ks.kernels...) + $M.$op(k::Kernel, ks::$T{<:AbstractVector{<:Kernel}}) = $T(vcat(k, ks.kernels)) + + $M.$op(ks::$T, k::Kernel) = $T(ks.kernels..., k) + $M.$op(ks::$T{<:AbstractVector{<:Kernel}}, k::Kernel) = $T(vcat(ks.kernels, k)) + end +end diff --git a/src/kernels/tensorproduct.jl b/src/kernels/tensorproduct.jl index c1b18212b..912a4d844 100644 --- a/src/kernels/tensorproduct.jl +++ b/src/kernels/tensorproduct.jl @@ -1,30 +1,58 @@ """ - TensorProduct(kernels...) + KernelTensorProduct <: Kernel -Create a tensor product kernel from kernels ``k_1, \\ldots, k_n``, i.e., -a kernel ``k`` that is given by +Tensor product kernel. + +## Definition + +For inputs ``x = (x_1, \\ldots, x_n)`` and ``x' = (x'_1, \\ldots, x'_n)``, a tensor +product kernel from kernels ``k_1, \\ldots, k_n`` is defined as ```math -k(x, y) = \\prod_{i=1}^n k_i(x_i, y_i). +k(x, x'; k_1, \\ldots, k_n) = \\prod_{i=1}^n k_i(x_i, y_i). +``` + +## Construction + +The simplest way to specify a `KernelTensorProduct` is to use the overloaded `tensor` +operator or its alias `⊗` (can be typed by `\\otimes`). +```jldoctest tensorproduct +julia> k1 = SqExponentialKernel(); k2 = LinearKernel(); X = [rand(2) for _ in 1:5]; + +julia> kernelmatrix(k1 ⊗ k2, X) == kernelmatrix(k1, first.(X)) .* kernelmatrix(k2, last.(X)) +true + +julia> kernelmatrix(k, X) == kernelmatrix(k1 + k2, X) +true ``` -The `kernels` can be specified as individual arguments, a tuple, or an iterable data -structure such as an array. Using a tuple or individual arguments guarantees that -`TensorProduct` is concretely typed but might lead to large compilation times if the -number of kernels is large. +You can also specify a `KernelTensorProduct` by providing kernels as individual arguments +or as an iterable data structure such as a `Tuple` or a `Vector`. Using a tuple or +individual arguments guarantees that `KernelTensorProduct` is concretely typed but might +lead to large compilation times if the number of kernels is large. +```jldoctest tensorproduct +julia> KernelTensorProduct(k1, k2) == k1 ⊗ k2 +true + +julia> KernelTensorProduct((k1, k2)) == k1 ⊗ k2 +true + +julia> KernelTensorProduct([k1, k2]) == k1 ⊗ k2 +true +``` """ -struct TensorProduct{K} <: Kernel +struct KernelTensorProduct{K} <: Kernel kernels::K end -function TensorProduct(kernel::Kernel, kernels::Kernel...) - return TensorProduct((kernel, kernels...)) +function KernelTensorProduct(kernel::Kernel, kernels::Kernel...) + return KernelTensorProduct((kernel, kernels...)) end -@functor TensorProduct +@functor KernelTensorProduct -Base.length(kernel::TensorProduct) = length(kernel.kernels) +Base.length(kernel::KernelTensorProduct) = length(kernel.kernels) -function (kernel::TensorProduct)(x, y) +function (kernel::KernelTensorProduct)(x, y) if !(length(x) == length(y) == length(kernel)) throw(DimensionMismatch("number of kernels and number of features are not consistent")) @@ -32,7 +60,7 @@ are not consistent")) return prod(k(xi, yi) for (k, xi, yi) in zip(kernel.kernels, x, y)) end -function validate_domain(k::TensorProduct, x::AbstractVector) +function validate_domain(k::KernelTensorProduct, x::AbstractVector) return dim(x) == length(k) || error("number of kernels and groups of features are not consistent") end @@ -42,7 +70,7 @@ slices(x::AbstractVector{<:Real}) = (x,) slices(x::ColVecs) = eachrow(x.X) slices(x::RowVecs) = eachcol(x.X) -function kernelmatrix!(K::AbstractMatrix, k::TensorProduct, x::AbstractVector) +function kernelmatrix!(K::AbstractMatrix, k::KernelTensorProduct, x::AbstractVector) validate_inplace_dims(K, x) validate_domain(k, x) @@ -56,7 +84,7 @@ function kernelmatrix!(K::AbstractMatrix, k::TensorProduct, x::AbstractVector) end function kernelmatrix!( - K::AbstractMatrix, k::TensorProduct, x::AbstractVector, y::AbstractVector + K::AbstractMatrix, k::KernelTensorProduct, x::AbstractVector, y::AbstractVector ) validate_inplace_dims(K, x, y) validate_domain(k, x) @@ -70,7 +98,7 @@ function kernelmatrix!( return K end -function kerneldiagmatrix!(K::AbstractVector, k::TensorProduct, x::AbstractVector) +function kerneldiagmatrix!(K::AbstractVector, k::KernelTensorProduct, x::AbstractVector) validate_inplace_dims(K, x) validate_domain(k, x) @@ -83,31 +111,31 @@ function kerneldiagmatrix!(K::AbstractVector, k::TensorProduct, x::AbstractVecto return K end -function kernelmatrix(k::TensorProduct, x::AbstractVector) +function kernelmatrix(k::KernelTensorProduct, x::AbstractVector) validate_domain(k, x) return mapreduce(kernelmatrix, hadamard, k.kernels, slices(x)) end -function kernelmatrix(k::TensorProduct, x::AbstractVector, y::AbstractVector) +function kernelmatrix(k::KernelTensorProduct, x::AbstractVector, y::AbstractVector) validate_domain(k, x) return mapreduce(kernelmatrix, hadamard, k.kernels, slices(x), slices(y)) end -function kerneldiagmatrix(k::TensorProduct, x::AbstractVector) +function kerneldiagmatrix(k::KernelTensorProduct, x::AbstractVector) validate_domain(k, x) return mapreduce(kerneldiagmatrix, hadamard, k.kernels, slices(x)) end -Base.show(io::IO, kernel::TensorProduct) = printshifted(io, kernel, 0) +Base.show(io::IO, kernel::KernelTensorProduct) = printshifted(io, kernel, 0) -function Base.:(==)(x::TensorProduct, y::TensorProduct) +function Base.:(==)(x::KernelTensorProduct, y::KernelTensorProduct) return ( length(x.kernels) == length(y.kernels) && all(kx == ky for (kx, ky) in zip(x.kernels, y.kernels)) ) end -function printshifted(io::IO, kernel::TensorProduct, shift::Int) +function printshifted(io::IO, kernel::KernelTensorProduct, shift::Int) print(io, "Tensor product of ", length(kernel), " kernels:") for k in kernel.kernels print(io, "\n") diff --git a/src/matrix/kernelkroneckermat.jl b/src/matrix/kernelkroneckermat.jl index 384bd1de7..d9d2ff075 100644 --- a/src/matrix/kernelkroneckermat.jl +++ b/src/matrix/kernelkroneckermat.jl @@ -1,11 +1,11 @@ -using .Kronecker +using .Kronecker: Kronecker export kernelkronmat function kernelkronmat(κ::Kernel, X::AbstractVector, dims::Int) @assert iskroncompatible(κ) "The chosen kernel is not compatible for kroenecker matrices (see [`iskroncompatible`](@ref))" k = kernelmatrix(κ, X) - return kronecker(k, dims) + return Kronecker.kronecker(k, dims) end function kernelkronmat( @@ -13,7 +13,7 @@ function kernelkronmat( ) @assert iskroncompatible(κ) "The chosen kernel is not compatible for Kronecker matrices" Ks = kernelmatrix.(κ, X) - return K = reduce(⊗, Ks) + return K = reduce(Kronecker.:⊗, Ks) end """ diff --git a/src/utils.jl b/src/utils.jl index 487405613..7e28c82ef 100644 --- a/src/utils.jl +++ b/src/utils.jl @@ -1,5 +1,3 @@ -hadamard(x, y) = x .* y - # Macro for checking arguments macro check_args(K, param, cond, desc=string(cond)) quote diff --git a/test/kernels/tensorproduct.jl b/test/kernels/tensorproduct.jl index d762112c0..27e64385a 100644 --- a/test/kernels/tensorproduct.jl +++ b/test/kernels/tensorproduct.jl @@ -8,11 +8,11 @@ # kernels k1 = SqExponentialKernel() k2 = ExponentialKernel() - kernel1 = TensorProduct(k1, k2) - kernel2 = TensorProduct([k1, k2]) + kernel1 = KernelTensorProduct(k1, k2) + kernel2 = KernelTensorProduct([k1, k2]) @test kernel1 == kernel2 - @test kernel1.kernels === (k1, k2) === TensorProduct((k1, k2)).kernels + @test kernel1.kernels === (k1, k2) === KernelTensorProduct((k1, k2)).kernels @test length(kernel1) == length(kernel2) == 2 @test_throws DimensionMismatch kernel1(rand(3), rand(3)) @@ -24,14 +24,31 @@ end end + # Overloaded operators. + k3 = LinearKernel() + @test (KernelTensorProduct([k1, k2]) ⊗ KernelTensorProduct([k2, k1])).kernels == + [k1, k2, k2, k1] + @test (KernelTensorProduct([k1, k2]) ⊗ k3).kernels == [k1, k2, k3] + @test (k3 ⊗ KernelTensorProduct([k1, k2])).kernels == [k3, k1, k2] + + @test (KernelTensorProduct((k1, k2)) ⊗ KernelTensorProduct((k2, k1))).kernels == + (k1, k2, k2, k1) + @test (KernelTensorProduct((k1, k2)) ⊗ k3).kernels == (k1, k2, k3) + @test (k3 ⊗ KernelTensorProduct((k1, k2))).kernels == (k3, k1, k2) + + # Deprecations + @test (@test_deprecated TensorProduct(k1, k2)) == k1 ⊗ k2 + @test (@test_deprecated TensorProduct((k1, k2))) == k1 ⊗ k2 + @test (@test_deprecated TensorProduct([k1, k2])) == k1 ⊗ k2 + # Standardised tests. TestUtils.test_interface(kernel1, ColVecs{Float64}) TestUtils.test_interface(kernel1, RowVecs{Float64}) - test_ADs(() -> TensorProduct(SqExponentialKernel(), LinearKernel()); dims=[2, 2]) # ADs = [:ForwardDiff, :ReverseDiff]) - test_params(TensorProduct(k1, k2), (k1, k2)) + test_ADs(() -> KernelTensorProduct(SqExponentialKernel(), LinearKernel()); dims=[2, 2]) + test_params(KernelTensorProduct(k1, k2), (k1, k2)) @testset "single kernel" begin - kernel = TensorProduct(k1) + kernel = KernelTensorProduct(k1) @test length(kernel) == 1 @testset "eval" begin From 3336f2966f87c2c1e12e7b7b3f66f4270f2b1d6f Mon Sep 17 00:00:00 2001 From: David Widmann Date: Sun, 17 Jan 2021 18:38:54 +0100 Subject: [PATCH 02/12] Reorganize tests --- src/KernelFunctions.jl | 2 +- ...ensorproduct.jl => kerneltensorproduct.jl} | 0 test/kernels/kernelproduct.jl | 25 +------------ test/kernels/kernelsum.jl | 25 +------------ ...ensorproduct.jl => kerneltensorproduct.jl} | 14 +------- test/kernels/overloads.jl | 36 +++++++++++++++++++ test/runtests.jl | 5 +-- 7 files changed, 43 insertions(+), 64 deletions(-) rename src/kernels/{tensorproduct.jl => kerneltensorproduct.jl} (100%) rename test/kernels/{tensorproduct.jl => kerneltensorproduct.jl} (71%) create mode 100644 test/kernels/overloads.jl diff --git a/src/KernelFunctions.jl b/src/KernelFunctions.jl index 07c89f2dd..080de52aa 100644 --- a/src/KernelFunctions.jl +++ b/src/KernelFunctions.jl @@ -107,7 +107,7 @@ include(joinpath("kernels", "scaledkernel.jl")) include(joinpath("matrix", "kernelmatrix.jl")) include(joinpath("kernels", "kernelsum.jl")) include(joinpath("kernels", "kernelproduct.jl")) -include(joinpath("kernels", "tensorproduct.jl")) +include(joinpath("kernels", "kerneltensorproduct.jl")) include(joinpath("kernels", "overloads.jl")) include(joinpath("approximations", "nystrom.jl")) include("generic.jl") diff --git a/src/kernels/tensorproduct.jl b/src/kernels/kerneltensorproduct.jl similarity index 100% rename from src/kernels/tensorproduct.jl rename to src/kernels/kerneltensorproduct.jl diff --git a/test/kernels/kernelproduct.jl b/test/kernels/kernelproduct.jl index d35ed3f56..01663b41b 100644 --- a/test/kernels/kernelproduct.jl +++ b/test/kernels/kernelproduct.jl @@ -1,36 +1,13 @@ @testset "kernelproduct" begin - rng = MersenneTwister(123456) - x = rand(rng) * 2 - v1 = rand(rng, 3) - v2 = rand(rng, 3) - k1 = LinearKernel() k2 = SqExponentialKernel() - k3 = RationalQuadraticKernel() - X = rand(rng, 2, 2) - k = KernelProduct(k1, k2) - ks1 = 2.0 * k1 - ks2 = 0.5 * k2 + @test k == KernelProduct([k1, k2]) == KernelProduct((k1, k2)) @test length(k) == 2 @test string(k) == ( "Product of 2 kernels:\n\tLinear Kernel (c = 0.0)\n\tSquared " * "Exponential Kernel" ) - @test k(v1, v2) == (k1 * k2)(v1, v2) - @test (k * k3)(v1, v2) ≈ (k3 * k)(v1, v2) - @test (k1 * k2)(v1, v2) == KernelProduct(k1, k2)(v1, v2) - @test (k * ks1)(v1, v2) ≈ (ks1 * k)(v1, v2) - @test (k * k)(v1, v2) == KernelProduct([k1, k2, k1, k2])(v1, v2) - @test KernelProduct([k1, k2]) == KernelProduct((k1, k2)) == k1 * k2 - - @test (KernelProduct([k1, k2]) * KernelProduct([k2, k1])).kernels == [k1, k2, k2, k1] - @test (KernelProduct([k1, k2]) * k3).kernels == [k1, k2, k3] - @test (k3 * KernelProduct([k1, k2])).kernels == [k3, k1, k2] - - @test (KernelProduct((k1, k2)) * KernelProduct((k2, k1))).kernels == (k1, k2, k2, k1) - @test (KernelProduct((k1, k2)) * k3).kernels == (k1, k2, k3) - @test (k3 * KernelProduct((k1, k2))).kernels == (k3, k1, k2) # Standardised tests. TestUtils.test_interface(k, Float64) diff --git a/test/kernels/kernelsum.jl b/test/kernels/kernelsum.jl index a76f8270b..1042c13f6 100644 --- a/test/kernels/kernelsum.jl +++ b/test/kernels/kernelsum.jl @@ -1,35 +1,12 @@ @testset "kernelsum" begin - rng = MersenneTwister(123456) - x = rand(rng) * 2 - v1 = rand(rng, 3) - v2 = rand(rng, 3) - k1 = LinearKernel() k2 = SqExponentialKernel() - k3 = RationalQuadraticKernel() - X = rand(rng, 2, 2) - k = KernelSum(k1, k2) - ks1 = 2.0 * k1 - ks2 = 0.5 * k2 + @test k == KernelSum([k1, k2]) == KernelSum((k1, k2)) @test length(k) == 2 @test string(k) == ( "Sum of 2 kernels:\n\tLinear Kernel (c = 0.0)\n\tSquared " * "Exponential Kernel" ) - @test k(v1, v2) == (k1 + k2)(v1, v2) - @test (k + k3)(v1, v2) ≈ (k3 + k)(v1, v2) - @test (k1 + k2)(v1, v2) == KernelSum(k1, k2)(v1, v2) - @test (k + ks1)(v1, v2) ≈ (ks1 + k)(v1, v2) - @test (k + k)(v1, v2) == KernelSum([k1, k2, k1, k2])(v1, v2) - @test KernelSum([k1, k2]) == KernelSum((k1, k2)) == k1 + k2 - - @test (KernelSum([k1, k2]) + KernelSum([k2, k1])).kernels == [k1, k2, k2, k1] - @test (KernelSum([k1, k2]) + k3).kernels == [k1, k2, k3] - @test (k3 + KernelSum([k1, k2])).kernels == [k3, k1, k2] - - @test (KernelSum((k1, k2)) + KernelSum((k2, k1))).kernels == (k1, k2, k2, k1) - @test (KernelSum((k1, k2)) + k3).kernels == (k1, k2, k3) - @test (k3 + KernelSum((k1, k2))).kernels == (k3, k1, k2) # Standardised tests. TestUtils.test_interface(k, Float64) diff --git a/test/kernels/tensorproduct.jl b/test/kernels/kerneltensorproduct.jl similarity index 71% rename from test/kernels/tensorproduct.jl rename to test/kernels/kerneltensorproduct.jl index 27e64385a..166f7967a 100644 --- a/test/kernels/tensorproduct.jl +++ b/test/kernels/kerneltensorproduct.jl @@ -1,4 +1,4 @@ -@testset "tensorproduct" begin +@testset "kerneltensorproduct" begin rng = MersenneTwister(123456) u1 = rand(rng, 10) u2 = rand(rng, 10) @@ -24,18 +24,6 @@ end end - # Overloaded operators. - k3 = LinearKernel() - @test (KernelTensorProduct([k1, k2]) ⊗ KernelTensorProduct([k2, k1])).kernels == - [k1, k2, k2, k1] - @test (KernelTensorProduct([k1, k2]) ⊗ k3).kernels == [k1, k2, k3] - @test (k3 ⊗ KernelTensorProduct([k1, k2])).kernels == [k3, k1, k2] - - @test (KernelTensorProduct((k1, k2)) ⊗ KernelTensorProduct((k2, k1))).kernels == - (k1, k2, k2, k1) - @test (KernelTensorProduct((k1, k2)) ⊗ k3).kernels == (k1, k2, k3) - @test (k3 ⊗ KernelTensorProduct((k1, k2))).kernels == (k3, k1, k2) - # Deprecations @test (@test_deprecated TensorProduct(k1, k2)) == k1 ⊗ k2 @test (@test_deprecated TensorProduct((k1, k2))) == k1 ⊗ k2 diff --git a/test/kernels/overloads.jl b/test/kernels/overloads.jl new file mode 100644 index 000000000..eb79d41f8 --- /dev/null +++ b/test/kernels/overloads.jl @@ -0,0 +1,36 @@ +@testset "overloads" begin + rng = MersenneTwister(123456) + + k1 = LinearKernel() + k2 = SqExponentialKernel() + k3 = RationalQuadraticKernel() + + for (op, T) in ((+, KernelSum), (*, KernelProduct), (⊗, KernelTensorProduct)) + if T === KernelTensorProduct + v2_1 = rand(rng, 2) + v2_2 = rand(rng, 2) + v3_1 = rand(rng, 3) + v3_2 = rand(rng, 3) + v4_1 = rand(rng, 4) + v4_2 = rand(rng, 4) + else + v2_1 = v3_1 = v4_1 = rand(rng, 3) + v2_2 = v3_2 = v4_2 = rand(rng, 3) + end + k = T(k1, k2) + + @test op(k1, k2)(v2_1, v2_2) == k(v2_1, v2_2) + @test op(k, k3)(v3_1, v3_2) == T((k1, k2, k3))(v3_1, v3_2) + @test op(k3, k)(v3_1, v3_2) == T((k3, k1, k2))(v3_1, v3_2) + @test op(k, k)(v4_1, v4_2) == T((k1, k2, k1, k2))(v4_1, v4_2) + @test op(k1, k2) == T([k1, k2]) == T((k1, k2)) + + @test op(T([k1, k2]), T([k2, k1])).kernels == [k1, k2, k2, k1] + @test op(T([k1, k2]), k3).kernels == [k1, k2, k3] + @test op(k3, T([k1, k2])).kernels == [k3, k1, k2] + + @test op(T((k1, k2)), T((k2, k1))).kernels == (k1, k2, k2, k1) + @test op(T((k1, k2)), k3).kernels == (k1, k2, k3) + @test op(k3, T((k1, k2))).kernels == (k3, k1, k2) + end +end diff --git a/test/runtests.jl b/test/runtests.jl index f04c934ac..aae8a6e75 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -1,7 +1,7 @@ using KernelFunctions using AxisArrays using Distances -using Kronecker +using Kronecker: Kronecker using LinearAlgebra using PDMats using Random @@ -119,8 +119,9 @@ include("test_utils.jl") @testset "kernels" begin include(joinpath("kernels", "kernelproduct.jl")) include(joinpath("kernels", "kernelsum.jl")) + include(joinpath("kernels", "kerneltensorproduct.jl")) + include(joinpath("kernels", "overloads.jl")) include(joinpath("kernels", "scaledkernel.jl")) - include(joinpath("kernels", "tensorproduct.jl")) include(joinpath("kernels", "transformedkernel.jl")) end @info "Ran tests on Kernel" From 1351d0c809d653d5e6fde25186afb8597d440349 Mon Sep 17 00:00:00 2001 From: David Widmann Date: Sun, 17 Jan 2021 18:57:35 +0100 Subject: [PATCH 03/12] Fix docs --- docs/src/api.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/src/api.md b/docs/src/api.md index 256e6bc8f..4bea296ed 100644 --- a/docs/src/api.md +++ b/docs/src/api.md @@ -50,7 +50,7 @@ TransformedKernel ScaledKernel KernelSum KernelProduct -TensorProduct +KernelTensorProduct ``` ## Transforms From 88e2d09dc7abba4b99ad9c0312f88647016f49f0 Mon Sep 17 00:00:00 2001 From: David Widmann Date: Sun, 17 Jan 2021 19:10:48 +0100 Subject: [PATCH 04/12] Fix docstring --- src/kernels/kerneltensorproduct.jl | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/kernels/kerneltensorproduct.jl b/src/kernels/kerneltensorproduct.jl index 912a4d844..799a690de 100644 --- a/src/kernels/kerneltensorproduct.jl +++ b/src/kernels/kerneltensorproduct.jl @@ -1,14 +1,14 @@ """ KernelTensorProduct <: Kernel -Tensor product kernel. +Tensor product of kernels. ## Definition -For inputs ``x = (x_1, \\ldots, x_n)`` and ``x' = (x'_1, \\ldots, x'_n)``, a tensor -product kernel from kernels ``k_1, \\ldots, k_n`` is defined as +For inputs ``x = (x_1, \\ldots, x_n)`` and ``x' = (x'_1, \\ldots, x'_n)``, the tensor +product of kernels ``k_1, \\ldots, k_n`` is defined as ```math -k(x, x'; k_1, \\ldots, k_n) = \\prod_{i=1}^n k_i(x_i, y_i). +k(x, x'; k_1, \\ldots, k_n) = \\Big(\\bigotimes_{i=1}^n k_i\\Big)(x, x') = \\prod_{i=1}^n k_i(x_i, x'_i). ``` ## Construction From 4863c24d103684a7be1659c86872c705da22866e Mon Sep 17 00:00:00 2001 From: David Widmann Date: Sun, 17 Jan 2021 21:06:07 +0100 Subject: [PATCH 05/12] Fix docstring --- src/kernels/kerneltensorproduct.jl | 3 --- 1 file changed, 3 deletions(-) diff --git a/src/kernels/kerneltensorproduct.jl b/src/kernels/kerneltensorproduct.jl index 799a690de..c43ef3cc7 100644 --- a/src/kernels/kerneltensorproduct.jl +++ b/src/kernels/kerneltensorproduct.jl @@ -20,9 +20,6 @@ julia> k1 = SqExponentialKernel(); k2 = LinearKernel(); X = [rand(2) for _ in 1: julia> kernelmatrix(k1 ⊗ k2, X) == kernelmatrix(k1, first.(X)) .* kernelmatrix(k2, last.(X)) true - -julia> kernelmatrix(k, X) == kernelmatrix(k1 + k2, X) -true ``` You can also specify a `KernelTensorProduct` by providing kernels as individual arguments From 28bbb6533edc0dbf807d5b26f21cad17c2a5e098 Mon Sep 17 00:00:00 2001 From: David Widmann Date: Sun, 17 Jan 2021 22:38:44 +0100 Subject: [PATCH 06/12] More fixes --- docs/make.jl | 1 + src/basekernels/sm.jl | 2 +- src/kernels/kerneltensorproduct.jl | 4 ++-- test/kernels/kernelproduct.jl | 4 +--- test/kernels/kernelsum.jl | 6 +----- test/kernels/kerneltensorproduct.jl | 6 +++++- 6 files changed, 11 insertions(+), 12 deletions(-) diff --git a/docs/make.jl b/docs/make.jl index ff16b17de..63813838d 100644 --- a/docs/make.jl +++ b/docs/make.jl @@ -23,6 +23,7 @@ makedocs(; "Custom Kernels" => "create_kernel.md", "API" => "api.md", ], + strict=true, ) deploydocs(; diff --git a/src/basekernels/sm.jl b/src/basekernels/sm.jl index b7d51f41a..833945295 100644 --- a/src/basekernels/sm.jl +++ b/src/basekernels/sm.jl @@ -92,7 +92,7 @@ function spectral_mixture_product_kernel( if !(size(αs) == size(γs) == size(ωs)) throw(DimensionMismatch("The dimensions of αs, γs, ans ωs do not match")) end - return TensorProduct( + return KernelTensorProduct( spectral_mixture_kernel(h, α, reshape(γ, 1, :), reshape(ω, 1, :)) for (α, γ, ω) in zip(eachrow(αs), eachrow(γs), eachrow(ωs)) ) diff --git a/src/kernels/kerneltensorproduct.jl b/src/kernels/kerneltensorproduct.jl index c43ef3cc7..d51cea78b 100644 --- a/src/kernels/kerneltensorproduct.jl +++ b/src/kernels/kerneltensorproduct.jl @@ -16,9 +16,9 @@ k(x, x'; k_1, \\ldots, k_n) = \\Big(\\bigotimes_{i=1}^n k_i\\Big)(x, x') = \\pro The simplest way to specify a `KernelTensorProduct` is to use the overloaded `tensor` operator or its alias `⊗` (can be typed by `\\otimes`). ```jldoctest tensorproduct -julia> k1 = SqExponentialKernel(); k2 = LinearKernel(); X = [rand(2) for _ in 1:5]; +julia> k1 = SqExponentialKernel(); k2 = LinearKernel(); X = rand(5, 2); -julia> kernelmatrix(k1 ⊗ k2, X) == kernelmatrix(k1, first.(X)) .* kernelmatrix(k2, last.(X)) +julia> kernelmatrix(k1 ⊗ k2, RowVecs(X)) == kernelmatrix(k1, X[:, 1]) .* kernelmatrix(k2, X[:, 2]) true ``` diff --git a/test/kernels/kernelproduct.jl b/test/kernels/kernelproduct.jl index 01663b41b..670a5bdf0 100644 --- a/test/kernels/kernelproduct.jl +++ b/test/kernels/kernelproduct.jl @@ -12,9 +12,7 @@ # Standardised tests. TestUtils.test_interface(k, Float64) test_ADs( - x -> SqExponentialKernel() * LinearKernel(; c=x[1]), - rand(1); - ADs=[:ForwardDiff, :ReverseDiff, :Zygote], + x -> KernelProduct(SqExponentialKernel(), LinearKernel(; c=exp(x[1]))), rand(1) ) test_params(k1 * k2, (k1, k2)) diff --git a/test/kernels/kernelsum.jl b/test/kernels/kernelsum.jl index 1042c13f6..da69fdcc7 100644 --- a/test/kernels/kernelsum.jl +++ b/test/kernels/kernelsum.jl @@ -10,11 +10,7 @@ # Standardised tests. TestUtils.test_interface(k, Float64) - test_ADs( - x -> KernelSum(SqExponentialKernel(), LinearKernel(; c=x[1])), - rand(1); - ADs=[:ForwardDiff, :ReverseDiff, :Zygote], - ) + test_ADs(x -> KernelSum(SqExponentialKernel(), LinearKernel(; c=exp(x[1]))), rand(1)) test_params(k1 + k2, (k1, k2)) end diff --git a/test/kernels/kerneltensorproduct.jl b/test/kernels/kerneltensorproduct.jl index 166f7967a..6b54f3dd0 100644 --- a/test/kernels/kerneltensorproduct.jl +++ b/test/kernels/kerneltensorproduct.jl @@ -32,7 +32,11 @@ # Standardised tests. TestUtils.test_interface(kernel1, ColVecs{Float64}) TestUtils.test_interface(kernel1, RowVecs{Float64}) - test_ADs(() -> KernelTensorProduct(SqExponentialKernel(), LinearKernel()); dims=[2, 2]) + test_ADs( + x -> KernelTensorProduct(SqExponentialKernel(), LinearKernel(; c=exp(x[1]))), + rand(1); + dims=[2, 2], + ) test_params(KernelTensorProduct(k1, k2), (k1, k2)) @testset "single kernel" begin From 96603ce69fc827ecee13b673ef0b8c159b547f2e Mon Sep 17 00:00:00 2001 From: David Widmann Date: Sun, 17 Jan 2021 22:49:47 +0100 Subject: [PATCH 07/12] Add documentation --- docs/make.jl | 1 + docs/src/api.md | 10 ++++++++++ 2 files changed, 11 insertions(+) diff --git a/docs/make.jl b/docs/make.jl index 63813838d..fab7cc2d4 100644 --- a/docs/make.jl +++ b/docs/make.jl @@ -24,6 +24,7 @@ makedocs(; "API" => "api.md", ], strict=true, + checkdocs=:exports, ) deploydocs(; diff --git a/docs/src/api.md b/docs/src/api.md index 4bea296ed..6cb250024 100644 --- a/docs/src/api.md +++ b/docs/src/api.md @@ -23,12 +23,17 @@ EyeKernel ZeroKernel CosineKernel SqExponentialKernel +GaussianKernel +RBFKernel +SEKernel ExponentialKernel +LaplacianKernel GammaExponentialKernel ExponentiatedKernel FBMKernel GaborKernel MaternKernel +Matern12Kernel Matern32Kernel Matern52Kernel NeuralNetworkKernel @@ -41,6 +46,9 @@ spectral_mixture_kernel spectral_mixture_product_kernel PeriodicKernel WienerKernel +MOKernel +IndependentMOKernel +LatentFactorMOKernel ``` ## Composite Kernels @@ -64,6 +72,7 @@ LinearTransform FunctionTransform SelectTransform ChainTransform +PeriodicTransform ``` ## Functions @@ -84,6 +93,7 @@ transform ```@docs ColVecs RowVecs +MOInput NystromFact ``` From 343fef0941fe53590bba4cbbc69b069e91e42622 Mon Sep 17 00:00:00 2001 From: David Widmann Date: Sun, 17 Jan 2021 23:01:36 +0100 Subject: [PATCH 08/12] Fix errors in docs --- docs/Project.toml | 4 ++++ docs/make.jl | 8 +++++++- docs/src/userguide.md | 2 +- 3 files changed, 12 insertions(+), 2 deletions(-) diff --git a/docs/Project.toml b/docs/Project.toml index 1b7a1fadc..41b8530da 100644 --- a/docs/Project.toml +++ b/docs/Project.toml @@ -1,5 +1,9 @@ [deps] Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4" +Kronecker = "2c470bb0-bcc8-11e8-3dad-c9649493f05e" +PDMats = "90014a1f-27ba-587c-ab20-58faa44d9150" [compat] Documenter = "0.23, 0.24, 0.25, 0.26" +Kronecker = "0.4" +PDMats = "0.10" diff --git a/docs/make.jl b/docs/make.jl index fab7cc2d4..836d46f57 100644 --- a/docs/make.jl +++ b/docs/make.jl @@ -4,7 +4,13 @@ using KernelFunctions DocMeta.setdocmeta!( KernelFunctions, :DocTestSetup, - :(using KernelFunctions, LinearAlgebra, Random); + quote + using KernelFunctions + using LinearAlgebra + using Random + using Kronecker: Kronecker + using PDMats: PDMats + end; recursive=true, ) diff --git a/docs/src/userguide.md b/docs/src/userguide.md index 622686043..240830624 100644 --- a/docs/src/userguide.md +++ b/docs/src/userguide.md @@ -52,7 +52,7 @@ To evaluate the kernel function on two vectors you simply call the kernel object ## Creating a kernel matrix Kernel matrices can be created via the `kernelmatrix` function or `kerneldiagmatrix` for only the diagonal. -An important argument to give is the data layout of the input `obsdim`. It specifies whether the number of observed data points is along the first dimension (`obsdim=1`, i.e. the matrix shape is number of samples times number of features) or along the second dimension (`obsdim=2`, i.e. the matrix shape is number of features times number of samples), similarly to [Distances.jl](https://github.com/JuliaStats/Distances.jl). If not given explicitly, `obsdim` defaults to [`defaultobs`](@ref). +An important argument to give is the data layout of the input `obsdim`. It specifies whether the number of observed data points is along the first dimension (`obsdim=1`, i.e. the matrix shape is number of samples times number of features) or along the second dimension (`obsdim=2`, i.e. the matrix shape is number of features times number of samples), similarly to [Distances.jl](https://github.com/JuliaStats/Distances.jl). If not given explicitly, `obsdim` defaults to `2`. For example: ```julia k = SqExponentialKernel() From 200016dd142cd3c8a4b6dfb82d1e95612fc2cc1a Mon Sep 17 00:00:00 2001 From: David Widmann Date: Sun, 17 Jan 2021 23:06:09 +0100 Subject: [PATCH 09/12] Remove `kernelkronmat` from docs --- docs/Project.toml | 2 -- docs/make.jl | 1 - docs/src/api.md | 1 - 3 files changed, 4 deletions(-) diff --git a/docs/Project.toml b/docs/Project.toml index 41b8530da..0e0cc849e 100644 --- a/docs/Project.toml +++ b/docs/Project.toml @@ -1,9 +1,7 @@ [deps] Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4" -Kronecker = "2c470bb0-bcc8-11e8-3dad-c9649493f05e" PDMats = "90014a1f-27ba-587c-ab20-58faa44d9150" [compat] Documenter = "0.23, 0.24, 0.25, 0.26" -Kronecker = "0.4" PDMats = "0.10" diff --git a/docs/make.jl b/docs/make.jl index 836d46f57..fb0726753 100644 --- a/docs/make.jl +++ b/docs/make.jl @@ -8,7 +8,6 @@ DocMeta.setdocmeta!( using KernelFunctions using LinearAlgebra using Random - using Kronecker: Kronecker using PDMats: PDMats end; recursive=true, diff --git a/docs/src/api.md b/docs/src/api.md index 6cb250024..38bc18cc8 100644 --- a/docs/src/api.md +++ b/docs/src/api.md @@ -83,7 +83,6 @@ kernelmatrix! kerneldiagmatrix kerneldiagmatrix! kernelpdmat -kernelkronmat nystrom transform ``` From a1132d2e594605981ef6ba768ad5c172a5f0a03d Mon Sep 17 00:00:00 2001 From: David Widmann Date: Mon, 18 Jan 2021 09:50:42 +0100 Subject: [PATCH 10/12] Change heading --- src/kernels/kerneltensorproduct.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/kernels/kerneltensorproduct.jl b/src/kernels/kerneltensorproduct.jl index d51cea78b..5ca88677b 100644 --- a/src/kernels/kerneltensorproduct.jl +++ b/src/kernels/kerneltensorproduct.jl @@ -3,7 +3,7 @@ Tensor product of kernels. -## Definition +# Definition For inputs ``x = (x_1, \\ldots, x_n)`` and ``x' = (x'_1, \\ldots, x'_n)``, the tensor product of kernels ``k_1, \\ldots, k_n`` is defined as @@ -11,7 +11,7 @@ product of kernels ``k_1, \\ldots, k_n`` is defined as k(x, x'; k_1, \\ldots, k_n) = \\Big(\\bigotimes_{i=1}^n k_i\\Big)(x, x') = \\prod_{i=1}^n k_i(x_i, x'_i). ``` -## Construction +# Construction The simplest way to specify a `KernelTensorProduct` is to use the overloaded `tensor` operator or its alias `⊗` (can be typed by `\\otimes`). From ea1803b7f2dc07bdd7eb47682bc7997e16e693e9 Mon Sep 17 00:00:00 2001 From: David Widmann Date: Mon, 18 Jan 2021 13:22:25 +0100 Subject: [PATCH 11/12] Add explanation for Kronecker import --- src/matrix/kernelkroneckermat.jl | 3 +++ 1 file changed, 3 insertions(+) diff --git a/src/matrix/kernelkroneckermat.jl b/src/matrix/kernelkroneckermat.jl index d9d2ff075..26f7d65a9 100644 --- a/src/matrix/kernelkroneckermat.jl +++ b/src/matrix/kernelkroneckermat.jl @@ -1,3 +1,6 @@ +# Since Kronecker does not implement `TensorCore.:⊗` but instead exports its own function +# `Kronecker.:⊗`, only the module is imported and Kronecker.:⊗ and Kronecker.kronecker are +# called explicitly. using .Kronecker: Kronecker export kernelkronmat From 4456f610217b600252cb51c8bf5ce464c444a575 Mon Sep 17 00:00:00 2001 From: David Widmann Date: Wed, 20 Jan 2021 13:12:56 +0100 Subject: [PATCH 12/12] Add test of representation --- src/kernels/kerneltensorproduct.jl | 3 +-- test/kernels/kerneltensorproduct.jl | 3 +++ 2 files changed, 4 insertions(+), 2 deletions(-) diff --git a/src/kernels/kerneltensorproduct.jl b/src/kernels/kerneltensorproduct.jl index 5ca88677b..16bed14f3 100644 --- a/src/kernels/kerneltensorproduct.jl +++ b/src/kernels/kerneltensorproduct.jl @@ -1,5 +1,5 @@ """ - KernelTensorProduct <: Kernel + KernelTensorProduct Tensor product of kernels. @@ -139,7 +139,6 @@ function printshifted(io::IO, kernel::KernelTensorProduct, shift::Int) for _ in 1:(shift + 1) print(io, "\t") end - print(io, "- ") printshifted(io, k, shift + 2) end end diff --git a/test/kernels/kerneltensorproduct.jl b/test/kernels/kerneltensorproduct.jl index 6b54f3dd0..e01fa3871 100644 --- a/test/kernels/kerneltensorproduct.jl +++ b/test/kernels/kerneltensorproduct.jl @@ -14,6 +14,9 @@ @test kernel1 == kernel2 @test kernel1.kernels === (k1, k2) === KernelTensorProduct((k1, k2)).kernels @test length(kernel1) == length(kernel2) == 2 + @test string(kernel1) == ( + "Tensor product of 2 kernels:\n\tSquared Exponential Kernel\n\tExponential Kernel" + ) @test_throws DimensionMismatch kernel1(rand(3), rand(3)) @testset "val" begin