From 08cb10bb49ba63f00a39100e061859a7bb967f73 Mon Sep 17 00:00:00 2001 From: Oliver Schulz Date: Wed, 13 Oct 2021 21:07:02 +0200 Subject: [PATCH 1/9] Add test_with_logabsdet_jacobian --- Project.toml | 4 +- src/ChangesOfVariables.jl | 2 + src/test.jl | 78 +++++++++++++++++++++++++++++++++++++ test/runtests.jl | 2 + test/test_test.jl | 36 +++++++++++++++++ test/test_with_ladj.jl | 82 ++++++++++++++++++--------------------- 6 files changed, 157 insertions(+), 47 deletions(-) create mode 100644 src/test.jl create mode 100644 test/test_test.jl diff --git a/Project.toml b/Project.toml index 8ca9df8..9e363ce 100644 --- a/Project.toml +++ b/Project.toml @@ -4,6 +4,7 @@ version = "0.1.0" [deps] LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" +Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" [compat] julia = "1" @@ -11,7 +12,6 @@ julia = "1" [extras] Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4" ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210" -Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" [targets] -test = ["Documenter", "ForwardDiff", "Test"] +test = ["Documenter", "ForwardDiff"] diff --git a/src/ChangesOfVariables.jl b/src/ChangesOfVariables.jl index 9315fb7..b2b034d 100644 --- a/src/ChangesOfVariables.jl +++ b/src/ChangesOfVariables.jl @@ -10,7 +10,9 @@ transformations). module ChangesOfVariables using LinearAlgebra +using Test include("with_ladj.jl") +include("test.jl") end # module diff --git a/src/test.jl b/src/test.jl new file mode 100644 index 0000000..d7fe106 --- /dev/null +++ b/src/test.jl @@ -0,0 +1,78 @@ +# This file is a part of ChangesOfVariables.jl, licensed under the MIT License (MIT). + + +_to_realvec_and_back(V::AbstractVector{<:Real}) = V, identity +_to_realvec_and_back(x::Real) = [x], V -> V[1] +_to_realvec_and_back(x::Complex) = [real(x), imag(x)], V -> Complex(V[1], V[2]) +_to_realvec_and_back(x::NTuple{N}) where N = [x...], V -> ntuple(i -> V[i], Val(N)) + +function _to_realvec_and_back(x::Ref) + xval = x[] + V, to_xval = _to_realvec_and_back(xval) + back_to_ref(V) = Ref(to_xval(V)) + return (V, back_to_ref) +end + +_to_realvec_and_back(A::AbstractArray{<:Real}) = vec(A), V -> reshape(V, size(A)) + +function _to_realvec_and_back(A::AbstractArray{Complex{T}, N}) where {T<:Real, N} + RA = cat(real.(A), imag.(A), dims = N+1) + V, to_array = _to_realvec_and_back(RA) + function back_to_complex(V) + RA = to_array(V) + Complex.(view(RA, map(_ -> :, size(A))..., 1), view(RA, map(_ -> :, size(A))..., 2)) + end + return (V, back_to_complex) +end + + +_to_realvec(x) = _to_realvec_and_back(x)[1] + + +function _auto_with_logabsdet_jacobian(f, x, getjacobian) + y = f(x) + V, to_x = _to_realvec_and_back(x) + vf(V) = _to_realvec(f(to_x(V))) + ladj = logabsdet(getjacobian(vf, V))[1] + return (y, ladj) +end + + +""" + ChangesOfVariables.test_with_logabsdet_jacobian( + f, x, getjacobian; + test_inferred::Bool = true, kwargs... + ) + +Test if [`with_logabsdet_jacobian(f, x)`](@ref) is implemented correctly. + +Checks if the result of `with_logabsdet_jacobian(f, x)` is approximately +equal to `(f(x), logabsdet(getjacobian(f, x)))` + +So the test uses `getjacobian(f, x)` to calculate a reference Jacobian for +`f` at `x`. Passing `ForwardDiff.jabobian`, `Zygote.jacobian` or similar as +the `getjacobian` function will do fine in most cases. + +If `x` or `f(x)` are real-valued scalars or complex-valued scalars or arrays, +the test will try to reshape them automatically, to account for limitations +of (e.g.) `ForwardDiff` and to ensure the result of `getjacobian` is a real +matrix. + +If `test_inferred == true` will test type inference on +`with_logabsdet_jacobian`. + +`kwargs...` are forwarded to `isapprox`. +""" +function test_with_logabsdet_jacobian(f, x, getjacobian; compare=isapprox, test_inferred::Bool = true, kwargs...) + @testset "test_with_logabsdet_jacobian: $f with input $x" begin + y, ladj = if test_inferred + @inferred with_logabsdet_jacobian(f, x) + else + with_logabsdet_jacobian(f, x) + end + ref_y, ref_ladj = _auto_with_logabsdet_jacobian(f, x, getjacobian) + @test compare(y, ref_y; kwargs...) + @test compare(ladj, ref_ladj; kwargs...) + end + return nothing +end diff --git a/test/runtests.jl b/test/runtests.jl index 3c63814..035b048 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -5,6 +5,7 @@ import ChangesOfVariables import Documenter Test.@testset "Package ChangesOfVariables" begin + include("test_test.jl") include("test_with_ladj.jl") # doctests @@ -16,3 +17,4 @@ Test.@testset "Package ChangesOfVariables" begin ) Documenter.doctest(ChangesOfVariables) end # testset + diff --git a/test/test_test.jl b/test/test_test.jl new file mode 100644 index 0000000..dee7698 --- /dev/null +++ b/test/test_test.jl @@ -0,0 +1,36 @@ +# This file is a part of ChangesOfVariables.jl, licensed under the MIT License (MIT). + +using ChangesOfVariables +using Test + +using LinearAlgebra +import ForwardDiff + +using ChangesOfVariables: test_with_logabsdet_jacobian + + +@testset "test_with_logabsdet_jacobian" begin + @testset "ChangesOfVariables._to_realvec_and_back" begin + for x in (rand(3), 0.5, Complex(0.2,0.7), (3,5,9), Ref(42), rand(3, 4, 5), Complex.(rand(3,5), rand(3,5))) + V, to_x = ChangesOfVariables._to_realvec_and_back(x) + @test V isa AbstractVector{<:Real} + @test V == ChangesOfVariables._to_realvec(x) + @test x isa Ref ? to_x(V)[] == x[] : to_x(V) == x + end + end + + x = Complex(0.2, -0.7) + y, ladj_y = ChangesOfVariables._auto_with_logabsdet_jacobian(inv, x, ForwardDiff.jacobian) + @test y == inv(x) + @test ladj_y ≈ -4 * log(abs(x)) + + X = Complex.(randn(3,3), randn(3,3)) + Y, ladj_Y = ChangesOfVariables._auto_with_logabsdet_jacobian(inv, X, ForwardDiff.jacobian) + @test Y == inv(X) + @test ladj_Y == -4 * 3 * logabsdet(X)[1] + + myisapprox(a, b; kwargs...) = isapprox(a, b; kwargs...) + test_with_logabsdet_jacobian(inv, X, ForwardDiff.jacobian, test_inferred = true, atol = 10^-6) + test_with_logabsdet_jacobian(inv, X, ForwardDiff.jacobian, test_inferred = false, atol = 10^-6) + test_with_logabsdet_jacobian(inv, X, ForwardDiff.jacobian, compare = myisapprox, atol = 10^-6) +end diff --git a/test/test_with_ladj.jl b/test/test_with_ladj.jl index 4c210ea..115d57f 100644 --- a/test/test_with_ladj.jl +++ b/test/test_with_ladj.jl @@ -4,66 +4,58 @@ using ChangesOfVariables using Test using LinearAlgebra -using ForwardDiff: derivative, jacobian +import ForwardDiff +using ChangesOfVariables: test_with_logabsdet_jacobian -fwddiff_ladj(f, x::Real) = log(abs(derivative(f, x))) -fwddiff_ladj(f, x::AbstractArray{<:Real}) = logabsdet(jacobian(f, x))[1] -fwddiff_with_ladj(f, x) = (f(x), fwddiff_ladj(f, x)) -ascomplex(A::AbstractArray{T}) where T = reinterpret(Complex{T}, A) -asreal(A::AbstractArray{Complex{T}}) where T = reinterpret(T, A) - -isaprx(a, b) = isapprox(a,b) -isaprx(a::NTuple{N,Any}, b::NTuple{N,Any}) where N = all(map(isaprx, a, b)) - - -foo(x) = inv(exp(-x) + 1) +@testset "with_logabsdet_jacobian" begin + foo(x) = inv(exp(-x) + 1) -function ChangesOfVariables.with_logabsdet_jacobian(::typeof(foo), x) - y = foo(x) - ladj = -x + 2 * log(y) - (y, ladj) -end + function ChangesOfVariables.with_logabsdet_jacobian(::typeof(foo), x) + y = foo(x) + ladj = -x + 2 * log(y) + (y, ladj) + end -@testset "with_logabsdet_jacobian" begin x = 4.2 X = rand(10) A = rand(5, 5) - CA = rand(10, 5) + CA = Complex.(rand(5, 5), rand(5, 5)) - @test isaprx(@inferred(with_logabsdet_jacobian(foo, x)), fwddiff_with_ladj(foo, x)) + + getjacobian = ForwardDiff.jacobian + + isaprx(a, b; kwargs...) = isapprox(a,b; kwargs...) + isaprx(a::NTuple{N,Any}, b::NTuple{N,Any}; kwargs...) where N = all(map((a,b) -> isaprx(a, b; kwargs...), a, b)) + + + test_with_logabsdet_jacobian(foo, x, getjacobian) @static if VERSION >= v"1.6" - log_foo = log ∘ foo - @test isaprx(@inferred(with_logabsdet_jacobian(log_foo, x)), fwddiff_with_ladj(log_foo, x)) + test_with_logabsdet_jacobian(log ∘ foo, x, getjacobian) + end + + @testset "getjacobian on mapped and broadcasted" begin + for f in (Base.Fix1(map, foo), Base.Fix1(broadcast, foo)) + for arg in (x, fill(x,), Ref(x), (x,), X) + test_with_logabsdet_jacobian(f, arg, getjacobian, compare = isaprx) + end + end end - mapped_foo = Base.Fix1(map, foo) - @test isaprx(@inferred(with_logabsdet_jacobian(mapped_foo, x)), fwddiff_with_ladj(mapped_foo, x)) - @test isaprx(@inferred(with_logabsdet_jacobian(mapped_foo, fill(x))), fwddiff_with_ladj(mapped_foo, fill(x))) - @test isaprx(@inferred(with_logabsdet_jacobian(mapped_foo, Ref(x))), fwddiff_with_ladj(mapped_foo, fill(x))) - @test isaprx(@inferred(with_logabsdet_jacobian(mapped_foo, (x,))), (mapped_foo((x,)), fwddiff_ladj(mapped_foo, x))) - @test isaprx(@inferred(with_logabsdet_jacobian(mapped_foo, X)), fwddiff_with_ladj(mapped_foo, X)) - - broadcasted_foo = Base.Fix1(broadcast, foo) - @test isaprx(@inferred(with_logabsdet_jacobian(broadcasted_foo, x)), fwddiff_with_ladj(broadcasted_foo, x)) - @test isaprx(@inferred(with_logabsdet_jacobian(broadcasted_foo, fill(x))), fwddiff_with_ladj(broadcasted_foo, x)) - @test isaprx(@inferred(with_logabsdet_jacobian(broadcasted_foo, Ref(x))), fwddiff_with_ladj(broadcasted_foo, x)) - @test isaprx(@inferred(with_logabsdet_jacobian(broadcasted_foo, (x,))), (mapped_foo((x,)), fwddiff_ladj(mapped_foo, x))) - @test isaprx(@inferred(with_logabsdet_jacobian(broadcasted_foo, X)), fwddiff_with_ladj(broadcasted_foo, X)) - - for f in (identity, adjoint, transpose) - @test isaprx(@inferred(with_logabsdet_jacobian(f, x)), fwddiff_with_ladj(f, x)) - @test isaprx(@inferred(with_logabsdet_jacobian(f, A)), fwddiff_with_ladj(f, A)) + @testset "getjacobian on identity, adjoint and transpose" begin + for f in (identity, adjoint, transpose) + for arg in (x, A) + test_with_logabsdet_jacobian(f, arg, getjacobian) + end + end end - - @test isaprx(@inferred(with_logabsdet_jacobian(inv, x)), fwddiff_with_ladj(inv, x)) - @test isaprx(@inferred(with_logabsdet_jacobian(inv, A)), fwddiff_with_ladj(inv, A)) - @test isaprx(@inferred(with_logabsdet_jacobian(inv, ascomplex(CA))), (inv(ascomplex(CA)), fwddiff_ladj(CA -> asreal(inv(ascomplex(CA))), CA))) - for f in (exp, log, exp2, log2, exp10, log10, expm1, log1p) - @test isaprx(@inferred(with_logabsdet_jacobian(f, x)), fwddiff_with_ladj(f, x)) + @testset "getjacobian on inv" begin + for arg in (x, A, CA) + test_with_logabsdet_jacobian(inv, arg, getjacobian) + end end end From 73aa3c331315f122a933fa79b30552cd4dd1c8de Mon Sep 17 00:00:00 2001 From: Oliver Schulz Date: Thu, 14 Oct 2021 00:54:11 +0200 Subject: [PATCH 2/9] Add rv_and_back argument to test_with_logabsdet_jacobian --- src/test.jl | 70 ++++++++++++++++-------------------------- test/rv_and_back.jl | 26 ++++++++++++++++ test/test_test.jl | 33 +++++++++++--------- test/test_with_ladj.jl | 12 +++++--- 4 files changed, 79 insertions(+), 62 deletions(-) create mode 100644 test/rv_and_back.jl diff --git a/src/test.jl b/src/test.jl index d7fe106..2ecf998 100644 --- a/src/test.jl +++ b/src/test.jl @@ -1,47 +1,22 @@ # This file is a part of ChangesOfVariables.jl, licensed under the MIT License (MIT). -_to_realvec_and_back(V::AbstractVector{<:Real}) = V, identity -_to_realvec_and_back(x::Real) = [x], V -> V[1] -_to_realvec_and_back(x::Complex) = [real(x), imag(x)], V -> Complex(V[1], V[2]) -_to_realvec_and_back(x::NTuple{N}) where N = [x...], V -> ntuple(i -> V[i], Val(N)) +_generalized_logabsdet(A) = logabsdet(A) +_generalized_logabsdet(x::Real) = log(abs(x)) -function _to_realvec_and_back(x::Ref) - xval = x[] - V, to_xval = _to_realvec_and_back(xval) - back_to_ref(V) = Ref(to_xval(V)) - return (V, back_to_ref) -end - -_to_realvec_and_back(A::AbstractArray{<:Real}) = vec(A), V -> reshape(V, size(A)) - -function _to_realvec_and_back(A::AbstractArray{Complex{T}, N}) where {T<:Real, N} - RA = cat(real.(A), imag.(A), dims = N+1) - V, to_array = _to_realvec_and_back(RA) - function back_to_complex(V) - RA = to_array(V) - Complex.(view(RA, map(_ -> :, size(A))..., 1), view(RA, map(_ -> :, size(A))..., 2)) - end - return (V, back_to_complex) -end - - -_to_realvec(x) = _to_realvec_and_back(x)[1] - - -function _auto_with_logabsdet_jacobian(f, x, getjacobian) +function _auto_with_logabsdet_jacobian(f, x, getjacobian, rv_and_back) y = f(x) - V, to_x = _to_realvec_and_back(x) - vf(V) = _to_realvec(f(to_x(V))) - ladj = logabsdet(getjacobian(vf, V))[1] + V, to_x = rv_and_back(x) + vf(V) = rv_and_back(f(to_x(V)))[1] + ladj = _generalized_logabsdet(getjacobian(vf, V))[1] return (y, ladj) end """ ChangesOfVariables.test_with_logabsdet_jacobian( - f, x, getjacobian; - test_inferred::Bool = true, kwargs... + f, x, getjacobian, rv_and_back = x -> (x, identity); + compare = isapprox, test_inferred::Bool = true, kwargs... ) Test if [`with_logabsdet_jacobian(f, x)`](@ref) is implemented correctly. @@ -51,26 +26,35 @@ equal to `(f(x), logabsdet(getjacobian(f, x)))` So the test uses `getjacobian(f, x)` to calculate a reference Jacobian for `f` at `x`. Passing `ForwardDiff.jabobian`, `Zygote.jacobian` or similar as -the `getjacobian` function will do fine in most cases. +the `getjacobian` function will do fine in most cases. If input and output +of `f` are real scalar values, use `ForwardDiff.derivative`. + +If `getjacobian(f, x)` can't handle the type of `x` of `f(x)` because they +are not real-valued vectors, use the `rv_and_back` argument to pass a +function with the following behavior -If `x` or `f(x)` are real-valued scalars or complex-valued scalars or arrays, -the test will try to reshape them automatically, to account for limitations -of (e.g.) `ForwardDiff` and to ensure the result of `getjacobian` is a real -matrix. +```julia +v, back = rv_and_back(x) +v isa AbstractVector{<:Real} +back(v) == x +``` -If `test_inferred == true` will test type inference on -`with_logabsdet_jacobian`. +If `test_inferred == true`, type inference on `with_logabsdet_jacobian` will +be tested. -`kwargs...` are forwarded to `isapprox`. +`kwargs...` are forwarded to `compare`. """ -function test_with_logabsdet_jacobian(f, x, getjacobian; compare=isapprox, test_inferred::Bool = true, kwargs...) +function test_with_logabsdet_jacobian( + f, x, getjacobian, rv_and_back = x -> (x, identity); + compare = isapprox, test_inferred::Bool = true, kwargs... +) @testset "test_with_logabsdet_jacobian: $f with input $x" begin y, ladj = if test_inferred @inferred with_logabsdet_jacobian(f, x) else with_logabsdet_jacobian(f, x) end - ref_y, ref_ladj = _auto_with_logabsdet_jacobian(f, x, getjacobian) + ref_y, ref_ladj = _auto_with_logabsdet_jacobian(f, x, getjacobian, rv_and_back) @test compare(y, ref_y; kwargs...) @test compare(ladj, ref_ladj; kwargs...) end diff --git a/test/rv_and_back.jl b/test/rv_and_back.jl new file mode 100644 index 0000000..31b33f9 --- /dev/null +++ b/test/rv_and_back.jl @@ -0,0 +1,26 @@ +# This file is a part of ChangesOfVariables.jl, licensed under the MIT License (MIT). + + +rv_and_back(V::AbstractVector{<:Real}) = V, identity +rv_and_back(x::Real) = [x], V -> V[1] +rv_and_back(x::Complex) = [real(x), imag(x)], V -> Complex(V[1], V[2]) +rv_and_back(x::NTuple{N}) where N = [x...], V -> ntuple(i -> V[i], Val(N)) + +function rv_and_back(x::Ref) + xval = x[] + V, to_xval = rv_and_back(xval) + back_to_ref(V) = Ref(to_xval(V)) + return (V, back_to_ref) +end + +rv_and_back(A::AbstractArray{<:Real}) = vec(A), V -> reshape(V, size(A)) + +function rv_and_back(A::AbstractArray{Complex{T}, N}) where {T<:Real, N} + RA = cat(real.(A), imag.(A), dims = N+1) + V, to_array = rv_and_back(RA) + function back_to_complex(V) + RA = to_array(V) + Complex.(view(RA, map(_ -> :, size(A))..., 1), view(RA, map(_ -> :, size(A))..., 2)) + end + return (V, back_to_complex) +end diff --git a/test/test_test.jl b/test/test_test.jl index dee7698..0450261 100644 --- a/test/test_test.jl +++ b/test/test_test.jl @@ -9,28 +9,33 @@ import ForwardDiff using ChangesOfVariables: test_with_logabsdet_jacobian -@testset "test_with_logabsdet_jacobian" begin - @testset "ChangesOfVariables._to_realvec_and_back" begin - for x in (rand(3), 0.5, Complex(0.2,0.7), (3,5,9), Ref(42), rand(3, 4, 5), Complex.(rand(3,5), rand(3,5))) - V, to_x = ChangesOfVariables._to_realvec_and_back(x) - @test V isa AbstractVector{<:Real} - @test V == ChangesOfVariables._to_realvec(x) - @test x isa Ref ? to_x(V)[] == x[] : to_x(V) == x - end +include("rv_and_back.jl") + +@testset "rv_and_back" begin + for x in (rand(3), 0.5, Complex(0.2,0.7), (3,5,9), Ref(42), rand(3, 4, 5), Complex.(rand(3,5), rand(3,5))) + V, to_x = rv_and_back(x) + @test V isa AbstractVector{<:Real} + @test V == rv_and_back(x)[1] + @test x isa Ref ? to_x(V)[] == x[] : to_x(V) == x end +end + +@testset "test_with_logabsdet_jacobian" begin x = Complex(0.2, -0.7) - y, ladj_y = ChangesOfVariables._auto_with_logabsdet_jacobian(inv, x, ForwardDiff.jacobian) + y, ladj_y = ChangesOfVariables._auto_with_logabsdet_jacobian(inv, x, ForwardDiff.jacobian, rv_and_back) @test y == inv(x) @test ladj_y ≈ -4 * log(abs(x)) X = Complex.(randn(3,3), randn(3,3)) - Y, ladj_Y = ChangesOfVariables._auto_with_logabsdet_jacobian(inv, X, ForwardDiff.jacobian) + Y, ladj_Y = ChangesOfVariables._auto_with_logabsdet_jacobian(inv, X, ForwardDiff.jacobian, rv_and_back) @test Y == inv(X) - @test ladj_Y == -4 * 3 * logabsdet(X)[1] + @test ladj_Y ≈ -4 * 3 * logabsdet(X)[1] myisapprox(a, b; kwargs...) = isapprox(a, b; kwargs...) - test_with_logabsdet_jacobian(inv, X, ForwardDiff.jacobian, test_inferred = true, atol = 10^-6) - test_with_logabsdet_jacobian(inv, X, ForwardDiff.jacobian, test_inferred = false, atol = 10^-6) - test_with_logabsdet_jacobian(inv, X, ForwardDiff.jacobian, compare = myisapprox, atol = 10^-6) + test_with_logabsdet_jacobian(inv, 0.5, ForwardDiff.derivative, test_inferred = true, atol = 10^-6) + test_with_logabsdet_jacobian(inv, rand(2,2), ForwardDiff.jacobian, test_inferred = true, atol = 10^-6) + test_with_logabsdet_jacobian(inv, X, ForwardDiff.jacobian, rv_and_back, test_inferred = true, atol = 10^-6) + test_with_logabsdet_jacobian(inv, X, ForwardDiff.jacobian, rv_and_back, test_inferred = false, atol = 10^-6) + test_with_logabsdet_jacobian(inv, X, ForwardDiff.jacobian, rv_and_back, compare = myisapprox, atol = 10^-6) end diff --git a/test/test_with_ladj.jl b/test/test_with_ladj.jl index 115d57f..50c6087 100644 --- a/test/test_with_ladj.jl +++ b/test/test_with_ladj.jl @@ -8,6 +8,8 @@ import ForwardDiff using ChangesOfVariables: test_with_logabsdet_jacobian +include("rv_and_back.jl") + @testset "with_logabsdet_jacobian" begin foo(x) = inv(exp(-x) + 1) @@ -31,16 +33,16 @@ using ChangesOfVariables: test_with_logabsdet_jacobian isaprx(a::NTuple{N,Any}, b::NTuple{N,Any}; kwargs...) where N = all(map((a,b) -> isaprx(a, b; kwargs...), a, b)) - test_with_logabsdet_jacobian(foo, x, getjacobian) + test_with_logabsdet_jacobian(foo, x, getjacobian, rv_and_back) @static if VERSION >= v"1.6" - test_with_logabsdet_jacobian(log ∘ foo, x, getjacobian) + test_with_logabsdet_jacobian(log ∘ foo, x, getjacobian, rv_and_back) end @testset "getjacobian on mapped and broadcasted" begin for f in (Base.Fix1(map, foo), Base.Fix1(broadcast, foo)) for arg in (x, fill(x,), Ref(x), (x,), X) - test_with_logabsdet_jacobian(f, arg, getjacobian, compare = isaprx) + test_with_logabsdet_jacobian(f, arg, getjacobian, rv_and_back, compare = isaprx) end end end @@ -48,14 +50,14 @@ using ChangesOfVariables: test_with_logabsdet_jacobian @testset "getjacobian on identity, adjoint and transpose" begin for f in (identity, adjoint, transpose) for arg in (x, A) - test_with_logabsdet_jacobian(f, arg, getjacobian) + test_with_logabsdet_jacobian(f, arg, getjacobian, rv_and_back) end end end @testset "getjacobian on inv" begin for arg in (x, A, CA) - test_with_logabsdet_jacobian(inv, arg, getjacobian) + test_with_logabsdet_jacobian(inv, arg, getjacobian, rv_and_back) end end end From 3ec06fa4b5ec945b66dc146a8f421b2f4d279d1c Mon Sep 17 00:00:00 2001 From: Oliver Schulz Date: Thu, 14 Oct 2021 01:21:12 +0200 Subject: [PATCH 3/9] Automatic type inference test in test_with_logabsdet_jacobian --- src/test.jl | 34 ++++++++++++++++++---------------- test/test_test.jl | 30 ++++++++++++++++-------------- 2 files changed, 34 insertions(+), 30 deletions(-) diff --git a/src/test.jl b/src/test.jl index 2ecf998..911265f 100644 --- a/src/test.jl +++ b/src/test.jl @@ -1,22 +1,10 @@ # This file is a part of ChangesOfVariables.jl, licensed under the MIT License (MIT). -_generalized_logabsdet(A) = logabsdet(A) -_generalized_logabsdet(x::Real) = log(abs(x)) - -function _auto_with_logabsdet_jacobian(f, x, getjacobian, rv_and_back) - y = f(x) - V, to_x = rv_and_back(x) - vf(V) = rv_and_back(f(to_x(V)))[1] - ladj = _generalized_logabsdet(getjacobian(vf, V))[1] - return (y, ladj) -end - - """ ChangesOfVariables.test_with_logabsdet_jacobian( f, x, getjacobian, rv_and_back = x -> (x, identity); - compare = isapprox, test_inferred::Bool = true, kwargs... + compare = isapprox, kwargs... ) Test if [`with_logabsdet_jacobian(f, x)`](@ref) is implemented correctly. @@ -46,17 +34,31 @@ be tested. """ function test_with_logabsdet_jacobian( f, x, getjacobian, rv_and_back = x -> (x, identity); - compare = isapprox, test_inferred::Bool = true, kwargs... + compare = isapprox, kwargs... ) @testset "test_with_logabsdet_jacobian: $f with input $x" begin - y, ladj = if test_inferred + ref_y, test_type_inference = try + @inferred(f(x)), true + catch err + f(x), false + end + + y, ladj = if test_type_inference @inferred with_logabsdet_jacobian(f, x) else with_logabsdet_jacobian(f, x) end - ref_y, ref_ladj = _auto_with_logabsdet_jacobian(f, x, getjacobian, rv_and_back) + + V, to_x = rv_and_back(x) + vf(V) = rv_and_back(f(to_x(V)))[1] + ref_ladj = _generalized_logabsdet(getjacobian(vf, V))[1] + @test compare(y, ref_y; kwargs...) @test compare(ladj, ref_ladj; kwargs...) end return nothing end + + +_generalized_logabsdet(A) = logabsdet(A) +_generalized_logabsdet(x::Real) = log(abs(x)) diff --git a/test/test_test.jl b/test/test_test.jl index 0450261..32b90d7 100644 --- a/test/test_test.jl +++ b/test/test_test.jl @@ -22,20 +22,22 @@ end @testset "test_with_logabsdet_jacobian" begin - x = Complex(0.2, -0.7) - y, ladj_y = ChangesOfVariables._auto_with_logabsdet_jacobian(inv, x, ForwardDiff.jacobian, rv_and_back) - @test y == inv(x) - @test ladj_y ≈ -4 * log(abs(x)) - - X = Complex.(randn(3,3), randn(3,3)) - Y, ladj_Y = ChangesOfVariables._auto_with_logabsdet_jacobian(inv, X, ForwardDiff.jacobian, rv_and_back) - @test Y == inv(X) - @test ladj_Y ≈ -4 * 3 * logabsdet(X)[1] + rx = 0.5 + cx = Complex(0.2, -0.7) + X = rand(3, 3) + CX = Complex.(randn(3,3), randn(3,3)) myisapprox(a, b; kwargs...) = isapprox(a, b; kwargs...) - test_with_logabsdet_jacobian(inv, 0.5, ForwardDiff.derivative, test_inferred = true, atol = 10^-6) - test_with_logabsdet_jacobian(inv, rand(2,2), ForwardDiff.jacobian, test_inferred = true, atol = 10^-6) - test_with_logabsdet_jacobian(inv, X, ForwardDiff.jacobian, rv_and_back, test_inferred = true, atol = 10^-6) - test_with_logabsdet_jacobian(inv, X, ForwardDiff.jacobian, rv_and_back, test_inferred = false, atol = 10^-6) - test_with_logabsdet_jacobian(inv, X, ForwardDiff.jacobian, rv_and_back, compare = myisapprox, atol = 10^-6) + + noninferrable_inv(x) = x!=rand(size(x)...) ? inv(x) : "" + ChangesOfVariables.with_logabsdet_jacobian(::typeof(noninferrable_inv), x) = noninferrable_inv(x), with_logabsdet_jacobian(inv, x)[2] + @test_throws ErrorException @inferred with_logabsdet_jacobian(noninferrable_inv, rand(2, 2)) + + test_with_logabsdet_jacobian(inv, rx, ForwardDiff.derivative, atol = 10^-6) + test_with_logabsdet_jacobian(inv, cx, ForwardDiff.jacobian, rv_and_back, atol = 10^-6) + test_with_logabsdet_jacobian(inv, X, ForwardDiff.jacobian, atol = 10^-6) + test_with_logabsdet_jacobian(inv, CX, ForwardDiff.jacobian, rv_and_back, atol = 10^-6) + test_with_logabsdet_jacobian(inv, CX, ForwardDiff.jacobian, rv_and_back, atol = 10^-6) + test_with_logabsdet_jacobian(inv, CX, ForwardDiff.jacobian, rv_and_back, compare = myisapprox, atol = 10^-6) + test_with_logabsdet_jacobian(noninferrable_inv, CX, ForwardDiff.jacobian, rv_and_back, atol = 10^-6) end From e81bf358b81b164d590510a8a65c01468ee3a7e2 Mon Sep 17 00:00:00 2001 From: Oliver Schulz Date: Thu, 14 Oct 2021 01:23:52 +0200 Subject: [PATCH 4/9] Add test_with_logabsdet_jacobian to docs --- docs/src/api.md | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/docs/src/api.md b/docs/src/api.md index 84413a6..46df42d 100644 --- a/docs/src/api.md +++ b/docs/src/api.md @@ -5,3 +5,9 @@ ```@docs with_logabsdet_jacobian ``` + +## Test utility + +```@docs +ChangesOfVariables.test_with_logabsdet_jacobian +``` From 98eee466cc385fb8fa30c887c168d599cbcf75ef Mon Sep 17 00:00:00 2001 From: Oliver Schulz Date: Thu, 14 Oct 2021 01:48:52 +0200 Subject: [PATCH 5/9] Fix typo in test_with_logabsdet_jacobian docstring --- src/test.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/test.jl b/src/test.jl index 911265f..d759508 100644 --- a/src/test.jl +++ b/src/test.jl @@ -17,7 +17,7 @@ So the test uses `getjacobian(f, x)` to calculate a reference Jacobian for the `getjacobian` function will do fine in most cases. If input and output of `f` are real scalar values, use `ForwardDiff.derivative`. -If `getjacobian(f, x)` can't handle the type of `x` of `f(x)` because they +If `getjacobian(f, x)` can't handle the type of `x` or `f(x)` because they are not real-valued vectors, use the `rv_and_back` argument to pass a function with the following behavior From 3abd5a289c9cb6fd58623414bc743e55eae57c3e Mon Sep 17 00:00:00 2001 From: Oliver Schulz Date: Fri, 15 Oct 2021 03:01:48 +0200 Subject: [PATCH 6/9] Simplify test_with_logabsdet_jacobian --- src/test.jl | 30 +++++++----------------------- test/getjacobian.jl | 36 ++++++++++++++++++++++++++++++++++++ test/rv_and_back.jl | 26 -------------------------- test/test_test.jl | 19 +++++++++---------- test/test_with_ladj.jl | 19 ++++++------------- 5 files changed, 58 insertions(+), 72 deletions(-) create mode 100644 test/getjacobian.jl delete mode 100644 test/rv_and_back.jl diff --git a/src/test.jl b/src/test.jl index d759508..e68183f 100644 --- a/src/test.jl +++ b/src/test.jl @@ -2,10 +2,7 @@ """ - ChangesOfVariables.test_with_logabsdet_jacobian( - f, x, getjacobian, rv_and_back = x -> (x, identity); - compare = isapprox, kwargs... - ) + ChangesOfVariables.test_with_logabsdet_jacobian(f, x, getjacobian; compare = isapprox, kwargs...) Test if [`with_logabsdet_jacobian(f, x)`](@ref) is implemented correctly. @@ -17,25 +14,14 @@ So the test uses `getjacobian(f, x)` to calculate a reference Jacobian for the `getjacobian` function will do fine in most cases. If input and output of `f` are real scalar values, use `ForwardDiff.derivative`. -If `getjacobian(f, x)` can't handle the type of `x` or `f(x)` because they -are not real-valued vectors, use the `rv_and_back` argument to pass a -function with the following behavior - -```julia -v, back = rv_and_back(x) -v isa AbstractVector{<:Real} -back(v) == x -``` - -If `test_inferred == true`, type inference on `with_logabsdet_jacobian` will -be tested. +Note that the result of `getjacobian(f, x)` must be a real-valued matrix +or a real scalar, so you may need to use a custom `getjacobian` function +that transforms the shape of `x` and `f(x)` internally, in conjunction +with automatic differentiation. `kwargs...` are forwarded to `compare`. """ -function test_with_logabsdet_jacobian( - f, x, getjacobian, rv_and_back = x -> (x, identity); - compare = isapprox, kwargs... -) +function test_with_logabsdet_jacobian(f, x, getjacobian; compare = isapprox, kwargs...) @testset "test_with_logabsdet_jacobian: $f with input $x" begin ref_y, test_type_inference = try @inferred(f(x)), true @@ -49,9 +35,7 @@ function test_with_logabsdet_jacobian( with_logabsdet_jacobian(f, x) end - V, to_x = rv_and_back(x) - vf(V) = rv_and_back(f(to_x(V)))[1] - ref_ladj = _generalized_logabsdet(getjacobian(vf, V))[1] + ref_ladj = _generalized_logabsdet(getjacobian(f, x))[1] @test compare(y, ref_y; kwargs...) @test compare(ladj, ref_ladj; kwargs...) diff --git a/test/getjacobian.jl b/test/getjacobian.jl new file mode 100644 index 0000000..4766d8a --- /dev/null +++ b/test/getjacobian.jl @@ -0,0 +1,36 @@ +# This file is a part of ChangesOfVariables.jl, licensed under the MIT License (MIT). + +import ForwardDiff + +torv_and_back(V::AbstractVector{<:Real}) = V, identity +torv_and_back(x::Real) = [x], V -> V[1] +torv_and_back(x::Complex) = [real(x), imag(x)], V -> Complex(V[1], V[2]) +torv_and_back(x::NTuple{N}) where N = [x...], V -> ntuple(i -> V[i], Val(N)) + +function torv_and_back(x::Ref) + xval = x[] + V, to_xval = torv_and_back(xval) + back_to_ref(V) = Ref(to_xval(V)) + return (V, back_to_ref) +end + +torv_and_back(A::AbstractArray{<:Real}) = vec(A), V -> reshape(V, size(A)) + +function torv_and_back(A::AbstractArray{Complex{T}, N}) where {T<:Real, N} + RA = cat(real.(A), imag.(A), dims = N+1) + V, to_array = torv_and_back(RA) + function back_to_complex(V) + RA = to_array(V) + Complex.(view(RA, map(_ -> :, size(A))..., 1), view(RA, map(_ -> :, size(A))..., 2)) + end + return (V, back_to_complex) +end + + +function getjacobian(f, x) + V, to_x = torv_and_back(x) + vf(V) = torv_and_back(f(to_x(V)))[1] + ForwardDiff.jacobian(vf, V) +end + +foo(x) = inv(exp(-x) + 1) diff --git a/test/rv_and_back.jl b/test/rv_and_back.jl deleted file mode 100644 index 31b33f9..0000000 --- a/test/rv_and_back.jl +++ /dev/null @@ -1,26 +0,0 @@ -# This file is a part of ChangesOfVariables.jl, licensed under the MIT License (MIT). - - -rv_and_back(V::AbstractVector{<:Real}) = V, identity -rv_and_back(x::Real) = [x], V -> V[1] -rv_and_back(x::Complex) = [real(x), imag(x)], V -> Complex(V[1], V[2]) -rv_and_back(x::NTuple{N}) where N = [x...], V -> ntuple(i -> V[i], Val(N)) - -function rv_and_back(x::Ref) - xval = x[] - V, to_xval = rv_and_back(xval) - back_to_ref(V) = Ref(to_xval(V)) - return (V, back_to_ref) -end - -rv_and_back(A::AbstractArray{<:Real}) = vec(A), V -> reshape(V, size(A)) - -function rv_and_back(A::AbstractArray{Complex{T}, N}) where {T<:Real, N} - RA = cat(real.(A), imag.(A), dims = N+1) - V, to_array = rv_and_back(RA) - function back_to_complex(V) - RA = to_array(V) - Complex.(view(RA, map(_ -> :, size(A))..., 1), view(RA, map(_ -> :, size(A))..., 2)) - end - return (V, back_to_complex) -end diff --git a/test/test_test.jl b/test/test_test.jl index 32b90d7..60dc722 100644 --- a/test/test_test.jl +++ b/test/test_test.jl @@ -4,18 +4,17 @@ using ChangesOfVariables using Test using LinearAlgebra -import ForwardDiff using ChangesOfVariables: test_with_logabsdet_jacobian +include("getjacobian.jl") -include("rv_and_back.jl") -@testset "rv_and_back" begin +@testset "torv_and_back" begin for x in (rand(3), 0.5, Complex(0.2,0.7), (3,5,9), Ref(42), rand(3, 4, 5), Complex.(rand(3,5), rand(3,5))) - V, to_x = rv_and_back(x) + V, to_x = torv_and_back(x) @test V isa AbstractVector{<:Real} - @test V == rv_and_back(x)[1] + @test V == torv_and_back(x)[1] @test x isa Ref ? to_x(V)[] == x[] : to_x(V) == x end end @@ -34,10 +33,10 @@ end @test_throws ErrorException @inferred with_logabsdet_jacobian(noninferrable_inv, rand(2, 2)) test_with_logabsdet_jacobian(inv, rx, ForwardDiff.derivative, atol = 10^-6) - test_with_logabsdet_jacobian(inv, cx, ForwardDiff.jacobian, rv_and_back, atol = 10^-6) + test_with_logabsdet_jacobian(inv, cx, getjacobian, atol = 10^-6) test_with_logabsdet_jacobian(inv, X, ForwardDiff.jacobian, atol = 10^-6) - test_with_logabsdet_jacobian(inv, CX, ForwardDiff.jacobian, rv_and_back, atol = 10^-6) - test_with_logabsdet_jacobian(inv, CX, ForwardDiff.jacobian, rv_and_back, atol = 10^-6) - test_with_logabsdet_jacobian(inv, CX, ForwardDiff.jacobian, rv_and_back, compare = myisapprox, atol = 10^-6) - test_with_logabsdet_jacobian(noninferrable_inv, CX, ForwardDiff.jacobian, rv_and_back, atol = 10^-6) + test_with_logabsdet_jacobian(inv, CX, getjacobian, atol = 10^-6) + test_with_logabsdet_jacobian(inv, CX, getjacobian, atol = 10^-6) + test_with_logabsdet_jacobian(inv, CX, getjacobian, compare = myisapprox, atol = 10^-6) + test_with_logabsdet_jacobian(noninferrable_inv, CX, getjacobian, atol = 10^-6) end diff --git a/test/test_with_ladj.jl b/test/test_with_ladj.jl index 50c6087..eff2997 100644 --- a/test/test_with_ladj.jl +++ b/test/test_with_ladj.jl @@ -4,45 +4,38 @@ using ChangesOfVariables using Test using LinearAlgebra -import ForwardDiff using ChangesOfVariables: test_with_logabsdet_jacobian -include("rv_and_back.jl") +include("getjacobian.jl") @testset "with_logabsdet_jacobian" begin - foo(x) = inv(exp(-x) + 1) - function ChangesOfVariables.with_logabsdet_jacobian(::typeof(foo), x) y = foo(x) ladj = -x + 2 * log(y) (y, ladj) end - x = 4.2 X = rand(10) A = rand(5, 5) CA = Complex.(rand(5, 5), rand(5, 5)) - - getjacobian = ForwardDiff.jacobian - isaprx(a, b; kwargs...) = isapprox(a,b; kwargs...) isaprx(a::NTuple{N,Any}, b::NTuple{N,Any}; kwargs...) where N = all(map((a,b) -> isaprx(a, b; kwargs...), a, b)) - test_with_logabsdet_jacobian(foo, x, getjacobian, rv_and_back) + test_with_logabsdet_jacobian(foo, x, getjacobian) @static if VERSION >= v"1.6" - test_with_logabsdet_jacobian(log ∘ foo, x, getjacobian, rv_and_back) + test_with_logabsdet_jacobian(log ∘ foo, x, getjacobian) end @testset "getjacobian on mapped and broadcasted" begin for f in (Base.Fix1(map, foo), Base.Fix1(broadcast, foo)) for arg in (x, fill(x,), Ref(x), (x,), X) - test_with_logabsdet_jacobian(f, arg, getjacobian, rv_and_back, compare = isaprx) + test_with_logabsdet_jacobian(f, arg, getjacobian, compare = isaprx) end end end @@ -50,14 +43,14 @@ include("rv_and_back.jl") @testset "getjacobian on identity, adjoint and transpose" begin for f in (identity, adjoint, transpose) for arg in (x, A) - test_with_logabsdet_jacobian(f, arg, getjacobian, rv_and_back) + test_with_logabsdet_jacobian(f, arg, getjacobian) end end end @testset "getjacobian on inv" begin for arg in (x, A, CA) - test_with_logabsdet_jacobian(inv, arg, getjacobian, rv_and_back) + test_with_logabsdet_jacobian(inv, arg, getjacobian) end end end From 6cbb1c3d32e9b927196c8a373cb17c19cbb9740b Mon Sep 17 00:00:00 2001 From: Oliver Schulz Date: Fri, 15 Oct 2021 10:01:25 +0200 Subject: [PATCH 7/9] Use jldoctest on with_logabsdet_jacobian example --- src/with_ladj.jl | 25 ++++++++++++++++++++----- 1 file changed, 20 insertions(+), 5 deletions(-) diff --git a/src/with_ladj.jl b/src/with_ladj.jl index ce38af6..4b5ddd6 100644 --- a/src/with_ladj.jl +++ b/src/with_ladj.jl @@ -14,9 +14,9 @@ For `(y, ladj) = with_logabsdet_jacobian(f, x)`, the following must hold true: `with_logabsdet_jacobian` comes with support for broadcasted/mapped functions (via `Base.Fix1`) and (Julia >=v1.6 only) `ComposedFunction`. -Example: +# Examples -```julia +```jldoctest a foo(x) = inv(exp(-x) + 1) function ChangesOfVariables.with_logabsdet_jacobian(::typeof(foo), x) @@ -28,14 +28,29 @@ end x = 4.2 y, ladj_y = with_logabsdet_jacobian(foo, x) -X = rand(10) +# output + +(0.9852259683067269, -4.229768509343836) +``` + +```jldoctest a +X = [3, 7, 5] broadcasted_foo = Base.Fix1(broadcast, foo) Y, ladj_Y = with_logabsdet_jacobian(broadcasted_foo, X) +# output + +([0.9525741268224334, 0.9990889488055994, 0.9933071490757153], -15.112428333033268) +``` + +```jldoctest a # Requires Julia >= v1.6: z, ladj_z = with_logabsdet_jacobian(log ∘ foo, x) -z == log(foo(x)) -ladj_z == ladj_y + with_logabsdet_jacobian(log, y)[2] +z == log(foo(x)) && ladj_z == ladj_y + with_logabsdet_jacobian(log, y)[2] + +# output + +true ``` """ function with_logabsdet_jacobian end From 69e26e818393dfb04bb6101f1e58d821bfd968dc Mon Sep 17 00:00:00 2001 From: Oliver Schulz Date: Fri, 15 Oct 2021 10:06:15 +0200 Subject: [PATCH 8/9] Refer to test routine in with_logabsdet_jacobian docstring. --- src/with_ladj.jl | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/src/with_ladj.jl b/src/with_ladj.jl index 4b5ddd6..1a7259d 100644 --- a/src/with_ladj.jl +++ b/src/with_ladj.jl @@ -52,6 +52,10 @@ z == log(foo(x)) && ladj_z == ladj_y + with_logabsdet_jacobian(log, y)[2] true ``` + +Implementations of with_logabsdet_jacobian can be tested (as a +`Test.@testset`) using +[`ChangesOfVariables.test_with_logabsdet_jacobian`](@ref). """ function with_logabsdet_jacobian end export with_logabsdet_jacobian From cabe46d17dc693cbc6362564b7eec20b9e307be4 Mon Sep 17 00:00:00 2001 From: Oliver Schulz Date: Fri, 15 Oct 2021 14:43:18 +0200 Subject: [PATCH 9/9] Stabilize doctests --- docs/Project.toml | 2 ++ src/with_ladj.jl | 19 +++++++++++++------ 2 files changed, 15 insertions(+), 6 deletions(-) diff --git a/docs/Project.toml b/docs/Project.toml index f471c80..f656746 100644 --- a/docs/Project.toml +++ b/docs/Project.toml @@ -1,5 +1,7 @@ [deps] Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4" +ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210" +LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" Markdown = "d6f4376e-aef5-505a-96c1-9c027394607a" [compat] diff --git a/src/with_ladj.jl b/src/with_ladj.jl index 1a7259d..7f5dc9a 100644 --- a/src/with_ladj.jl +++ b/src/with_ladj.jl @@ -17,6 +17,8 @@ For `(y, ladj) = with_logabsdet_jacobian(f, x)`, the following must hold true: # Examples ```jldoctest a +using ChangesOfVariables + foo(x) = inv(exp(-x) + 1) function ChangesOfVariables.with_logabsdet_jacobian(::typeof(foo), x) @@ -28,25 +30,30 @@ end x = 4.2 y, ladj_y = with_logabsdet_jacobian(foo, x) +using LinearAlgebra, ForwardDiff +y == foo(x) && ladj_y ≈ log(abs(ForwardDiff.derivative(foo, x))) + # output -(0.9852259683067269, -4.229768509343836) +true ``` ```jldoctest a -X = [3, 7, 5] +X = rand(10) broadcasted_foo = Base.Fix1(broadcast, foo) Y, ladj_Y = with_logabsdet_jacobian(broadcasted_foo, X) +Y == broadcasted_foo(X) && ladj_Y ≈ logabsdet(ForwardDiff.jacobian(broadcasted_foo, X))[1] # output -([0.9525741268224334, 0.9990889488055994, 0.9933071490757153], -15.112428333033268) +true ``` ```jldoctest a -# Requires Julia >= v1.6: -z, ladj_z = with_logabsdet_jacobian(log ∘ foo, x) -z == log(foo(x)) && ladj_z == ladj_y + with_logabsdet_jacobian(log, y)[2] +VERSION < v"1.6" || begin # Support for ∘ requires Julia >= v1.6 + z, ladj_z = with_logabsdet_jacobian(log ∘ foo, x) + z == log(foo(x)) && ladj_z == ladj_y + with_logabsdet_jacobian(log, y)[2] +end # output