Homogeneous Spaces

Homogeneous spaces are very important in GeometricMachineLearning as we can generalize existing neural network optimizers from vector spaces to such homogenous spaces. They are intricately linked to the notion of a Lie Group and its Lie Algebra^[1].

Now fix a distinct element $E\in\mathcal{M}$; we will refer to this as the canonical element. We can also establish an isomorphism between $\mathcal{M}$ and the quotient space $G/\sim$ with the equivalence relation:

\[A_1 \sim A_2 \iff A_1E = A_2E.\]

Note that this is independent of the chosen $E$.

The tangent spaces of $\mathcal{M}$ are of the form $T_Y\mathcal{M} = \mathfrak{g}\cdot{}Y$, i.e. can be fully described through its Lie algebra. Based on this we can perform a splitting of $\mathfrak{g}$ into two parts:

We will refer to the mapping from $T_Y\mathcal{M}$ to $\mathfrak{g}^{\mathrm{hor}, Y}$ by $\Omega$. We will give explicit examples of $\Omega$ below. If we have now defined a metric $\langle\cdot,\cdot\rangle$ on $\mathfrak{g}$, then this induces a Riemannian metric on $\mathcal{M}$:

\[g_Y(\Delta_1, \Delta_2) = \langle\Omega(Y,\Delta_1),\Omega(Y,\Delta_2)\rangle\text{ for $\Delta_1,\Delta_2\in{}T_Y\mathcal{M}$.}\]

Two examples of homogeneous spaces implemented in GeometricMachineLearning are the Stiefel and the Grassmann manifold. The Lie group $SO(N)$ acts transitively on both of these manifolds, i.e. turns them into homogeneous spaces. The Lie algebra of $SO(N)$ are the skew-symmetric matrices $\mathfrak{so}(N):=\{V\in\mathbb{R}^{N\times{}N}:V^T + V = 0\}$ and the canonical metric associated with it is simply $(V_1,V_2)\mapsto\frac{1}{2}\mathrm{Tr}(V_1^TV_2)$.

The Stiefel Manifold

The Stiefel manifold $St(n, N)$ is the space of all orthonormal frames in $\mathbb{R}^{N\times{}n}$, i.e. matrices $Y\in\mathbb{R}^{N\times{}n}$ s.t. $Y^TY = \mathbb{I}_n$. It can also be seen as $SO(N)$ modulo an equivalence relation: $A\sim{}B\iff{}AE = BE$ for

\[E = \begin{bmatrix} \mathbb{I}_n \\ \mathbb{O} \end{bmatrix}\in{}St(n, N),\]

which is the canonical element of the Stiefel manifold. In words: the first $n$ columns of $A$ and $B$ are the same. We also use this principle to draw random elements from the Stiefel manifold.

The tangent space to the element $Y\in{}St(n,N)$ can be determined by considering $C^\infty$ curves on $SO(N)$ through $\mathbb{I}$ which we write $t\mapsto{}A(t)$. Because $SO(N)$ acts transitively on $St(n, N)$ each $C^\infty$ curve on $St(n, N)$ through $Y$ can be written as $A(t)Y$ and we get:

\[T_YSt(n,N)=\{BY : B\in\mathfrak{g}\} = \{\Delta\in\mathbb{R}^{N\times{}n}: \Delta^TY + Y^T\Delta = \mathbb{O}\},\]

where the last equality can be established through the isomorphism

\[\Omega: T_YSt(n, N) \to \mathfrak{g}^{\mathrm{vec}, Y}, \Delta \mapsto (\mathbb{I} - \frac{1}{2}YY^T)\Delta{}Y^T - Y\Delta^T(\mathbb{I} - \frac{1}{2}YY^T).\]

That this is an isomorphism can be easily checked:

\[ \Omega(\Delta)Y = (\mathbb{I} - \frac{1}{2}YY^T)\Delta - \frac{1}{2}Y\Delta^TY = \Delta.\]

The isomorphism is also implemented in GeometricMachineLearning:

using GeometricMachineLearning

Y = rand(StiefelManifold{Float32}, 5, 3)
Δ = rgrad(Y, rand(Float32, 5, 3))
GeometricMachineLearning.Ω(Y, Δ) * Y.A ≈ Δ

true

The function rgrad is introduced below.

The Riemannian Gradient for the Stiefel Manifold

We defined the Riemannian gradient to be a vector field $\mathrm{grad}^gL$ such that it is compatible with the Riemannian metric in some sense; the definition we gave relied on an explicit coordinate chart. We can also express the Riemannian gradient for matrix manifolds by not relying on an explicit coordinate representation (which would be computationally expensive) [6].

For the Stiefel manifold the Riemannian gradient is given by:

\[ \mathrm{grad}_YL = \nabla_YL - Y(\nabla_YL)^TY =: \mathtt{rgrad}(Y, \nabla_YL),\]

where $\nabla_YL$ refers to the Euclidean gradient, i.e.

\[ [\nabla_YL]_{ij} = \frac{\partial{}L}{\partial{}y_{ij}}.\]

The Euclidean gradient $\nabla{}L$ can in practice be obtained with an AD routine. We then use the function rgrad to map $\nabla_YL$ from $\mathbb{R}^{N\times{}n}$ to $T_YSt(n,N)$. We can check that this mapping indeed maps to the Riemannian gradient

using GeometricMachineLearning
using LinearAlgebra: tr

Y = rand(StiefelManifold{Float32}, 5, 3)
∇L = rand(Float32, 5, 3)
gradL = rgrad(Y, ∇L)
Δ = rgrad(Y, rand(Float32, 5, 3))

metric(Y, gradL, Δ) ≈ tr(∇L' * Δ)

true

The Grassmann Manifold

The Grassmann manifold is closely related to the Stiefel manifold, and an element of the Grassmann manifold can be represented through an element of the Stiefel manifold (but not vice-versa). An element of the Grassmann manifold $G(n,N)$ is a vector subspace $\subset\mathbb{R}^N$ of dimension $n$. Each such subspace (i.e. element of the Grassmann manifold) can be represented by a full-rank matrix $A\in\mathbb{R}^{N\times{}n}$ and we identify two elements with the following equivalence relation:

\[ A_1 \sim A_2 \iff \exists{}C\in\mathbb{R}^{n\times{}n}\text{ s.t. }A_1C = A_2.\]

The resulting manifold is of dimension $n(N-n)$. One can find a parametrization of the manifold the following way: Because the matrix $Y$ has full rank, there have to be $n$ independent columns in it: $i_1, \ldots, i_n$. For simplicity assume that $i_1 = 1, i_2=2, \ldots, i_n=n$ and call the matrix made up of these columns $C$. Then the mapping to the coordinate chart is: $YC^{-1}$ and the last $N-n$ columns are the coordinates.

We can also define the Grassmann manifold based on the Stiefel manifold since elements of the Stiefel manifold are already full-rank matrices. In this case we have the following equivalence relation (for $Y_1, Y_2\in{}St(n,N)$):

\[ Y_1 \sim Y_2 \iff \exists{}C\in{}SO(n)\text{ s.t. }Y_1C = Y_2.\]

In GeometricMachineLearning elements of the Grassmann manifold are drawn the same way as elements of the Stiefel manifold:

using GeometricMachineLearning

rand(GrassmannManifold{Float32}, 5, 3)

5×3 GrassmannManifold{Float32, Matrix{Float32}}:
 -0.157512    0.185112  -0.216194
  0.0595288  -0.523402  -0.567974
 -0.379467    0.669637  -0.545434
  0.908048    0.323416  -0.259806
 -0.0556817  -0.372523  -0.51543

The Riemannian Gradient of the Grassmann Manifold

Obtaining the Riemannian Gradient for the Grassmann manifold is slightly more difficult than it is in the case of the Stiefel manifold [6]. Since the Grassmann manifold can be obtained from the Stiefel manifold through an equivalence relation, we can however use this as a starting point.

Library Functions

rand(manifold_type::Type{MT}, N::Integer, n::Integer) where MT <: Manifold

using GeometricMachineLearning
using GeometricMachineLearning: _round # hide
import Random
Random.seed!(123)

N, n = 5, 3
Y = rand(StiefelManifold{Float32}, N, n)
_round(Y; digits = 5) # hide

# output

5×3 StiefelManifold{Float32, Matrix{Float32}}:
 -0.27575   0.32991   0.77275
 -0.62485  -0.33224  -0.0686
 -0.69333   0.36724  -0.18988
 -0.09295  -0.73145   0.46064
  0.2102    0.33301   0.38717

rgrad(Y::StiefelManifold, e_grad::AbstractMatrix)

using GeometricMachineLearning

Y = StiefelManifold([1 0 ; 0 1 ; 0 0; 0 0])
Δ = [1 2; 3 4; 5 6; 7 8]
rgrad(Y, Δ)

# output

4×2 Matrix{Int64}:
 0  -1
 1   0
 5   6
 7   8

rgrad(Y::GrassmannManifold, e_grad::AbstractMatrix)

using GeometricMachineLearning

Y = GrassmannManifold([1 0 ; 0 1 ; 0 0; 0 0])
Δ = [1 2; 3 4; 5 6; 7 8]
rgrad(Y, Δ)

# output

4×2 Matrix{Int64}:
 0  0
 0  0
 5  6
 7  8

metric(Y::GrassmannManifold, Δ₁::AbstractMatrix, Δ₂::AbstractMatrix)

Ω(Y::StiefelManifold{T}, Δ::AbstractMatrix{T}) where T

SkewSymMatrix(2 * (I(n) - .5 * Y * Y') * Δ * Y')

using GeometricMachineLearning
E = StiefelManifold(StiefelProjection(5, 2))
Δ = [0. -1.; 1. 0.; 2. 3.; 4. 5.; 6. 7.]
GeometricMachineLearning.Ω(E, Δ)

# output

5×5 SkewSymMatrix{Float64, Vector{Float64}}:
 0.0  -1.0  -2.0  -4.0  -6.0
 1.0   0.0  -3.0  -5.0  -7.0
 2.0   3.0   0.0  -0.0  -0.0
 4.0   5.0   0.0   0.0  -0.0
 6.0   7.0   0.0   0.0   0.0

Ω(Y::GrassmannManifold{T}, Δ::AbstractMatrix{T}) where T

using GeometricMachineLearning
E = GrassmannManifold(StiefelProjection(5, 2))
Δ = [0. 0.; 0. 0.; 2. 3.; 4. 5.; 6. 7.]
GeometricMachineLearning.Ω(E, Δ)

# output

5×5 SkewSymMatrix{Float64, Vector{Float64}}:
 0.0  -0.0  -2.0  -4.0  -6.0
 0.0   0.0  -3.0  -5.0  -7.0
 2.0   3.0   0.0  -0.0  -0.0
 4.0   5.0   0.0   0.0  -0.0
 6.0   7.0   0.0   0.0   0.0

References