Merge branch 'master' into tgw/heterotopic

Author: thomasgudjonwright

.github/workflows/ci.yml

@@ -58,20 +58,8 @@ jobs:
           GROUP: ${{ matrix.group }}
       - uses: julia-actions/julia-processcoverage@v1
         if: matrix.version == '1' && matrix.os == 'ubuntu-latest'
-      - name: Coveralls parallel
+      - name: Send coverage to CodeCov
         if: matrix.version == '1' && matrix.os == 'ubuntu-latest'
-        uses: coverallsapp/github-action@master
+        uses: codecov/codecov-action@v2
         with:
-          github-token: ${{ secrets.GITHUB_TOKEN }}
-          path-to-lcov: ./lcov.info
-          flag-name: run-${{ matrix.group }}
-          parallel: true
-  finish:
-    needs: test
-    runs-on: ubuntu-latest
-    steps:
-      - name: Send coverage
-        uses: coverallsapp/github-action@master
-        with:
-          github-token: ${{ secrets.github_token }}
-          parallel-finished: true
+            file: lcov.info

Project.toml

@@ -1,6 +1,6 @@
 name = "KernelFunctions"
 uuid = "ec8451be-7e33-11e9-00cf-bbf324bd1392"
-version = "0.10.24"
+version = "0.10.27"
 
 [deps]
 ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4"

README.md

@@ -1,7 +1,7 @@
 # KernelFunctions.jl
 
 ![CI](https://github.com/JuliaGaussianProcesses/KernelFunctions.jl/workflows/CI/badge.svg?branch=master)
-[![Coverage Status](https://coveralls.io/repos/github/JuliaGaussianProcesses/KernelFunctions.jl/badge.svg?branch=master)](https://coveralls.io/github/JuliaGaussianProcesses/KernelFunctions.jl?branch=master)
+[![codecov](https://codecov.io/gh/JuliaGaussianProcesses/KernelFunctions.jl/branch/master/graph/badge.svg?token=rmDh3gb7hN)](https://codecov.io/gh/JuliaGaussianProcesses/KernelFunctions.jl)
 [![Documentation (stable)](https://img.shields.io/badge/docs-stable-blue.svg)](https://juliagaussianprocesses.github.io/KernelFunctions.jl/stable)
 [![Documentation (latest)](https://img.shields.io/badge/docs-dev-blue.svg)](https://juliagaussianprocesses.github.io/KernelFunctions.jl/dev)
 [![ColPrac: Contributor's Guide on Collaborative Practices for Community Packages](https://img.shields.io/badge/ColPrac-Contributor's%20Guide-blueviolet)](https://github.com/SciML/ColPrac)

docs/src/api.md

@@ -84,7 +84,6 @@ To find out more about the background, read this [review of kernels for vector-v
 
 KernelFunctions also provides miscellaneous utility functions.
 ```@docs
-kernelpdmat
 nystrom
 NystromFact
 ```
@@ -96,4 +95,11 @@ To keep the dependencies of KernelFunctions lean, some functionality is only ava
 [*https://github.com/MichielStock/Kronecker.jl*](https://github.com/MichielStock/Kronecker.jl)
 ```@docs
 kronecker_kernelmatrix
+kernelkronmat
 ```
+
+### PDMats.jl
+[*https://github.com/JuliaStats/PDMats.jl*](https://github.com/JuliaStats/PDMats.jl)
+```@docs
+kernelpdmat
+```

examples/gaussian-process-priors/script.jl

@@ -120,6 +120,7 @@ kernels = [
     LinearKernel(),
     compose(PeriodicKernel(), ScaleTransform(0.2)),
     NeuralNetworkKernel(),
+    GibbsKernel(; lengthscale=x -> sum(exp ∘ sin, x)),
 ]
 plot(
     [visualize(k) for k in kernels]...;

src/basekernels/rational.jl

@@ -5,11 +5,11 @@ Rational kernel with shape parameter `α` and given `metric`.
 
 # Definition
 
-For inputs ``x, x'``, the rational kernel with shape parameter
-``\\alpha > 0`` is defined as
+For inputs ``x, x'`` and metric ``d(\\cdot, \\cdot)``, the rational kernel with shape parameter ``\\alpha > 0`` is defined as
 ```math
-k(x, x'; \\alpha) = \\bigg(1 + \\frac{\\|x - x'\\|}{\\alpha}\\bigg)^{-\\alpha}.
+k(x, x'; \\alpha) = \\bigg(1 + \\frac{d(x, x')}{\\alpha}\\bigg)^{-\\alpha}.
 ```
+By default, ``d`` is the Euclidean metric ``d(x, x') = \\|x - x'\\|_2``.
 
 The [`ExponentialKernel`](@ref) is recovered in the limit as ``\\alpha \\to \\infty``.
 

src/distances/delta.jl

@@ -1,10 +1,18 @@
 # Delta is not following the PreMetric rules since d(x, x) == 1
 struct Delta <: Distances.UnionPreMetric end
 
-@inline Distances.eval_op(::Delta, a::Real, b::Real) = a == b
-@inline Distances.eval_reduce(::Delta, a, b) = a && b
-@inline Distances.eval_start(::Delta, a, b) = true
-@inline (dist::Delta)(a::AbstractArray, b::AbstractArray) = Distances._evaluate(dist, a, b)
-@inline (dist::Delta)(a::Number, b::Number) = a == b
+(dist::Delta)(a::Number, b::Number) = a == b
+Base.@propagate_inbounds function (dist::Delta)(
+    a::AbstractArray{<:Number}, b::AbstractArray{<:Number}
+)
+    @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
 
 Distances.result_type(::Delta, Ta::Type, Tb::Type) = Bool

src/distances/pairwise.jl

@@ -7,7 +7,7 @@ end
 pairwise(d::PreMetric, X::AbstractVector) = pairwise(d, X, X)
 
 function pairwise!(out::AbstractMatrix, d::PreMetric, X::AbstractVector, Y::AbstractVector)
-    return broadcast!(d, out, X, Y')
+    return broadcast!(d, out, X, permutedims(Y))
 end
 
 pairwise!(out::AbstractMatrix, d::PreMetric, X::AbstractVector) = pairwise!(out, d, X, X)

src/kernels/kerneltensorproduct.jl

@@ -71,10 +71,17 @@ function kernelmatrix!(K::AbstractMatrix, k::KernelTensorProduct, x::AbstractVec
     validate_inplace_dims(K, x)
     validate_domain(k, x)
 
-    kernels_and_inputs = zip(k.kernels, slices(x))
-    kernelmatrix!(K, first(kernels_and_inputs)...)
-    for (k, xi) in Iterators.drop(kernels_and_inputs, 1)
-        K .*= kernelmatrix(k, xi)
+    first_kernels, tail_kernels = Iterators.peel(k.kernels)
+    first_x, tail_x = Iterators.peel(slices(x))
+
+    # handle first kernel and input
+    kernelmatrix!(K, first_kernels, first_x)
+
+    # handle remaining kernels and inputs
+    Ktmp = similar(K)
+    for (ki, xi) in zip(tail_kernels, tail_x)
+        kernelmatrix!(Ktmp, ki, xi)
+        hadamard!(K, K, Ktmp)
     end
 
     return K
@@ -86,10 +93,18 @@ function kernelmatrix!(
     validate_inplace_dims(K, x, y)
     validate_domain(k, x)
 
-    kernels_and_inputs = zip(k.kernels, slices(x), slices(y))
-    kernelmatrix!(K, first(kernels_and_inputs)...)
-    for (k, xi, yi) in Iterators.drop(kernels_and_inputs, 1)
-        K .*= kernelmatrix(k, xi, yi)
+    first_kernels, tail_kernels = Iterators.peel(k.kernels)
+    first_x, tail_x = Iterators.peel(slices(x))
+    first_y, tail_y = Iterators.peel(slices(y))
+
+    # handle first kernel and inputs
+    kernelmatrix!(K, first_kernels, first_x, first_y)
+
+    # handle remaining kernels and inputs
+    Ktmp = similar(K)
+    for (ki, xi, yi) in zip(tail_kernels, tail_x, tail_y)
+        kernelmatrix!(Ktmp, ki, xi, yi)
+        hadamard!(K, K, Ktmp)
     end
 
     return K
@@ -99,10 +114,40 @@ function kernelmatrix_diag!(K::AbstractVector, k::KernelTensorProduct, x::Abstra
     validate_inplace_dims(K, x)
     validate_domain(k, x)
 
-    kernels_and_inputs = zip(k.kernels, slices(x))
-    kernelmatrix_diag!(K, first(kernels_and_inputs)...)
-    for (k, xi) in Iterators.drop(kernels_and_inputs, 1)
-        K .*= kernelmatrix_diag(k, xi)
+    first_kernels, tail_kernels = Iterators.peel(k.kernels)
+    first_x, tail_x = Iterators.peel(slices(x))
+
+    # handle first kernel and input
+    kernelmatrix_diag!(K, first_kernels, first_x)
+
+    # handle remaining kernels and inputs
+    Ktmp = similar(K)
+    for (ki, xi) in zip(tail_kernels, tail_x)
+        kernelmatrix_diag!(Ktmp, ki, xi)
+        hadamard!(K, K, Ktmp)
+    end
+
+    return K
+end
+
+function kernelmatrix_diag!(
+    K::AbstractVector, k::KernelTensorProduct, x::AbstractVector, y::AbstractVector
+)
+    validate_inplace_dims(K, x, y)
+    validate_domain(k, x)
+
+    first_kernels, tail_kernels = Iterators.peel(k.kernels)
+    first_x, tail_x = Iterators.peel(slices(x))
+    first_y, tail_y = Iterators.peel(slices(y))
+
+    # handle first kernel and inputs
+    kernelmatrix_diag!(K, first_kernels, first_x, first_y)
+
+    # handle remaining kernels and inputs
+    Ktmp = similar(K)
+    for (ki, xi, yi) in zip(tail_kernels, tail_x, tail_y)
+        kernelmatrix_diag!(Ktmp, ki, xi, yi)
+        hadamard!(K, K, Ktmp)
     end
 
     return K

src/kernels/scaledkernel.jl

@@ -61,6 +61,14 @@ function kernelmatrix_diag!(K::AbstractVector, κ::ScaledKernel, x::AbstractVect
     return K
 end
 
+function kernelmatrix_diag!(
+    K::AbstractVector, κ::ScaledKernel, x::AbstractVector, y::AbstractVector
+)
+    kernelmatrix_diag!(K, κ.kernel, x, y)
+    K .*= κ.σ²
+    return K
+end
+
 Base.:*(w::Real, k::Kernel) = ScaledKernel(k, w)
 
 Base.show(io::IO, κ::ScaledKernel) = printshifted(io, κ, 0)

src/matrix/kernelkroneckermat.jl

@@ -6,27 +6,58 @@ using .Kronecker: Kronecker
 export kernelkronmat
 export kronecker_kernelmatrix
 
-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.kronecker(k, dims)
+@doc raw"""
+    kernelkronmat(κ::Kernel, X::AbstractVector{<:Real}, dims::Int) -> KroneckerPower
+
+Return a `KroneckerPower` matrix on the `D`-dimensional input grid constructed by ``\otimes_{i=1}^D X``,
+where `D` is given by `dims`.
+
+!!! warning
+
+    Require `Kronecker.jl` and for `iskroncompatible(κ)` to return `true`.
+"""
+function kernelkronmat(κ::Kernel, X::AbstractVector{<:Real}, dims::Int)
+    checkkroncompatible(κ)
+    K = kernelmatrix(κ, X)
+    return Kronecker.kronecker(K, dims)
 end
 
-function kernelkronmat(
-    κ::Kernel, X::AbstractVector{<:AbstractVector}; obsdim::Int=defaultobs
-)
-    @assert iskroncompatible(κ) "The chosen kernel is not compatible for Kronecker matrices"
+@doc raw"""
+
+    kernelkronmat(κ::Kernel, X::AbstractVector{<:AbstractVector}) -> KroneckerProduct
+
+Returns a `KroneckerProduct` matrix on the grid built with the collection of vectors ``\{X_i\}_{i=1}^D``: ``\otimes_{i=1}^D X_i``.
+
+!!! warning
+
+    Requires `Kronecker.jl` and for `iskroncompatible(κ)` to return `true`.
+"""
+function kernelkronmat(κ::Kernel, X::AbstractVector{<:AbstractVector})
+    checkkroncompatible(κ)
     Ks = kernelmatrix.(κ, X)
-    return K = reduce(Kronecker.:⊗, Ks)
+    return reduce(Kronecker.:⊗, Ks)
 end
 
-"""
-    To be compatible with kroenecker constructions the kernel must satisfy
-    the property : for x,x' ∈ ℜᴰ
-    k(x,x') = ∏ᵢᴰ k(xᵢ,x'ᵢ)
+@doc raw"""
+    iskroncompatible(k::Kernel)
+
+Determine whether kernel `k` is compatible with Kronecker constructions such as [`kernelkronmat`](@ref)
+
+The function returns `false` by default. If `k` is compatible it must satisfy for all ``x, x' \in \mathbb{R}^D`:
+```math
+k(x, x') = \prod_{i=1}^D k(x_i, x'_i).
+```
 """
 @inline iskroncompatible(κ::Kernel) = false # Default return for kernels
 
+function checkkroncompatible(κ::Kernel)
+    return iskroncompatible(κ) || throw(
+        ArgumentError(
+            "The chosen kernel is not compatible for Kronecker matrices (see [`iskroncompatible`](@ref))",
+        ),
+    )
+end
+
 function _kernelmatrix_kroneckerjl_helper(
     ::Type{<:MOInputIsotopicByFeatures}, Kfeatures, Koutputs
 )

src/test_utils.jl

@@ -87,6 +87,7 @@ function test_interface(
 
     tmp_diag = Vector{Float64}(undef, length(x0))
     @test kernelmatrix_diag!(tmp_diag, k, x0) ≈ kernelmatrix_diag(k, x0)
+    @test kernelmatrix_diag!(tmp_diag, k, x0, x1) ≈ kernelmatrix_diag(k, x0, x1)
 end
 
 function test_interface(

src/utils.jl

@@ -197,17 +197,27 @@ function validate_inplace_dims(K::AbstractMatrix, x::AbstractVector, y::Abstract
     end
 end
 
-function validate_inplace_dims(K::AbstractMatrix, x::AbstractVector)
-    return validate_inplace_dims(K, x, x)
-end
-
-function validate_inplace_dims(K::AbstractVector, x::AbstractVector)
-    if length(K) != length(x)
+function validate_inplace_dims(K::AbstractVector, x::AbstractVector, y::AbstractVector)
+    validate_inputs(x, y)
+    n = length(x)
+    if length(y) != n
+        throw(
+            DimensionMismatch(
+                "Length of input x ($n) not consistent with length of input y " *
+                "($(length(y))",
+            ),
+        )
+    end
+    if length(K) != n
         throw(
             DimensionMismatch(
-                "Length of target vector K ($(length(K))) not consistent with length of input" *
-                "vector x ($(length(x))",
+                "Length of target vector K ($(length(K))) not consistent with length of " *
+                "inputs ($n)",
             ),
         )
     end
 end
+
+function validate_inplace_dims(K::AbstractVecOrMat, x::AbstractVector)
+    return validate_inplace_dims(K, x, x)
+end

test/Project.toml

@@ -17,7 +17,7 @@ Zygote = "e88e6eb3-aa80-5325-afca-941959d7151f"
 
 [compat]
 AxisArrays = "0.4.3"
-Distances = "0.9, 0.10"
+Distances = "= 0.10.0, = 0.10.1, = 0.10.2, = 0.10.3, = 0.10.4"
 Documenter = "0.25, 0.26, 0.27"
 FiniteDifferences = "0.10.8, 0.11, 0.12"
 ForwardDiff = "0.10"

test/matrix/kernelkroneckermat.jl

@@ -6,7 +6,7 @@
 
     @test all(collect(kernelkronmat(k, collect(x), 2)) .≈ kernelmatrix(k, X; obsdim=1))
     @test all(collect(kernelkronmat(k, [x, x])) .≈ kernelmatrix(k, X; obsdim=1))
-    @test_throws AssertionError kernelkronmat(LinearKernel(), collect(x), 2)
+    @test_throws ArgumentError kernelkronmat(LinearKernel(), collect(x), 2)
 
     @testset "lazy kernelmatrix" begin
         rng = MersenneTwister(123)