Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Introduce Kendall's shape space #550

Merged
merged 17 commits into from
Nov 12, 2022
Merged
Show file tree
Hide file tree
Changes from 1 commit
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
25 changes: 21 additions & 4 deletions docs/src/manifolds/shapespace.md
Original file line number Diff line number Diff line change
@@ -1,6 +1,9 @@
# Shape spaces

Kendall's pre-shape and shape spaces.
Shape spaces are spaces of ``k`` points in ``\mathbb{R}^n`` up to simultaneous action of a group on all points.
The most commonly encountered are Kendall's pre-shape and shape spaces.
In the case of the Kendall's pre-shape spaces the action is translation and scaling.
In the case of the Kendall's shape spaces the action is translation, scaling and rotation.

```@example
kellertuer marked this conversation as resolved.
Show resolved Hide resolved
using Manifolds, Plots
Expand All @@ -26,16 +29,30 @@ rot_q = apply(A, a, q)
scatter!(fig, rot_q[1,:], rot_q[2,:], label="q aligned to p")
```

A more extensive usage example is available in the `hand_gestures.jl` tutorial.

```@autodocs
Modules = [Manifolds, ManifoldsBase]
Pages = ["manifolds/KendallsPreShapeSpace.jl"]
Order = [:type]
```

```@autodocs
Modules = [Manifolds, ManifoldsBase]
Pages = ["manifolds/ShapeSpace.jl"]
Pages = ["manifolds/KendallsShapeSpace.jl"]
Order = [:type]
```

## Provided functions

```@autodocs
Modules = [Manifolds, ManifoldsBase]
Pages = ["manifolds/ShapeSpace.jl"]
Pages = ["manifolds/KendallsPreShapeSpace.jl"]
Order = [:function]
```
```

```@autodocs
Modules = [Manifolds, ManifoldsBase]
Pages = ["manifolds/KendallsShapeSpace.jl"]
Order = [:function]
```
3 changes: 2 additions & 1 deletion src/Manifolds.jl
Original file line number Diff line number Diff line change
Expand Up @@ -368,7 +368,8 @@ include("manifolds/Rotations.jl")
include("manifolds/Orthogonal.jl")

# shape spaces require Sphere
include("manifolds/ShapeSpace.jl")
include("manifolds/KendallsPreShapeSpace.jl")
include("manifolds/KendallsShapeSpace.jl")

# Introduce the quotient, Grassmann, only after Stiefel
include("manifolds/Grassmann.jl")
Expand Down
12 changes: 12 additions & 0 deletions src/groups/rotation_action.jl
Original file line number Diff line number Diff line change
Expand Up @@ -276,6 +276,18 @@ function inverse_apply(::LeftColumnwiseMultiplicationAction, a, p)
return a \ p
end

"""
optimal_alignment(A::LeftColumnwiseMultiplicationAction, p, q)

Compute optimal alignment for the left [`ColumnwiseMultiplicationAction`](@ref), i.e. the
group element that, when it acts on `p`, returns the point closest to `q`. Details of
computation are described in Section 2.2.1 of [^Srivastava2016].

# References

[^Srivastava2016]:
> A. Srivastava and E. P. Klassen, Functional and Shape Data Analysis. Springer New York, 2016.
"""
function optimal_alignment(A::LeftColumnwiseMultiplicationAction, p, q)
kellertuer marked this conversation as resolved.
Show resolved Hide resolved
is_point(A.manifold, p, true)
is_point(A.manifold, q, true)
Expand Down
149 changes: 149 additions & 0 deletions src/manifolds/KendallsPreShapeSpace.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,149 @@

@doc raw"""
KendallsPreShapeSpace{n,k} <: AbstractSphere{ℝ}

Kendall's pre-shape space of ``k`` landmarks in ``ℝ^n`` represented by n×k matrices.
In each row the sum of elements of a matrix is equal to 0. The Frobenius norm of the matrix
is equal to 1 [^Kendall1984][^Kendall1989].

The space can be interpreted as tuples of ``k`` points in ``ℝ^n`` up to simultaneous
translation and scaling of all points, so this can be thought of as a quotient manifold.

# Constructor

KendallsPreShapeSpace(n::Int, k::Int)

# References

[^Kendall1989]:
> D. G. Kendall, “A Survey of the Statistical Theory of Shape,” Statist. Sci., vol. 4,
> no. 2, pp. 87–99, May 1989, doi: 10.1214/ss/1177012582.
mateuszbaran marked this conversation as resolved.
Show resolved Hide resolved
[^Kendall1984]:
> D. G. Kendall, “Shape Manifolds, Procrustean Metrics, and Complex Projective Spaces,”
> Bull. London Math. Soc., vol. 16, no. 2, pp. 81–121, Mar. 1984, doi: 10.1112/blms/16.2.81.
"""
struct KendallsPreShapeSpace{n,k} <: AbstractSphere{ℝ} end

KendallsPreShapeSpace(n::Int, k::Int) = KendallsPreShapeSpace{n,k}()

function active_traits(f, ::KendallsPreShapeSpace, args...)
return merge_traits(IsEmbeddedSubmanifold())
end

representation_size(::KendallsPreShapeSpace{n,k}) where {n,k} = (n, k)

"""
check_point(M::KendallsPreShapeSpace, p; atol=sqrt(max_eps(X, Y)), kwargs...)

Check whether `p` is a valid point on [`KendallsPreShapeSpace`](@ref), i.e. whether
each row has zero mean. Other conditions are checked via embedding in [`ArraySphere`](@ref).
"""
function check_point(M::KendallsPreShapeSpace, p; atol=sqrt(eps(eltype(p))), kwargs...)
for p_row in eachrow(p)
if !isapprox(mean(p_row), 0; atol, kwargs...)
return DomainError(
mean(p_row),
"The point $(p) does not lie on the $(M) since one of the rows does not have zero mean.",
)
end
end
return nothing
end

"""
check_vector(M::KendallsPreShapeSpace, p, X; kwargs... )

Check whether `X` is a valid tangent vector on [`KendallsPreShapeSpace`](@ref), i.e. whether
each row has zero mean. Other conditions are checked via embedding in [`ArraySphere`](@ref).
"""
function check_vector(M::KendallsPreShapeSpace, p, X; atol=sqrt(max_eps(X, Y)), kwargs...)
for X_row in eachrow(X)
if !isapprox(mean(X_row), 0; atol, kwargs...)
return DomainError(
mean(X_row),
"The vector $(X) is not a tangent vector to $(p) on $(M), since one of the rows does not have zero mean.",
)
end
end
return nothing
end

embed(::KendallsPreShapeSpace, p) = p
embed(::KendallsPreShapeSpace, p, X) = X

function get_embedding(::KendallsPreShapeSpace{N,K}) where {N,K}
return ArraySphere(N, K)
end

@doc raw"""
manifold_dimension(M::KendallsPreShapeSpace)

Return the dimension of the [`KendallsPreShapeSpace`](@ref) manifold `M`. The dimension is
given by ``n(k - 1) - 1``.
"""
manifold_dimension(::KendallsPreShapeSpace{n,k}) where {n,k} = n * (k - 1) - 1

"""
project(M::KendallsPreShapeSpace, p)

Project point `p` from the embedding to [`KendallsPreShapeSpace`](@ref) by selecting
the right element from the orthogonal section representing the quotient manifold `M`.
See Section 3.7 of [^Srivastava2016] for details.

# References

[^Srivastava2016]:
> A. Srivastava and E. P. Klassen, Functional and Shape Data Analysis. Springer New York, 2016.
mateuszbaran marked this conversation as resolved.
Show resolved Hide resolved
"""
project(::KendallsPreShapeSpace, p)

function project!(::KendallsPreShapeSpace, q, p)
q .= p .- mean(p, dims=2)
q ./= norm(q)
return q
end

"""
project(M::KendallsPreShapeSpace, p, X)

Project tangent vector `X` at point `p` from the embedding to [`KendallsPreShapeSpace`](@ref)
by selecting the right element from the tangent space to orthogonal section representing the
quotient manifold `M`. See Section 3.7 of [^Srivastava2016] for details.

# References

[^Srivastava2016]:
kellertuer marked this conversation as resolved.
Show resolved Hide resolved
> A. Srivastava and E. P. Klassen, Functional and Shape Data Analysis. Springer New York, 2016.
"""
project(::KendallsPreShapeSpace, p, X)

function project!(::KendallsPreShapeSpace, Y, p, X)
Y .= X .- mean(X, dims=2)
Y .-= dot(p, Y) .* p
return Y
end

function Random.rand!(M::KendallsPreShapeSpace, pX; vector_at=nothing, σ=one(eltype(pX)))
if vector_at === nothing
project!(M, pX, randn(representation_size(M)))
else
n = σ * randn(size(pX)) # Gaussian in embedding
project!(M, pX, vector_at, n) # project to TpM (keeps Gaussianness)
end
return pX
end
function Random.rand!(
rng::AbstractRNG,
M::KendallsPreShapeSpace,
pX;
vector_at=nothing,
σ=one(eltype(pX)),
)
if vector_at === nothing
project!(M, pX, randn(rng, representation_size(M)))
else
n = σ * randn(rng, size(pX)) # Gaussian in embedding
project!(M, pX, vector_at, n) #project to TpM (keeps Gaussianness)
end
return pX
end
129 changes: 20 additions & 109 deletions src/manifolds/ShapeSpace.jl → src/manifolds/KendallsShapeSpace.jl
Original file line number Diff line number Diff line change
@@ -1,122 +1,27 @@

@doc raw"""
KendallsPreShapeSpace{n,k} <: AbstractSphere{ℝ}

Kendall's pre-shape space of ``k`` landmarks in $ℝ^n$ represented by n×k matrices.

# Constructor

KendallsPreShapeSpace(n::Int, k::Int)
"""
struct KendallsPreShapeSpace{n,k} <: AbstractSphere{ℝ} end

KendallsPreShapeSpace(n::Int, k::Int) = KendallsPreShapeSpace{n,k}()

function active_traits(f, ::KendallsPreShapeSpace, args...)
return merge_traits(IsEmbeddedSubmanifold())
end

representation_size(::KendallsPreShapeSpace{n,k}) where {n,k} = (n, k)

"""
check_point(M::KendallsPreShapeSpace, p; atol=sqrt(max_eps(X, Y)), kwargs...)

Check whether `p` is a valid point on [`KendallsPreShapeSpace`](@ref), i.e. whether
each row has zero mean. Other conditions are checked via embedding in [`ArraySphere`](@ref).
"""
function check_point(M::KendallsPreShapeSpace, p; atol=sqrt(eps(eltype(p))), kwargs...)
for p_row in eachrow(p)
if !isapprox(mean(p_row), 0; atol, kwargs...)
return DomainError(
mean(p_row),
"The point $(p) does not lie on the $(M) since one of the rows does not have zero mean.",
)
end
end
return nothing
end

"""
check_vector(M::KendallsPreShapeSpace, p, X; kwargs... )

Check whether `X` is a valid tangent vector on [`KendallsPreShapeSpace`](@ref), i.e. whether
each row has zero mean. Other conditions are checked via embedding in [`ArraySphere`](@ref).
"""
function check_vector(M::KendallsPreShapeSpace, p, X; atol=sqrt(max_eps(X, Y)), kwargs...)
for X_row in eachrow(X)
if !isapprox(mean(X_row), 0; atol, kwargs...)
return DomainError(
mean(X_row),
"The vector $(X) is not a tangent vector to $(p) on $(M), since one of the rows does not have zero mean.",
)
end
end
return nothing
end

embed(::KendallsPreShapeSpace, p) = p
embed(::KendallsPreShapeSpace, p, X) = X

function get_embedding(::KendallsPreShapeSpace{N,K}) where {N,K}
return ArraySphere(N, K)
end

@doc raw"""
manifold_dimension(M::KendallsPreShapeSpace)

Return the dimension of the [`KendallsPreShapeSpace`](@ref) manifold `M`. The dimension is
given by ``n(k - 1) - 1``.
"""
manifold_dimension(::KendallsPreShapeSpace{n,k}) where {n,k} = n * (k - 1) - 1

function project!(::KendallsPreShapeSpace, q, p)
q .= p .- mean(p, dims=2)
q ./= norm(q)
return q
end

function project!(::KendallsPreShapeSpace, Y, p, X)
Y .= X .- mean(X, dims=2)
Y .-= dot(p, Y) .* p
return Y
end

function Random.rand!(M::KendallsPreShapeSpace, pX; vector_at=nothing, σ=one(eltype(pX)))
if vector_at === nothing
project!(M, pX, randn(representation_size(M)))
else
n = σ * randn(size(pX)) # Gaussian in embedding
project!(M, pX, vector_at, n) # project to TpM (keeps Gaussianness)
end
return pX
end
function Random.rand!(
rng::AbstractRNG,
M::KendallsPreShapeSpace,
pX;
vector_at=nothing,
σ=one(eltype(pX)),
)
if vector_at === nothing
project!(M, pX, randn(rng, representation_size(M)))
else
n = σ * randn(rng, size(pX)) # Gaussian in embedding
project!(M, pX, vector_at, n) #project to TpM (keeps Gaussianness)
end
return pX
end

###

@doc raw"""
KendallsShapeSpace{n,k} <: AbstractDecoratorManifold{ℝ}

Kendall's shape space, defined as quotient of a [`KendallsPreShapeSpace`](@ref)
(represented by n×k matrices) by the action [`ColumnwiseMultiplicationAction`](@ref).

The space can be interpreted as tuples of ``k`` points in ``ℝ^n`` up to simultaneous
translation and scaling and rotation of all points [^Kendall1984][^Kendall1989].

This manifold possesses the [`IsQuotientManifold`](@ref) trait.

# Constructor

KendallsShapeSpace(n::Int, k::Int)

# References

[^Kendall1989]:
> D. G. Kendall, “A Survey of the Statistical Theory of Shape,” Statist. Sci., vol. 4,
> no. 2, pp. 87–99, May 1989, doi: 10.1214/ss/1177012582.
[^Kendall1984]:
> D. G. Kendall, “Shape Manifolds, Procrustean Metrics, and Complex Projective Spaces,”
> Bull. London Math. Soc., vol. 16, no. 2, pp. 81–121, Mar. 1984, doi: 10.1112/blms/16.2.81.
"""
struct KendallsShapeSpace{n,k} <: AbstractDecoratorManifold{ℝ} end

Expand Down Expand Up @@ -165,6 +70,12 @@ end
embed(::KendallsShapeSpace, p) = p
embed(::KendallsShapeSpace, p, X) = X

"""
get_embedding(M::KendallsShapeSpace)

Get the manifold in which [`KendallsShapeSpace`](@ref) `M` is embedded, i.e.
[`KendallsPreShapeSpace`](@ref) of matrices of the same shape.
"""
function get_embedding(::KendallsShapeSpace{N,K}) where {N,K}
return KendallsPreShapeSpace(N, K)
end
Expand Down