Deprecate implicit obsdim keyword argument (#429)

* Deprecate implicit `obsdim` keyword argument

* Fix format

* Fix typo

* Relax type annotation

* Fix tests

* Update src/utils.jl

Co-authored-by: st-- <[email protected]>

* Fix format

Co-authored-by: github-actions[bot] <41898282+github-actions[bot]@users.noreply.github.com>

* Throw `ArgumentError` instead of assertion

* Update Project.toml

* Improve docstring of `kernelpdmat`

* Update src/utils.jl

Co-authored-by: st-- <[email protected]>

* Update src/utils.jl

Co-authored-by: st-- <[email protected]>

* Fix utils.jl

* Additional test of `vec_of_vecs`

Co-authored-by: st-- <[email protected]>
Co-authored-by: github-actions[bot] <41898282+github-actions[bot]@users.noreply.github.com>

Author: 3 people

Project.toml

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

docs/src/design.md

@@ -139,11 +139,9 @@ In short: many people like matrices, and are familiar with `obsdim`-style keywor
 arguments.
 
 All internals are implemented using `AbstractVector`s though, and the `obsdim` interface
-is just a thin layer of utility functionality which sits on top of this.
-
-
-
-
+is just a thin layer of utility functionality which sits on top of this. To avoid
+confusion and silent errors, we do not favour a specific convention (rows or columns)
+but instead it is necessary to specify the `obsdim` keyword argument explicitly.
 
 ## [Kernels for Multiple-Outputs](@id inputs_for_multiple_outputs)
 

src/approximations/nystrom.jl

@@ -9,9 +9,9 @@ function sampleindex(X::AbstractVector, r::Real)
     return S
 end
 
-@deprecate sampleindex(X::AbstractMatrix, r::Real; obsdim::Integer=defaultobs) sampleindex(
-    vec_of_vecs(X; obsdim=obsdim), r
-) false
+@deprecate sampleindex(
+    X::AbstractMatrix, r::Real; obsdim::Union{Integer,Nothing}=defaultobs
+) sampleindex(vec_of_vecs(X; obsdim=obsdim === nothing ? nothing : Int(obsdim)), r) false
 
 function nystrom_sample(k::Kernel, X::AbstractVector, S::AbstractVector{<:Integer})
     Xₘ = @view X[S]
@@ -21,8 +21,11 @@ function nystrom_sample(k::Kernel, X::AbstractVector, S::AbstractVector{<:Intege
 end
 
 @deprecate nystrom_sample(
-    k::Kernel, X::AbstractMatrix, S::Vector{<:Integer}; obsdim::Integer=defaultobs
-) nystrom_sample(k, vec_of_vecs(X; obsdim=obsdim), S) false
+    k::Kernel,
+    X::AbstractMatrix,
+    S::Vector{<:Integer};
+    obsdim::Union{Integer,Nothing}=defaultobs,
+) nystrom_sample(k, vec_of_vecs(X; obsdim=obsdim === nothing ? nothing : Int(obsdim)), S) false
 
 function nystrom_pinv!(Cs::Matrix{T}, tol::T=eps(T) * size(Cs, 1)) where {T<:Real}
     # Compute eigendecomposition of sampled component of K
@@ -97,13 +100,32 @@ function nystrom(k::Kernel, X::AbstractVector, r::Real)
     return nystrom(k, X, S)
 end
 
+"""
+    nystrom(k::Kernel, X::AbstractMatrix, S::AbstractVector{<:Integer}; obsdim)
+
+If `obsdim=1`, equivalent to `nystrom(k, RowVecs(X), S)`.
+If `obsdim=2`, equivalent to `nystrom(k, ColVecs(X), S)`.
+
+See also: [`ColVecs`](@ref), [`RowVecs`](@ref)
+"""
 function nystrom(
-    k::Kernel, X::AbstractMatrix, S::AbstractVector{<:Integer}; obsdim::Int=defaultobs
+    k::Kernel,
+    X::AbstractMatrix,
+    S::AbstractVector{<:Integer};
+    obsdim::Union{Int,Nothing}=defaultobs,
 )
     return nystrom(k, vec_of_vecs(X; obsdim=obsdim), S)
 end
 
-function nystrom(k::Kernel, X::AbstractMatrix, r::Real; obsdim::Int=defaultobs)
+"""
+    nystrom(k::Kernel, X::AbstractMatrix, r::Real; obsdim)
+
+If `obsdim=1`, equivalent to `nystrom(k, RowVecs(X), r)`.
+If `obsdim=2`, equivalent to `nystrom(k, ColVecs(X), r)`.
+
+See also: [`ColVecs`](@ref), [`RowVecs`](@ref)
+"""
+function nystrom(k::Kernel, X::AbstractMatrix, r::Real; obsdim::Union{Int,Nothing}=nothing)
     return nystrom(k, vec_of_vecs(X; obsdim=obsdim), r)
 end
 

src/matrix/kernelmatrix.jl

@@ -5,18 +5,19 @@
 In-place version of [`kernelmatrix`](@ref) where pre-allocated matrix `K` will be
 overwritten with the kernel matrix.
 
-    kernelmatrix!(K::AbstractMatrix, κ::Kernel, X::AbstractMatrix; obsdim::Integer=2)
+    kernelmatrix!(K::AbstractMatrix, κ::Kernel, X::AbstractMatrix; obsdim)
     kernelmatrix!(
         K::AbstractMatrix,
         κ::Kernel,
         X::AbstractMatrix,
         Y::AbstractMatrix;
-        obsdim::Integer=2,
+        obsdim,
     )
 
-Equivalent to `kernelmatrix!(K, κ, ColVecs(X))` and
-`kernelmatrix(K, κ, ColVecs(X), ColVecs(Y))` respectively.
-Set `obsdim=1` to get `RowVecs`.
+If `obsdim=1`, equivalent to `kernelmatrix!(K, κ, RowVecs(X))` and
+`kernelmatrix(K, κ, RowVecs(X), RowVecs(Y))`, respectively.
+If `obsdim=2`, equivalent to `kernelmatrix!(K, κ, ColVecs(X))` and
+`kernelmatrix(K, κ, ColVecs(X), ColVecs(Y))`, respectively.
 
 See also: [`ColVecs`](@ref), [`RowVecs`](@ref)
 """
@@ -35,12 +36,13 @@ Compute the kernel `κ` for each pair of inputs in `x` and `y`.
 Returns a matrix of size `(length(x), length(y))` satisfying
 `kernelmatrix(κ, x, y)[p, q] == κ(x[p], y[q])`.
 
-    kernelmatrix(κ::Kernel, X::AbstractMatrix; obsdim::Int=2)
-    kernelmatrix(κ::Kernel, X::AbstractMatrix, Y::AbstractMatrix; obsdim::Int=2)
+    kernelmatrix(κ::Kernel, X::AbstractMatrix; obsdim)
+    kernelmatrix(κ::Kernel, X::AbstractMatrix, Y::AbstractMatrix; obsdim)
 
-Equivalent to `kernelmatrix(κ, ColVecs(X))` and `kernelmatrix(κ, ColVecs(X), ColVecs(Y))`
-respectively.
-Set `obsdim=1` to get `RowVecs`.
+If `obsdim=1`, equivalent to `kernelmatrix(κ, RowVecs(X))` and
+`kernelmatrix(κ, RowVecs(X), RowVecs(Y))`, respectively.
+If `obsdim=2`, equivalent to `kernelmatrix(κ, ColVecs(X))` and
+`kernelmatrix(κ, ColVecs(X), ColVecs(Y))`, respectively.
 
 See also: [`ColVecs`](@ref), [`RowVecs`](@ref)
 """
@@ -52,18 +54,19 @@ kernelmatrix
 
 In place version of [`kernelmatrix_diag`](@ref).
 
-    kernelmatrix_diag!(K::AbstractVector, κ::Kernel, X::AbstractMatrix; obsdim::Int=2)
+    kernelmatrix_diag!(K::AbstractVector, κ::Kernel, X::AbstractMatrix; obsdim)
     kernelmatrix_diag!(
         K::AbstractVector,
         κ::Kernel,
         X::AbstractMatrix,
         Y::AbstractMatrix;
-        obsdim::Int=2,
+        obsdim
     )
 
-Equivalent to `kernelmatrix_diag!(K, κ, ColVecs(X))` and
-`kernelmatrix_diag!(K, κ, ColVecs(X), ColVecs(Y))` respectively.
-Set `obsdim=1` to get `RowVecs`.
+If `obsdim=1`, equivalent to `kernelmatrix_diag!(K, κ, RowVecs(X))` and
+`kernelmatrix_diag!(K, κ, RowVecs(X), RowVecs(Y))`, respectively.
+If `obsdim=2`, equivalent to `kernelmatrix_diag!(K, κ, ColVecs(X))` and
+`kernelmatrix_diag!(K, κ, ColVecs(X), ColVecs(Y))`, respectively.
 
 See also: [`ColVecs`](@ref), [`RowVecs`](@ref)
 """
@@ -79,11 +82,13 @@ Compute the diagonal of `kernelmatrix(κ, x)` efficiently.
 Compute the diagonal of `kernelmatrix(κ, x, y)` efficiently.
 Requires that `x` and `y` are the same length.
 
-    kernelmatrix_diag(κ::Kernel, X::AbstractMatrix; obsdim::Int=2)
-    kernelmatrix_diag(κ::Kernel, X::AbstractMatrix, Y::AbstractMatrix; obsdim::Int=2)
+    kernelmatrix_diag(κ::Kernel, X::AbstractMatrix; obsdim)
+    kernelmatrix_diag(κ::Kernel, X::AbstractMatrix, Y::AbstractMatrix; obsdim)
 
-Equivalent to `kernelmatrix_diag(κ, ColVecs(X))` and
-`kernelmatrix_diag(κ, ColVecs(X), ColVecs(Y))` respectively.
+If `obsdim=1`, equivalent to `kernelmatrix_diag(κ, RowVecs(X))` and
+`kernelmatrix_diag(κ, RowVecs(X), RowVecs(Y))`, respectively.
+If `obsdim=2`, equivalent to `kernelmatrix_diag(κ, ColVecs(X))` and
+`kernelmatrix_diag(κ, ColVecs(X), ColVecs(Y))`, respectively.
 
 See also: [`ColVecs`](@ref), [`RowVecs`](@ref)
 """
@@ -162,10 +167,10 @@ end
 # Wrapper methods for AbstractMatrix inputs to maintain obsdim interface.
 #
 
-const defaultobs = 2
+const defaultobs = nothing
 
 function kernelmatrix!(
-    K::AbstractMatrix, κ::Kernel, X::AbstractMatrix; obsdim::Int=defaultobs
+    K::AbstractMatrix, κ::Kernel, X::AbstractMatrix; obsdim::Union{Int,Nothing}=defaultobs
 )
     return kernelmatrix!(K, κ, vec_of_vecs(X; obsdim=obsdim))
 end
@@ -175,12 +180,12 @@ function kernelmatrix!(
     κ::Kernel,
     X::AbstractMatrix,
     Y::AbstractMatrix;
-    obsdim::Int=defaultobs,
+    obsdim::Union{Int,Nothing}=defaultobs,
 )
     return kernelmatrix!(K, κ, vec_of_vecs(X; obsdim=obsdim), vec_of_vecs(Y; obsdim=obsdim))
 end
 
-function kernelmatrix(κ::Kernel, X::AbstractMatrix; obsdim::Int=defaultobs)
+function kernelmatrix(κ::Kernel, X::AbstractMatrix; obsdim::Union{Int,Nothing}=defaultobs)
     return kernelmatrix(κ, vec_of_vecs(X; obsdim=obsdim))
 end
 
@@ -189,7 +194,7 @@ function kernelmatrix(κ::Kernel, X::AbstractMatrix, Y::AbstractMatrix; obsdim=d
 end
 
 function kernelmatrix_diag!(
-    K::AbstractVector, κ::Kernel, X::AbstractMatrix; obsdim::Int=defaultobs
+    K::AbstractVector, κ::Kernel, X::AbstractMatrix; obsdim::Union{Int,Nothing}=defaultobs
 )
     return kernelmatrix_diag!(K, κ, vec_of_vecs(X; obsdim=obsdim))
 end
@@ -199,19 +204,21 @@ function kernelmatrix_diag!(
     κ::Kernel,
     X::AbstractMatrix,
     Y::AbstractMatrix;
-    obsdim::Int=defaultobs,
+    obsdim::Union{Int,Nothing}=defaultobs,
 )
     return kernelmatrix_diag!(
         K, κ, vec_of_vecs(X; obsdim=obsdim), vec_of_vecs(Y; obsdim=obsdim)
     )
 end
 
-function kernelmatrix_diag(κ::Kernel, X::AbstractMatrix; obsdim::Int=defaultobs)
+function kernelmatrix_diag(
+    κ::Kernel, X::AbstractMatrix; obsdim::Union{Int,Nothing}=defaultobs
+)
     return kernelmatrix_diag(κ, vec_of_vecs(X; obsdim=obsdim))
 end
 
 function kernelmatrix_diag(
-    κ::Kernel, X::AbstractMatrix, Y::AbstractMatrix; obsdim::Int=defaultobs
+    κ::Kernel, X::AbstractMatrix, Y::AbstractMatrix; obsdim::Union{Int,Nothing}=defaultobs
 )
     return kernelmatrix_diag(
         κ, vec_of_vecs(X; obsdim=obsdim), vec_of_vecs(Y; obsdim=obsdim)

src/matrix/kernelpdmat.jl

@@ -3,18 +3,13 @@ using .PDMats: PDMat
 export kernelpdmat
 
 """
-    kernelpdmat(k::Kernel, X::AbstractMatrix; obsdim::Int=2)
     kernelpdmat(k::Kernel, X::AbstractVector)
 
 Compute a positive-definite matrix in the form of a `PDMat` matrix (see [PDMats.jl](https://github.com/JuliaStats/PDMats.jl)),
 with the Cholesky decomposition precomputed.
 The algorithm adds a diagonal "nugget" term to the kernel matrix which is increased until positive
 definiteness is achieved. The algorithm gives up with an error if the nugget becomes larger than 1% of the largest value in the kernel matrix.
 """
-function kernelpdmat(κ::Kernel, X::AbstractMatrix; obsdim::Int=defaultobs)
-    return kernelpdmat(κ, vec_of_vecs(X; obsdim=obsdim))
-end
-
 function kernelpdmat(κ::Kernel, X::AbstractVector)
     K = kernelmatrix(κ, X)
     Kmax = maximum(K)
@@ -31,3 +26,15 @@ function kernelpdmat(κ::Kernel, X::AbstractVector)
     end
     return PDMat(K + α * I)
 end
+
+"""
+    kernelpdmat(k::Kernel, X::AbstractMatrix; obsdim)
+
+If `obsdim=1`, equivalent to `kernelpdmat(k, RowVecs(X))`.
+If `obsdim=2`, equivalent to `kernelpdmat(k, ColVecs(X))`.
+
+See also: [`ColVecs`](@ref), [`RowVecs`](@ref)
+"""
+function kernelpdmat(κ::Kernel, X::AbstractMatrix; obsdim::Union{Int,Nothing}=defaultobs)
+    return kernelpdmat(κ, vec_of_vecs(X; obsdim=obsdim))
+end

src/utils.jl

@@ -21,12 +21,28 @@ macro check_args(K, param, cond, desc=string(cond))
     end
 end
 
-function vec_of_vecs(X::AbstractMatrix; obsdim::Int=2)
-    @assert obsdim ∈ (1, 2) "obsdim should be 1 or 2, see docs of kernelmatrix"
-    if obsdim == 1
-        RowVecs(X)
+function deprecated_obsdim(obsdim::Union{Int,Nothing})
+    _obsdim = if obsdim === nothing
+        Base.depwarn(
+            "implicit `obsdim=2` argument is deprecated and now has to be passed " *
+            "explicitly to specify that each column corresponds to one observation",
+            :vec_of_vecs,
+        )
+        2
+    else
+        obsdim
+    end
+    return _obsdim
+end
+
+function vec_of_vecs(X::AbstractMatrix; obsdim::Union{Int,Nothing}=nothing)
+    _obsdim = deprecated_obsdim(obsdim)
+    if _obsdim == 1
+        return RowVecs(X)
+    elseif _obsdim == 2
+        return ColVecs(X)
     else
-        ColVecs(X)
+        throw(ArgumentError("`obsdim` keyword argument should be 1 or 2"))
     end
 end
 

test/approximations/nystrom.jl

@@ -14,4 +14,8 @@
         @test kernelmatrix(k, X; obsdim=obsdim) ≈
             kernelmatrix(nystrom(k, X, collect(1:dims[obsdim]); obsdim=obsdim))
     end
+    @test kernelmatrix(@test_deprecated(nystrom(k, X, 1.0))) ≈
+        kernelmatrix(nystrom(k, ColVecs(X), 1.0))
+    @test kernelmatrix(@test_deprecated(nystrom(k, X, 1:5))) ≈
+        kernelmatrix(nystrom(k, ColVecs(X), 1:5))
 end

test/basekernels/fbm.jl

@@ -5,7 +5,7 @@
     v1 = rand(rng, 3)
     v2 = rand(rng, 3)
     @test k(v1, v2) ≈
-          (
+        (
         sqeuclidean(v1, zero(v1))^h + sqeuclidean(v2, zero(v2))^h -
         sqeuclidean(v1 - v2, zero(v1 - v2))^h
     ) / 2 atol = 1e-5

test/basekernels/matern.jl

@@ -45,7 +45,7 @@
         k = Matern52Kernel()
         @test kappa(k, x) ≈ (1 + sqrt(5) * x + 5 / 3 * x^2)exp(-sqrt(5) * x)
         @test k(v1, v2) ≈
-              (
+            (
             1 + sqrt(5) * norm(v1 - v2) + 5 / 3 * norm(v1 - v2)^2
         )exp(-sqrt(5) * norm(v1 - v2))
         @test kappa(Matern52Kernel(), x) == kappa(k, x)

test/basekernels/piecewisepolynomial.jl

@@ -20,7 +20,7 @@
     @test PiecewisePolynomialKernel(; dim=D) isa PiecewisePolynomialKernel{0}
 
     @test repr(k) ==
-          "Piecewise Polynomial Kernel (degree = $(degree), ⌊dim/2⌋ = $(div(D, 2)), metric = Euclidean(0.0))"
+        "Piecewise Polynomial Kernel (degree = $(degree), ⌊dim/2⌋ = $(div(D, 2)), metric = Euclidean(0.0))"
 
     k3 = PiecewisePolynomialKernel(;
         degree=degree, dim=D, metric=WeightedEuclidean(ones(D))

test/basekernels/wiener.jl

@@ -27,9 +27,9 @@
     @test k0(v1, v2) ≈ minXY
     @test k1(v1, v2) ≈ 1 / 3 * minXY^3 + 1 / 2 * minXY^2 * euclidean(v1, v2)
     @test k2(v1, v2) ≈
-          1 / 20 * minXY^5 + 1 / 12 * minXY^3 * euclidean(v1, v2) * (X + Y - 1 / 2 * minXY)
+        1 / 20 * minXY^5 + 1 / 12 * minXY^3 * euclidean(v1, v2) * (X + Y - 1 / 2 * minXY)
     @test k3(v1, v2) ≈
-          1 / 252 * minXY^7 +
+        1 / 252 * minXY^7 +
           1 / 720 *
           minXY^4 *
           euclidean(v1, v2) *

test/kernels/gibbskernel.jl

@@ -8,6 +8,6 @@
     k_gibbs = GibbsKernel(ell)
 
     @test k_gibbs(x, y) ≈
-          sqrt((2 * ell(x) * ell(y)) / (ell(x)^2 + ell(y)^2)) *
+        sqrt((2 * ell(x) * ell(y)) / (ell(x)^2 + ell(y)^2)) *
           exp(-(x - y)^2 / (ell(x)^2 + ell(y)^2))
 end

test/kernels/neuralkernelnetwork.jl

@@ -43,12 +43,12 @@ using KernelFunctions: NeuralKernelNetwork, LinearLayer, product, Primitive
             # Vector input.
             @test kernelmatrix_diag(nkn_add_kernel, x0) ≈ kernelmatrix_diag(sum_k, x0)
             @test kernelmatrix_diag(nkn_add_kernel, x0, x1) ≈
-                  kernelmatrix_diag(sum_k, x0, x1)
+                kernelmatrix_diag(sum_k, x0, x1)
 
             # ColVecs input.
             @test kernelmatrix_diag(nkn_add_kernel, X0) ≈ kernelmatrix_diag(sum_k, X0)
             @test kernelmatrix_diag(nkn_add_kernel, X0, X1) ≈
-                  kernelmatrix_diag(sum_k, X0, X1)
+                kernelmatrix_diag(sum_k, X0, X1)
         end
         @testset "product" begin
             nkn_prod_kernel = NeuralKernelNetwork(primitives, product)

test/kernels/transformedkernel.jl

@@ -19,8 +19,8 @@
     @test ktard(v1, v2) == (k ∘ ARDTransform(v))(v1, v2)
     @test ktard(v1, v2) == k(v .* v1, v .* v2)
     @test (k ∘ LinearTransform(P') ∘ ScaleTransform(s))(v1, v2) ==
-          ((k ∘ LinearTransform(P')) ∘ ScaleTransform(s))(v1, v2) ==
-          (k ∘ (LinearTransform(P') ∘ ScaleTransform(s)))(v1, v2)
+        ((k ∘ LinearTransform(P')) ∘ ScaleTransform(s))(v1, v2) ==
+        (k ∘ (LinearTransform(P') ∘ ScaleTransform(s)))(v1, v2)
 
     @test repr(kt) == repr(k) * "\n\t- " * repr(ScaleTransform(s))