point_extrinsic (array-like, shape=[…, *embedding_space.point_shape]) – Point in the embedded manifold in extrinsic coordinates,
i. e. in the coordinates of the embedding manifold.
Returns:
point_intrinsic (array-lie, shape=[…, *point_shape]) – Point in the embedded manifold in intrinsic coordinates.
Compute the coefficients of matrices in the given basis.
This takes a matrix (the matrix representation of a point) and
transforms it into its corresponding vector representation
(the coefficients wrt a given basis).
Previously, this method was called to_vector. basis_representation
makes it more clear that the vector representation depends on the chosen
basis.
Points are sampled from the embedding space using the distribution set
for that manifold and then projected to the manifold. As a result, this
is not a uniform distribution on the manifold itself.
Parameters:
n_samples (int) – Number of samples.
Optional, default: 1.
bound (float) – Bound of the interval in which to sample for the embedding space.
Optional, default: 1.
Returns:
samples (array-like, shape=[…, *point_shape]) – Points sampled on the hypersphere.
This method is not recommended for statistical purposes, as the
tangent vectors generated are not drawn from a distribution related
to the Riemannian metric.
Parameters:
n_samples (int) – Number of samples.
Optional, default: 1.
The Poincaré complex disk is a representation of the Hyperbolic space of dimension 2.
Lead author: Yann Cabanes.
References
[Cabanes2022]
Yann Cabanes. Multidimensional complex stationary
centered Gaussian autoregressive time series machine learning
in Poincaré and Siegel disks: application for audio and radar
clutter classification, PhD thesis, 2022
Marc Arnaudon, Frédéric Barbaresco and Le Yang.
Riemannian medians and means with applications to radar signal processing,
IEEE Journal of Selected Topics in Signal Processing, vol. 7, no. 4, pp. 595-604,
Aug. 2013, doi: 10.1109/JSTSP.2013.2261798.
https://ieeexplore.ieee.org/document/6514112
Generate a random unit tangent vector at a given point.
Parameters:
base_point (array-like, shape=[…, dim]) – Point.
n_vectors (float) – Number of vectors to be generated at base_point.
For vectorization purposes n_vectors can be greater than 1 iff
base_point consists of a single point.
Returns:
normalized_vector (array-like, shape=[…, n_vectors, dim]) – Random unit tangent vector at base_point.
Compute the covariant derivative of the directional curvature.
For tangent vector fields at base_point \(P\):
\(X|_P = tangent\_vec\_a\),
\(Y|_P = tangent\_vec\_b\),
the covariant derivative (in the direction X)
\((\nabla_X R_Y)(X) |_P = (\nabla_X R)(Y, X) Y |_P\) of the
directional curvature (in the direction Y)
\(R_Y(X) = R(Y, X) Y\)
is a quadratic tensor in X and Y that
plays an important role in the computation of the moments of the
empirical Fréchet mean.
References
[Pennec]
Pennec, Xavier. Curvature effects on the empirical mean in
Riemannian and affine Manifolds: a non-asymptotic high
concentration expansion in the small-sample regime. Preprint. 2019.
https://arxiv.org/abs/1906.07418
Generate parameterized function for the geodesic curve.
Geodesic curve defined by either:
an initial point and an initial tangent vector,
an initial point and an end point.
Parameters:
initial_point (array-like, shape=[…, dim]) – Point on the manifold, initial point of the geodesic.
end_point (array-like, shape=[…, dim], optional) – Point on the manifold, end point of the geodesic. If None,
an initial tangent vector must be given.
initial_tangent_vec (array-like, shape=[…, dim],) – Tangent vector at base point, the initial speed of the geodesics.
Optional, default: None.
If None, an end point must be given and a logarithm is computed.
Returns:
path (callable) – Time parameterized geodesic curve. If a batch of initial
conditions is passed, the output array’s first dimension
represents the different initial conditions, and the second
corresponds to time.
This is is the supremum of radii r for which the exponential map is a
diffeomorphism from the open ball of radius r centered at the base
point onto its image.
Parameters:
base_point (array-like, shape=[…, {dim, [n, m]}]) – Point on the manifold.
Approximate parallel transport using the pole ladder scheme.
Approximate Parallel transport using either the pole ladder or the
Schild’s ladder scheme [LP2013b]. Pole ladder is exact in symmetric
spaces and of order two in general while Schild’s ladder is a first
order approximation [GP2020]. Both schemes are available on any affine
connection manifolds whose exponential and logarithm maps are
implemented. tangent_vec is transported along the geodesic starting
at the base_point with initial tangent vector direction.
Parameters:
tangent_vec (array-like, shape=[…, dim]) – Tangent vector at base point to transport.
direction (array-like, shape=[…, dim]) – Tangent vector at base point, initial speed of the geodesic along
which to transport.
base_point (array-like, shape=[…, dim]) – Point on the manifold, initial position of the geodesic along
which to transport.
n_rungs (int) – Number of steps of the ladder.
Optional, default: 1.
scheme (str, {‘pole’, ‘schild’}) – The scheme to use for the construction of the ladder at each step.
Optional, default: ‘pole’.
alpha (float) – Exponent for the scaling of the vector to transport. Must be
greater or equal to 1, 2 is optimal. See [GP2020].
Optional, default: 2
Returns:
ladder (dict of array-like and callable with following keys) –
Approximation of the parallel transport of tangent vector a.
trajectorylist of list of callable, len=n_steps
List of lists containing the geodesics of the
construction, only if return_geodesics=True in the step
function. The geodesics are methods of the class connection.
Lorenzi, Marco, and Xavier Pennec. “Efficient Parallel
Transport of Deformations in Time Series of Images: From Schild to
Pole Ladder.” Journal of Mathematical Imaging and Vision 50, no. 1
(September 1, 2014): 5–17.
https://doi.org/10.1007/s10851-013-0470-3.
Guigui, Nicolas, and Xavier Pennec. “Numerical Accuracy
of Ladder Schemes for Parallel Transport on Manifolds.”
Foundations of Computational Mathematics, June 18, 2021.
https://doi.org/10.1007/s10208-021-09515-x.
Compute the parallel transport of a tangent vector.
Closed-form solution for the parallel transport of a tangent vector
along the geodesic between two points base_point and end_point
or alternatively defined by \(t \mapsto exp_{(base\_point)}(
t*direction)\).
Parameters:
tangent_vec (array-like, shape=[…, {dim, [n, m]}]) – Tangent vector at base point to be transported.
base_point (array-like, shape=[…, {dim, [n, m]}]) – Point on the manifold. Point to transport from.
direction (array-like, shape=[…, {dim, [n, m]}]) – Tangent vector at base point, along which the parallel transport
is computed.
Optional, default: None.
end_point (array-like, shape=[…, {dim, [n, m]}]) – Point on the manifold. Point to transport to.
Optional, default: None.
Note that geomstats puts the contravariant index on
the first dimension of the Christoffel symbols.
Parameters:
base_point (array-like, shape=[…, dim]) – Point on the manifold.
Returns:
riemann_curvature (array-like, shape=[…, dim, dim, dim, dim]) – riemann_tensor[…,i,j,k,l] = R_{ijk}^l
Riemannian tensor curvature,
with the contravariant index on the last dimension.
image_tangent_vec (array-like, shape=[…, *image_shape]) – Image tangent vector at image of the base point.
Notes
Signature choice (i.e. having the possibility to pass both base and image
points) comes from performance considerations.
In several methods (e.g. PullbackDiffeoMetric.inner_product) the need
to call tangent_diffeomorphism comes after having called
diffeomorphism, which means we have both base_point and
image_point available.
In some cases, tangent_diffeomorphism needs access to base_point
(e.g. SRVTransform), while in others, it needs access to image_point
(e.g. CholeskyMap).
By passing both, we give the corresponding implementation the possibility
to choose the one that is more convenient (performancewise).
If we pass only one of the two, it has the mechanims to perform the
necessary transformations (e.g. calling diffeomorphism or
inverse_diffeomorphism on it).
Each individual curve is represented by a 2d-array of shape [
k_sampling_points - 1, ambient_dim].
This space corresponds to the space of immersions defined below, i.e. the
space of smooth functions from an interval I into the ambient Euclidean
space M, with non-vanishing derivative.
\[Imm(I, M)=\{c \in C^{\infty}(I, M) \|c'(s)\|\neq 0 \forall s \in I \},\]
where the interval of parameters I is taken to be I = [0, 1]
without loss of generality.
Parameters:
ambient_dim (int) – Dimension of the ambient Euclidean space in which curves take values.
k_sampling_points (int) – Number of sampling points.
equip (bool) – If True, equip space with default metric.
Find the reparametrization gamma of end_curve that minimizes the distance
between initial_curve and end_curve reparametrized by gamma, and output
the corresponding distance, using a dynamic programming algorithm.
The objective can be expressed in terms of square root velocity (SRV)
representations: it is equivalent to finding the gamma that maximizes
the L2 scalar product between \(initial_{srv}\) and \(end_{srv}@\gamma\),
where \(initial_{srv}\) is the SRV representation of the initial curve
and \(end_{srv}@\gamma\) is the SRV representation of the end curve
reparametrized by \(\gamma\), i.e
The dynamic programming algorithm assumes that for every subinterval
\(\left[\frac{i}{n},\frac{i+1}{n}\right]\) of \(\left[0,1\right]\),
gamma is linear.
Parameters:
total_space (Manifold) – Total space with reparametrizations fiber bundle structure.
n_space_grid (int) – Number of sampling points of the unit interval.
max_slope (int) – Maximum slope allowed for a reparametrization.
References
[WAJ2007] M. Washington, S. Anuj & H. Joshi,
“On Shape of Plane Elastic Curves”, in International Journal of Computer
Vision. 73(3):307-324, 2007.
S. Kurtek and T. Needham,
“Simplifying transforms for general elastic metrics on the space of
plane curves”, arXiv:1803.10894 [math.DG], 29 Mar 2018.
[BCKKNP2022]
Martin Bauer, Nicolas Charon, Eric Klassen, Sebastian Kurtek,
Tom Needham, and Thomas Pierron.
“Elastic Metrics on Spaces of Euclidean Curves: Theory and Algorithms.”
arXiv, September 20, 2022. https://doi.org/10.48550/arXiv.2209.09862.
The f_transform is defined on the space of curves
quotiented by translations, which is identified with the space
of curves with their first sampling point located at 0:
Align two curves through iterative horizontal geodesic algorithm.
This algorithm computes the horizontal geodesic between two curves in the shape
bundle of curves modulo reparametrizations, and at the same time, aligns the
end curve with respect to the initial curve. This is done through an iterative
procedure where the initial curve stays fixed and the sampling points are moved
on the end curve to obtain its optimal parametrization with respect to the initial
curve. This procedure is based on the decomposition of any path of curves into a
horizontal path of curves composed with a path of reparametrizations:
\(c(t, s) = c_h(t, phi(t, s))\) where \(d/dt c_h(t, .)\) is horizontal.
Here t is the time parameter of the path and s the space parameter of the curves.
The algorithm sets current_end_curve to be the end curve and iterates three steps:
1) compute the geodesic between the initial curve and current_end_curve
2) compute the path of reparametrizations such that the path of its inverses
transforms this geodesic into a horizontal path of curves
3) invert this path of reparametrizations to find the horizontal path and update
current_end_curve to be its end point.
The algorithm stops when the new current_end_curve is sufficiently
close to the former current_end_curve.
Parameters:
total_space (Manifold) – Total space with reparametrizations fiber bundle structure.
n_time_grid (int) – Number of times in which compute the geodesic.
threshold (float) – When the difference between the new end curve and the current end
curve becomes lower than this threshold, the algorithm stops.
Optional, default: 1e-3.
max_iter (int) – Maximum number of iterations.
Optional, default: 20.
tol (float) – Minimal spacing between time samples in the unit segment when
reparametrizing the end curve.
Optional, default: 1e-3.
verbose (boolean) – Optional, default: False.
save_history (bool) – If True, history is saved in a self.history.
References
[LAB2017]
A. Le Brigant, “Optimal matching between curves in a manifold”,
in Geometric Science of Information. Springer Lecture Notes in Computer
Science 10589 (2017), 57 - 64. https://hal.science/hal-04374199.
L2 metric on the space of regularly sampled discrete curves
defined on the unit interval. The inner product between two tangent vectors
is given by the integral of the ambient space inner product, approximated by
a left Riemann sum.
Compute the left Riemann sum approximation of the integral.
Compute the left Riemann sum approximation of the integral of a
function func defined on the unit interval, given by sample points
at regularly spaced times
\[\begin{split}t_i = i / k, \\
i = 0, ..., k - 1\end{split}\]
where \(k\) is the number of landmarks (last time is missing).
Parameters:
func (array-like, shape=[…, k_landmarks]) – Sample points of a function at regularly spaced times.
Returns:
riemann_sum (array-like, shape=[…, ]) – Left Riemann sum.
Principal bundle of curves modulo reparametrizations with an elastic metric.
The space of parameterized curves is the total space of a principal bundle
where the group action is given by reparametrization and the base space is
the shape space of curves modulo reparametrization, i.e.unparametrized
curves. In the discrete case, reparametrization corresponds to resampling.
Each tangent vector to the space of parameterized curves can be split into a
vertical part (tangent to the fibers of the principal bundle) and a
horizontal part (orthogonal to the vertical part with respect to the SRV
metric). The geodesic between the shapes of two curves is computed by
aligning (i.e. reparametrizing) one of the two curves with respect to the
other, and computing the geodesic between the aligned curves. This geodesic
will be horizontal, and will project to a geodesic on the shape space.
Two different aligners are available:
- IterativeHorizontalGeodesicAligner (default)
- DynamicProgrammingAligner.
Parameters:
total_space (DiscreteCurvesStartingAtOrigin) – Space of discrete curves starting at the origin
Principal bundle of curves modulo rotations with an elastic metric.
This is the fiber bundle where the total space is the space of parameterized
curves equipped with an elastic metric, the action is given by rotations, and
the base space is the shape space of curves modulo rotations.
Parameters:
total_space (DiscreteCurvesStartingAtOrigin) – Space of discrete curves starting at the origin
A. Srivastava, E. Klassen, S. H. Joshi and I. H. Jermyn,
“Shape Analysis of Elastic Curves in Euclidean Spaces,”
in IEEE Transactions on Pattern Analysis and Machine Intelligence,
vol. 33, no. 7, pp. 1415-1428, July 2011.
tangent_vec (array-like, shape=[…, ambient_dim]) – Inverse of the differential of the square root velocity transform at
curve evaluated at tangent_vec.
d_srv_vec (array-like, shape=[…, k_sampling_points - 1, ambient_dim,]) – Differential of the square root velocity transform at curve
evaluated at tangent_vec.
Emmanuel Hartman, Yashil Sukurdeep, Eric Klassen,
Nicolas Charon, and Martin Bauer.
“Elastic shape analysis of surfaces with second-order Sobolev metrics:
a comprehensive numerical framework”. arXiv:2204.04238 [cs.CV], 25 Sep 2022
[HPBDC2023]
Emmanuel Hartman, Emery Pierson, Martin Bauer,
Mohamed Daoudi, and Nicolas Charon.
“Basis Restricted Elastic Shape Analysis on the Space of Unregistered Surfaces.”
arXiv, November 7, 2023. https://doi.org/10.48550/arXiv.2311.04382.
Each surface is sampled with fixed n_vertices vertices and n_faces
faces in \(\mathbb{R}^3\).
Each individual surface is represented by a 2d-array of shape
[n_vertices, 3]. This space corresponds to the space of immersions
defined below, i.e. the
space of smooth functions from a template to manifold \(M\)
into \(\mathbb{R}^3\), with non-vanishing Jacobian.
\[Imm(M,\mathbb{R}^3)=\{ f \in C^{\infty}(M, \mathbb{R}^3)
\|Df(x)\|\neq 0 \forall x \in M \}.\]
Parameters:
faces (integer array-like, shape=[n_faces, 3]) – Triangulation of the surface.
Each face is given by 3 indices that indicate its vertices.
Compute the areas for each face of a triangulated surface.
The corresponds to the volume area for the surface metric, that is
the volume area of the pullback metric of the immersion defining the
surface metric.
Parameters:
point (array-like, shape=[…, n_vertices, 3]) – Surface, as the 3D coordinates of the vertices of its triangulation.
Returns:
_ (array-like, shape=[…, n_faces,]) – Area computed at each face of the triangulated surface.
Compute the mesh Laplacian operator of a triangulated surface.
Denoting q the surface, i.e. the point in the DiscreteSurfaces manifold,
the laplacian at \(q\) is defined as the operator:
\(\Delta_q = - Tr(g_q^{-1} \nabla^2)\)
where \(g_q\) is the surface metric matrix of \(q\).
The area of the triangles is computed using Heron’s formula.
Parameters:
point (array-like, shape=[…, n_vertices, 3]) – Surface, as the 3D coordinates of the vertices of its triangulation.
Returns:
_laplacian (callable) – Function that evaluates the mesh Laplacian operator at a
tangent vector field to the surface.
The matrices are evaluated at the faces of a triangulated surface.
The surface metric is the pullback metric of the immersion q
defining the surface, i.e. of
the map q: M -> R3, where M is the parameterization manifold.
Parameters:
point (array like, shape=[…, n_vertices, 3]) – Surface, as the 3D coordinates of the vertices of its triangulation.
Returns:
metric_mats (array-like, shape=[…, n_faces, 2, 2]) – Surface metric matrices evaluated at each face of
the triangulated surface.
Vertex area is the area of all of the triangles who are in contact (incident)
with a specific vertex, according to the formula:
vertex_areas = 2 * sum_incident_areas / 3.0
Parameters:
point (array-like, shape=[…, n_vertices, 3]) – Surface, as the 3D coordinates of the vertices of its triangulation.
Returns:
vertex_areas (array-like, shape=[…, n_vertices]) – Vertex area for each vertex.
Elastic metric defined by a family of second order Sobolev metrics.
Each individual discrete surface is represented by a 2D-array of shape
[n_vertices, 3]. See [HSKCB2022] and [HPBDC2023] (appendix) for details.
The parameters a0, a1, b1, c1, d1, a2 (detailed below) are non-negative weighting
coefficients for the different terms in the metric.
Parameters:
space (DiscreteSurfaces) – Instantiated DiscreteSurfaces manifold.
a0 (float) – First order parameter.
Default: 1.
a1 (float) – Stretching parameter.
Default: 1.
b1 (float) – Shearing parameter.
Default: 1.
c1 (float) – Bending parameter.
Default: 1.
d1 (float) – Additonal first order parameter.
Default: 1.
a2 (float) – Second order parameter.
Default: 1.
References
[HSKCB2022]
Emmanuel Hartman, Yashil Sukurdeep, Eric Klassen,
Nicolas Charon, and Martin Bauer.
“Elastic shape analysis of surfaces with second-order Sobolev metrics:
a comprehensive numerical framework”. arXiv:2204.04238 [cs.CV], 25 Sep 2022
[HPBDC2023]
Emmanuel Hartman, Emery Pierson, Martin Bauer,
Mohamed Daoudi, and Nicolas Charon.
“Basis Restricted Elastic Shape Analysis on the Space of Unregistered Surfaces.”
arXiv, November 7, 2023. https://doi.org/10.48550/arXiv.2311.04382.
Adds base_manifold (manifold that is repeated n_copies
times) and n_copies to space as attribute if it does
not have them as their expected by NFoldMetric.
Class to solve the geodesic boundary value problem with path-straightening.
Parameters:
space (Manifold) – Equipped manifold.
n_nodes (int) – Number of path discretization points.
lambda_ (float) – Discrepancy loss weight.
discrepancy_loss (callable) – A generic discrepancy term. Receives point and outputs a callable
which receives another point and outputs a scalar measuring
discrepancy between point and another point.
Scalar sums to energy of a path.
path_energy (callable) – Method to compute Riemannian path energy.
optimizer (ScipyMinimize) – An optimizer to solve path energy minimization problem.
initialization (callable or array-like, shape=[n_nodes - 2, n_vertices, 3]) – A method to get initial guess for optimization or an initial path.
References
[HSKCB2022]
“Elastic shape analysis of surfaces with second-order
Sobolev metrics: a comprehensive numerical framework”.
arXiv:2204.04238 [cs.CV], 25 Sep 2022
Generate parameterized function for the geodesic curve.
Geodesic curve defined by either:
an initial point and an initial tangent vector,
an initial point and an end point.
Parameters:
initial_point (array-like, shape=[…, dim]) – Point on the manifold, initial point of the geodesic.
end_point (array-like, shape=[…, dim], optional) – Point on the manifold, end point of the geodesic. If None,
an initial tangent vector must be given.
initial_tangent_vec (array-like, shape=[…, dim],) – Tangent vector at base point, the initial speed of the geodesics.
Optional, default: None.
If None, an end point must be given and a logarithm is computed.
Returns:
path (callable) – Time parameterized geodesic curve. If a batch of initial
conditions is passed, the output array’s first dimension
represents the different initial conditions, and the second
corresponds to time.
This is is the supremum of radii r for which the exponential map is a
diffeomorphism from the open ball of radius r centered at the base
point onto its image.
Parameters:
base_point (array-like, shape=[…, {dim, [n, m]}]) – Point on the manifold.
Compute derivative of the inner prod matrix at base point.
Writing \(g_{ij}\) the inner-product matrix at base point,
this computes \(mat_{ijk} = \partial_k g_{ij}\), where the
index k of the derivation is put last.
Parameters:
base_point (array-like, shape=[…, dim]) – Base point.
Returns:
metric_derivative (array-like, shape=[…, dim, dim, dim]) – Derivative of the inner-product matrix, where the index
k of the derivation is last: \(mat_{ijk} = \partial_k g_{ij}\).
Assumes total space is equipped with several group actions.
Aligns points wrt these group actions by alternate minimization
wrt each of them (see e.g. [SK2016] for more details).
This alignment results in an approximate quotient of the product
group action.
Parameters:
total_space (Manifold) – Manifold equipped with a quotient structure.
threshold (float) – Distance between consecutive aligned points for which
convergence is considered reached.
max_iter (int) – Maximum number of iterations.
verbose (boolean) – If log number of iterations need for convergence.
Lift a tangent vector to a horizontal vector in the total space.
It means that horizontal lift is the inverse of the restriction of the
tangent submersion to the horizontal space at point, where point must
be in the fiber above the base point. By default, the base manifold
is not explicit but is identified with a horizontal section of the
fiber bundle, so the horizontal lift is the horizontal projection.
fiber_point (array-like, shape=[…, {ambient_dim, [n, m]}]) – Point of the total space.
Optional, default : None. The lift method is used to compute a
point at which to compute a tangent vector.
base_point (array-like, shape=[…, {base_dim, [n, m]}]) – Point of the base space.
Optional, default : None. In this case, point must be given,
and submersion is used to compute the base_point if needed.
Returns:
horizontal_lift (array-like, shape=[…, {total_space.dim, [n, m]}]) – Horizontal tangent vector to the total space at point.
Compute the fundamental tensor A of the submersion.
The fundamental integrability tensor A is defined for tangent vectors
\(X = tangent\_vec\_a\) and \(Y = tangent\_vec\_b\) of the
total space by [ONeill] as
\(A_X Y = ver\nabla_{hor X} (hor Y) + hor \nabla_{hor X}( ver Y)\)
where \(hor, ver\) are the horizontal and vertical projections
and \(\nabla\) is the connection of the total space.
base_point (array-like, shape=[…, {total_space.dim, [n, m]}]) – Point of the total space.
Returns:
vector (array-like, shape=[…, {total_space.dim, [n, m]}]) – Tangent vector at base_point, result of the A tensor applied to
tangent_vec_a and tangent_vec_b.
References
[ONeill]
O’Neill, Barrett. The Fundamental Equations of a
Submersion, Michigan Mathematical Journal 13, no. 4
(December 1966): 459–69. https://doi.org/10.1307/mmj/1028999604.
Compute the covariant derivative of the integrability tensor A.
The covariant derivative \(\nabla_X (A_Y E)\) in total space is
necessary to compute the covariant derivative of the directional
curvature in a submersion. The components \(\nabla_X (A_Y E)\)
and \(A_Y E\) are computed at base-point \(P = base\_point\)
for horizontal vector fields \(X, Y\) extending the values
given in argument \(X|_P = horizontal\_vec\_x\),
\(Y|_P = horizontal\_vec\_y\) and a general vector field
\(E\) extending \(E|_x = tangent\_vec\_e\)
in a neighborhood of x with covariant derivatives
\(\nabla_X Y |_P = nabla_x y\) and
\(\nabla_X E |_P = nabla_x e\).
This is a section of the fiber bundle, defined on the base manifold,
with values in the total space. It means that submersion applied after
lift results in the identity map. By default, the base manifold
is not explicit but is identified with a section of the fiber bundle,
so the lift is the identity map.
Parameters:
point (array-like, shape=[…, {base_dim, [n, m]}]) – Point of the base manifold.
Returns:
lift (array-like, shape=[…, {total_space.dim, [n, m]}]) – Point of the total space.
This is the projection of the fiber bundle, defined on the total
space, with values in the base manifold. This map is surjective.
By default, the base manifold is not explicit but is identified with a
local section of the fiber bundle, so the submersion is the identity
map.
Parameters:
point (array-like, shape=[…, {total_space.dim, [n, m]}]) – Point of the total space.
Returns:
projection (array-like, shape=[…, {base_dim, [n, m]}]) – Point of the base manifold.
This is the differential of the projection of the fiber bundle,
defined on the tangent space of a point of the total space,
with values in the tangent space of the projection of this point in the
base manifold. This map is surjective. By default, the base manifold
is not explicit but is identified with a horizontal section of the
fiber bundle, so the tangent submersion is the horizontal projection.
Parameters:
tangent_vec (array-like, shape=[…, {total_space.dim, [n, m]}]) – Tangent vector to the total space at base_point.
base_point (array-like, shape=[…, {total_space.dim, [n, m]}]) – Point of the total space.
Returns:
projection (array-like, shape=[…, {base_dim, [n, m]}]) – Tangent vector to the base manifold.
Fiber bundle to construct the quotient metric on correlation matrices.
Correlation matrices are obtained as the quotient of the space of SPD
matrices by the action by congruence of diagonal matrices.
References
[TP21]
Thanwerdas, Yann, and Xavier Pennec. “Geodesics and Curvature of
the Quotient-Affine Metrics on Full-Rank CorrelationMatrices.”
In Proceedings of Geometric Science of Information.
Paris, France, 2021.
https://hal.archives-ouvertes.fr/hal-03157992.
Class for the quotient of the affine-invariant metric.
The affine-invariant metric on SPD matrices is invariant under the
action of diagonal matrices, thus it induces a quotient metric on the
manifold of full-rank correlation matrices.
Compute the correlation matrix associated to an SPD matrix.
The correlation matrix associated to an SPD matrix (the covariance)
\(\Sigma\) is given by \(D \Sigma D\) where \(D\) is
the inverse square-root of the diagonal of \(\Sigma\).
Parameters:
point (array-like, shape=[…, n, n]) – Symmetric Positive definite matrix.
Returns:
corr (array-like, shape=[…, n, n]) – Correlation matrix obtained by dividing all elements by the
diagonal entries.
Diffeormorphism between full-rank correlation matrices and
the space of lower triangular matrices with positive diagonal
and unit normed rows.
Since this image space is also diffeomorphic to another space, the
product space of successively increasing factor-dimension open hemispheres,
we take advantage of ComposedDiffeo to avoid explicitly representing
the intermediate space.
Find unique positive diagonal matrix corresponding to an SPD matrix.
That is, for every symmetric positive-definite matrix \(\Sigma\),
there exists a unique positive diagonal matrix \(\Delta\) such that
\(\Delta \Sigma \Delta\) is a symmetric positive-definite matrix
with unit row sums.
This result is known as the existence and uniqueness of the scaling of
SPD matrices ([T2024], [MO1968], [JR2009]).
It finds the roots of the gradient of the strictly convex map
\[F(D) = \frac{1}{2} \mathbb{1}^{\top} D^{\top}
\Sigma D \mathbb{1}-\operatorname{tr}(\log (D))\]
To ensure convergence, it is better to use a damped Newton method, which
leverages function being standard self-concordant and strictly convex
(check out section 5.1.5 of [N2018]).
Parameters:
root_finder (RootFinder)
References
[T2024]
Thanwerdas, Yann. “Permutation-Invariant Log-Euclidean Geometries
on Full-Rank Correlation Matrices.”
SIAM Journal on Matrix Analysis and Applications, 2024, 930–53.
https://doi.org/10.1137/22M1538144.
Marshall, Albert W., and Ingram Olkin.
“Scaling of Matrices to Achieve Specified Row and Column Sums.”
Numerische Mathematik 12, no. 1 (August 1, 1968): 83–90.
https://doi.org/10.1007/BF02170999.
Johnson, Charles R., and Robert Reams.
“Scaling of Symmetric Matrices by Positive Diagonal Congruence.”
Linear and Multilinear Algebra 57, no. 2 (March 1, 2009): 123–40.
https://doi.org/10.1080/03081080600872327.
Project a matrix to the set of full rank matrices.
As the space of full rank matrices is dense in the space of matrices,
this is not a projection per se, but a regularization if the matrix input X
is not already full rank: \(X + \epsilon [I_{rank}, 0]\) is returned
where \(\epsilon=gs.atol\)
Parameters:
point (array-like, shape=[…, n, k]) – Point in embedding manifold.
Returns:
projected (array-like, shape=[…, n, k]) – Projected point.
Class for the general linear group GL(n) and its identity component.
If positive_det=True, this is the connected component of the identity,
i.e. the space of matrices with positive determinant.
Parameters:
n (int) – Integer representing the shape of the matrices: n x n.
positive_det (bool) – Whether to restrict to the identity connected component of the
general linear group, i.e. matrices with positive determinant.
Optional, default: False.
As GL(n) is dense in the space of matrices, this is not a projection
per se, but a regularization if the matrix is not already invertible:
\(X + \epsilon I_n\) is returned where \(\epsilon=gs.atol\)
is returned for an input X.
Parameters:
point (array-like, shape=[…, dim_embedding]) – Point in embedding manifold.
The Grassmannian \(Gr(n, k)\) is the manifold of k-dimensional
subspaces in n-dimensional Euclidean space.
Lead author: Olivier Peltre.
\(Gr(n, p)\) is represented by
\(n \times n\) matrices
of rank \(p\) satisfying \(P^2 = P\) and \(P^T = P\).
Each \(P \in Gr(n, p)\) is identified with the unique
orthogonal projector onto \({\rm Im}(P)\).
\(Gr(n, p)\) is a homogoneous space, quotient of the special orthogonal
group by the subgroup of rotations stabilising a \(p\)-dimensional subspace:
It is therefore customary to represent the Grassmannian
by equivalence classes of orthogonal \(p\)-frames in \({\mathbb R}^n\).
For such a representation, work in the Stiefel manifold instead.
\[Gr(n, p) \simeq St(n, p) / SO(p)\]
References
[Batzies15]
Batzies, E., K. Hüper, L. Machado, and F. Silva Leite.
“Geometric Mean and Geodesic Regression on Grassmannians.”
Linear Algebra and Its Applications 466 (February 1, 2015):
83–101. https://doi.org/10.1016/j.laa.2014.10.003.
Following [Chikuse03], \(n\_samples * n * p\) scalars are sampled
from a standard normal distribution and reshaped to matrices,
the projectors on their first \(p\) columns follow a uniform
distribution.
Parameters:
n_samples (int) – The number of points to sample
Optional. default: 1.
Returns:
projectors (array-like, shape=[…, n, n]) – Points following a uniform distribution.
References
[Chikuse03]
Yasuko Chikuse, Statistics on special manifolds,
New York: Springer-Verlag. 2003, 10.1007/978-0-387-21540-2
Following [Chikuse03], \(n\_samples * n * p\) scalars are sampled
from a standard normal distribution and reshaped to matrices,
the projectors on their first \(p\) columns follow a uniform
distribution.
Parameters:
n_samples (int) – The number of points to sample
Optional. default: 1.
Returns:
projectors (array-like, shape=[…, n, n]) – Points following a uniform distribution.
References
[Chikuse03]
Yasuko Chikuse, Statistics on special manifolds,
New York: Springer-Verlag. 2003, 10.1007/978-0-387-21540-2
The Grassmann manifold is defined here as embedded in the set of
symmetric matrices, as the pre-image of the function defined around the
projector on the space spanned by the first \(p\) columns of the identity
matrix by (see Exercise E.25 in [Pau07]).
\[\begin{split}\begin{pmatrix} I_p + A & B^T \\ B & D \end{pmatrix} \mapsto
(D - B(I_p + A)^{-1}B^T, A + A^2 + B^TB)\end{split}\]
This map is a submersion and its zero space is the set of orthogonal
rank-\(p\) projectors.
This is is the supremum of radii r for which the exponential map is a
diffeomorphism from the open ball of radius r centered at the base
point onto its image.
In this case it is Pi / 2 everywhere.
Parameters:
base_point (array-like, shape=[…, n, n]) – Point on the manifold.
Compute the parallel transport of a tangent vector.
Closed-form solution for the parallel transport of a tangent vector
along the geodesic between two points base_point and end_point
or alternatively defined by \(t \mapsto exp_{(base\_point)}(
t*direction)\).
Parameters:
tangent_vec (array-like, shape=[…, n, n]) – Tangent vector at base point to be transported.
base_point (array-like, shape=[…, n, n]) – Point on the Grassmann manifold. Point to transport from.
direction (array-like, shape=[…, n, n]) – Tangent vector at base point, along which the parallel transport
is computed.
Optional, default: None
end_point (array-like, shape=[…, n, n]) – Point on the Grassmann manifold to transport to. Unused if
direction is given.
Optional, default: None
Returns:
transported_tangent_vec (array-like, shape=[…, n, n]) – Transported tangent vector at exp_(base_point)(direction).
References
[BZA20]
Bendokat, Thomas, Ralf Zimmermann, and P.-A. Absil.
“A Grassmann Manifold Handbook: Basic Geometry and Computational
Aspects.” ArXiv:2011.13699 [Cs, Math], November 27, 2020.
https://arxiv.org/abs/2011.13699.
lie_algebra_repr (bool) – If True, group elements are represented by the vector
representation of Lie algebra elements.
This is amenable to perform gradient-based optimizations
on the Lie algebra.
n (int) – Integer representing the shapes of the matrices : n x n.
lie_algebra_repr (bool) – If True, group elements are represented by the vector
representation of Lie algebra elements.
This is amenable to perform gradient-based optimizations
on the Lie algebra.
Class for the 3D Heisenberg group in the vector representation.
The 3D Heisenberg group represented as R^3. It is a step-2 Carnot Lie
group. It can be equipped with a natural sub-Riemannian structure, and it
is a fundamental example in sub-Riemannian geometry.
Compute the Jacobian matrix of left/right translation by a point.
This calculates the differential of the left translation L_(point)
evaluated at point. Note that it only depends on the point we are
left-translating by, not on the point where the differential is
evaluated.
Parameters:
point (array-like, shape=[…, 3]) – Point.
left (bool) – Indicate whether to calculate the differential of left or right
translations.
Optional, default: True.
Returns:
_ (array-like, shape=[…, 3, 3]) – Jacobian of the left/right translation by point.
The manifold of Hermitian positive definite (HPD) matrices.
Lead author: Yann Cabanes.
References
[Cabanes2022]
Yann Cabanes. Multidimensional complex stationary
centered Gaussian autoregressive time series machine learning
in Poincaré and Siegel disks: application for audio and radar
clutter classification, PhD thesis, 2022
Compute the Riemannian exponential at point base_point
of tangent vector tangent_vec wrt the metric defined in inner_product.
This gives a Hermitian positive definite matrix.
Parameters:
tangent_vec (array-like, shape=[…, n, n]) – Tangent vector at base point.
base_point (array-like, shape=[…, n, n]) – Base point.
Returns:
exp (array-like, shape=[…, n, n]) – Riemannian exponential.
Compute the Riemannian logarithm at point base_point,
of point wrt the metric defined in inner_product.
This gives a tangent vector at point base_point.
Parameters:
point (array-like, shape=[…, n, n]) – Point.
base_point (array-like, shape=[…, n, n]) – Base point.
Returns:
log (array-like, shape=[…, n, n]) – Riemannian logarithm of point at base_point.
Closed-form solution for the parallel transport of a tangent vector
along the geodesic between two points base_point and end_point
or alternatively defined by \(t \mapsto exp_{(base\_point)}(
t*direction)\).
Denoting tangent_vec_a by S, base_point by A, and end_point
by B or B = Exp_A(tangent_vec_b) and \(E = (BA^{- 1})^{( 1
/ 2)}\). Then the parallel transport to B is:
\[S' = ESE^T\]
Parameters:
tangent_vec (array-like, shape=[…, n, n]) – Tangent vector at base point to be transported.
base_point (array-like, shape=[…, n, n]) – Point on the manifold of HPD matrices. Point to transport from
direction (array-like, shape=[…, n, n]) – Tangent vector at base point, initial speed of the geodesic along
which the parallel transport is computed. Unused if end_point is given.
Optional, default: None.
end_point (array-like, shape=[…, n, n]) – Point on the manifold of HPD matrices. Point to transport to.
Optional, default: None.
Returns:
transported_tangent_vec (array-like, shape=[…, n, n]) – Transported tangent vector at exp_(base_point)(tangent_vec_b).
Compute the inner-product of tangent_vec_a \(A\) and tangent_vec_b
\(B\) at point base_point \(S=PDP^\top\) using the
Bures-Wasserstein Riemannian metric:
Compute the parallel transport of a tangent vec along a geodesic.
Approximation of the solution of the parallel transport of a tangent
vector a along the geodesic defined by \(t \mapsto exp_{(
base\_point)}(t* tangent\_vec\_b)\). The parallel transport equation is
formulated in this case in [TP2021].
Parameters:
tangent_vec (array-like, shape=[…, n, n]) – Tangent vector at base_point to transport.
base_point (array-like, shape=[…, n, n]) – Initial point of the geodesic.
direction (array-like, shape=[…, n, n]) – Tangent vector at base_point, initial velocity of the geodesic to
transport along.
end_point (array-like, shape=[…, n, n]) – Point to transport to.
Optional, default: None.
n_steps (int) – Number of steps to use to approximate the solution of the
ordinary differential equation.
Optional, default: 100
step (str, {‘euler’, ‘rk2’, ‘rk4’}) – Scheme to use in the integration scheme.
Optional, default: ‘rk4’.
Returns:
transported (array-like, shape=[…, n, n]) – Transported tangent vector at exp_(base_point)(tangent_vec_b).
This class is a common interface to the different models of hyperbolic
geometry:
the hyperboloid, embedded in Minkowski space of dimension dim + 1 as the set of
points whose squared norm is equal to -1. This representation is called
extrinsic here.
the Poincare ball, the open ball of the Euclidean space of dimension dim.
the Poincare half-space, the open space of points of the Euclidean
space of dimension dim, whose last coordinate is positive.
Parameters:
dim (int) – Dimension of the hyperbolic space.
coords_type (str, {‘extrinsic’, ‘ball’, ‘half-space’}) – Default coordinates to represent points in hyperbolic space.
Optional, default: ‘extrinsic’.
Class for the n-dimensional hyperboloid space as embedded in (n+1)-dimensional
Minkowski space as the set of points with squared norm equal to -1, i.e.
\[- x_0^2 + x_1^2 + ... + x_n^2 = - 1.\]
For other
representations of hyperbolic spaces see the Hyperbolic class.
point_extrinsic (array-like, shape=[…, dim + 1]) – Point in the embedded manifold in extrinsic coordinates,
i. e. in the coordinates of the embedding manifold.
Returns:
point_intrinsic (array-like, shape=[…, dim]) – Point in intrinsic coordinates.
Project point onto geodesic going through base point in direction
of tangent vector. See reference below.
Parameters:
point (array-like, shape=[…, dim + 1]) – Point in hyperbolic space.
base_point (array-like, shape=[…, dim + 1]) – Point through which the geodesic passes.
tangent_vec (array-like, shape=[…, dim + 1]) – Tangent vector in Minkowski space, direction of the geodesic.
Returns:
proj (array-like, shape=[…, dim + 1]) – Projected point on the geodesic.
References
[CSV2016]
R. Chakraborty, D. Seo, and B. C. Vemuri,
“An efficient exact-pga algorithm for constant curvature manifolds.”
Proceedings of the IEEE conference on computer vision and pattern
recognition. 2016.
This is is the supremum of radii r for which the exponential map is a
diffeomorphism from the open ball of radius r centered at the base
point onto its image.
In the case of the hyperbolic space, it does not depend on the base
point and is infinite everywhere, because of the negative curvature.
Parameters:
base_point (array-like, shape=[…, dim+1]) – Point on the manifold.
Compute Riemannian logarithm of a point wrt a base point.
If coords_type is poincare then base_point belongs
to the Poincare ball and point is a vector in the Euclidean
space of the same dimension as the ball.
Parameters:
point (array-like, shape=[…, dim + 1]) – Point in hyperbolic space.
base_point (array-like, shape=[…, dim + 1]) – Point in hyperbolic space.
Returns:
log (array-like, shape=[…, dim + 1]) – Tangent vector at the base point equal to the Riemannian logarithm
of point at the base point.
Compute the parallel transport of a tangent vector.
Closed-form solution for the parallel transport of a tangent vector
along the geodesic between two points base_point and end_point
or alternatively defined by \(t \mapsto exp_{(base\_point)}(
t*direction)\).
Parameters:
tangent_vec (array-like, shape=[…, dim + 1]) – Tangent vector at base point to be transported.
base_point (array-like, shape=[…, dim + 1]) – Point on the hyperboloid.
direction (array-like, shape=[…, dim + 1]) – Tangent vector at base point, along which the parallel transport
is computed.
Optional, default : None.
end_point (array-like, shape=[…, dim + 1]) – Point on the hyperboloid. Point to transport to. Unused if tangent_vec_b
is given.
Optional, default : None.
Returns:
transported_tangent_vec (array-like, shape=[…, dim + 1]) – Transported tangent vector at exp_(base_point)(tangent_vec_b).
Class for the n-dimensional hypersphere embedded in the
(n+1)-dimensional Euclidean space.
By default, points are parameterized by their extrinsic
(n+1)-coordinates. For dimensions 1 and 2, this can be changed with the
intrinsic parameter. For dimensions 1 (the circle),
the intrinsic coordinates correspond angles in radians, with 0. mapping
to point [1., 0.]. For dimension 2, the intrinsic coordinates are the
spherical coordinates from the north pole, i.e. where angles [0.,
0.] correspond to point [0., 0., 1.].
Parameters:
dim (int) – Dimension of the hypersphere.
intrinsic (bool) – Type of representation for dimensions 1 and 2.
For three tangent vectors at a base point \(x,y,z\),
the curvature is defined by
\(R(x, y)z = \nabla_{[x,y]}z
- \nabla_z\nabla_y z + \nabla_y\nabla_x z\), where \(\nabla\)
is the Levi-Civita connection. In the case of the hypersphere,
we have the closed formula
\(R(x,y)z = \langle y,z \rangle x - \langle x, z \rangle y\).
Parameters:
tangent_vec_a (array-like, shape=[…, dim]) – Tangent vector at base_point.
tangent_vec_b (array-like, shape=[…, dim]) – Tangent vector at base_point.
tangent_vec_c (array-like, shape=[…, dim]) – Tangent vector at base_point.
base_point (array-like, shape=[…, dim]) – Point on the hypersphere.
Returns:
curvature (array-like, shape=[…, dim]) – Tangent vector at base_point.
This is is the supremum of radii r for which the exponential map is a
diffeomorphism from the open ball of radius r centered at the base
point onto its image.
In the case of the sphere, it does not depend on the base point and is
Pi everywhere.
Parameters:
base_point (array-like, shape=[…, dim+1]) – Point on the manifold.
Compute the parallel transport of a tangent vector.
Closed-form solution for the parallel transport of a tangent vector
along the geodesic between two points base_point and end_point
or alternatively defined by \(t \mapsto exp_{(base\_point)}(
t*direction)\).
Parameters:
tangent_vec (array-like, shape=[…, dim + 1]) – Tangent vector at base point to be transported.
base_point (array-like, shape=[…, dim + 1]) – Point on the hypersphere.
direction (array-like, shape=[…, dim + 1]) – Tangent vector at base point, along which the parallel transport
is computed.
Optional, default : None.
end_point (array-like, shape=[…, dim + 1]) – Point on the hypersphere. Point to transport to. Unused if
tangent_vec_b is given.
Optional, default : None.
Returns:
transported_tangent_vec (array-like, shape=[…, dim + 1]) – Transported tangent vector at
end_point=exp_(base_point)(tangent_vec_b).
Class for bi-invariant metrics on compact Lie groups.
Compact Lie groups and direct products of compact Lie groups with vector
spaces admit bi-invariant metrics [Gallier]. Products Lie groups are not
implemented. Other groups such as SE(3) admit bi-invariant pseudo-metrics.
Parameters:
space (LieGroup) – The group to equip with the bi-invariant metric
References
[Gallier]
Gallier, Jean, and Jocelyn Quaintance. Differential Geometry
and Lie Groups: A Computational Perspective.
Geonger International Publishing, 2020.
https://doi.org/10.1007/978-3-030-46040-2.
Compute Riemannian exponential of tangent vector from the identity.
For a bi-invariant metric, this corresponds to the group exponential.
Parameters:
tangent_vec – Tangent vector at identity.
base_point (array-like, shape=[…, {dim, [n, n]}]) – Point in the group.
Optional, default : identity.
Returns:
exp (array-like, shape=[…, {dim, [n, n]}]) – Point in the group.
References
[Gallier]
Gallier, Jean, and Jocelyn Quaintance. Differential
Geometry and Lie Groups: A Computational Perspective.
Geonger International Publishing, 2020.
https://doi.org/10.1007/978-3-030-46040-2.
This is is the supremum of radii r for which the exponential map is a
diffeomorphism from the open ball of radius r centered at the base
point onto its image.
In the case of a bi-invariant metric, it does not depend on the base
point.
Parameters:
base_point (array-like, shape=[…, n, n]) – Point on the manifold.
Compute Riemannian logarithm of a point wrt the identity.
For a bi-invariant metric this corresponds to the group logarithm.
Parameters:
point (array-like, shape=[…, {dim, [n, n]}]) – Point in the group.
base_point (array-like, shape=[…, {dim, [n, n]}]) – Point in the group.
Optional, default : identity.
Returns:
log (array-like, shape=[…, {dim, [n, n]}]) – Tangent vector at the identity equal to the Riemannian logarithm
of point at the identity.
References
[Gallier]
Gallier, Jean, and Jocelyn Quaintance. Differential
Geometry and Lie Groups: A Computational Perspective.
Geonger International Publishing, 2020.
https://doi.org/10.1007/978-3-030-46040-2.
Compute the parallel transport of a tangent vec along a geodesic.
Closed-form solution for the parallel transport of a tangent vector a
along the geodesic between the base point and an end point, or alternatively
defined by \(t \mapsto exp_{(base\_point)}(t*direction)\).
As a compact Lie group endowed with its canonical bi-invariant metric is a
symmetric space, parallel transport is achieved by a geodesic symmetry, or
equivalently, one step of the pole ladder scheme.
Parameters:
tangent_vec (array-like, shape=[…, n, n]) – Tangent vector at base point to be transported.
base_point (array-like, shape=[…, n, n]) – Point on the manifold.
direction (array-like, shape=[…, n, n]) – Tangent vector at base point, along which the parallel transport
is computed.
Optional, default: None.
end_point (array-like, shape=[…, n, n]) – Point on the manifold. Point to transport to.
Optional, default: None.
Returns:
transported_tangent_vec (array-like, shape=[…, n, n]) – Transported tangent vector at
end_point=exp_(base_point)(tangent_vec_b).
where \(ad^*\) is the dual adjoint map with respect to the
metric. For a right-invariant metric, \(dR\) is used instead of
\(dL\) and \(ad^*\) is replaced by \(-ad^*\). The
exponential map is approximated by numerical integration
of this equation, with initial conditions \(\dot{\gamma}(0)\)
given by the argument tangent_vec and \(\gamma(0)\) by
base_point.
Parameters:
tangent_vec (array-like, shape=[…, n, n]) – Tangent vector at a base point.
base_point (array-like, shape=[…, n, n]) – Point in the group.
Optional, defaults to identity if None.
Returns:
exp (array-like, shape=[…, n, n]) – Point in the group equal to the Riemannian exponential
of tangent_vec at the base point.
Kolev, Boris. “Lie Groups and Mechanics: An Introduction.”
Journal of Nonlinear Mathematical Physics 11, no. 4, 2004:
480–98. https://doi.org/10.2991/jnmp.2004.11.4.5.
Class for the Klein Bottle, a two dimensional manifold which is not orientable.
Points are understood to be on the Klein Bottle if they are in the interval
\([0,1]^2\).
Each point in \(\mathbb R^2\) can be understood as a point on the Klein Bottle
by considering the equivalence relation
\begin{align}
& (x_1,y_1) \sim (x_2, y_2) \\
\Leftrightarrow \quad & y_1-y_2 \in \mathbb Z
\text{ and } x_1-x_2 \in 2\mathbb Z \\
\text{ or } \quad &
y_1 + y_2 \in \mathbb Z \text{ and } x_1-x_2 \in \mathbb Z\setminus 2\mathbb Z.
\end{align}
Evaluate whether two points represent equivalent points on the Klein Bottle.
This method uses the equivalence stated in the class description.
This means points are equivalent if one walks an even number of squares in
x-direction and any number of squares in y-direction or an uneven number of
squares in x-direction and any number of squares in y-direction while
“mirroring” the y coordinate.
Parameters:
point_a (array-like, shape=[…, 2]) – Representation of point on Klein Bottle
point_b (array-like, shape=[…, 2]) – Representation of point on Klein Bottle
atol (Absolute tolerance to test for belonging to mathbb Z.)
Returns:
is_equivalent (array-like, shape=[…, 2]) – Boolean evaluating if points are equivalent
Convert point to coordinates in R^3 parametrizing the Klein bagel.
Convert from the intrinsic coordinates in the Klein bottle (2 parameters),
to the coordinates of the Klein bagel parametrization in 3d Euclidean space
(3 parameters).
For intrinsic parameters (theta,v) the Klein bagel parametrization is
[https://en.wikipedia.org/wiki/Klein_bottle#The_figure_8_immersion]:
x = left(r + cos(theta/2)sin(v) - sin(theta/2)sin(2v)right)cos(theta)
y = left(r + cos(theta/2)sin(v) - sin(theta/2)sin(2v)right)sin(theta)
z = sin(theta/2)sin(v) + cos(theta/2)sin(2v)
for 0leqtheta<2pi and 0leq v<2pi. r is a constant to determine
the aspect ratio.
Parameters:
point (array-like, shape=[…, 2]) – Point on the Klein bottle, in intrinsic coordinates.
radius (float) – Radius of the circle.
Returns:
point_extrinsic (array-like, shape=[…, 4]) – Point on the Klein bottle, in the Klein bagel parametrization.
Convert point to coordinates in R^3 parametrizing the Klein bottle.
Convert from the intrinsic coordinates in the Klein bottle (2 parameters),
to the coordinates of the Klein bottle parametrization in 3d Euclidean space
(3 parameters).
For intrinsic parameters (theta,v) the Klein bottle parametrization is
[https://en.wikipedia.org/wiki/Klein_bottle#Bottle_shape]
Parameters:
point (array-like, shape=[…, 2]) – Point on the Klein bottle, in intrinsic coordinates.
Returns:
point_extrinsic (array-like, shape=[…, 4]) – Point on the Klein bottle, in the Klein bagel parametrization.
Convert from the intrinsic coordinates in the Klein bottle (2 parameters),
to the extrinsic coordinates in Euclidean space (4 parameters).
For intrinsic parameters (theta,v) the extrinsic Euclidean parametrization is
[https://en.wikipedia.org/wiki/Klein_bottle#4-D_non-intersecting]:
x = R left(cos(theta/2)cos(v) - sin(theta/2)sin(2v)right)
y = R left(sin(theta/2)cos(v) + cos(theta/2)sin(2v)right)
z = P cos(theta)(1+epsilon sin(v))
w = P sin(theta)(1+epsilon sin(v))
for 0leqtheta<2pi and 0leq v<2pi. P and R are constants to determine
the aspect ratio.
ε is any small constant .
Parameters:
point (array-like, shape=[…, 2]) – Point on the Klein bottle, in intrinsic coordinates.
Returns:
point_extrinsic (array-like, shape=[…, 4]) – Point on the Klein bottle, in extrinsic coordinates in
Euclidean space.
This is is the supremum of radii r for which the exponential map is a
diffeomorphism from the open ball of radius r centered at the base
point onto its image.
Parameters:
base_point (array-like, shape=[…, 2]) – Point on the manifold. Unused.
The landmark space is a product manifold where all manifolds in the
product are the same. The default metric is the product metric and
is often referred to as the L2 metric.
Parameters:
ambient_manifold (Manifold) – Manifold to which landmarks belong.
Module providing an implementation of MatrixLieAlgebras.
There are two main forms of representation for elements of a MatrixLieAlgebra
implemented here. The first one is as a matrix, as elements of R^(n x n).
The second is by choosing a base and remembering the coefficients of an element
in that base. This base will be provided in child classes
(e.g. SkewSymmetricMatrices).
Calculate the Baker-Campbell-Hausdorff approximation of given order.
The implementation is based on [CM2009a] with the pre-computed
constants taken from [CM2009b]. Our coefficients are truncated to
enable us to calculate BCH up to order 15.
This represents Z = log(exp(X)exp(Y)) as an infinite linear combination
of the form Z = sum z_i e_i where z_i are rational numbers and e_i are
iterated Lie brackets starting with e_1 = X, e_2 = Y, each e_i is given
by some i’,i’’: e_i = [e_i’, e_i’’].
order (int) – The order to which the approximation is calculated. Note that this
is NOT the same as using only e_i with i < order.
Optional, default 2.
F. Casas and A. Murua. An efficient algorithm for
computing the Baker–Campbell–Hausdorff series and some of its
applications. Journal of Mathematical Physics 50, 2009
In this class, point_type (‘vector’ or ‘matrix’) will be used to describe
the format of the points on the Lie group.
If point_type is ‘vector’, the format of the inputs is
[…, dimension], where dimension is the dimension of the Lie group.
If point_type is ‘matrix’ the format of the inputs is
[…, n, n] where n is the parameter of GL(n) e.g. the amount of rows
and columns of the matrix.
For matrix Lie groups with tangent vectors A,B at the same base point P
this is given by (translate to identity, compute commutator, go back)
\([A,B] = A_P^{-1}B - B_P^{-1}A\)
Parameters:
tangent_vec_a (array-like, shape=[…, n, n]) – Tangent vector at base point.
tangent_vec_b (array-like, shape=[…, n, n]) – Tangent vector at base point.
base_point (array-like, shape=[…, n, n]) – Base point.
Returns:
bracket (array-like, shape=[…, n, n]) – Lie bracket.
Compute the push-forward map by the left/right translation.
Compute the push-forward map, of the left/right translation by the
point. It corresponds to the tangent map, or differential of the
group multiplication by the point or its inverse. For groups with a
vector representation, it is only implemented at identity, but it can
be used at other points by passing inverse=True. This method wraps
the jacobian translation which actually computes the matrix
representation of the map.
Parameters:
point (array-like, shape=[…, {dim, [n, n]]) – Point.
left (bool) – Whether to calculate the differential of left or right
translations.
Optional, default: True.
inverse (bool,) – Whether to inverse the jacobian matrix. If True, the push forward
by the translation by the inverse of point is returned.
Optional, default: False.
Returns:
tangent_map (callable) – Tangent map of the left/right translation by point. It can be
applied to tangent vectors.
Exponentiate a left-invariant vector field from a base point.
The vector input is not an element of the Lie algebra, but of the
tangent space at base_point: if \(g\) denotes base_point,
\(v\) the tangent vector, and \(V = g^{-1} v\) the associated
Lie algebra vector, then
\[\exp(v, g) = mul(g, \exp(V))\]
Therefore, the Lie exponential is obtained when base_point is None, or
the identity.
Parameters:
tangent_vec (array-like, shape=[…, n, n]) – Tangent vector at base point.
base_point (array-like, shape=[…, n, n]) – Base point.
Optional, defaults to identity if None.
Returns:
point (array-like, shape=[…, n, n]) – Left multiplication of exp(algebra_mat) with base_point.
For matrix Lie groups with tangent vectors A,B at the same base point P
this is given by (translate to identity, compute commutator, go back)
\([A,B] = A_P^{-1}B - B_P^{-1}A\)
Parameters:
tangent_vec_a (array-like, shape=[…, n, n]) – Tangent vector at base point.
tangent_vec_b (array-like, shape=[…, n, n]) – Tangent vector at base point.
base_point (array-like, shape=[…, n, n]) – Base point.
Returns:
bracket (array-like, shape=[…, n, n]) – Lie bracket.
Compute the push-forward map by the left/right translation.
Compute the push-forward map, of the left/right translation by the
point. It corresponds to the tangent map, or differential of the
group multiplication by the point or its inverse. For groups with a
vector representation, it is only implemented at identity, but it can
be used at other points by passing inverse=True. This method wraps
the jacobian translation which actually computes the matrix
representation of the map.
Parameters:
point (array-like, shape=[…, {dim, [n, n]]) – Point.
left (bool) – Whether to calculate the differential of left or right
translations.
Optional, default: True
inverse (bool,) – Whether to inverse the jacobian matrix. If True, the push forward
by the translation by the inverse of point is returned.
Optional, default: False.
Returns:
tangent_map (callable) – Tangent map of the left/right translation by point. It can be
applied to tangent vectors.
This method is not recommended for statistical purposes,
as the tangent vectors generated are not drawn from a
distribution related to the Riemannian metric.
Parameters:
n_samples (int) – Number of samples.
Optional, default: 1.
Flatten a matrix (compatible with vectorization on data axis 0).
The reverse operation is reshape. These operations are often called
matrix vectorization / matricization in mathematics
(https://en.wikipedia.org/wiki/Tensor_reshaping).
The names flatten / reshape were chosen to avoid confusion with
vectorization on data axis 0.
Parameters:
mat (array-like, shape=[…, m, n]) – Matrix.
Returns:
vec (array-like, shape=[…, m * n]) – Flatten copy of mat.
Matricize a vector (compatible with vectorization on data axis 0).
The reverse operation is matrices.flatten. These operations are often
called matrix vectorization / matricization in mathematics
(https://en.wikipedia.org/wiki/Tensor_reshaping).
The names flatten / reshape were chosen to avoid confusion with
vectorization on data axis 0.
Parameters:
vec (array-like, shape=[…, m * n]) – Vector.
Returns:
mat (array-like, shape=[…, m, n]) – Matricized copy of vec.
This is is the supremum of radii r for which the exponential map is a
diffeomorphism from the open ball of radius r centered at the base
point onto its image.
In the case of the hyperbolic space, it does not depend on the base
point and is infinite everywhere, because of the negative curvature.
Parameters:
base_point (array-like, shape=[…, dim]) – Point on the manifold.
This is is the supremum of radii r for which the exponential map is a
diffeomorphism from the open ball of radius r centered at the base
point onto its image.
In the case of the hyperbolic space, it does not depend on the base
point and is infinite everywhere, because of the negative curvature.
Parameters:
base_point (array-like, shape=[…, dim]) – Point on the manifold.
The Poincare polydisk is defined as a product manifold of the Hyperbolic space
of dimension 2. The Poincare polydisk has a product metric. The metric on each
space is the natural Poincare metric multiplied by a constant.
Compute the Riemannian exponential at point base_point
of tangent vector tangent_vec wrt the Cholesky metric.
This gives a lower triangular matrix with positive elements.
Parameters:
tangent_vec (array-like, shape=[…, n, n]) – Tangent vector at base point.
base_point (array-like, shape=[…, n, n]) – Base point.
Returns:
exp (array-like, shape=[…, n, n]) – Riemannian exponential.
Manifold of lower triangular matrices with >0 diagonal.
This manifold is also called the cholesky space.
Parameters:
n (int) – Integer representing the shape of the matrices: n x n.
References
[TP2019]
Z Lin. “Riemannian Geometry of Symmetric
Positive Definite Matrices Via Cholesky Decomposition”
SIAM journal on Matrix Analysis and Applications , 2019.
https://arxiv.org/abs/1908.09326
A diffeomorphism from UnitNormedRowsPLTDMatrices to OpenHemispheresProduct.
A diffeomorphism between the space of lower triangular matrices with
positive diagonal and unit normed rows and the product space of
successively increasing factor-dimension open hemispheres.
Given the way OpenHemisphere is implemented, i.e. the first component
is positive, rows need to be flipped.
This manifold is a particular case in dimension 1
of the SPD matrices endowed with the SPD affine-invariant metric
with a power affine coefficient equal to one.
This manifold is used as part of a product manifold involving
complex Poincare disks to classify complex stationary centered
Gaussian autoregressive time series in [Cabanes2022], [Cabanes_CESAR_2019],
[Cabanes_GSI_2019] and [Cabanes_RADAR_2019].
This product manifold is called ProductPositiveRealsAndComplexPoincareDisks
in geomstats. It is a product manifold of complex type, therefore the manifold
PositiveReals is compatible with complex input points and vectors even if the
input values should be reals.
Lead author: Yann Cabanes.
References
[Cabanes2022]
Yann Cabanes. Multidimensional complex stationary
centered Gaussian autoregressive time series machine learning
in Poincaré and Siegel disks: application for audio and radar
clutter classification, PhD thesis, 2022
[Cabanes_CESAR_2019]
Yann Cabanes, Frédéric Barbaresco, Marc Arnaudon,
Jérémie Bigot. Unsupervised Machine Learning for Pathological Radar
Clutter Clustering: the P-Mean-Shift Algorithm, C&ESAR 2019, Nov 2019,
Rennes, France. hal-02875430
https://hal.archives-ouvertes.fr/hal-02875430/document
[Cabanes_GSI_2019]
Yann Cabanes, Frédéric Barbaresco, Marc Arnaudon,
Jérémie Bigot. Toeplitz Hermitian Positive Definite Matrix Machine Learning
based on Fisher Metric. Geometric Science of Information, 2019,
10.1007/978-3-030-26980-7_27. hal-02875403
https://hal.archives-ouvertes.fr/hal-02875403/document
[Cabanes_RADAR_2019]
Yann Cabanes, Frédéric Barbaresco, Marc Arnaudon,
Jérémie Bigot. Non-Supervised High Resolution Doppler Machine Learning
for Pathological Radar Clutter. RADAR 2019, Sep 2020, Toulon, France.
10.1109/RADAR41533.2019.171295. hal-02875415
https://hal.archives-ouvertes.fr/hal-02875415/document
Compute the covariant derivative of the directional curvature.
For two vectors fields \(X|_P =\) tangent_vec_a,
\(Y|_P =\) tangent_vec_b with tangent vector value specified in argument
at the base_point \(P\),
the covariant derivative (in the direction \(X\))
\((\nabla_X R_Y)(X) |_P = (\nabla_X R)(Y, X) Y |_P\) of the
directional curvature (in the direction \(Y\))
\(R_Y(X) = R(Y, X) Y\) is a quadratic tensor in \(X\) and \(Y\)
that plays an important role in the computation of the moments of the
empirical Fréchet mean [Pennec].
In more details, let \(X, Y\) be the horizontal lift of parallel
vector fields extending the tangent vectors given in argument by
parallel transport in a neighborhood of the`base_point` \(P\) in the
base-space. Such vector fields verify \(\nabla^T_X X=0\) and
\(\nabla^T_X Y = A_X Y\) using the connection \(\nabla^T\)
of the total space. Then the covariant derivative of the
directional curvature tensor is given by
\(\nabla_X (R_Y(X)) = hor \nabla^T_X (R^T_Y(X)) - A_X( ver R^T_Y(X))
- 3 (\nabla_X^T A_Y A_X Y - A_X A_Y A_X Y )\), where \(R^T_Y(X)\)
is the directional curvature tensor of the total space.
Parameters:
tangent_vec_a (array-like, shape=[…, k_landmarks, ambient_dim]) – Tangent vector at base_point.
tangent_vec_b (array-like, shape=[…, k_landmarks, ambient_dim]) – Tangent vector at base_point.
base_point (array-like, shape=[…, k_landmarks, ambient_dim]) – Point on the group.
Returns:
curvature_derivative (array-like, shape=[…, k_landmarks, ambient_dim]) – Tangent vector at base_point.
Compute the parallel transport of a tangent vec along a geodesic.
Approximation of the solution of the parallel transport of a tangent
vector a along the geodesic between two points base_point and
end_point or alternatively defined by
\(t\mapsto exp_{(base\_point)}(t*direction)\)
Parameters:
tangent_vec (array-like, shape=[…, k_landmarks, ambient_dim]) – Tangent vector at base_point to transport.
base_point (array-like, shape=[…, k_landmarks, ambient_dim]) – Initial point of the geodesic to transport along.
direction (array-like, shape=[…, k_landmarks, ambient_dim]) – Tangent vector ar base_point, initial velocity of the geodesic to
transport along.
Optional, default: None.
end_point (array-like, shape=[…, k_landmarks, ambient_dim]) – Point to transport to. Unused if tangent_vec_b is given.
Optional, default: None.
n_steps (int) – Number of steps to use to approximate the solution of the
ordinary differential equation.
Optional, default: 100.
step (str, {‘euler’, ‘rk2’, ‘rk4’}) – Scheme to use in the integration scheme.
Optional, default: ‘rk4’.
Guigui, Nicolas, Elodie Maignant, Alain Trouvé, and Xavier
Pennec. “Parallel Transport on Kendall Shape Spaces.”
5th conference on Geometric Science of Information,
Paris 2021. Lecture Notes in Computer Science.
Springer, 2021. https://hal.inria.fr/hal-03160677.
Compute the fundamental tensor A of the submersion.
The fundamental tensor A is defined for tangent vectors of the total
space by [ONeill]\(A_X Y = ver\nabla^M_{hor X}(hor Y) + hor \nabla^M_{hor X}(ver Y)\)
where \(hor, ver\) are the horizontal and vertical projections.
For the Kendall shape space, we have the closed-form expression at
base-point P [Pennec]:
\(A_X E = P Sylv_P(E^\top hor(X)) + F + <F,P> P\) where
\(F = hor(X) Sylv_P(P^\top E)\) and \(Sylv_P(B)\) is the
unique skew-symmetric matrix \(\Omega\) solution of
\(P^\top P \Omega + \Omega P^\top P = B - B^\top\).
Parameters:
tangent_vec_x (array-like, shape=[…, k_landmarks, ambient_dim]) – Tangent vector at base_point.
tangent_vec_e (array-like, shape=[…, k_landmarks, ambient_dim]) – Tangent vector at base_point.
base_point (array-like, shape=[…, k_landmarks, ambient_dim]) – Point of the total space.
Returns:
vector (array-like, shape=[…, k_landmarks, ambient_dim]) – Tangent vector at base_point, result of the A tensor applied to
tangent_vec_x and tangent_vec_e.
References
[ONeill]
O’Neill, Barrett. The Fundamental Equations of a
Submersion, Michigan Mathematical Journal 13, no. 4
(December 1966): 459–69. https://doi.org/10.1307/mmj/1028999604.
[Pennec]
Pennec, Xavier. Computing the curvature and its gradient
in Kendall shape spaces. Unpublished.
Compute the covariant derivative of the integrability tensor A.
The horizontal covariant derivative \(\nabla_X (A_Y E)\) is
necessary to compute the covariant derivative of the curvature in a
submersion.
The components \(\nabla_X (A_Y E)\) and \(A_Y E\) are
computed here for the Kendall shape space at base-point
\(P = base\_point\) for horizontal vector fields fields :math:
X, Y extending the values \(X|_P = horizontal\_vec\_x\),
\(Y|_P = horizontal\_vec\_y\) and a general vector field
\(E\) extending \(E|_P = tangent\_vec\_e\) in a neighborhood
of the base-point P with covariant derivatives
\(\nabla_X Y |_P = nabla_x y\) and
\(\nabla_X E |_P = nabla_x e\).
Compute derivative of the integrability tensor A (special case).
The horizontal covariant derivative \(\nabla_X (A_Y Z)\) of the
integrability tensor A may be computed more efficiently in the case of
parallel vector fields in the quotient space.
\(\nabla_X (A_Y Z)\) and \(A_Y Z\) are computed here for the
Kendall shape space with quotient-parallel vector fields \(X,
Y, Z\) extending the values horizontal_vec_x, horizontal_vec_y and
horizontal_vec_z by parallel transport in a neighborhood of the
base-space. Such vector fields verify \(\nabla_X^X = A_X X =
0\), \(\nabla_X^Y = A_X Y\) and similarly for Z.
Parameters:
horizontal_vec_x (array-like, shape=[…, k_landmarks, ambient_dim]) – Tangent vector at base_point.
horizontal_vec_y (array-like, shape=[…, k_landmarks, ambient_dim]) – Tangent vector at base_point.
horizontal_vec_z (array-like, shape=[…, k_landmarks, ambient_dim]) – Tangent vector at base_point.
base_point (array-like, shape=[…, k_landmarks, ambient_dim]) – Point of the total space.
Returns:
nabla_x_a_y_z (array-like, shape=[…, k_landmarks, ambient_dim]) – Tangent vector at base_point, result of
\(\nabla_X (A_Y Z)\) with X = horizontal_vec_x,
Y = horizontal_vec_y and Z = horizontal_vec_z.
a_y_z (array-like, shape=[…, k_landmarks, ambient_dim]) – Tangent vector at base_point, result of \(A_Y Z\)
with Y = horizontal_vec_y and Z = horizontal_vec_z.
References
[Pennec]
Pennec, Xavier. Computing the curvature and its gradient
in Kendall shape spaces. Unpublished.
Compute iterated derivatives of the integrability tensor A.
The iterated horizontal covariant derivative
\(\nabla_X (A_Y A_X Y)\) (where \(X\) and \(Y\) are
horizontal vector fields) is a key ingredient in the computation of
the covariant derivative of the directional curvature in a submersion.
The components \(\nabla_X (A_Y A_X Y)\), \(A_X A_Y A_X Y\),
\(\nabla_X (A_X Y)\), and intermediate computations
\(A_Y A_X Y\) and \(A_X Y\) are computed here for the
Kendall shape space in the special case of quotient-parallel vector
fields \(X, Y\) extending the values horizontal_vec_x and
horizontal_vec_y by parallel transport in a neighborhood.
Such vector fields verify \(\nabla_X^X = A_X X\) and
\(\nabla_X^Y = A_X Y\).
Parameters:
horizontal_vec_x (array-like, shape=[…, k_landmarks, ambient_dim]) – Tangent vector at base_point.
horizontal_vec_y (array-like, shape=[…, k_landmarks, ambient_dim]) – Tangent vector at base_point.
base_point (array-like, shape=[…, k_landmarks, ambient_dim]) – Point of the total space.
Returns:
nabla_x_a_y_a_x_y (array-like, shape=[…, k_landmarks, ambient_dim]) – Tangent vector at base_point, result of
\(\nabla_X^S (A_Y A_X Y)\) with
X = horizontal_vec_x and Y = horizontal_vec_y.
a_x_a_y_a_x_y (array-like, shape=[…, k_landmarks, ambient_dim]) – Tangent vector at base_point, result of
\(A_X A_Y A_X Y\) with
X = horizontal_vec_x and Y = horizontal_vec_y.
nabla_x_a_x_y (array-like, shape=[…, k_landmarks, ambient_dim]) – Tangent vector at base_point, result of
\(\nabla_X^S (A_X Y)\) with
X = horizontal_vec_x and Y = horizontal_vec_y.
a_y_a_x_y (array-like, shape=[…, k_landmarks, ambient_dim]) – Tangent vector at base_point, result of \(A_Y A_X Y\) with
X = horizontal_vec_x and Y = horizontal_vec_y.
a_x_y (array-like, shape=[…, k_landmarks, ambient_dim]) – Tangent vector at base_point, result of \(A_X Y\) with
X = horizontal_vec_x and Y = horizontal_vec_y.
References
[Pennec]
Pennec, Xavier. Computing the curvature and its gradient
in Kendall shape spaces. Unpublished.
Compute the covariant derivative of the curvature.
For four vectors fields \(H|_P =\) tangent_vec_a,
\(X|_P =\) tangent_vec_b, \(Y|_P =\) tangent_vec_c,
\(Z|_P =\) tangent_vec_d with
tangent vector value specified in argument at the base_point \(P\),
the covariant derivative of the curvature
\((\nabla_H R)(X, Y) Z |_P\) is computed at the base_point \(P\).
Since the sphere is a constant curvature space this vanishes identically.
Parameters:
tangent_vec_a (array-like, shape=[…, k_landmarks, ambient_dim]) – Tangent vector at base_point along which the curvature is
derived.
This is is the supremum of radii r for which the exponential map is a
diffeomorphism from the open ball of radius r centered at the base
point onto its image.
In the case of the sphere, it does not depend on the base point and is
Pi everywhere.
Parameters:
base_point (array-like, shape=[…, k_landmarks, ambient_dim]) – Point on the manifold.
The pre-shape space is the sphere of the space of centered k-ad of
landmarks in \(R^m\) (for the Frobenius norm). It is endowed with the
spherical Procrustes metric d(x, y):= arccos(tr(xy^t)).
Points are represented by \(k \times m\) centred matrices as in
[Nava]. Beware that this is not the usual convention from the literature.
Parameters:
k_landmarks (int) – Number of landmarks
ambient_dim (int) – Number of coordinates of each landmark.
Nava-Yazdani, E., H.-C. Hege, T. J.Sullivan, and C. von Tycowicz.
“Geodesic Analysis in Kendall’s Shape Space with Epidemiological
Applications.”
Journal of Mathematical Imaging and Vision 62, no. 4 549–59.
https://doi.org/10.1007/s10851-020-00945-w.
Project a vector in the embedding matrix space
to the tangent space of the pre-shape space at a base point.
Parameters:
vector (array-like, shape=[…, k_landmarks, ambient_dim]) – Vector in Matrix space.
base_point (array-like, shape=[…, k_landmarks, ambient_dim]) – Point on the pre-shape space defining the tangent space,
where the vector will be projected.
Returns:
tangent_vec (array-like, shape=[…, k_landmarks, ambient_dim]) – Tangent vector in the tangent space of the pre-shape space
at the base point.
The ProductHPDMatricesAndSiegelDisks manifold is defined as a
product manifold of the HPD manifold and (n-1) Siegel disks.
The HPD Siegel disks product has a product metric.
The product metric on the HPD Siegel disks product space is the usual
HPD matrices affine-invariant metric (with power affine parameter equal 1)
and Siegel metrics multiplied by constants.
This product manifold can be used to represent Block-Toeplitz HPD matrices.
Lead author: Yann Cabanes.
References
[Cabanes2022]
Yann Cabanes. Multidimensional complex stationary
centered Gaussian autoregressive time series machine learning
in Poincaré and Siegel disks: application for audio and radar
clutter classification, PhD thesis, 2022
Class for the ProductHPDMatricesAndSiegelDisks manifold.
The HPD and Siegel product manifold is a direct product of the HPD manifold
and (n-1) Siegel disks. Each manifold of the product is a square matrix
manifold of the same dimension.
Parameters:
n_manifolds (int) – Number of manifolds of the product.
Class defining the ProductHPDMatricesAndSiegelDisks metric.
The HPD Siegel disks product metric is a product of the HPD metric
and (n-1) Siegel metrics, each of them being multiplied by a specific
constant factor (see [Cabanes2022], [Cabanes2021] and [JV2016]).
This metric comes from a model used to represent multidimensional complex
stationary centered Gaussian autoregressive times series.
A multidimensional times series can be seen as a realization of
a multidimensional complex Gaussian distributions with zero mean,
a block-Toeplitz HPD covariance matrix and a zero relation matrix.
The ProductHPDMatricesAndSiegelDisks metric is inspired by information geometry
on this specific set of Gaussian distributions.
References
[Cabanes2022]
Yann Cabanes. Multidimensional complex stationary
centered Gaussian autoregressive time series machine learning
in Poincaré and Siegel disks: application for audio and radar
clutter classification, PhD thesis, 2022
In contrast to the classes NFoldManifold, Landmarks, or DiscretizedCurves,
the manifolds M_1, …, M_n need not be the same, nor of
same dimension, but the list of manifolds needs to be provided.
Parameters:
factors (tuple) – Collection of manifolds in the product.
point_ndim (int or None) – If None, defaults to 1, unless all factors have the same shape.
Sample on the tangent space from the product distribution.
The distribution used is the product of the distributions used by the
random_tangent_vec methods of each individual factor manifold.
Parameters:
base_point (array-like, shape=[…, n, n]) – Base point of the tangent space.
Optional, default: None.
n_samples (int) – Number of samples.
Optional, default: 1.
Returns:
samples (array-like, shape=[…, {dim, embedding_space.dim, [n_manifolds, dim_each]}]) – Points sampled in the tangent space of the product manifold at base_point.
Generate parameterized function for the geodesic curve.
Geodesic curve defined by either:
an initial point and an initial tangent vector,
an initial point and an end point.
Parameters:
initial_point (array-like, shape=[…, dim]) – Point on the manifold, initial point of the geodesic.
end_point (array-like, shape=[…, dim], optional) – Point on the manifold, end point of the geodesic. If None,
an initial tangent vector must be given.
initial_tangent_vec (array-like, shape=[…, dim]) – Tangent vector at base point, the initial speed of the geodesics.
Optional, default: None.
If None, an end point must be given and a logarithm is computed.
Returns:
path (callable) – Time parameterized geodesic curve. If a batch of initial
conditions is passed, the output array’s first dimension
represents the different initial conditions, and the second
corresponds to time.
Matrix of the inner-product defined by the Riemmanian metric
at point base_point of the manifold.
Parameters:
base_point (array-like, shape=[…, self.shape]) – Point on the manifold at which to compute the inner-product matrix.
Optional, default: None.
Returns:
matrix (array-like, shape as described below) – Matrix of the inner-product at the base point.
The matrix is in block diagonal form with a block for each factor.
Each block is the same size as the metric_matrix for that factor.
The ProductPositiveRealsAndComplexPoincareDisks manifold.
The ProductPositiveRealsAndComplexPoincareDisks manifold is defined as a
product manifold of the positive reals manifold and (n-1) complex Poincaré disks.
The positive reals and complex Poincaré disks product has a product metric.
The product metric on the positive reals and complex Poincaré disks product space is
the positive reals metric and (n - 1) complex Poincaré metrics multiplied by constants.
This product manifold can be used to represent Toeplitz HPD matrices.
The ProductPositiveRealsAndComplexPoincareDisks corresponds to
the one-dimensional case of ProductHPDMatricesAndSiegelDisks.
Indeed, PositiveReals is the one-dimensional case of HPDMatrices and
ComplexPoincareDisk is the one-dimensional case of Siegel.
In these one-dimensional manifolds, many simplifications occur compared with
the multidimensional manifolds since matrices commute in dimension 1.
Lead author: Yann Cabanes.
References
[Cabanes_2022]
Yann Cabanes. Multidimensional complex stationary
centered Gaussian autoregressive time series machine learning
in Poincaré and Siegel disks: application for audio and radar
clutter classification, PhD thesis, tel-03708515, 2022.
https://theses.hal.science/tel-03708515
[Cabanes_CESAR_2019]
Yann Cabanes, Frédéric Barbaresco, Marc Arnaudon et
Jérémie Bigot. Unsupervised Machine Learning for Pathological Radar Clutter
Clustering: the P-Mean-Shift Algorithm, IEEE, C&ESAR 2019, Rennes, France, 2019.
https://hal.archives-ouvertes.fr/hal-02875430
[Cabanes_RADAR_2019]
Yann Cabanes, Frédéric Barbaresco, Marc Arnaudon et
Jérémie Bigot. Non-Supervised High Resolution Doppler Machine Learning for
Pathological Radar Clutter, IEEE, RADAR 2019, Toulon, France, 2019.
https://hal.archives-ouvertes.fr/hal-02875415
[Cabanes_GSI_2019]
Yann Cabanes, Frédéric Barbaresco, Marc Arnaudon et
Jérémie Bigot. Toeplitz Hermitian Positive Definite Matrix Machine Learning
based on Fisher Metric, IEEE, GSI 2019, Toulouse, France, 2019.
https://hal.archives-ouvertes.fr/hal-02875403
[Le_Brigant_2017]
Alice Le Brigant. Probability on the spaces of curves and
the associated metric spaces using information geometry; radar applications,
PhD thesis, tel-01635258, 2017.
https://theses.hal.science/tel-01635258
Class for the ProductPositiveRealsAndComplexPoincareDisks manifold.
The positive reals and complex Poincaré disks product manifold is a
direct product of the positive reals manifold and (n-1) complex Poincaré disks.
Each manifold of the product is a one-dimensional manifold.
Parameters:
n_manifolds (int) – Number of manifolds of the product.
Class defining the ProductPositiveRealsAndComplexPoincareDisks metric.
The positive reals and complex Poincaré disks product metric is a product
of the positive reals metric and (n-1) complex Poincaré metrics, each of them
being multiplied by a specific constant factor (see [Cabanes_2022],
[Cabanes_CESAR_2019], [Cabanes_RADAR_2019], [Cabanes_GSI_2019],
[Le_Brigant_2017], [Jeuris_2016] and [Yang_2013]).
This metric comes from a model used to represent one-dimensional complex
stationary centered Gaussian autoregressive times series.
Such a times series can be seen as a realization of
a multidimensional complex Gaussian distributions with zero mean,
a Toeplitz HPD covariance matrix and a zero relation matrix.
The ProductPositiveRealsAndComplexPoincareDisks metric is inspired by
information geometry on this specific set of Gaussian distributions.
References
[Cabanes_2022]
Yann Cabanes. Multidimensional complex stationary
centered Gaussian autoregressive time series machine learning
in Poincaré and Siegel disks: application for audio and radar
clutter classification, PhD thesis, tel-03708515, 2022.
https://theses.hal.science/tel-03708515
[Cabanes_CESAR_2019]
Yann Cabanes, Frédéric Barbaresco, Marc Arnaudon et
Jérémie Bigot. Unsupervised Machine Learning for Pathological Radar Clutter
Clustering: the P-Mean-Shift Algorithm, IEEE, C&ESAR 2019, Rennes, France, 2019.
https://hal.archives-ouvertes.fr/hal-02875430
[Cabanes_RADAR_2019]
Yann Cabanes, Frédéric Barbaresco, Marc Arnaudon et
Jérémie Bigot. Non-Supervised High Resolution Doppler Machine Learning for
Pathological Radar Clutter, IEEE, RADAR 2019, Toulon, France, 2019.
https://hal.archives-ouvertes.fr/hal-02875415
[Cabanes_GSI_2019]
Yann Cabanes, Frédéric Barbaresco, Marc Arnaudon et
Jérémie Bigot. Toeplitz Hermitian Positive Definite Matrix Machine Learning
based on Fisher Metric, IEEE, GSI 2019, Toulouse, France, 2019.
https://hal.archives-ouvertes.fr/hal-02875403
[Le_Brigant_2017]
Alice Le Brigant. Probability on the spaces of curves and
the associated metric spaces using information geometry; radar applications,
PhD thesis, tel-01635258, 2017.
https://theses.hal.science/tel-01635258
Let \(f\) be an immersion \(f: M \rightarrow N\)
of one manifold \(M\) into the Riemannian manifold \(N\)
with metric \(g\).
The pull-back metric \(f^*g\) is defined on \(M\) for a
base point \(p\) as:
\((f^*g)_p(u, v) = g_{f(p)}(df_p u , df_p v)
\quad \forall u, v \in T_pM\)
Note
The pull-back metric is currently only implemented for an
immersion into the Euclidean space, i.e. for
The derivative of the metrix matrix is given by
\(\partial_k g_{ij}(p)\)
where \(p\) is the base_point.
The index k of the derivation is last.
Parameters:
base_point (array-like, shape=[…, *shape]) – Base point.
Optional, default: None.
Returns:
inner_prod_deriv_mat (array-like, shape=[…, dim, dim, dim]) – Inner-product derivative matrix, where the index of the derivation
is last: \(mat_{ijk} = \partial_k g_{ij}\).
The mean curvature vector is defined at base point \(p\) by
\(H_p^\alpha= \frac{1}{d} (f^{*}g)_{p}^{ij} (\partial_{i j}^2 f^\alpha(p)\)\(-\Gamma_{i j}^k(p) \partial_k f^\alpha(p))\)
where \(f^{*}g\) is the pullback of the metric \(g\) by the
immersion \(f\).
Parameters:
base_point (array-like, shape=[…, dim]) – Base point.
Returns:
mean_curvature_vector (array-like, shape=[…, embedding_dim]) – Mean curvature vector.
Metric matrix at the tangent space at a base point.
Let \(f\) be the immersion
\(f: M \rightarrow \mathbb{R}^n\) of the manifold
\(M\) into the Euclidean space \(\mathbb{R}^n\).
The elements of the metric matrix at a base point \(p\)
are defined as:
\((f*g)_{ij}(p) = <df_p e_i , df_p e_j>\),
for \(e_i, e_j\) basis elements of \(M\).
Parameters:
base_point (array-like, shape=[…, dim]) – Base point.
Optional, default: None.
Returns:
mat (array-like, shape=[…, dim, dim]) – Inner-product matrix.
In the case of an immersion f, the second fundamental form is
given by the formula:
\(\RN{2}(p)_{ij}^\alpha = \partial_{i j}^2 f^\alpha(p)\)\(-\Gamma_{i j}^k(p) \partial_k f^\alpha(p)\)
at base_point \(p\).
Parameters:
base_point (array-like, shape=[…, dim]) – Base point.
Returns:
second_fundamental_form (array-like, shape=[…, embedding_dim, dim, dim]) – Second fundamental form \(\RN{2}(p)_{ij}^\alpha\) where the
\(\alpha\) index is first.
Given a (principal) fiber bundle, or more generally a manifold with a
Lie group acting on it by the right, the quotient space is the space of
orbits under this action.
In general, the quotient space is not a manifold, as it can exhibit singularities.
Instead, the quotient space is a stratified space, whose principal stratum is a manifold.
We restrict computations on the quotient space to computations on its principal stratum,
that we equip with the quotient metric.
The quotient metric is defined such that the
canonical projection is a Riemannian submersion, i.e., it is isometric to
the restriction of the metric of the total space to horizontal subspaces.
Parameters:
space (Manifold) – Base.
total_space (Manifold) – Total space equipped with a fiber bundle structure.
For three vectors fields \(X|_P = tangent\_vec\_a,
Y|_P = tangent\_vec\_b, Z|_P = tangent\_vec\_c\) with tangent vector
specified in argument at the base point \(P\),
the curvature is defined by \(R(X,Y)Z = \nabla_{[X,Y]}Z
- \nabla_X\nabla_Y Z + \nabla_Y\nabla_X Z\).
In the case of quotient metrics, the fundamental equations of a
Riemannian submersion allow to compute the curvature of the base
manifold from the one of the total space and a correction term that
uses the integrability tensor A [ONeill].
In more details, let \(X, Y, Z\) be the horizontal lift of
vector fields extending the tangent vectors given in argument in a
neighborhood of the base-point P in the base-space. Then the
curvature of the base-space at the base-points is
\(R(X,Y) Z = hor( R^T(X,Y) Z) - 2 A_Z A_X Y + A_X A_Y Z - A_Y
A_X Z\), where \(R^T(X,Y)Z\) is the curvature tensor of the
total space.
O’Neill, Barrett. The Fundamental Equations of a
Submersion, Michigan Mathematical Journal 13, no. 4
(December 1966): 459–69. https://doi.org/10.1307/mmj/1028999604.
Compute the covariant derivative of the curvature.
For four vectors fields \(H|_P =\) tangent_vec_a,
\(X|_P =\) tangent_vec_b, \(Y|_P =\) tangent_vec_c,
\(Z|_P =\) tangent_vec_d with
tangent vector value specified in argument at the base_point \(P\),
the covariant derivative of the curvature
\((\nabla_H R)(X, Y)Z |_P\) is computed at the base_point \(P\).
In the case of quotient metrics, the fundamental equations of a
Riemannian submersion allow to compute this tensor on the base manifold
from the one of the total space T and its covariant derivative with
additional correction terms involving the integrability tensor A and
its covariant derivatives [Pennec].
In more details, let \(H, X, Y, Z\) be the horizontal lift of
parallel vector fields extending the tangent vectors given in argument
by parallel transport in a neighborhood of the base_point \(P\) in the
base-space. Such vector fields verify \(\nabla^T_H H=0\),
\(\nabla^T_H X = A_H X\) (and similarly for Y and Z) using the
connection \(\nabla^T\) of the total space. Then the covariant
derivative of the curvature tensor is given by
Compute the covariant derivative of the directional curvature.
For two vectors fields \(X|_P =\) tangent_vec_a,
\(Y|_P =\) tangent_vec_b with tangent vector value specified in
argument at the base_point \(P\),
the covariant derivative (in the direction \(X\))
\((\nabla_X R_Y)(X) |_P = (\nabla_X R)(Y, X) Y |_P\) of the
directional curvature (in the direction \(Y\))
\(R_Y(X) = R(Y, X) Y\) is a quadratic tensor in \(X\) and \(Y\)
that plays an important role in the computation of the moments of the
empirical Fréchet mean.
This tensor can be computed from the covariant derivative of the
curvature tensor as in done generically the Connection class.
However, in the case of quotient metrics, a simplified expression can
be implemented based on the directional curvature of the total space T
and its covariant derivative with additional correction terms
involving the integrability tensor A and its covariant derivatives
[Pennec].
In more details, let \(X, Y\) be the horizontal lift of parallel
vector fields extending the tangent vectors given in argument by
parallel transport in a neighborhood of the base_point \(P\) in the
base-space. Such vector fields verify \(\nabla^T_X X=0\) and
\(\nabla^T_X^Y = A_X Y\) using the connection \(\nabla^T\)
of the total space. Then the covariant derivative of the
directional curvature tensor is given by
\[\nabla_X (R_Y(X)) = hor \nabla^T_X (R^T_Y(X)) - A_X( ver R^T_Y(X))
+ 3 A_X A_Y A_X Y - 3 \nabla_X^T A_Y A_X Y\]
where \(R^T_Y(X)\) is the directional curvature tensor of the total space.
Parameters:
tangent_vec_a (array-like, shape=[…, n, n]) – Tangent vector at base_point.
tangent_vec_b (array-like, shape=[…, n, n]) – Tangent vector at base_point.
base_point (array-like, shape=[…, n, n]) – Point on the base manifold.
Returns:
curvature_derivative (array-like, shape=[…, n, n]) – Tangent vector at base_point.
References
[Pennec]
Pennec, Xavier. Computing the curvature and its gradient
in Kendall shape spaces. Unpublished.
Generate parameterized function for the geodesic curve.
Geodesic curve defined by either:
an initial point and an initial tangent vector,
an initial point and an end point.
Parameters:
initial_point (array-like, shape=[…, dim]) – Point on the manifold, initial point of the geodesic.
end_point (array-like, shape=[…, dim], optional) – Point on the manifold, end point of the geodesic. If None,
an initial tangent vector must be given.
initial_tangent_vec (array-like, shape=[…, dim],) – Tangent vector at base point, the initial speed of the geodesics.
Optional, default: None.
If None, an end point must be given and a logarithm is computed.
Returns:
path (callable) – Time parameterized geodesic curve. If a batch of initial
conditions is passed, the output array’s first dimension
represents the different initial conditions, and the second
corresponds to time.
Compute the inner-product of two tangent vectors at a base point.
Parameters:
tangent_vec_a (array-like, shape=[…, {dim, [n, n]}]) – Tangent vector to the quotient manifold.
tangent_vec_b (array-like, shape=[…, {dim, [n, n]}]) – Tangent vector to the quotient manifold.
base_point (array-like, shape=[…, {dim, [n, n]}]) – Point on the quotient manifold.
Optional, default: None.
fiber_point (array-like, shape=[…, {dim, [n, n]}]) – Point on the total space, lift of base_point, i.e. such that
submersion applied to point results in base_point.
Optional, default: None. In this case, it is computed using the
method lift.
The class is redirecting to the correct embedding manifold.
The stratum PSD rank k if the matrix is not full rank.
The top stratum SPD if the matrix is full rank.
The whole stratified space of PSD if no rank is specified.
Parameters:
n (int) – Integer representing the shapes of the matrices : n x n.
k (int) – Integer representing the rank of the matrices.
Project a matrix to the space of PSD matrices of rank k.
The nearest symmetric positive semidefinite matrix in the
Frobenius norm to an arbitrary real matrix A is shown to be (B + H)/2,
where H is the symmetric polar factor of B=(A + A’)/2.
As [Higham1988] is turning the matrix into a PSD, the rank
is then forced to be k.
Parameters:
point (array-like, shape=[…, n, n]) – Matrix to project.
Returns:
projected (array-like, shape=[…, n, n]) – PSD matrix rank k.
Highamm, N. J.
“Computing a nearest symmetric positive semidefinite matrix.”
Linear Algebra and Its Applications 103 (May 1, 1988):
103-118. https://doi.org/10.1016/0024-3795(88)90223-6
Compute Christoffel symbols of the Levi-Civita connection.
The Koszul formula defining the Levi-Civita connection gives the
expression of the Christoffel symbols with respect to the metric:
\(\Gamma^k_{ij}(p) = \frac{1}{2} g^{lk}(
\partial_i g_{jl} + \partial_j g_{li} - \partial_l g_{ij})\),
where:
\(p\) represents the base point, and
\(g\) represents the Riemannian metric tensor.
Note that the function computing the derivative of the metric matrix
puts the index of the derivation last.
Parameters:
base_point (array-like, shape=[…, dim]) – Base point.
Returns:
christoffels (array-like, shape=[…, dim, dim, dim]) – Christoffel symbols, where the contravariant index is first.
Compute purely covariant version of Riemannian tensor at base_point.
In the literature the covariant riemannian tensor is noted \(R_{ijkl}\).
Convention used in the literature (tensor index notation, ref. Wikipedia) is:
\(R_{ijkl} = <R(x_k, x_l)x_j, x_i>\)
where \(x_i\) is the i-th basis vector of the tangent space at base point.
In other words:
[cov_riemann_tensor]_{ijkl} = [metric_matrix]_{im} [riemann_tensor]_{jkl}^m
Parameters:
base_point (array-like, shape=[…, dim]) – Point on the manifold.
Returns:
covariant_tensor (array-like, shape=[…, dim, dim, dim, dim]) – covariant_riemann_tensor[…, i, j, k, l] = R_{ijkl}
Covariant version of Riemannian curvature tensor.
If n_samples_a == n_samples_b then dist is the element-wise
distance result of a point in points_a with the point from
points_b of the same index. If n_samples_a not equal to
n_samples_b then dist is the result of applying geodesic
distance for each point from points_a to all points from
points_b.
Parameters:
point_a (array-like, shape=[n_samples_a, dim]) – Set of points in the Poincare ball.
point_b (array-like, shape=[n_samples_b, dim]) – Second set of points in the Poincare ball.
Returns:
dist (array-like, shape=[n_samples_a, dim] or [n_samples_a, n_samples_b, dim]) – Geodesic distance between the two points.
points (array-like, shape=[n_samples, dim]) – Set of points in the manifold.
n_jobs (int) – Number of jobs to run in parallel, using joblib. Note that a
higher number of jobs may not be beneficial when one computation
of a geodesic distance is cheap.
Optional. Default: 1.
**joblib_kwargs (dict) – Keyword arguments to joblib.Parallel
Returns:
dist (array-like, shape=[n_samples, n_samples]) – Pairwise distance matrix between all the points.
Compute derivative of the inner prod matrix at base point.
Writing \(g_{ij}\) the inner-product matrix at base point,
this computes \(mat_{ijk} = \partial_k g_{ij}\), where the
index k of the derivation is put last.
Parameters:
base_point (array-like, shape=[…, dim]) – Base point.
Returns:
metric_derivative (array-like, shape=[…, dim, dim, dim]) – Derivative of the inner-product matrix, where the index
k of the derivation is last: \(mat_{ijk} = \partial_k g_{ij}\).
Generate a random unit tangent vector at a given point.
Parameters:
base_point (array-like, shape=[…, dim]) – Point.
n_vectors (float) – Number of vectors to be generated at base_point.
For vectorization purposes n_vectors can be greater than 1 iff
base_point consists of a single point.
Returns:
normalized_vector (array-like, shape=[…, n_vectors, dim]) – Random unit tangent vector at base_point.
In the literature scalar_curvature is noted S and writes:
\(S = g^{ij} Ric_{ij}\),
with Einstein notation, where we recognize the trace of the Ricci
tensor according to the Riemannian metric \(g\).
Parameters:
base_point (array-like, shape=[…, dim]) – Point on the manifold.
Implements of the Sasaki metric on the tangent bundle TM of a Riem. manifold M.
The Sasaki metric is characterized by the following three properties:
the canonical projection of TM becomes a Riemannian submersion,
parallel vector fields along curves are orthogonal to their fibres, and
its restriction to any tangent space is Euclidean.
Geodesic computations are realized via a discrete formulation of the
geodesic equation on TM that involve geodesics, parallel translation, and
the curvature tensor on the base manifold M (see [1]_ for details).
However, as the implemented energy in the discrete-geodesics-optimization
as well as the approximations of its gradient slightly differ from those
proposed in [1]_, we also refer to [2] for additional details.
Parameters:
space (Manifold) – Tangent bundle.
n_jobs (int) – Number of jobs for parallel computing.
Optional, default: 1.
n_steps (int) – Number of discrete time steps.
Optional, default: 3.
Class for scalar products of Riemannian and pseudo-Riemannian metrics.
This class multiplies the (0,2) metric tensor ‘space.metric’ by a
scalar ‘scaling_factor’. Note that this does not scale distances by
‘scaling_factor’. That would require multiplication by the square of the
scalar.
The space is not automatically equipped with the ScalarProductMetric.
An object of this type can also be instantiated by the expression
scaling_factor * space.metric.
This class acts as a wrapper for the underlying Riemannian metric. All
public attributes apart from ‘underlying_metric’ and ‘scaling_factor’ are
loaded from the underlying metric at initialization and rescaled by the
appropriate factor. Changes to the underlying metric at runtime will not
affect the attributes of this object.
One exception to this is when the ‘underlying_metric’ is itself of type
ScalarProductMetric. In this case, rather than wrapping the wrapper, the
‘underlying_metric’ of the first ScalarProductMetric object is wrapped a
second time with a new ‘scaling_factor’.
Parameters:
space (Manifold or ComplexManifold) – A manifold equipped with a metric which is being scaled.
scale (float) – The value by which to scale the metric. Note that this rescales the
(0,2) metric tensor, so distances are rescaled by the square root of
this.
Register the scaling factor of a method of a RiemannianMetric.
The ScalarProductMetric class rescales various methods of a
RiemannianMetric by the correct factor. The default behaviour is to
rescale linearly. This method allows the user to add a new method to be
rescaled according to a different rule.
Note that this method must be called before the ScalarProductMetric is
instantiated. It does not affect objects which already exist.
Parameters:
func_name (str) – The name of a method from a RiemannianMetric object which must be
rescaled.
scaling_type (str, {‘sqrt’,) – ‘linear’,
‘quadratic’,
‘inverse’,
‘inverse_sqrt’}
How the method should be rescaled as a function of
ScalarProductMetric.scale.
The Siegel manifold is a generalization of the Poincare disk to complex matrices
of singular values strictly lower than one.
It is defined as the set of complex matrices M such that:
I - M @ M.conj().T is a positive definite matrix.
Warning: another more restrictive definition of the Siegel disk also exists
which add a symmetry condition on the matrices.
It has been proven in [Cabanes2022] that the sub-manifold of symmetric Siegel
matrices is a totally geodesic sub-manifold of the Siegel space.
The sub-manifold of real Siegel matrices is also a totally geodesic sub-manifold
of the Siegel space.
Lead author: Yann Cabanes.
References
[Cabanes2022]
Yann Cabanes. Multidimensional complex stationary
centered Gaussian autoregressive time series machine learning
in Poincaré and Siegel disks: application for audio and radar
clutter classification, PhD thesis, 2022
The Siegel space is the set of complex matrices of singular values
strictly lower than one.
The Frobenius norm of a matrix is greater than or equal to the spectral norm
which corresponds to the largest singular value of a matrix.
Then, simulating a matrix with Frobenius norm strictly lower than one,
its singular values are also strictly lower than one,
therefore this matrix belongs to the Siegel disk.
Parameters:
n_samples (int) – Number of samples.
Optional, default: 1.
bound (float) – Bound of the interval in which to sample in the tangent space.
Optional, default: 1.
Returns:
samples (array-like, shape=[…, n, n]) – Points sampled in the Siegel space.
Compute the inner-product of tangent_vec_a and tangent_vec_b
at point base_point using the Siegel Riemannian metric.
The expression of the inner product between the vectors v and w
at point O is \(<v, w>_{O}
= 1/2 * trace((I - O O^{H})^{-1} v (I - O^{H} O)^{-1} w^{H})
+ 1/2 * trace((I - O O^{H})^{-1} w (I - O^{H} O)^{-1} v^{H})
= Re(trace((I - O O^{H})^{-1} v (I - O^{H} O)^{-1} w^{H}))\)
where H denotes the conjugate transpose operator.
Parameters:
tangent_vec_a (array-like, shape=[…, n, n]) – Tangent vector at base point.
tangent_vec_b (array-like, shape=[…, n, n]) – Tangent vector at base point.
base_point (array-like, shape=[…, n, n]) – Base point.
For two orthonormal tangent vectors at a base point \(x,y\),
the sectional curvature is defined by \(<R(x, y)x,
y>\). Non-orthonormal vectors can be given.
Parameters:
tangent_vec_a (array-like, shape=[…, n, n]) – Tangent vector at base point.
tangent_vec_b (array-like, shape=[…, n, n]) – Tangent vector at base point.
base_point (array-like, shape=[…, n, n]) – Base point. Optional, default is zero
This is the Lie algebra of the Special Orthogonal Group.
As basis we choose the matrices with a single 1 on the upper triangular part
of the matrices (and a -1 in its lower triangular part), except in dim 2 and
3 to match usual conventions.
Class for the affine-invariant metric on the SPD manifold.
Parameters:
power_affine (int) – Power transformation of the classical SPD metric.
Optional, default: 1.
References
[TP2019]
Thanwerdas, Pennec. “Is affine-invariance well defined on
SPD matrices? A principled continuum of metrics” Proc. of GSI,
2019. https://arxiv.org/abs/1906.01349
Compute the Riemannian exponential at point base_point
of tangent vector tangent_vec wrt the metric defined in inner_product.
This gives a symmetric positive definite matrix.
Parameters:
tangent_vec (array-like, shape=[…, n, n]) – Tangent vector at base point.
base_point (array-like, shape=[…, n, n]) – Base point.
Returns:
exp (array-like, shape=[…, n, n]) – Riemannian exponential.
Compute the Riemannian logarithm at point base_point,
of point wrt the metric defined in inner_product.
This gives a tangent vector at point base_point.
Parameters:
point (array-like, shape=[…, n, n]) – Point.
base_point (array-like, shape=[…, n, n]) – Base point.
Returns:
log (array-like, shape=[…, n, n]) – Riemannian logarithm of point at base_point.
Closed-form solution for the parallel transport of a tangent vector
along the geodesic between two points base_point and end_point
or alternatively defined by \(t \mapsto exp_{(base\_point)}(
t*direction)\).
Denoting tangent_vec_a by S, base_point by A, and end_point
by B or B = Exp_A(tangent_vec_b) and \(E = (BA^{- 1})^{( 1
/ 2)}\). Then the parallel transport to B is:
\[S' = ESE^T\]
Parameters:
tangent_vec (array-like, shape=[…, n, n]) – Tangent vector at base point to be transported.
base_point (array-like, shape=[…, n, n]) – Point on the manifold of SPD matrices. Point to transport from
direction (array-like, shape=[…, n, n]) – Tangent vector at base point, initial speed of the geodesic along
which the parallel transport is computed. Unused if end_point is given.
Optional, default: None.
end_point (array-like, shape=[…, n, n]) – Point on the manifold of SPD matrices. Point to transport to.
Optional, default: None.
Returns:
transported_tangent_vec (array-like, shape=[…, n, n]) – Transported tangent vector at exp_(base_point)(tangent_vec_b).
Compute the inner-product of tangent_vec_a \(A\) and tangent_vec_b
\(B\) at point base_point \(S=PDP^\top\) using the
Bures-Wasserstein Riemannian metric:
Compute the parallel transport of a tangent vec along a geodesic.
Approximation of the solution of the parallel transport of a tangent
vector a along the geodesic defined by \(t \mapsto exp_{(
base\_point)}(t* tangent\_vec\_b)\). The parallel transport equation is
formulated in this case in [TP2021].
Parameters:
tangent_vec (array-like, shape=[…, n, n]) – Tangent vector at base_point to transport.
base_point (array-like, shape=[…, n, n]) – Initial point of the geodesic.
direction (array-like, shape=[…, n, n]) – Tangent vector at base_point, initial velocity of the geodesic to
transport along.
end_point (array-like, shape=[…, n, n]) – Point to transport to.
Optional, default: None.
n_steps (int) – Number of steps to use to approximate the solution of the
ordinary differential equation.
Optional, default: 100
step (str, {‘euler’, ‘rk2’, ‘rk4’}) – Scheme to use in the integration scheme.
Optional, default: ‘rk4’.
Returns:
transported (array-like, shape=[…, n, n]) – Transported tangent vector at exp_(base_point)(tangent_vec_b).
Thanwerdas, Y. (2022) Riemannian and stratified
geometries of covariance and correlation matrices. Theses.
Université Côte d’Azur. Available at:
https://hal.archives-ouvertes.fr/tel-03698752 (Accessed: 29 September 2022).
Compute the generalized eigenvalues of SPD matrix pair.
Steps (check section 7.2 of [GKC2023]):
1. compute eigendecomposition of point_b
2. get matrix turning point_b into identity by congruence
3. apply congruence to point_a and get generalized eigenvalues
Parameters:
point_a (array_like, shape=[…, n, n]) – Symmetric positive definite matrix.
point_b (array_like, shape=[…, n, n]) – Symmetric positive definite matrix.
Class for the canonical left-invariant metric on SE(n).
The canonical left-invariant metric is defined by endowing the tangent
space at the identity with the Frobenius inned-product, and to define the
metric at any point by left-translation. This results in a direct product
metric between rotations and translations, whose geodesics are therefore
easily computable with the matrix exponential and straight lines.
Parameters:
space (SpecialEuclidean) – Instance of the class SpecialEuclidean with point_type=’matrix’.
Exponential map associated to the cannonical metric.
Exponential map at base_point of tangent_vec. The geodesics of this
metric correspond to a direct product metric between rotation and
translation: the translation part is a straight line, while the
rotation part has constant angular velocity (which corresponds to one-
parameter subgroups of the rotation group).
Parameters:
tangent_vec (array-like, shape=[…, n + 1, n + 1]) – Tangent vector at the base point.
base_point (array-like, shape=[…, n + 1, n + 1]) – Point on the manifold.
Returns:
exp (array-like, shape=[…, n + 1, n + 1]) – Point on the manifold.
Generate parameterized function for the geodesic curve.
Geodesic curve defined by either:
an initial point and an initial tangent vector,
an initial point and an end point.
Parameters:
initial_point (array-like, shape=[…, dim]) – Point on the manifold, initial point of the geodesic.
end_point (array-like, shape=[…, dim], optional) – Point on the manifold, end point of the geodesic. If None,
an initial tangent vector must be given.
initial_tangent_vec (array-like, shape=[…, dim],) – Tangent vector at base point, the initial speed of the geodesics.
Optional, default: None.
If None, an end point must be given and a logarithm is computed.
Returns:
path (callable) – Time parameterized geodesic curve. If a batch of initial
conditions is passed, the output array’s first dimension
represents the different initial conditions, and the second
corresponds to time.
This is is the supremum of radii r for which the exponential map is a
diffeomorphism from the open ball of radius r centered at the base point onto
its image.
In this case, it does not depend on the base point. If the rotation part is
null, then the radius is infinite, otherwise it is the same as the special
orthonormal group.
Parameters:
base_point (array-like, shape=[…, n + 1, n + 1]) – Point on the manifold.
Compute logarithm map associated to the canonical metric.
Log map at base_point of point. The geodesics of this
metric correspond to a direct product metric between rotation and
translation: the translation part is a straight line, while the
rotation part has constant angular velocity (which corresponds to one-
parameter subgroups of the rotation group).
Parameters:
point (array-like, shape=[…, n + 1, n + 1]) – Point on the manifold.
base_point (array-like, shape=[…, n + 1, n + 1]) – Point on the manifold.
Returns:
tangent_vec (array-like, shape=[…, n + 1, n + 1]) – Tangent vector at the base point.
References
[Zefran98]
Zefran, M., V. Kumar, and C.B. Croke.
“On the Generation of Smooth Three-Dimensional Rigid Body Motions.”
IEEE Transactions on Robotics and Automation 14,
no. 4 (August 1998): 576–89.
https://doi.org/10.1109/70.704225.
Compute the parallel transport of a tangent vector.
Closed-form solution for the parallel transport of a tangent vector a
along the geodesic between two points base_point and end_point
or alternatively defined by \(t \mapsto exp_{(base\_point)}(
t*direction)\). As the special Euclidean group endowed with its
canonical left-invariant metric is a symmetric space, parallel
transport is achieved by a geodesic symmetry, or equivalently, one step
of the pole ladder scheme.
Parameters:
tangent_vec (array-like, shape=[…, n + 1, n + 1]) – Tangent vector at base point to be transported.
base_point (array-like, shape=[…, n + 1, n + 1]) – Point on the hypersphere.
direction (array-like, shape=[…, n + 1, n + 1]) – Tangent vector at base point, along which the parallel transport
is computed.
Optional, default: None
end_point (array-like, shape=[…, n + 1, n + 1]) – Point on the Grassmann manifold to transport to. Unused if
tangent_vec_b is given.
Optional, default: None
Returns:
transported_tangent_vec (array-like, shape=[…, n + 1, n + 1]) – Transported tangent vector at exp_(base_point)(tangent_vec_b).
Embed rotation, translation couples into n+1 square matrices.
Construct a block matrix of size \(n + 1 \times n + 1\) of the form
\[\begin{split}\begin{pmatrix} R & t \\ 0 & c \end{pmatrix}\end{split}\]
where \(R\) is a square matrix, \(t\) a vector of size
\(n\), and \(c\) a constant (either 0 or 1 should be used).
Parameters:
rotation (array-like, shape=[…, n, n]) – Square Matrix.
translation (array-like, shape=[…, n]) – Vector.
constant (float or array-like of shape […]) – Constant to use at the last line and column of the square matrix.
Optional, default: 1.
Returns:
mat (array-like, shape=[…, n + 1, n + 1]) – Square Matrix of size n + 1. It can represent an element of the
special euclidean group or its Lie algebra.
This is is the supremum of radii r for which the exponential map is a
diffeomorphism from the open ball of radius r centered at the base
point onto its image.
In this case the exact injectivity radius is not known, and we use here
a lower bound given by [Rentmeesters2015].
Parameters:
base_point (array-like, shape=[…, n, p]) – Point on the manifold.
R Zimmermann. A matrix-algebraic algorithm for the
Riemannian logarithm on the Stiefel manifold under the canonical
metric. SIAM Journal on Matrix Analysis and Applications 38 (2),
322-342, 2017. https://epubs.siam.org/doi/pdf/10.1137/16M1074485
Exponential map at base_point of cotangent_vec computed by integration
of the Hamiltonian equation (initial value problem), using the cometric.
In the Riemannian case this yields the exponential associated to the
Levi-Civita connection of the metric (the inverse of the cometric).
Parameters:
cotangent_vec (array-like, shape=[…, dim]) – Tangent vector at the base point.
base_point (array-like, shape=[…, dim]) – Point on the manifold.
n_steps (int) – Number of discrete time steps to take in the integration.
Optional, default: N_STEPS.
Returns:
exp (array-like, shape=[…, dim]) – Point on the manifold.
Generate parameterized function for the normal geodesic curve.
Normal geodesic curve defined by an initial point and an initial
cotangent vector.
Parameters:
initial_point (array-like, shape=[…, dim]) – Point on the manifold, initial point of the geodesic.
initial_cotangent_vec (array-like, shape=[…, dim]) – Cotangent vector at base point, the initial speed of the geodesics.
Returns:
path (callable) – Time parameterized normal geodesic curve. If a batch of initial
conditions is passed, the output array’s first dimension
represents the different initial conditions, and the second
corresponds to time.
This is the sub-Riemannian sharp map, mapping a covector at base_point to a
tangent vector in the distribution subspace at base_point. For an orthonormal
frame (F_i)_{i=1..dist_dim}, the sharp map is given by