<figure><img src="image_1.png" /></figure>

HAL
open science

# A geometric framework for asymptotic inference of principal subspaces in PCA

Dimbihery Rabenoro, Xavier Pennec

► To cite this version:

Dimbihery Rabenoro, Xavier Pennec. A geometric framework for asymptotic inference of principal subspaces in PCA. 2022. hal-03842125

HAL Id: hal-03842125

https://inria.hal.science/hal-03842125

Preprint submitted on 7 Nov 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.
# A geometric framework for asymptotic inference of principal subspaces in PCA

Dimbihery Rabenoro, Xavier Pennec

<figure><img src="image_1.png" /></figure>

## To cite this version:

Dimbihery Rabenoro, Xavier Pennec. A geometric framework for asymptotic inference of principal subspaces in PCA. 2024. hal-03842125v2

HAL Id: hal-03842125

https://inria.hal.science/hal-03842125v2

Preprint submitted on 12 Dec 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# A geometric framework for asymptotic inference of principal subspaces
in PCA

Dimbihery Rabenoro† and Xavier Pennec†

†Université Côte d'Azur and INRIA, Epione team

September 5, 2022

## Abstract

In this article, we develop an asymptotic method for testing hypothesis on the set of all linear subspaces arising from PCA and for constructing confidence regions for this set. This procedure is derived from intrinsic estimation in each Grassmannian, endowed with a structure of Riemannian manifold, to which each of these subspaces belong.

## 1 Introduction

Let $X$ be Gaussian a random vector valued in $\mathbb{R}^d$, $d \ge 1$, with covariance matrix $\Sigma$. Denote by $\lambda_1 > \dots > \lambda_r$ the different eigenvalues of $\Sigma$ and by $I = (q_i)_{1 \le i \le r}$ their respective multiplicities. The sequence $I$ induces a partition of the indices $(1, \dots, d)$ into $r$ parts $(B_i)_{1 \le i \le r}$, of respective sizes $(q_i)_{1 \le i \le r}$. For example, $B_1 = (1, \dots, q_1)$ and $B_r = (q_{r-1} + 1, \dots, d)$. In PCA, the principal subspaces (PS's) are defined as the eigenspaces of $\Sigma$, denoted by $(P_i(\Sigma))_{1 \le i \le r}$. Now, $\Sigma$ is estimated by the sample covariance matrices $(\hat{\Sigma}_n)_{n>1}$ from iid samples of $X$. Recall that almost surely, for all $n \ge 1$, $\hat{\Sigma}_n$ has $d$ distinct (random) eigenvalues $\mu_1(\hat{\Sigma}_n) > \dots > \mu_d(\hat{\Sigma}_n)$. To estimate the eigenspaces of $\Sigma$, these eigenvalues need to be regrouped according to the sequence $I$: the PS's are estimated by the eigenspaces $(P_i(\hat{\Sigma}_n))_{1 \le i \le r}$ which are the linear subspaces spanned by the eigenvectors of $\hat{\Sigma}_n$ associated to its eigenvalues in $M_n^i := \{\mu_k(\hat{\Sigma}_n) : k \in B_i\}$.

The landmark papers (Anderson, 1963) and (Tyler, 1981) investigated the inference of a single PS from the asymptotic distribution of statistics involving the eigenvectors of $\hat{\Sigma}_n$. The inference obtained for a PS of dimension $q = 1$ in (Anderson, 1963) was generalized to the case of one PS of any dimension $q \ge 1$ in (Tyler, 1981) using linear algebraic methods. In contrast, we introduce in this paper a geometric setting which allow us to derive the inference of all PS's together, for any sequence of dimensions. Namely, a given PS belongs to a Grassmannian, while the collection $F^I(\Sigma) := (P_i(\Sigma))_{1 \le i \le r}$ of them lies in a Flag manifold. Here, a flag is viewed as a collection $\mathcal{P} = (P_i)_{1 \le i \le r}$ of mutually orthogonal subspaces spanning $\mathbb{R}^d$, indexed by its type, that is the sequence $I = (q_i)_{1 \le i \le r}$ of dimensions of the $P_i$'s. This incremental subspaces representation is equivalent to the more classical one as properly nested sequence of linear subspaces. We denote by $\mathcal{F}^I$ the set of all flags of type I. In order to estimate the PS's, the idea is then to endow Grassmannians and flags of given type with a structure of Riemannian manifold, allowing to perform intrinsic estimation, based on Central Limit Theorems (CLT's), as in (Bhattacharya and Patrangenaru, 2005). It is well-known that $\mathcal{F}^I$ is diffeomorphic to $O(d)/O(I)$, where $O(d)$ is the orthogonal group of $\mathbb{R}^d$ and $O(I) := O(q_1) \times \dots \times O(q_r)$. While in (Tyler, 1981) each PS is estimated separately, our geometric approach enables us to estimate the whole flag $F^I(\Sigma)$ from the limiting distribution of the sample flag $F^I(\hat{\Sigma}_n) := (P_i(\hat{\Sigma}_n))_{1 \le i \le r}$.

In order to obtain Confidence Regions for $F^I(\Sigma)$ in $\mathcal{F}^I$, a first method could be to establish a CLT in $\mathcal{F}^I$ in normal coordinates, following (Bhattacharya and Patrangenaru, 2005). However, the Riemannian Logarithm and even the geodesic distance are unknown in closed form in $\mathcal{F}^I$. The representation of a flag by its incremental subspaces is a key point of our method: it allows embedding the flag space $\mathcal{F}^I$ in the product of Grassmannians $G^I$, for which the geometry is much simpler and the Logarithm is available in closed form. We establish CLT's
A geometric framework for asymptotic inference of
principal subspaces in PCA

Dimbihery Rabenoro and Xavier Pennec

Université Côte d'Azur and Inria, Epione Project Team, Sophia Antipolis
2004 route des Lucioles, 06902 Valbonne Cedex, France

Abstract

Consider data assumed to be iid samples from a multivariate Gaussian distribu-
tion whose covariance matrix has repeated eigenvalues. In our model, these eigen-
values are unknown but their multiplicities are supposed to be known. In this paper,
we develop an asymptotic method to infer the collection of all principal subspaces
together, i.e. the eigenspaces of this covariance matrix. Our approach is based on
the geometry of the flag manifold to which the collection of all principal subspaces
and our estimators of it belong.

# 1 Introduction

This article deals with the asymptotic study of principal component analysis (PCA), in the light of methods from Riemannian geometry. Throughout this paper, we assume that the data are iid samples from a Gaussian random vector X defined on a probability space (Ω, A, Pr) and valued in R^d, d ≥ 1, whose covariance matrix Σ has possibly repeated eigenvalues. In PCA, the principal subspaces (PS's) are the eigenspaces of Σ associated respectively to eigenvalues in decreasing order. The sequence of multiplicities of the latter is called the type of Σ. Then, a classical problem is to estimate the PS's by deriving Central Limit Theorems (CLT's) involving the eigenvectors of sample covariance matrices. This question is addressed notably in the Anderson's celebrated paper [1]. The matrix Σ may have repeated eigenvalues, although the set of symmetric matrices having multiple eigenvalues is negligible. In fact, it is implicit in [1] that Σ is obtained by equalizing beforehand nearly equal but distinct eigenvalues of a sample covariance matrix. In [13], a principled procedure is developed to equalize such eigenvalues and thus to select the type of Σ.

The method of [1] yields only an estimation of each PS separately and when its
dimension is equal to 1. This result is improved in [14] where other CLT's, issued from
Perturbation Theory [8], provide an estimation of each PS of any dimension. However,
each PS remains to be addressed separately. In contrast, our geometric method allows us
notably to estimate the collection of all PS together, regardless of their dimensions. Thus,
the main novelty of our paper is to build statistics based on the Riemannian geometry of
manifolds of linear subspaces of R^d. Namely, each PS belongs to a Grassmannian, while
the collection of them lies in a flag manifold. Here, a flag is viewed as a collection of
mutually orthogonal subspaces spanning R^d, whose sequence of dimensions is also called
the type of the flag. This incremental subspaces representation of a flag is equivalent to
the more usual one as nested sequence of linear subspaces. In our PCA setting, we denote
by λ₁ > ... > λr the distinct eigenvalues of Σ, of respective multiplicities q₁, ..., qr. Thus,
the collection of all PS's forms a flag of type I = (qᵢ)₁≤i≤r, denoted by F¹(Σ). Let F be the set of all flags of type I = (qᵢ)ᵢ, so that F¹(Σ) ∈ F. It is well-known that F is endowed with a structure of homogeneous space for the orthogonal group O(d): see [12].
The introduction of such manifolds in a PCA setting is notably inspired by [11], where
it is proved that PCA can be rephrased as an optimization on a flag manifold. We also
refer to [9] for the use of flag manifolds for studying PCA.
in Grassmannians for each $P_i(\hat{\Sigma}_n)$ and we prove that their rate of convergence to a Gaussian distribution are of order $\frac{1}{\sqrt{n}}$. Thanks to the fact that the tangent space to all Grassmannians can be embedded in the same symmetric matrix space with the projector representation, we then derive the estimation of the whole flag $F^1(\Sigma)$. Namely, we prove a result of the form

$$n\Delta_{\text{Spec}(\Sigma)}(F^1(\Sigma), F^1(\hat{\Sigma}_n)) \xrightarrow[n \to \infty]{d} \chi_D^2, \qquad (1.1)$$

for some $D > 0$. Here, $\Delta_{\text{Spec}(\Sigma)}(\cdot, \cdot)$ can be viewed as a discrepancy function on $\mathcal{F}^1$, indexed by the spectrum of $\Sigma$. The convergence in (1.1) is proved rigorously and checked numerically on synthetic experiments. In an estimation setting, the unknown eigenvalues of $\Sigma$ are replaced by the block-means of eigenvalues of the empirical matrix $\hat{\Sigma}_n$. Thanks to the consistency of this estimator, the convergence to the $\chi_D^2$ distribution is preserved, and this provides confidence regions for $F^1(\Sigma)$.

This article is organized as follows. In Section 2, we develop the background of geometry of Grassmannians and Flag manifolds needed for the asymptotic inference methods derived in Section 3.

# 2 Grassmannian and Flags

For Riemannian manifolds, we refer to (Lee, 1963). Our reference for the geometry of Grassmannians is (Bendokat Zimmermann and Absil, 2020) and for that of Flag manifolds is (Ye Wong and Lim, 2022). The key formula for the closed form for the Riemannian Logarithm for Grassmannians (Theorem 1 below) is due to (Batzies Huper Machado and Silva Leite, 2015).

## 2.1 The Grassmannian manifold

For $1 \le q < d$, the Grassmannian $G(q, d)$ is by definition the set of all linear subspaces of dimension $q$ of $\mathbb{R}^d$. In order to perform calculations, it is identified with the set of orthogonal projectors of rank $q$:

$$G(q, d) \simeq \{ P \in \mathbb{R}^{d \times d} : P' = P, P^2 = P, \text{rank}(P) = q \}. \qquad (2.1)$$

Then, we introduce the Stiefel manifold $\text{St}(q, d)$. As a set, it is the set of all $d \times q$ matrices whose columns are orthonormal vectors. Denoting by $I_d$ the identity matrix of order $d$,

$$\text{St}(q, d) = \{ Y \in \mathbb{R}^{d \times q} : Y'Y = I_d \}. $$

The sets $G(q, d)$ and $\text{St}(q, d)$ are linked by the map $\pi_q^{\text{SG}}: \text{St}(q, d) \to G(q, d)$, which associates to any set of $q$ orthonormal vectors, the projector whose range is the subspace which they span, i.e.:

$$\pi_q^{\text{SG}}(Y) = YY'.$$

### 2.1.1 Action of the orthogonal group

The orthogonal group $O(d)$ acts on $G(q, d)$ by

$$(Q, P) \mapsto QPQ', \quad \text{for } Q \in O(d) \text{ and } P \in G(q, d). \qquad (2.2)$$

Let $P_0 \in G(q, d)$ whose range is spanned by the $q$ first vectors of the standard basis of $\mathbb{R}^d$. Let $\pi_q^{\text{OG}}$ be the map defined by

$$\pi_q^{\text{OG}}: O(d) \to G(q, d), \quad \pi_q^{\text{OG}}(Q) = QP_0Q'.$$

Fix $P \in G(q, d)$. Then, for any $Q \in O(d)$ whose $q$ first columns span the range of $P$, we have that $\pi_q^{\text{OG}}(Q) = P$. So, the orbit of $P_0$ is the whole $G(q, d)$. On the other hand, the stabilizer of $P_0$ is identified with $O(q) \times O(d-q)$, so that $\pi_q^{\text{OG}}$ induces a bijection

$$\tilde{\pi}_q^{\text{OG}}: O(d)/(O(q) \times O(d-q)) \simeq G(q, d).$$

Now, $O(d)/(O(q) \times O(d-q))$ is endowed with the manifold structure such that the canonical quotient map from $O(d)$ to $O(d)/(O(q) \times O(d-q))$ is a smooth submersion, i.e. with surjective differential at every point. Then, the manifold structure on $G(q, d)$ is the one for which $\tilde{\pi}_q^{\text{OG}}$ is a diffeomorphism.
A first idea to estimate $F^1(\Sigma)$ could be to exploit the structure of Riemannian manifold on $\mathcal{F}^1$ introduced in [12] and, following [3], to establish a CLT in $\mathcal{F}^1$ in normal coordinates i.e. defined from the Riemannian Logarithm on $\mathcal{F}^1$. However, the latter is *unknown* in closed form. In fact, there is no known metric on $\mathcal{F}^1$ for which an explicit expression of the Riemannian Logarithm is available. Here, the incremental subspaces representation of a flag allows us to overcome this lack, since it provides an embedding of $\mathcal{F}^1$ within a product of Grassmannians. Indeed, the Riemannian Logarithm is available in closed form in any Grassmannian.

In fact, our method is a wide geometric extension of that of [1]. For $n \ge 1$, let $\hat{\Sigma}_n$ be the sample covariance matrix from an iid sample of size $n$ of $X$. Then, in [1], a matrix $E_n \in O(d)$ of the form $E_n = \Gamma'C_n$ is considered, where $\Gamma$ and $C_n$ are respectively matrices of eigenvectors of $\Sigma$ and $\hat{\Sigma}_n$, such that all diagonal entries of $E_n$ are positive (See Section 2 for precise definitions). So, $E_n$ measures the deviation between the PS's and the eigenvectors of $\hat{\Sigma}_n$. Then, the asymptotic distribution of $E_n$ is described in [1]. Thus, the diagonal blocks of size $q_i \times q_i$ of $E_n$ converge in distribution to a Haar measure on the set of orthogonal $q_i \times q_i$ matrices whose diagonal entries are positive. The other blocks of $E_n$ converge to Gaussian distributions, which provide Central Limit Theorems (CLT's) that allow to infer the PS's.

Instead of considering the usual coordinates of $E_n$, we introduce a geometric splitting of $E_n$ into two parts, for which we recover respectively Gaussian distributions and Haar measures as $n \to \infty$. The first part yields CLT's from which we deduce our estimation of $F^1(\Sigma)$. The second part provides statistics that are valued, for all $n \ge 1$, in the orthogonal groups $(O(q_i))_i$, contrarily to the diagonal blocks of $E_n$. We derive such a geometric splitting of any orthogonal matrix from the structure of Riemannian submersion of the map which associates to a $q$-frame the $q$-linear subspace spanned by it, for any $1 \le q \le d$.

We describe briefly hereafter our method to infer $F^1(\Sigma)$. Thanks to the explicit expression of the Riemannian Logarithm in any Grassmannian, we establish a CLT in a Grassmannian for each PS separately. Then, contrarily to the CLT's for each PS in [14], some geometric properties of Grassmannians allow us to concatenate our CLT's to obtain a pivotal statistic for $F^1(\Sigma)$, whose limiting distribution is a $\chi^2$ one. Namely, recall that a.s., for all $n \ge 1$, the eigenvalues of $\hat{\Sigma}_n$ are all distinct. Then, one can construct a flag of type I from the eigenvectors of $\hat{\Sigma}_n$ as follows: those associated to the $q_1$ largest eigenvalues are grouped together, then those to the next $q_2$ ones, ..., and finally those to the $q_r$ smallest ones. Then, the linear subspaces generated by each group form a flag $F^1(\hat{\Sigma}_n) \in \mathcal{F}^1$ which may be compared to $F^1(\Sigma)$. Thus, we prove a result of the form

$$ \frac{n}{4} \left[ \tilde{\mathcal{D}}_{\hat{\mathbb{K}}_n} \left( F^1(\Sigma), F^1(\hat{\Sigma}_n) \right) \right]^2 \xrightarrow{n \to \infty} \chi_{\mathrm{D}^1}^2, \quad (1.1) $$

where $\tilde{\mathcal{D}}_{\hat{\mathbb{K}}_n}(\cdot, \cdot)$ is a discrepancy function on $\mathcal{F}^1$, indexed by a statistic $\hat{\mathbb{K}}_n$ which is a function of the spectrum of $\hat{\Sigma}_n$, and $\mathrm{D}^1 = \frac{1}{2}(d^2 - \sum_i q_i^2)$. The convergence in Eq. (1.1) is illustrated numerically on synthetic experiments. Finally, we derive, from Eq. (1.1), confidence regions for $F^1(\Sigma)$ that are interpreted as confidence ellipsoids in $\mathcal{F}^1$, centered at $F^1(\hat{\Sigma}_n)$.

Thus, our confidence regions are defined directly in the manifold $\mathcal{F}^1$. In contrast, the confidence regions usually met in the literature on estimation on manifolds are defined in a chart, as in [3]. Furthermore, in [14], for fixed $1 \le i \le r$, the confidence regions for the principal subspace $P_i(\Sigma)$ are based on a $\chi^2_{q_i(d-q_i)}$ distribution. Then, the sum of their degrees of freedom is $\sum_i q_i(d-q_i) = 2\mathrm{D}^1$, so that we obtain a gain of a factor of 2.

This article is organized as follows. In Section 2, we present Anderson's results and a brief sketch of their proofs. In Section 3, we develop the background of geometry needed for the rest of the paper and we present notably the aforementioned geometric splitting of $E_n$. Then, in Section 4, we state our main result in Theorem 5, which provides the
### 2.1.2 Tangent space

By definition, the Lie algebra of $O(d)$ is the tangent space $T_{I_d}O(d)$. It is the set $\mathfrak{so}(d)$ of skew-symmetric matrices. Now, for any $Q \in O(d)$,

$$T_Q O(d) = \{Q \tilde{\Omega} : \tilde{\Omega} \in \mathfrak{so}(d)\}.$$

Then, the tangent space at $P \in G(q, d)$ is the linear subspace of $\mathrm{Sym}_d$ defined by

$$T_P G(q, d) = \{\Delta \in \mathrm{Sym}_d : \Delta P + P \Delta = \Delta\}.$$

## 2.2 Riemannian geometry of the Grassmannian

### 2.2.1 Metric

A *metric* on a smooth manifold $M$ is a collection $g = (g_p)_{p \in M}$ of inner products on the tangent spaces $T_p M$, varying smoothly wrt $p$. We say that $(M, g)$ is a *Riemannian manifold* and we assume in the sequel that $M$ is *connected*. For $v \in T_p M$, its length is $\|v\|_p := (g_p(v, v))^{1/2}$. For $p, q \in M$, a $C^1$ curve between $p$ and $q$ is a $C^1$ map $\gamma : [a, b] \to M$ such that $\gamma(a) = p$ and $\gamma(b) = q$. Its length is $\int_a^b \|\frac{d}{dt}\gamma(t)\|_{\gamma(t)} dt$. Then, $g$ defines the geodesic distance $d_g$, where for any $p, q \in M$, $d_g(p, q)$ is the infimum of lengths of all $C^1$ curves between $p$ and $q$.

Let $\pi : M \to N$ be a smooth submersion. A metric $g$ on $M$ allows to split any tangent space $T_p M$ into a vertical part $\mathrm{Ver}_p^\pi M := \ker(d\pi_p)$ and a horizontal part $\mathrm{Hor}_p^\pi M$, which is the orthogonal complement of $\mathrm{Ver}_p^\pi M$ in $T_p M$ wrt $g_p$. We say that $\pi$ is a *Riemannian submersion* when, for all $p \in M$, the restriction of $d\pi_p$ to $\mathrm{Hor}_p^\pi M$ is an isometry.

As a compact Lie group, $O(d)$ is endowed with its canonical metric $g^O$ defined by

$$g_Q^O(\Omega_1, \Omega_2) = \frac{1}{2} \mathrm{tr}(\Omega'_1 \Omega_2), \quad \text{for } Q \in O(d), \Omega_i \in T_Q O(d).$$

We endow $G(q, d)$ with the unique metric $g^G$ for which $\pi^{OG}$ is a *Riemannian submersion*. Thus, for $P \in G(q, d)$ and $\Delta_1, \Delta_2 \in T_P G(q, d)$,

$$g_P^G(\Delta_1, \Delta_2) = g_Q^O(\Delta_{1,Q}^h, \Delta_{2,Q}^h)$$

where $\Delta_{i,Q}^h$ is the horizontal lift at $Q$ of $\Delta_i$, i.e. $\Delta_{i,Q}^h \in \mathrm{Hor}_Q O(d)$ and $d\pi_Q^{OG}(\Delta_{i,Q}^h) = \Delta_i$.

**Lemma 1.** The action of $O(d)$ on $G(q, d)$ defined above is an isometry. This means that, for any $P, R \in G(q, d)$ and $Q \in O(d)$,

$$d_g^G(QPQ', QRQ') = d_g^G(P, R),$$

where $d_g^G$ denotes the geodesic distance on $G(q, d)$ associated to $g^G$.

### 2.2.2 Geodesics

Let $\gamma$ be a smooth curve in $M$. If $M$ is a submanifold of $\mathbb{R}^n$, then the acceleration $\ddot{\gamma}(t)$ is not always a tangent vector at $\gamma(t)$. On a Riemannian manifold $(M, g)$, one builds from the metric an operation of derivation of vector fields whose outputs are tangent vectors, called the *covariant derivative* and denoted by $\nabla$. Thus, if $\gamma$ is in $(M, g)$, the acceleration vector of $\gamma$ at $\gamma(t)$ is the vector of $T_{\gamma(t)} M$ equals to the covariant derivative of the velocity $\dot{\gamma}$ at $\gamma(t)$ in the direction of $\dot{\gamma}(t)$, denoted by $(\nabla_{\dot{\gamma}}\dot{\gamma})_{\gamma(t)}$. Then, a *geodesic* is a smooth curve $\gamma$ on $M$ of zero acceleration, i.e. such that

$$\nabla_{\dot{\gamma}} \dot{\gamma} = 0.$$

This equation is locally an ODE, so that for any $p \in M$ and $v \in T_p M$, there exists a unique geodesic $\gamma_{p,v}$ such that $\gamma_{p,v}(0) = p$ and $\dot{\gamma}_{p,v}(0) = v$. By the *Hopf-Rinow theorem*, the metric space $(M, d_g)$ is complete if and only if all geodesics are defined on $\mathbb{R}$. In that case, for any $p, q \in M$, there exists at least one geodesic of minimal length between $p$ and $q$. Given a Riemannian submersion $\pi : M \to N$, the geodesics in $N$ are the images by $\pi$ of the geodesics with horizontal tangent vectors in $M$. Since $\pi^{OG}$ is a Riemannian submersion, the geodesics in $G(q, d)$ are given by those in $O(d)$. The geodesic in $O(d)$ starting at $Q$ with initial velocity $Q\tilde{\Omega}$ is

$$t \mapsto Q \exp_m(t\tilde{\Omega}),$$
limiting distribution of the statistics issued from our splitting of $E_n$ and we deduce the estimation of the flag $F^1(\Sigma)$. Finally, in Section 5, we summarize our results and propose some open questions raised by them. The proof of Theorem 5 is deferred to the Appendix.

## 2 Review of Anderson's results

We keep the notations of the Introduction. We fix a matrix $\Gamma \in O(d)$ of eigenvectors of $\Sigma$, associated respectively to $\lambda_1, ..., \lambda_r$. Then, $\Gamma'\Sigma\Gamma = \Delta$ where $\Delta$ is the following $d \times d$ block-diagonal matrix: $\Delta := \text{Diag}(\lambda_1 I_{q_1}, ..., \lambda_r I_{q_r})$.

### 2.1 Preliminaries

For $1 \le i \le r - 1$, set $\bar{q}_i := \sum_{j=1}^i q_j$. The sequence $I = (q_i)_{1 \le i \le r}$ induces a partition of the indices $\{1, ..., d\}$ into $r$ blocks $(\beta_i)_{1 \le i \le r}$ as follows:

$$ \beta_1 = \{1, ..., q_1\}, \dots, \beta_i = \{\bar{q}_i + 1, ..., \bar{q}_{i+1}\}, \dots, \beta_r = \{\bar{q}_{r-1} + 1, ..., d\}. $$

The sequence $(\beta_i)_{1 \le i \le r}$ is called the partition of $\{1, ..., d\}$ wrt I. Let $A = (a_{kl})_{1 \le k \le l \le d}$ be a $d \times d$ matrix. For fixed $1 \le i \le j \le r$, consider the submatrix $A^{(i,j)} := (a_{kl})_{k \in \beta_i, l \in \beta_j}$. Let $A^{(i)}$ be the $i$-th block of columns of $A$, i.e. $A^{(i)} = (A^{(1,i)} \dots A^{(r,i)})'$. Thus, setting $I_d^{(i)} := (I_d)^{(i)}$,

$$ A^{(i)} = A I_d^{(i)} \quad \text{and} \quad A^{(i,j)} = (I_d^{(i)})' A I_d^{(j)}. \qquad (2.1) $$

**Definition 1.** Let $\text{Sym}_d$ be the set of all $d \times d$ symmetric matrices and $\text{Sym}_d^{\pm}$ the subset of $\text{Sym}_d$ of matrices whose eigenvalues are all distinct. For $S \in \text{Sym}_d^{\pm}$, its spectrum is denoted by $\{\mu_1(S) > \dots > \mu_d(S)\}$. Let $\psi: \text{Sym}_d^{\pm} \to O(d)$ be the map such that for any $1 \le k \le d$, the $k$-th column of $\psi(S)$ is the eigenvector of norm 1 associated to $\mu_k(S)$ whose $k$-th entry is non-negative. Then, the map $\psi$ is measurable and is called the eigenvector map.

### 2.2 Statement of Anderson's results

Assume that $\hat{\Sigma}_n \in \text{Sym}_d^{\pm}$, which holds a.s. Set

$$ E_n := \psi(T_n) \quad \text{where} \quad T_n := \Gamma' \hat{\Sigma}_n \Gamma, \quad n \ge 1. $$

$E_n$ is well-defined, since $T_n \in \text{Sym}_d^{\pm}$. Then, $C_n := \Gamma E_n$ is an orthogonal matrix of eigenvectors of $\hat{\Sigma}_n$, also associated to eigenvalues in decreasing order. Thus, $E_n = \Gamma' C_n$ is the matrix mentioned in the Introduction, which measures the deviation between eigenvectors of $\Sigma$ (the PS's) and those of $\hat{\Sigma}_n$. The asymptotic distribution of $E_n = \Gamma' C_n$ is the main result of [1], stated in Theorem 1 below, for which we need the following definition.

**Definition 2.** The Haar measure is the unique probability distribution on the compact group $O(d)$ that is invariant under right-multiplications. It is the uniform measure on $O(d)$. The conditional Haar invariant distribution is defined in [1] as $2^d$ times the restriction of the Haar measure to the space of orthogonal matrices whose diagonal entries are non-negative.

**Theorem 1.** (i) For $1 \le i \le r$, $E_n^{(i,i)} \xrightarrow{n \to \infty} E^{i,i}$, where the distribution of $E^{i,i}$ is the conditional Haar invariant distribution.

(ii) For $n \ge 1$, set $F_n := \sqrt{n}E_n$. Then, for all $i \ne j$, $F_n^{(i,j)} \xrightarrow{n \to \infty} F^{i,j}$, where the entries of $F^{i,j}$ are iid rv's $N(0, \sigma_{i,j}^2)$, of standard deviation $\sigma_{i,j} := \frac{\sqrt{\lambda_i \lambda_j}}{|\lambda_i - \lambda_j|}$.

(iii) The blocks $\{E^{i,i}, F^{i,j} : 1 \le i \le r, 1 \le j \le r, i \ne j\}$ are mutually independent.
where $\exp_m$ is the matrix exponential. So, the geodesic in $G(q, d)$ starting at $P$ with initial velocity $\Delta$ is

$$t \mapsto \pi^{OG}(Q \exp_m(t\tilde{\Omega})) = Q \exp_m(t\tilde{\Omega})P_0 \exp_m(-t\tilde{\Omega})Q', \quad (2.3)$$

for any $Q \in (\pi^{OG})^{-1}(P)$ and where $Q\tilde{\Omega} = \Delta_Q^h$. Therefore, $(G(q, d), g^O)$ is complete.

### 2.2.3 Exponential map

In the sequel, $(M, g)$ is a complete Riemannian manifold. For any $p \in M$, the Riemannian exponential at $p$ is the map $\text{Exp}_p^M$ from $T_pM$ to $M$ defined by

$$\text{Exp}_p^M(v) = \gamma_{p,v}(1).$$

**Lemma 2.** By (2.3), for any $P \in G(q, d)$ and $\Delta \in T_pG(q, d)$,

$$\text{Exp}_p^{G(q,d)}(\Delta) = \exp_m([\Delta, P])P \exp_m(-[\Delta, P]).$$

### 2.2.4 Cut locus

It is known that, for all $p \in M$, any geodesic starting at $p$ is locally length-minimizing. The cut locus of $p \in M$ is the set $\text{Cut}(p)$ of points $q$ such that the geodesics starting at $p$ cease to be length-minimizing beyond $q$. It is established that if $q \notin \text{Cut}(p)$, then there exists a unique length-minimizing geodesic between $p$ and $q$. For example, on a sphere, the geodesics are the great circles, so that two antipodal points are each other's cut locus.

**Lemma 3.** Let $P \in G(q, d)$. Then, writing $P = YY'$,

$$\text{Cut}(P) = \{R = ZZ' \in G(q, d) : \text{rank}(Y'Z) < q\}.$$

This means that this characterization is independent of $Y \in \text{St}(q, d)$ such that $P = YY'$.

### 2.2.5 Riemannian Logarithm

Let $p \in M$ and $q \in M \setminus \text{Cut}(p)$. Then, there exists a unique shortest geodesic between $p$ and $q$, whose initial velocity is thus the smallest tangent vector $v \in T_pM$ such that $\text{Exp}_p^M(v) = q$. It is called the Riemannian Logarithm of $q$ at $p$, denoted by $\text{Log}_p^M(q)$. We have that

$$d_g(p, q) = \|\text{Log}_p^M(q)\|_p. \quad (2.4)$$

For $G(q, d)$, combining (2.4), Lemmas 3 and 1 yield the following result.

**Proposition 1.** For all $P, R \in G(q, d)$ and $Q \in O(d)$,

$$QRQ' \notin \text{Cut}(QPQ') \iff R \notin \text{Cut}(P), \quad (2.5)$$

and in that case,

$$\left\| \text{Log}_{QPQ'}^{G(q,d)} (QRQ') \right\|_{QPQ'} = \left\| \text{Log}_P^{G(q,d)} (R) \right\|_P. \quad (2.6)$$

Finally, the next result provides the closed form for the Riemannian Logarithm for $G(q, d)$.

**Theorem 1.** Let $P \in G(q, d)$ and $R \in G(q, d) \setminus \text{Cut}(P)$. Then, $\text{Log}_P^{G(q,d)}(R) = [\Omega, P]$, where $\Omega$ is determined by

$$\Omega = \frac{1}{2} \log_m ((I_d - 2R)(I_d - 2P)),$$

and $\log_m$ denotes the matrix logarithm.

## 2.3 The set of flags of fixed type

Let $I = (q_i)_{1 \le i \le r}$ be a sequence of positive integers such that $\sum_{i=1}^r q_i = d$.
**Remark 1.** As pointed out in the Introduction, Theorem 1 allows to estimate separately each PS's, but only when its dimension is equal to 1: See Appendix B in [1].

## 2.3 Method for the proof of Theorem 1

### 2.3.1 Preliminary results

The method consists of propagating the uncertainty of the estimation of $\hat{\Sigma}_n$ by $\Sigma$, provided by Theorem 2 hereafter which is obtained in [1], to derive the limiting distribution of the deviation between their respective eigenvectors measured by $E_n = \Gamma'C_n$. This propagation is realized by a generalized $\delta$-method stated in Theorem 3 below which can be found in Section 7 of [1]. A more actual reference is section 18.11 in [15]. In fact, our version extends slightly the latter by adding Condition (2.2) in the statement, because of constraints imposed by our geometric method.

**Theorem 2.** Set $T_n := \Gamma'\hat{\Sigma}_n\Gamma$. Then, the following CLT in $\text{Sym}_d$ holds.

$$U_n := \sqrt{n}(T_n - \Delta) \xrightarrow{n \to \infty} U,$$

where $U$ is a random matrix whose distribution is characterized as follows: $U \in \text{Sym}_d$, the blocks $\{U^{(i,j)}: 1 \le i < j \le r\}$ are mutually independent and for $i \neq j$, the entries of $U^{(i,j)}$ are iid rv's $N(0, s_{i,j}^2)$, of standard deviation $s_{i,j} := \sqrt{\lambda_i \lambda_j}$.

**Theorem 3.** Let $\mathcal{S}, \mathcal{T}$ be metric spaces endowed with their Borel $\sigma$-algebras. Let $(\mathcal{S}_n)_{n \ge 0}$ be a sequence of measurable subsets of $\mathcal{S}$ and for all $n \ge 0$, $g_n: \mathcal{S}_n \to \mathcal{T}$ a measurable map. Let $(V_n)_{n \ge 1}$ be a sequence of rv's valued in $(\mathcal{S}_n)_{n \ge 1}$ which converges in distribution to a rv $V$ valued in $\mathcal{S}$, with $\text{Pr}(V \in \mathcal{S}_0) = 1$. Assume that for all $n \ge 1$, there exists a measurable subset $\mathcal{S}_n^*$ of $\mathcal{S}_n$ such that

$$\text{Pr}(V_n \in \mathcal{S}_n^*) \xrightarrow{n \to \infty} 1 \qquad (2.2)$$

and that, for any $v \in \mathcal{S}_0$ and all sequence $(v_n)_{n \ge 1}$ valued in $(\mathcal{S}_n^*)_{n \ge 1}$,

$$v_n \xrightarrow{n \to \infty} v \implies g_n(v_n) \xrightarrow{n \to \infty} g_0(v). \qquad (2.3)$$

Then, $g_n(V_n)$ converges in distribution to $g_0(V)$.

*Proof.* We prove our extended version in subsection 6.5 of the Appendix. $\square$

# 3 Preliminaries of geometry

In subsection 3.1, we introduce the manifolds of linear subspaces in which the (collection of) PS's and their estimators lie. In subsection 3.2, we define a geometric parametrization of $E_n$ that is an alternative to that by its coordinates in the standard basis of $\mathbb{R}^{d \times d}$, used in Anderson's paper [1] to describe (Theorem 1) the limiting distribution of $E_n$. In contrast, our parametrization of $E_n$ provides notably confidence regions for the collection of PS's, for which geometric preliminaries are presented in subsection 3.3.

## 3.1 Grassmannians, Stiefel and flag manifolds

### 3.1.1 Sets of linear subspaces

Let $1 \le q \le d$. The Stiefel manifold is the set $\text{St}_q$ of all $q$-frames of $\mathbb{R}^d$ i.e. of all $q$-tuples of orthonormal vectors of $\mathbb{R}^d$. Thus,

$$\text{St}_q = \{U \in \mathbb{R}^{d \times q} : U'U = I_q\}.$$
### 2.3.1 Notations

The sequence I induces a partition of the indices {1, ..., d} into r parts $(B_i)_{1\le i\le r}$, i.e. $B_1 = (1, ..., q_1)$, $B_r = (q_{r-1} + 1, ..., d)$ and for $1 < i < r$, $B_i = (\bar{q}_{i-1} + 1, ..., \bar{q}_i)$, where $\bar{q}_i := \sum_{j=1}^{i} q_j$. Let A be a $d \times d$ matrix. If A is defined by its columns, i.e. $A = [A_1, ..., A_d]$, then, for all $1 \le i \le r$, set

$$A^{(i)} := [A_k : k \in B_i].$$

If A is defined by its entries, i.e. $A = (a_{k\ell})_{1\le k\le\ell\le d}$, then, for all $1 \le i \le j \le r$, set

$$A^{(ij)} := (a_{k\ell})_{\substack{k \in B_i \\ \ell \in B_j}}.$$

Then, setting $I_d^{(i)} := (I_d)^{(i)}$, we have that $A^{(i)} = AI_d^{(i)}$ and $A^{(ij)} = (I_d^{(i)})' AI_d^{(j)}$.

### 2.3.2 Definitions

We recall that a flag of type I is a collection $\mathcal{V} = (V_i)_{1\le i\le r}$ of r mutually orthogonal linear subspaces of $\mathbb{R}^d$ such that for all $1 \le i \le r$, $\dim(V_i) = q_i$. Thus, the set $\mathcal{F}^I$ of flags of type I is a subset of

$$G^I := \prod_{i=1}^{r} G^i \quad \text{where} \quad G^i := G(q_i, d).$$

Identifying Grassmannians with projectors as in (2.1),

$$\mathcal{F}^I \simeq \{ \mathcal{P} = (P_i)_{1\le i\le r} \in G^I : P_i P_j = 0, i \neq j \}.$$

Introduce the standard flag $\mathcal{P}_0^I = (P_0^i)_{1\le i\le r}$, where $P_0^i \in G^i$ is the block-diagonal matrix defined by

$$P_0^i := \text{Diag}[0_{q_1}, \dots, I_{q_i}, \dots, 0_{q_r}],$$

and $0_{q_i}$ is the null matrix of order $q_i$. Then,

$$P_0^i = I_d^{(i)} (I_d^{(i)})'. \qquad (2.7)$$

### 2.3.3 Action of the orthogonal group

The action of $O(d)$ on Grassmannians induces an action of $O(d)$ on $\mathcal{F}^I$, defined by

$$(Q, \mathcal{P}) \mapsto (QP_iQ')_{1\le i\le r} \in \mathcal{F}^I,$$

for $Q \in O(d)$ and $\mathcal{P} = (P_i)_{1\le i\le r} \in \mathcal{F}^I$. Indeed, $O(d)$ acts naturally on $G^I$ and this action preserves the mutual orthogonality of linear subspaces. Then, define the map $\pi^I : O(d) \to \mathcal{F}^I$ by

$$\pi^I(Q) = (QP_0^iQ')_{1\le i\le r}.$$

**Proposition 2.** The orbit of $\mathcal{P}_0^I$ under this action is the whole $\mathcal{F}^I$ and the stabilizer of $\mathcal{P}_0^I$ is the group $O(I)$, where

$$O(I) := \{ \text{Diag}[H_1, \dots, H_r] : H_i \in O(q_i), 1 \le i \le r \} = \prod_{i=1}^{r} O(q_i).$$

In other words, the map $\pi^I$ is surjective and for all $Q \in O(d)$ and $H \in O(I)$, $\pi^I(QH) = \pi^I(Q)$. Therefore, $\pi^I$ induces a bijection

$$\tilde{\pi}^I : O(d)/O(I) \simeq \mathcal{F}^I.$$

such that $\pi^I = \tilde{\pi}^I \circ \bar{p}^I$, where $\bar{p}^I$ is the canonical quotient map from $O(d)$ to $O(d)/O(I)$.
The Grassmannian $G_q$ is the set of all $q$-linear subspaces of $\mathbb{R}^d$, identified with the orthogonal projectors of rank $q$ i.e. such a projector $P$ is identified with its range $rg(P)$:

$$G_q = \{P \in \mathbb{R}^{d \times d} : P' = P, P^2 = P, \mathrm{rank}(P) = q\}. \quad (3.1)$$

Thus, $G_q$ is a subset of $\mathrm{Sym}_d$. Then, $\mathrm{St}_q$ and $G_q$ are linked by the map $\pi_q$ which associates to $\mathcal{U} \in \mathrm{St}_q$, the (projector onto the) linear subspace spanned by $\mathcal{U}$, defined by

$$\pi_q : \mathrm{St}_q \to G_q, \quad \pi_q(\mathcal{U}) = \mathcal{U}\mathcal{U}'.$$

**Remark 2.** Let $P \in G_q$ and $\mathcal{U}_P$ a fixed $q$-frame spanning $rg(P)$, i.e. $\mathcal{U}_P \in \pi_q^{-1}(P)$. Then, any $\mathcal{V} \in \pi_q^{-1}(P)$ is obtained from $\mathcal{U}_P$ by a unique orthonormal base change in $rg(P)$. Thus, there exists a unique $K \in O(q)$ such that $\mathcal{V} = \mathcal{U}_P K$ i.e. $K = \mathcal{U}_P' \mathcal{V}$.

**Definition 3.** (i) A flag of $\mathbb{R}^d$ is a sequence of mutually orthogonal linear subspaces spanning $\mathbb{R}^d$, whose sequence of respective dimensions is called the type of the flag.

(ii) A flag whose type is $(1, 1, ..., 1)$ is called a full flag.

(iii) For $\Sigma \in \mathrm{Sym}_d$, its flag of eigenspaces is the collection of its eigenspaces associated to eigenvalues in decreasing order, whose type is the sequence of their dimensions.

(iv) Given a sequence $I = (q_i)_{1 \le i \le r}$ with $\sum q_i = d$, consider the product manifold

$$G^I := \prod G^i \quad \text{where} \quad G^i := G_{q_i}, \quad 1 \le i \le r. \quad (3.2)$$

The flag manifold of type I is the set $\mathcal{F}^I$ of all flags of type I, identified with a subset of $G^I$:

$$\mathcal{F}^I = \{ \mathcal{P} = (P_i)_{i \in G^I} : \sum P_i = \mathbb{R}^d, P_i P_j = 0, i \neq j \}. \quad (3.3)$$

**Remark 3.** This definition of a flag is more suitable for computing with flags of eigenspaces than the usual equivalent one as a nested sequence of linear subspaces.

In the sequel, for $1 \le i \le r$, $\mathrm{St}_{q_i}$ and $\pi_{q_i}$ are respectively denoted by $\mathrm{St}^i$ and $\pi^i$.

**Definition 4.** (i) Let $\mathrm{St}_\perp$ be the set of all collections $\bar{\mathcal{U}} = (\mathcal{U}_i)_{1 \le i \le r}$ of frames such that for $1 \le i \le r$, $\mathcal{U}_i \in \mathrm{St}^i$ and the collection $(\pi^i(\mathcal{U}_i))_i$ lies in $\mathcal{F}^I$. Define the map $(\pi^i)_i$ by

$$(\pi^i)_i : \mathrm{St}_\perp^I \to \mathcal{F}^I \text{ with } (\pi^i)_i(\bar{\mathcal{U}}) := (\pi^i(\mathcal{U}_i))_i.$$

(ii) We say that $\Gamma \in O(d)$ spans a flag $\mathcal{P} = (P_i)_{i \in \mathcal{F}^I}$ when for all $1 \le i \le r$, $\pi^i(\Gamma^{(i)}) = P_i$.

(iii) For $Q \in O(d)$, the collection $\zeta(Q) := (Q^{(i)})_{1 \le i \le r}$ lies in $\mathrm{St}_\perp^I$, which defines a diffeomorphism $\zeta : O(d) \to \mathrm{St}_\perp^I$.

### 3.1.2 Actions of the orthogonal group

The group $O(d)$ acts on vectors of $\mathbb{R}^d$ by linear isometries, which induces an action on $\mathrm{St}_q$ by left-multiplication:

$$(Q, \mathcal{U}) \mapsto Q\mathcal{U}, \quad Q \in O(d), \mathcal{U} \in \mathrm{St}_q. \quad (3.4)$$

Let $P \in G_q$ and $\mathcal{U}_P$ a $q$-frame spanning $rg(P)$. Any $Q \in O(d)$ sends $\mathcal{U}_P$ to $Q\mathcal{U}_P$, which spans $rg(\pi_q(Q\mathcal{U}_P)) = rg(QPQ')$. The latter is independent of $\mathcal{U}_P$. So, $O(d)$ acts on $G_q$ by

$$(Q, P) \mapsto Q \cdot P := QPQ', \quad Q \in O(d), P \in G_q. \quad (3.5)$$

The action on $G_q$ extends to the product $G^I$: for $Q \in O(d)$ and $\mathcal{P} = (P_i)_{1 \le i \le r} \in G^I$, set

$$Q * \mathcal{P} := (Q \cdot P_i)_{1 \le i \le r} \in G^I. \quad (3.6)$$

Now, the action of $Q \in O(d)$ preserves the mutual orthogonality of linear subspaces. Thus, $\mathcal{P} \in \mathcal{F}^I \implies Q * \mathcal{P} \in \mathcal{F}^I$. So, the action defined in (3.6) induces an action of $O(d)$ on $\mathcal{F}^I$.
Now, let $\pi_I^{SG}$ be the map from $\mathrm{St}^I := \prod_{i=1}^r \mathrm{St}(q_i, d)$ to $\mathbb{G}^I$ defined by

$$ \mathcal{Y} = (Y_i)_{1 \le i \le r} \mapsto (\pi_{q_i}^{SG}(Y_i))_{1 \le i \le r}. $$

Set $\mathrm{St}_{\perp}^I = \{\mathcal{Y} = (Y_i)_{1 \le i \le r} \in \mathrm{St}^I : (Y_i)'Y_j = 0, i \ne j\}$. Then, the restriction of $\pi_I^{SG}$ to $\mathrm{St}_{\perp}^I$ is valued in $\mathcal{F}^I$.

**Lemma 4.** Let $\pi_I^{OS}$ be the bijection between $O(d)$ and $\mathrm{St}_{\perp}^I$ defined by $\pi_I^{OS}(Q) = (Q^{(i)})_{1 \le i \le r}$. Then, for all $Q \in O(d)$,

$$ \pi_I^{SG}(\pi_I^{OS}(Q)) = \pi^I(Q). $$

*Proof.* For $Q \in O(d)$, $\pi_I^{SG}(\pi_I^{OS}(Q)) = (\pi_{q_i}^{SG}(Q^{(i)}))_{1 \le i \le r}$ and for all $1 \le i \le r$, by (2.7),

$$ \pi_{q_i}^{SG}(Q^{(i)}) = Q^{(i)}(Q^{(i)})' = (QI_d^{(i)}) (QI_d^{(i)})' = QP_0^i Q'. $$

□

**Corollary 1.** The triangles of the following diagram commute.

<figure><img src="image_2.png" /></figure>

## 2.4 Flag of eigenspaces

For $S \in \mathrm{Sym}_d$, its spectrum is denoted by

$$ \mathrm{Spec}(S) = \{\mu_1(S) \ge \dots \ge \mu_d(S)\}. $$

We denote by $\mathrm{Sym}_d^\ne$ the set of $S \in \mathrm{Sym}_d$ such that $\mu_1(S) > \dots > \mu_d(S)$. For all $1 \le i \le r$, set $\bar{q}_i := \sum_{j=1}^i q_j$. Let $\overline{\mathrm{Sym}}_d^{(I)}$ be the set of $S \in \mathrm{Sym}_d$ such that for all $1 \le i \le r-1$, $\mu_{\bar{q}_i}(S) > \mu_{\bar{q}_i+1}(S)$. Denote by $\mathrm{Sym}_d^I$ the set of $S \in \overline{\mathrm{Sym}}_d^{(I)}$ such that for all $1 \le i \le r$ and $k, \ell$ in $B^i$, $\mu_k(S) = \mu_\ell(S)$. Then,

$$ \mathrm{Sym}_d^\ne \cup \mathrm{Sym}_d^I \subset \overline{\mathrm{Sym}}_d^{(I)}. $$

Let $S \in \overline{\mathrm{Sym}}_d^{(I)}$. For $1 \le i \le r$, let $\{u_k : 1 \le k \le d\}$ be eigenvectors respectively associated to the eigenvalues $\{\mu_k(S) : k \in B_i\}$. Then, the linear subspace $\mathcal{V}_i$ spanned by $\{u_k : 1 \le k \le d\}$ is independent of the choice of these eigenvectors. We define the *eigenprojection* associated to $\{\mu_k(S) : k \in B_i\}$ as the map $P_i : \overline{\mathrm{Sym}}_d^{(I)} \to \mathbb{G}^i$ such that $P_i(S)$ is the projector onto $\mathcal{V}_i$. Then, define the map $F^I : \overline{\mathrm{Sym}}_d^{(I)} \to \mathcal{F}^I$ by

$$ F^I(S) = (P_i(S))_{1 \le i \le r}. $$

$F^I(S)$ is called the *flag of eigenspaces* of $S$.

**Lemma 5.** Let $S \in S \in \overline{\mathrm{Sym}}_d^{(I)}$ and $Q \in O(d)$. Then, $\pi^I(Q) = F^I(S)$ if and only if $Q$ is a matrix of eigenvectors of $S$ such that for all $1 \le i \le r$, the columns of $Q^{(i)}$ are eigenvectors associated to $\{\mu_k(S) : k \in B_i\}$ if and only if

$$ QP_0^i Q' = P_i(S), \quad 1 \le i \le r. \qquad (2.8) $$
**Definition 5.** For $1 \le i \le r$, the standard $i$-th Grassmannian is $P_0^i \in G^i$ defined by

$$P_0^i := \pi^i(\{\epsilon_j : j \in \beta_i\}) \quad \text{i.e.} \quad P_0^i = \text{Diag}[0_{q_1}, \dots, I_{q_i}, \dots, 0_{q_r}].$$

The standard flag of type I is the flag $\mathcal{P}_0^I := (P_0^i)_i$. The main orbital map is the map

$$\pi^I : O(d) \to \mathcal{F}^I, \quad \pi^I(Q) = Q * \mathcal{P}_0^I.$$

**Remark 4.** (i) It is easily obtained from (2.1) (see [12] for details) that

$$(\pi^i)_i \circ \zeta = \pi^I \quad \text{i.e.} \quad \left(\pi^i(Q^{(i)})\right)_i = \pi^I(Q) = (QP_0^i Q')_i, \quad Q \in O(d). \qquad (3.7)$$

(ii) Let $Q \in O(d)$ which spans a flag $\mathcal{P} \in \mathcal{F}^I$. By definition, $\mathcal{P} = (\pi^i(Q^{(i)}))_i$. So, by (3.7),

$$\mathcal{P} = \pi^I(Q) := Q * \mathcal{P}_0^I. \qquad (3.8)$$

### 3.1.3 Manifolds of linear subspaces

These actions of $O(d)$ provide the manifold structures of $\text{St}_q$, $G_q$ and $\mathcal{F}^I$, through Proposition 1 hereafter (see [4, Proposition A.2]).

**Proposition 1.** Consider a compact Lie group $\mathbb{G}$ which acts smoothly on a smooth manifold $\mathbb{M}$. For $b_0 \in \mathbb{M}$, let $\mathbb{B}$ be its orbit and $\mathbb{K}$ its isotropy group. Let $\psi : \mathbb{G} \to \mathbb{G}/\mathbb{K}$ be the canonical quotient map. Define the map $\pi_0 : \mathbb{G} \to \mathbb{B}$ by $Q \mapsto Q \cdot b_0$. Then,

(i) $\mathbb{B}$ is an embedded submanifold of $\mathbb{M}$.

(ii) The map $\tilde{\pi}_0 : \mathbb{G}/\mathbb{K} \to \mathbb{B}$ such that $\tilde{\pi}_0 \circ \psi = \pi_0$ is a diffeomorphism.

$O(d)$ is a compact Lie group and its actions defined in subsection 3.1.2 are smooth. In [4], Proposition 1 is applied to prove that $\text{St}_q$ (resp. $G_q$) is an embedded submanifold of $\mathbb{R}^{d \times q}$ (resp. $\text{Sym}_d$) that is diffeomorphic to $O(d)/O(d-q)$ (resp. $O(d)/(O(q) \times O(d-q))$). The manifold structure of $\mathcal{F}^I$ is provided by the following Lemma, proved in [12].

**Lemma 1.** The main orbital map $\pi^I : O(d) \to \mathcal{F}^I$ with $\pi^I(Q) = Q * \mathcal{P}_0^I$ is surjective.

By Lemma 1, $\mathcal{F}^I$ is the orbit of $\mathcal{P}_0^I$ under the smooth action of $O(d)$ on $G^I$ defined in (3.6). Now, the isotropy group of $\mathcal{P}_0^I$ under this action is the group

$$O(I) := \{\text{Diag}(H_1, \dots, H_r) : H_i \in O(q_i), 1 \le i \le r\} \simeq \prod O(q_i).$$

By Proposition 1, $\mathcal{F}^I$ is an embedded submanifold of $G^I$ and the map $\tilde{\pi}^I : O(d)/O(I) \to \mathcal{F}^I$ s.t. $\pi^I = \tilde{\pi}^I \circ \psi^I$ is a diffeomorphism, where $\psi^I : O(d) \to O(d)/O(I)$ is the canonical map.

### 3.1.4 Background of Riemannian geometry

We introduce Riemannian structures on these manifolds, in order to measure intrinsic distances and perform estimation on them.

More generally, let $\mathbb{M}$ be a connected smooth manifold. A *metric* on $\mathbb{M}$ is a collection $g = (g_x)_{x \in \mathbb{M}}$ of inner products on the tangent spaces $T_x\mathbb{M}$, varying smoothly wrt $x$. Then, $g$ defines the *geodesic distance* $\rho_g$ on $\mathbb{M}$, where for $x, y \in \mathbb{M}$, $\rho_g(x, y)$ is the infimum of lengths (wrt $g$) of all $C^1$ curves between $x$ and $y$. A *geodesic* is a smooth curve on $\mathbb{M}$ that realizes locally this infimum. For $x \in \mathbb{M}$ and $v \in T_x\mathbb{M}$, there exists a unique geodesic $\gamma_{x,v}$ such that $\gamma_{x,v}(0) = x$ and $\dot{\gamma}_{x,v}(0) = v$. By the Hopf-Rinow theorem, the metric space $(\mathbb{M}, \rho_g)$ is complete iff all geodesics are defined on $\mathbb{R}$, which we assume in the sequel. In that case, for all $x, y \in \mathbb{M}$, there exists at least one geodesic of minimal length between $x$ and $y$. For $x \in \mathbb{M}$, the *exponential map* at $x$ is the map $\text{Exp}_x^M : T_x\mathbb{M} \to \mathbb{M}$ such that $\text{Exp}_x^M(v) = \gamma_{x,v}(1)$.
# 3 Estimation of $F^I(\Sigma)$

Recall that $X$ is a random vector valued in $\mathbb{R}^d$, $d \ge 1$, with covariance matrix $\Sigma$ whose different eigenvalues $\lambda_1 > \dots > \lambda_r$ have respective multiplicities $I = (q_i)_{1 \le i \le r}$. We aim at estimating $F^I(\Sigma)$, through $F^I(\hat{\Sigma}_n)$, where $\hat{\Sigma}_n$ is a sample covariance matrix from an iid sample of size $n$ of $X$. In the preceding section, we have not described the Riemannian geometry of $\mathcal{F}^I$, since the Riemannian logarithm is not available in closed form. Instead, when $X$ is normally distributed in $\mathbb{R}^d$, we derive the estimation of $F^I(\Sigma)$ from the collection of CLT's for all $P_i(\hat{\Sigma}_n)$, $1 \le i \le r$, established below. Throughout the sequel, $\Gamma$ is a matrix of eigenvectors of $\Sigma$ such that $\pi^I(\Gamma) = F^I(S)$. Such a $\Gamma$ is described by Lemma 5. So, for all $1 \le i \le r$, $P_i(\Sigma) = \Gamma P_0^i \Gamma'$.

## 3.1 Review of Anderson's results

The following CLT provides the uncertainty of the estimation of $\hat{\Sigma}_n$ by $\Sigma$, which one wishes to propagate, in order to derive that of flags of eigenspaces.

**Theorem 2.** Let $\Delta$ be the diagonal matrix of eigenvalues of $\Sigma$ defined by $\Delta = \Gamma'\Sigma\Delta$. Then,

$$U_n := \sqrt{n} (\Gamma' \hat{\Sigma}_n \Gamma - \Delta) \xrightarrow[n \to \infty]{d} U,$$

where $U$ is a random matrix whose distribution is characterized as follows: $U \in \text{Sym}_d$, the blocks $\{U^{(i,j)}: 1 \le i < j \le r\}$ are mutually independent and for $i \ne j$, the entries of $U^{(i,j)}$ are iid rv's $\mathcal{N}(0, s_{i,j}^2)$, of standard deviation $s_{i,j} := \sqrt{\lambda_i \lambda_j}$.

Let $\psi$ be the map from $\text{Sym}_d^\pm$ to $O(d)$ such that for $S \in \text{Sym}_d^\pm$ and $1 \le k \le d$, the $k$-th column of $\psi(S)$ is the eigenvector associated to $\mu_k(S)$ whose $k$-th entry is non-negative. Then, for all $n \ge 1$, set

$$C_n := \psi(\hat{\Sigma}_n) \quad \text{and} \quad E_n := \Gamma' C_n.$$

Thus, $C_n$ is a matrix of eigenvectors of $\hat{\Sigma}_n$ such that $\pi^I(C_n) = F^I(\hat{\Sigma}_n)$. The following Theorem is the main result of (Anderson, 1963).

**Theorem 3.** For $1 \le i \le r$,

$$E_n^{(i,i)} \xrightarrow[n \to \infty]{d} E^{(i,i)},$$

where $E^{(i,i)}$ is uniformly distributed over the set of orthogonal matrices with non-negative diagonal entries. For all $i \ne j$ and $n \ge 1$, set $F_n^{(i,j)} := \sqrt{n} E_n^{(i,j)}$. Then,

$$F_n^{(i,j)} \xrightarrow[n \to \infty]{d} F^{(i,j)},$$

where $F^{(i,j)}$ is a random matrix whose entries are iid rv's $\mathcal{N}(0, \sigma_{i,j}^2)$, of standard deviation

$$\sigma_{i,j} := \frac{\sqrt{\lambda_i \lambda_j}}{|\lambda_i - \lambda_j|}.$$

Furthermore, the blocks $\{E^{(i,i)}, F^{(i,j)}: 1 \le i \le r, 1 \le j \le r, i \ne j\}$ are mutually independent.

The strategy for the proof of Theorem 3 is to express $E_n$ in function of $U_n$ and then to apply Theorem 4 below. Recall that $U_n = \sqrt{n} (\Gamma' \hat{\Sigma}_n \Gamma - \Delta)$, so that $\hat{\Sigma}_n = \phi_n(U_n)$ where $\phi_n$ is the map from $\text{Sym}_d$ to $\text{Sym}_d$ defined by

$$\phi_n(u) = \Gamma \left( \Delta + \frac{u}{\sqrt{n}} \right) \Gamma.$$

Then, $E_n = f_n(U_n)$, where $f_n$ is the map from $(\phi_n)^{-1}(\text{Sym}_d^\pm)$ to $O(d)$ defined by

$$f_n: u \mapsto \Gamma' \psi(\phi_n(u)) \quad \text{that is,} \quad f_n(u) = \psi \left( \Delta + \frac{u}{\sqrt{n}} \right).$$
We recall that, for all $x \in \mathbb{M}$, any geodesic starting at $x$ is locally length-minimizing. The cut locus of $x \in \mathbb{M}$ is the set $\text{Cut}(x)$ of points $y$ such that the geodesics starting at $x$ cease to be length-minimizing beyond $y$. Thus, on a sphere, the geodesics are the great circles, so that two antipodal points are each other's cut locus. Let $x \in \mathbb{M}$ and $y \in \mathbb{M} \setminus \text{Cut}(x)$. Then, there exists a unique geodesic of minimal length between $x$ and $y$, whose initial velocity is thus the smallest tangent vector $v \in T_x \mathbb{M}$ such that $\text{Exp}_x^M(v) = y$. This vector is called the Riemannian Logarithm of $y$ at $x$, denoted by $\text{Log}_x^M(y)$, related to the geodesic distance by:

$$ \rho_g(x, y) = \|\text{Log}_x^M(y)\|_x. \qquad (3.9) $$

For any $x \in \mathbb{M}$, there exists an open set $\text{ID}_x \subset T_x \mathbb{M}$, called the injectivity domain at $x$, such that the restriction of $\text{Exp}_x^M$ to $\text{ID}_x$ is a diffeomorphism onto $\mathbb{M} \setminus \text{Cut}(x)$, whose inverse is a restriction of the map $\text{Log}_x^M(\cdot)$. Finally, we present results on Riemannian submersions.

**Definition 6.** Let $(\mathbb{E}, g^\mathbb{E})$ and $(\mathbb{B}, g^\mathbb{B})$ be Riemannian manifolds. Let $\pi : \mathbb{E} \to \mathbb{B}$ be a smooth submersion, i.e. for all $x \in \mathbb{E}$, $d_x\pi$ is a surjective linear map.

(i) For $x \in \mathbb{E}$, the vertical space at $x$ is $V_x\mathbb{E} := \ker(d_x\pi) \subset T_x\mathbb{E}$. The horizontal space at $x$ is its orthogonal complement in $T_x\mathbb{E}$ wrt the inner product $g_x^\mathbb{E}$, denoted by $H_x\mathbb{E} := (V_x\mathbb{E})^\perp$.

(ii) Let $y \in \mathbb{B}$ and $v \in T_y\mathbb{B}$. By (i), for all $x \in \mathbb{E}$, there exists a unique vector $v_x^\sharp \in H_x$ such that $(d_x\pi)(v_x^\sharp) = v$. The vector $v_x^\sharp$ is called the horizontal lift wrt $\pi$ of $v$ at $x$.

(iii) The map $\pi$ is a Riemannian submersion when for all $x \in \mathbb{E}$, the restriction of $d_x\pi$ to $H_x\mathbb{E}$ is a linear isometry between $(H_x\mathbb{E}, g_x^\mathbb{E})$ and $(T_{\pi(x)}\mathbb{B}, g_{\pi(x)}^\mathbb{B})$.

**Proposition 2.** Let $\pi : \mathbb{E} \to \mathbb{B}$ be a Riemannian submersion.

(i) Let $\gamma$ be a geodesic in $\mathbb{B}$. Then, for all $x \in \mathbb{E}$, there exists a unique geodesic $\hat{\gamma}^x$ in $\mathbb{E}$ through $x$ which is horizontal, i.e. all its tangent vectors are horizontal, and which projects to $\gamma$, i.e. $\pi(\hat{\gamma}^x) = \gamma$. The geodesic $\hat{\gamma}^x$ is called the horizontal lift wrt $\pi$ of $\gamma$ through $x$.

(ii) The geodesics in $\mathbb{B}$ are the images by $\pi$ of the horizontal geodesics in $\mathbb{E}$.

### 3.1.5 Metrics on Stiefel manifolds and Grassmannians

We recall hereafter the description of the tangent spaces of $\text{St}_q$ and $G_q$: See [4]. For any $U \in \text{St}_q$ and $P \in G_q$,

$$ T_U\text{St}_q = \{\mathcal{D} \in \mathbb{R}^{d \times q} : U'\mathcal{D} = -\mathcal{D}'U\} \quad \text{and} \quad T_P G_q = \{\Delta \in \text{Sym}_d : \Delta P + P\Delta = \Delta\}. \qquad (3.10) $$

In the sequel, we endow $\text{St}_q$ with the Euclidean metric $g^\text{St}$ defined by

$$ g_{\mathcal{U}}^{\text{St}}(\mathcal{D}_1, \mathcal{D}_2) = \text{tr}(\mathcal{D}_1^T \mathcal{D}_2), \quad U \in \text{St}_q \text{ and } \mathcal{D}_1, \mathcal{D}_2 \in T_{\mathcal{U}}\text{St}_q. $$

Then, $G_q$ is endowed with the metric $g^G$ (see [4] for a discussion on it) defined by

$$ g_P^G(\Delta_1, \Delta_2) = \frac{1}{2}\text{tr}(\Delta_1 \Delta_2), \quad P \in G_q \text{ and } \Delta_1, \Delta_2 \in T_P G_q. $$

The Riemannian manifold $(G_q, g^G)$ has a rich structure, notably of symmetric space. So, many Riemannian operations in $G_q$ are available in closed form: see paragraph 3.2.3. Moreover, the following result holds: see [12].

**Lemma 2.** The map $\pi_q$ is a Riemannian submersion from $(\text{St}_q, g^{\text{St}})$ onto $(G_q, g^G)$.
**Theorem 4.** Let $\mathcal{F}, \mathcal{T}$ be metric spaces and $(\mathcal{F}_n)_{n \ge 1}$ a sequence of subsets of $\mathcal{F}$. Let $(V_n)_{n \ge 1}$ be a sequence of rv's valued in $(\mathcal{F}_n)_{n \ge 1}$. Assume that

$$V_n \xrightarrow[n \to \infty]{d} V,$$

where $V$ is a rv valued in $\mathcal{F}$. For all $n \ge 1$, let $g_n$ be a map from $\mathcal{F}_n$ to $\mathcal{T}$ and $g$ a map from $\mathcal{F}$ to $\mathcal{T}$. Let $\mathcal{F}_0$ be a subset of $\mathcal{F}$ such that $g$ is continuous on $\mathcal{F}_0$ and

$$\mathrm{Pr}(V \in \mathcal{F}_0) = 1.$$

Assume that for all $u \in \mathcal{F}_0$, and all sequence $(u_n)_{n \ge 1}$ valued in $(\mathcal{F}_n)_{n \ge 1}$,

$$u_n \xrightarrow[n \to \infty]{\phantom{=}d} u \implies g_n(u_n) \xrightarrow[n \to \infty]{\phantom{=}d} g(u).$$

Then,

$$g_n(V_n) \xrightarrow[n \to \infty]{d} g(V).$$

*Proof.* See APPENDIX D in (Anderson, 1963) or section 18.11 in (Van der Vaart, 2000). $\square$

## 3.2 CLT for $P_i(\hat{\Sigma}_n)$

For any $u \in \mathrm{Sym}_d$, set $M^i_<(u) = [\frac{1}{\lambda_i-\lambda_j}u^{(i,j)} : i < j]$ and $M^i_>(u) = [\frac{1}{\lambda_i-\lambda_j}u^{(i,j)} : i > j]$.

For all $1 \le i \le r$, define the map $g^i : \mathrm{Sym}_d \to T_{P_0^i} G^i$ by

$$g^1(u) = \begin{pmatrix} 0_{q_1} & M^1_>(u) \\ (M^1_>(u))' & 0 \end{pmatrix}, \quad g^r(u) = \begin{pmatrix} 0 & (M^r_<(u))' \\ M^r_<(u) & 0_{q_r} \end{pmatrix}$$

$$g^i(u) = \begin{pmatrix} 0 & (M^i_<(u))' & 0 \\ M^i_<(u) & 0_{q_i} & M^i_>(u) \\ 0 & (M^i_>(u))' & 0 \end{pmatrix}, \quad i \ne 1, r.$$

Hereabove, when we write 0 as a block matrix, it is implicit that its size is the suitable one. In the sequel, given rv's $X : \Omega \to V$ and $Y : \Omega \to (E, \mathcal{B})$, where $V$ is a vector space and $(E, \mathcal{B})$ a measurable space, for $B \in \mathcal{B}$, we denote by $X1_{\{Y \in B\}}$ the rv valued in $V$ equals to $X$ if $Y \in B$ and to 0 else.

**Theorem 5.** For all $1 \le i \le r$,

$$\sqrt{n}\mathrm{Log}_{P_0^i}^{G^i}(\Gamma'P_i(\hat{\Sigma}_n)\Gamma) \mathbf{1}_{\mathcal{T}_n^i} = g^i(U) + O_{\mathrm{Pr}}(1/\sqrt{n}), \quad (3.1)$$

$$\text{where } \mathcal{T}_n^i := \{\Gamma' P_i(\hat{\Sigma}_n) \Gamma \notin \mathrm{Cut}(P_0^i)\}.$$

The proof of Theorem 5 is based on Theorem 4. We introduce below some sequences of different natures satisfying the assumptions of Theorem 4. For all $1 \le i \le r$ and $n \ge 1$, set

$$\mathcal{F}_n^i := \{u \in \phi_n^{-1}(\mathrm{Sym}_d^\nolimits) : \Gamma' P_i(\phi_n(u))\Gamma \notin \mathrm{Cut}(P_0^i)\}.$$

Then, consider the rv $V_n^i$ valued in $\mathcal{F}_n^i$, defined by

$$V_n^i := U_n \mathbf{1}_{\{U_n \in \mathcal{F}_n^i\}}$$

Define the map $g_n^i : \mathcal{F}_n^i \to T_{P_0^i} G^i$ by

$$g_n^i(u) = \sqrt{n} \log_{P_0^i}^{G^i}(\Gamma' P_i(\hat{\Sigma}_n) \Gamma).$$
## 3.2 Geodesic parametrization of a full flag

With the notations of Section 2, let $F^1(\Sigma)$ be the flag of eigenspaces of $\Sigma$ (I denotes its type) and $\hat{F}_n$ that of $\hat{\Sigma}_n$. Thus, $F^1(\Sigma)$ is the flag of PS's and $\hat{F}_n$ is a.s. a full flag. As in Anderson's paper [1], $\hat{F}_n$ is represented by the matrix $C_n \in O(d)$ of eigenvectors of $\hat{\Sigma}_n$. When $F^1(\Sigma)$ is estimated by $\hat{F}_n$, the main difficulty to measure the deviation between them is that their types are different in general.

To solve this problem, we define a parametrization of (almost) any full flag, represented by an orthogonal matrix, which determines its relative position wrt the pair formed by a given flag and an orthogonal matrix that spans the latter. Then, we prove that such a parametrization of $C_n$ wrt $(F^1(\Sigma), \Gamma)$ is equivalent to that of $E_n := \Gamma'C_n$ wrt $(\mathcal{P}_0^I, I_d)$.

### 3.2.1 Geodesic projection of a frame

Fix $R \in G_q$ and set

$$V_R := \pi_q^{-1}(G_q \setminus \text{Cut}(R)) = \{\mathcal{U} \in \text{St}(q, d) : \pi_q(\mathcal{U}) \notin \text{Cut}(R)\}.$$

For $\mathcal{U} \in V_R$, set $P := \pi_q(\mathcal{U})$. Let $\gamma_{(P,R)}$ be the minimal geodesic between $P$ and $R$ in $G_q$ and $\tilde{\gamma}_{(P,R)}^U$ its horizontal lift wrt $\pi_q$ through $\mathcal{U}$. Then, define the $q$-frame $\mathcal{H}_{(P,R)}^q(\mathcal{U})$ as the endpoint of $\tilde{\gamma}_{(P,R)}^U$. By (ii) of Proposition 2, $\pi_q(\mathcal{H}_{(P,R)}^q(\mathcal{U})) = R$.

<figure><img src="image_3.png" /></figure>

**Remark 5.** Among all $q$-frames generating $\text{rg}(R)$, $\mathcal{H}_{(P,R)}^q(\mathcal{U})$ minimizes the geodesic distance in $\text{St}_q$ to $\mathcal{U}$, denoted by $d_g^{\text{St}}$: See Lemma 26.11 in [10]. Formally,

$$d_g^{\text{St}}(\mathcal{U}, \mathcal{H}_{(P,R)}^q(\mathcal{U})) = \min \{d_g^{\text{St}}(\mathcal{U}, \mathcal{W}) : \pi^{\text{SG}}(\mathcal{W}) = R\}. \quad (3.11)$$

So, $\mathcal{H}_{(P,R)}^q(\mathcal{U})$ is interpreted as the orthogonal projection of $\mathcal{U}$ on $\text{rg}(R)$ wrt the metric $g^{\text{St}}$.

**Definition 7.** The $q$-frame $\mathcal{H}_{(P,R)}^q(\mathcal{U})$ is called the geodesic projection of $\mathcal{U}$ onto $\text{rg}(R)$. For $v \in T_pG_q$, its horizontal lift wrt $\pi_q$ at $\mathcal{U}$ is denoted by $v_{\mathcal{U}}^{\sharp}$. By definition,

$$\mathcal{H}_{(P,R)}^q(\mathcal{U}) = \text{Exp}_{\mathcal{U}}^{\text{St}_q} \left( \left( \text{Log}_P^{G_q}(R) \right)_{\mathcal{U}}^{\sharp} \right). \quad (3.12)$$

### 3.2.2 Geodesic decomposition of a frame

Fix $R \in G_q$ and $W_R \in \pi_q^{-1}(R)$.

**Lemma 3.** (i) Define the map $\theta_R: V_R \to \pi_q(V_R) \times \pi_q^{-1}(R)$ by

$$\theta_R(\mathcal{U}) = (\pi_q(\mathcal{U}), \mathcal{H}_{(\pi_q(\mathcal{U}),R)}^q(\mathcal{U})) . \quad (3.13)$$

Then, the map $\theta_R$ is a diffeomorphism, whose inverse is defined by

$$\theta_R^{-1} : \pi_q(V_R) \times \pi_q^{-1}(R) \to V_R \quad \text{with} \quad (P, V) \mapsto \mathcal{H}_{(R,P)}^q(V).$$

(ii) Let $\text{ID}_R$ be the injectivity domain at $R$. Define the map $\xi_R^{W_R}: V_R \to \text{ID}_R \times O(q)$ by

$$\xi_R^{W_R}(\mathcal{U}) = (\text{Log}_R^{G_q}(\pi_q(\mathcal{U})), (W_R)' \mathcal{H}_{(\pi_q(\mathcal{U}),R)}(\mathcal{U})) . \quad (3.14)$$

Then, the map $\xi_R^{W_R}$ is a diffeomorphism.
Since $\phi_n(0) = \Gamma\Delta\Gamma' = \Sigma$ and $\Gamma^t P_i(\Sigma)\Gamma = P_0^i$, we have that $0 \in \mathcal{F}_n^i$ and $g_n^i(0) = 0$. Therefore,

$$g_n^i(V_n^i) = \sqrt{n} \log_{P_0^i} \left( \Gamma' P_i \left( \hat{\Sigma}_n \right) \Gamma \right) \mathbf{1}_{\tau_n}.$$

The following series of results imply that the assumptions of Theorem 4 are fullfilled.

**Proposition 3.** For all $1 \le i \le r$,

$$V_n^i \xrightarrow{n \to \infty} U.$$

*Proof.* By the CLT for $\hat{\Sigma}_n$, $U_n \xrightarrow{n \to \infty} U$, which combined to Lemma 6 below, concludes the proof. $\square$

**Lemma 6.** For all $1 \le i \le r$,

$$\mathrm{Pr} (U_n \in \mathcal{F}_n^i) \xrightarrow{n \to \infty} 0.$$

*Proof.* For all $1 \le i \le r$ and $n \ge 1$,

$$U_n \in \mathcal{F}_n^i \iff E_n^{(i)} (E_n^{(i)})' \notin \mathrm{Cut} (P_0^i) \iff \mathrm{rk} (E_n^{(ii)}) = q_i.$$

By Theorem 3, $E_n^{ii} \xrightarrow{n \to \infty} E^{ii}$ and $E^{ii}$ is uniformly distributed on the orthogonal group $O(q_i)$, so that, a.s., $\mathrm{rk}(E^{ii}) = q_i$. Now, for all $n \ge 1$, $\mathrm{rk}(E_n^{ii}) \le q_i$. We conclude by Lemma 7 below. $\square$

**Lemma 7.** Let $B_n$ and $B$ be random matrices of same size. Assume that $B_n \xrightarrow{n \to \infty} B$, that, a.s., $\mathrm{rk}(B) = b$ and that $\mathrm{Pr}(\mathrm{rk}(B_n) \le b) \xrightarrow{n \to \infty} 1$. Then, $\mathrm{Pr}(\mathrm{rk}(B_n) = b) \xrightarrow{n \to \infty} 1$.

*Proof.* See Lemma 2.6. in (Tyler, 1981) $\square$

**Remark 1.** For all $1 \le i \le r$, $g^i$ is continuous on $\mathrm{Sym}_d^\ddagger$ and

$$\mathrm{Pr} (U \in \mathrm{Sym}_d^\ddagger) = 1$$

**Proposition 4.** Fix $1 \le i \le r$. Let $v \in \mathrm{Sym}_d^\ddagger$. Let $(v_n)_{n \ge 1}$ be a sequence in $(\mathcal{F}_n^i)_{n \ge 1}$ such that $v_n \xrightarrow{n \to \infty} v$. Then, $g_n^i(v_n)$ is of the form

$$g_n^i(v_n) = [\Lambda_n^i(v_n), P_0^i] + (\sqrt{n})^{-1} [\Phi_n^i(v_n), P_0^i] + o((\sqrt{n})^{-1}),$$

where for all $n \ge 1$, $\Delta_n^i(v_n)$ is a matrix such that

$$[\Lambda_n^i(v_n), P_0^i] \in T_{P_0^i} G^i \quad \text{with} \quad [\Lambda_n^i(v_n), P_0^i] \xrightarrow{n \to \infty} g^i(v),$$

and the sequence $(\Phi_n^i(v_n))_{n \ge 1}$ converges to a matrix $\Phi^i(v)$, as $n \to \infty$.

*Proof.* We use the proof of Theorem 3 and expand the matrix logarithm involved in the closed form of the Riemannian Logarithm for Grassmannians. $\square$

### 3.3 A pivotal statistic

We aim at deriving a pivotal statistic of $F^I(\Sigma)$, i.e. which depends only on $F^I(\Sigma)$ and the sample, and whose asymptotic distribution does not depend on any unknown parameter.

In the sequel, we denote by $\mathcal{D}(I)$ the set of $d \times d$ diagonal matrices $D$ of the form

$$D = \mathrm{diag} [\beta_i I_{q_i} : 1 \le i \le r].$$

**Remark 2.** For all $H \in \mathcal{O}(I)$ and $K \in \mathcal{D}(I)$, $[H, K] = 0$.
**Proof.** One checks that the map $\theta_R^{-1}$ is indeed the inverse of $\theta_R$: see the diagram in Definition 7. Their smoothness follows from properties of Riemannian submersions. This proves (i). Then, (ii) follows readily from (i) and Remark 2. $\square$

**Definition 8.** For $U \in \mathbb{V}_R$, the pair $\theta_R(U)$ is called the geodesic decomposition of $U$ wrt $R$ and $\xi_R^{W_R}(U)$ its effective geodesic decomposition wrt $(R, W_R)$. By Lemma 3, any $U \in \mathbb{V}_R$ is uniquely determined by these decompositions.

**Remark 6.** The effective geodesic decomposition lies in a product of a vector space and a group. Thus, it provides a handleable measure of the relative position of any q-frame $U \in \mathbb{V}_R$ wrt a given q-linear subspace $R$ and a q-frame spanning $rg(R)$.

### 3.2.3 Explicit formulas

The explicit formulas in $G_q$ provide those for the (effective) geodesic decomposition of a frame. First, we describe the cut locus in $G_q$: see [4].

**Lemma 4.** (i) Let $P \in G_q$. Then, writing $P = YY'$,

$$ \mathrm{Cut}(P) = \{R = ZZ' \in G_q : \mathrm{rank}(Y'Z) < q\}. $$

(ii) Thus, the following symmetry holds: $R \in \mathrm{Cut}(P) \iff P \in \mathrm{Cut}(R)$.

(iii) The cut locus is compatible with the action of $O(d)$ on $G_q$, i.e. for all $Q \in O(d)$,

$$ R \in \mathrm{Cut}(P) \iff Q \cdot R \in \mathrm{Cut}(Q \cdot P). $$

Now, a closed form for the Riemannian Logarithm in $G_q$ is presented hereafter: see [2].

**Theorem 4.** Let $P \in G_q$ and $R \in G_q \setminus \mathrm{Cut}(P)$. Then,

$$ \mathrm{Log}_P^{G_q}(R) = [\Omega, P] \quad \text{where} \quad \Omega = \frac{1}{2} \log_m((I_d - 2R)(I_d - 2P)), \qquad (3.15) $$

and $\log_m$ denotes the matrix logarithm.

**Corollary 1.** Let $P, R \in G_q$ and $Q \in O(d)$. If $R \notin \mathrm{Cut}(P)$, then

$$ \mathrm{Log}_{Q \cdot P}^{G_q}(Q \cdot R) = Q(\mathrm{Log}_P^{G_q}(R))Q'. $$

*Proof.* By the properties of the matrix logarithm,

$$ \begin{aligned} \log_m((I_d - 2QRQ')(I_d - 2QPQ')) &= \log_m(Q(I_d - 2R)(I_d - 2P)Q') \\ &= Q \log_m((I_d - 2R)(I_d - 2P))Q'. \end{aligned} $$

So, with the notations of (3.15), $\mathrm{Log}_{Q \cdot P}^{G_q}(Q \cdot R) = [Q\Omega Q', QPQ'] = Q[\Omega, P]Q'$. $\square$

Then, we establish from (3.12) an explicit expression of $\mathcal{H}_{(P,R)}^q(\mathcal{U})$, proved in [12].

**Proposition 3.** For $U \in \mathbb{V}_R$, set $P := \pi_q(U)$ and $\Delta := \mathrm{Log}_P^{G_q}(R)$. Then, setting also $V := (U \ \Delta U) \in \mathbb{R}^{n \times 2q}$ and $C := U'\Delta^2 U$,

$$ \mathcal{H}_{(P,R)}^q(\mathcal{U}) = V \begin{pmatrix} \exp_m \begin{pmatrix} 0 & -C \\ I_q & 0 \end{pmatrix} \\ I_q \\ 0_q \end{pmatrix}. \qquad (3.16) $$

We deduce hereafter the effect on the geodesic decomposition of the action of $O(d)$.
For any $1 \le i \le r$, we aim at normalizing the entries of $g^i(U)$. By Theorem 2, for $i \ne j$, the entries of $\frac{1}{\lambda_i - \lambda_j} U^{(i,j)}$ are real i.i.d rv's $\mathcal{N}(0, \sigma_{i,j}^2)$, of standard deviation $\sigma_{i,j} = \frac{\sqrt{\lambda_i \lambda_j}}{|\lambda_i - \lambda_j|}$. For any $1 \le i \le r$, define $K^i \in \mathcal{D}(I)$ by

$$ K^i := \operatorname{diag} \left[ \frac{1}{\sigma_{i,1}} I_{q_1}, \dots, \frac{1}{\sigma_{i,r}} I_{q_r} \right]. \qquad (3.2) $$

Now, consider the rv $\bar{g}^i(U)$ defined by

$$ \bar{g}^i(U) := K^i g^i(U) K^i. $$

Then, $\bar{g}^i(U)$ is a rv valued in $T_{P_0^i} G^i$, with the same stochastic dependence between and within blocks as $\bar{g}^i(U)$. Now, left and right multiplications of $g^i(U)$ by $K^i$ imply that $\frac{1}{\lambda_i - \lambda_j} U^{(i,j)}$ is replaced by

$$ W^{(i,j)} := \frac{1}{\sigma_{i,j}} \left( \frac{1}{\lambda_i - \lambda_j} U^{(i,j)} \right), $$

whose entries are real iid rv's $\mathcal{N}(0, 1)$. Since $U$ is a symmetric random matrix,

$$ \sum_{i=1}^{r} \| \bar{g}^i(U) \|_F^2 = 4 \sum_{1 \le i < j \le r} \| W^{(i,j)} \|_F^2. $$

Now, the blocks $\{W^{(i,j)} : 1 \le i < j \le r\}$ are mutually independent. So,

$$ \frac{1}{4} \sum_{i=1}^{r} \| \bar{g}^i(U) \|_F^2 = \sum_{1 \le i < j \le r} \| W^{(i,j)} \|_F^2 \sim \chi_{\mathrm{D}^1}^2, $$

where

$$ \mathrm{D}^I := \frac{1}{2} \left( d^2 - \sum_{i=1}^{r} q_i^2 \right). \qquad (3.3) $$

By Theorem 5, for all $1 \le i \le r$,

$$ \sqrt{n} \left[ K^i \left( \operatorname{Log}_{\Pi_0^i} \left( \Gamma' P^j \left( \hat{\Sigma}_n \right) \Gamma \right) \right) K^i \right] \mathbf{1}_{\mathcal{T}_n^i} \xrightarrow{n \to \infty} \bar{g}^i(U). $$

Therefore, we have proved the following result.

**Proposition 5.** For $(K^i)_{1 \le i \le r}$ defined in (3.2), and $\mathrm{D}^I$ defined in (3.3),

$$ \frac{n}{4} \sum_{i=1}^{r} \| K^i \operatorname{Log}_{P_0^i}^{G^i} (\Gamma' P^i (\hat{\Sigma}_n) \Gamma) K^i \|_F^2 \mathbf{1}_{\mathcal{T}_n^i} \xrightarrow{n \to \infty} \chi_{\mathrm{D}^I}^2. $$

**Corollary 2.** The statistic $\hat{T}_n$ defined below is a a pivotal statistic of $F^I(\Sigma)$:

$$ \hat{T}_n := \frac{n}{4} \sum_{i=1}^{r} \| \hat{K}_n^i \operatorname{Log}_{P_0^i}^{G^i} (\Gamma' P^i (\hat{\Sigma}_n) \Gamma) \hat{K}_n^i \|_F^2 \mathbf{1}_{\mathcal{T}_n^i} \xrightarrow{n \to \infty} \chi_{\mathrm{D}^I}^2, \qquad (3.4) $$

where, for $1 \le i \le r$, $\hat{K}_n^i$ is derived from $K^i$ by replacing $\lambda_i$ by any consistent estimator $\hat{\lambda}_n^i$, for example

$$ \hat{\lambda}_n^i = \frac{1}{q_i} \sum_{k \in B_i} \mu_k (\hat{\Sigma}_n). $$

## 3.4 Simulation

Figure 1 illustrates the convergence in distribution of $\hat{T}_n$ to a $\chi^2$ distribution. The parameters for this simulation are $d = 4$ and $I = (1, 1, 1, 1)$, so that $\Sigma \in \mathrm{Sym}_d^\pm$. Then $\mathrm{D}^I = d(d-1)/2 = 6$. For the sample size, we take $n = 10000$. The histogram in blue represents the probability distribution of $\hat{T}_n$ and the curve in red is that of the probability distribution function of the $\chi_6^2$ distribution. We see that the empirical distribution, i.e. that of $\hat{T}_n$ is indeed very close to that of the $\chi_6^2$ one.
**Corollary 2.** (i) Let $R \in G_q$ and $U \in V_R$. Then, for all $Q \in O(d)$,

$$QU \in V_{Q \cdot R} \quad \text{and} \quad \mathrm{Log}_{Q \cdot R}^{G_q}(\pi_q(QU)) = Q(\mathrm{Log}_R^{G_q}(\pi_q(U)))' \quad (3.17)$$

(ii) For the geodesic projections, we have that

$$\mathcal{H}_{(\pi_q(QU), Q \cdot R)}(QU) = Q\mathcal{H}_{(\pi_q(U), R)}(U).$$

(iii) For any $W_R \in \pi_q^{-1}(R)$, the effective geodesic decomposition of $U$ wrt $(R, W_R)$ is fully determined by that of $QU$ wrt $(Q \cdot R, QW_R)$.

*Proof.* (i) and (ii) follow from Lemma 4, Corollary 1 and Proposition 3. Then, (iii) is deduced from (i) and (ii). □

### 3.2.4 Geodesic parametrization of an orthogonal matrix

In this paragraph, we consider a fixed flag $\mathcal{R} = (R_i)_i \in \mathcal{F}^I$ and a matrix $\Gamma \in O(d)$ which spans $\mathcal{R}$. Now, set

$$W_{\mathcal{R}} := \{ Q \in O(d) : \forall 1 \le i \le r, \pi^i(Q^{(i)}) \notin \mathrm{Cut}(R_i) \text{ i.e. } QP_0^i Q' \notin \mathrm{Cut}(R_i) \}. \quad (3.18)$$

For $C \in W_R$ identified with $\zeta(C) := (C^{(i)})_i$, consider the collection of geodesic decompositions of $C^{(i)}$ wrt $R_i$, illustrated hereafter, where (3.7) justifies that $\pi^I(C) = (\pi^i(C^{(i)}))_i$.

<figure><img src="image_4.png" /></figure>

**Definition 9.** The geodesic parametrization of $C \in W_R$ (or of the associated full flag) wrt $(\mathcal{R}, \Gamma)$ is defined as the collection of effective geodesic decompositions of $C^{(i)}$ wrt $(R_i, \Gamma^{(i)})$, for $1 \le i \le r$. This parametrization determines $C \in W_R$ uniquely.

**Lemma 5.** For $C \in W_R$, set $E := \Gamma'C$. Then, $E \in W_{P_0^I}$ and the geodesic parametrization of $C$ wrt $(\mathcal{R}, \Gamma)$ is fully determined by that of $E$ wrt $(P_0^I, I_d)$.

*Proof.* We apply Corollary 2 with $Q = \Gamma'$. Thus, the geodesic parametrization of $C$ wrt $(\mathcal{R}, \Gamma)$ is fully determined by that of $E := \Gamma'C$ wrt $(\Gamma' * \mathcal{R}, \Gamma'\Gamma) = (\Gamma^{-1} * \mathcal{R}, I_d)$. Now, $\Gamma \in O(d)$ spans the flag $\mathcal{R}$. So, by (3.8), $\mathcal{R} = \Gamma * P_0^I$, which implies that $\Gamma^{-1} * \mathcal{R} = P_0^I$. □

### 3.2.5 Tangent spaces at standard Grassmannians

A geodesic parametrization wrt $(P_0^I, I_d)$ involves the maps $\mathrm{Log}_{P_0^I}^{G_i}$, for $1 \le i \le r$, each valued in the tangent space $T_{P_0^I} G^i$. Thus, we describe hereafter the structure of the latter.

**Lemma 6.** For $1 \le i \le r$, set $\mathbb{T}_0^i := T_{P_0^i} G^i$. Then,

$$\mathbb{T}_0^1 = \left\{ \begin{pmatrix} 0_{q_1} & A_1^1 \\ (A_1^1)' & 0 \end{pmatrix} : A_1^1 \in \mathbb{R}^{q_1 \times (q_2 + \dots + q_r)} \right\}, \quad \mathbb{T}_0^r = \left\{ \begin{pmatrix} 0 & (A_r^r)' \\ A_r^r & 0_{q_r} \end{pmatrix} : A_r^r \in \mathbb{R}^{q_r \times (q_1 + \dots + q_{r-1})} \right\}$$

$$\mathbb{T}_0^i = \left\{ \begin{pmatrix} 0 & (A_i^i)' & 0 \\ A_i^i & 0_{q_i} & A_i^i \\ 0 & (A_i^i)' & 0 \end{pmatrix} : A_i^i \in \mathbb{R}^{q_i \times (q_1 + \dots + q_{i-1})}, \quad A_i^i \in \mathbb{R}^{q_i \times (q_{i+1} + \dots + q_r)} \right\}, \quad i \notin \{1, r\}.$$
<figure><img src="image_3.png" /><figcaption>Figure 1: Illustration of the convergence of $\hat{T}_n$ to the $\chi_0^2$ distribution</figcaption></figure>

## 3.5 Confidence regions and Tests

Let $\mathbb{K} = (K^i)_{1\le i\le r} \in (\mathcal{D}(\mathbb{I}))^r$. For $Q \in \mathcal{O}(d)$ and $\mathcal{R} = (R_i)_{1\le i\le r} \in \mathcal{F}^1$, set

$$ \Delta_{\mathbb{K}}(Q, \mathcal{R}) = \frac{1}{4} \sum_{i=1}^{r} \left\| K^i \mathrm{Log}_{P_0^i} (Q'R_iQ) K^i \right\|_F^2 \mathbf{1}_{\{Q'R_iQ \notin \mathrm{Cut}(P_0^i)\}}. $$

**Proposition 6.** For all $H \in \mathcal{O}(\mathbb{I})$,

$$ \Delta_{\mathbb{K}}(QH, \mathcal{R}) = \Delta_{\mathbb{K}}(Q, \mathcal{R}). \qquad (3.5) $$

*Proof.* First, we prove that for all $1 \le i \le r$ and $H \in \mathcal{O}(\mathbb{I})$,

$$ (QH)'R_i(QH) \notin \mathrm{Cut}(P_0^i) \iff Q'R_iQ \notin \mathrm{Cut}(P_0^i) \qquad (3.6) $$

Indeed, writing $R_i = C_i(C_i)'$ with $C_i \in \mathrm{St}(q_i, d)$,

$$ \begin{align*} (QH)'R_i(QH) \notin \mathrm{Cut}(P_0^i) &\iff (H'Q'C_i)(H'Q'C_i)' \notin \mathrm{Cut}(P_0^i) \\ &\iff \mathrm{rank}((I_d^{(i)})'(H'Q'C_i)) = q_i \\ &\iff \mathrm{rank}((HI_d^{(i)})'(Q'C_i)) = q_i \end{align*} $$

Now, for all $1 \le i \le r$, $(HI_d^{(i)}) (HI_d^{(i)})' = H P_0^i H' = P_0^i$. Therefore,

$$ \mathrm{rank}((HI_d^{(i)})'(Q'C_i)) = q_i \iff \mathrm{rank}((I_d^{(i)})'(Q'C_i)) = q_i \iff Q'R_iQ \notin \mathrm{Cut}(P_0^i). $$

This proves (3.6). Now, when $(QH)'R_i(QH) \notin \mathrm{Cut}(P_0^i)$, the closed formula for $\mathrm{Log}_{P_0^i}(\cdot)$ yields that

$$ \mathrm{Log}_{P_0^i}((QH)'R_i(QH)) = \mathrm{Log}_{P_0^i}(H'(Q'R_iQ)H) = H'\mathrm{Log}_{P_0^i}(Q'R_iQ)H. $$

Now, Remark 2 combined to the properties of invariance of the Frobenius norm imply that for all $1 \le i \le r$,

$$ \|K^i(H'\mathrm{Log}_{P_0^i} (Q'R_iQ) H) K^i\|_F = \|H'(K^i\mathrm{Log}_{P_0^i} (Q'R_iQ) K^i) H\|_F = \|K^i\mathrm{Log}_{P_0^i} (Q'R_iQ) K^i\|_F $$

□

For fixed $n \ge 1$, set $\hat{\mathbb{K}}_n := (\hat{K}_n^i)_{1\le i\le r} \in (\mathcal{D}(\mathbb{I}))^r$ and define the map $\hat{T}_n : \mathcal{O}(d) \to [0, \infty)$ by

$$ \hat{T}_n(Q) := n\Delta_{\hat{\mathbb{K}}_n} (Q, F^I (\hat{\Sigma}_n)). $$

Then, for any matrix $\Gamma$ of eigenvectors of $\Sigma$, $\hat{T}_n(\Gamma)$ is the pivotal statistic $\hat{T}_n$ defined in (3.4). Furthermore, Proposition 6 implies that the value of $\hat{T}_n(Q)$ is independent of the class of $Q$ modulo $\mathcal{O}(\mathbb{I})$.
Proof. By (3.10), for any $P \in G(q, d)$, $T_P G(q, d) = \{\Delta \in \text{Sym}_d : \Delta P + P\Delta = \Delta\}$. So,

$$ \Delta \in \mathbb{T}_0^i \iff \Delta \in \text{Sym}_d \text{ and } \Delta P_0^i + P_0^i \Delta = \Delta. \quad (3.19) $$

We readily check that the sets described in the rhs' in the statement of this Lemma satisfy (3.19). Thus, each of these sets is a linear subspace of $\mathbb{T}_0^i$, and clearly, its dimension is equal to that of $\mathbb{T}_0^i$. This concludes the proof. $\square$

**Corollary 3.** For any $d \times d$ matrix $A = (a_{k\ell})_{1\le k\le\ell\le d}$, its Frobenius norm is

$$ \|A\|_F := \text{tr}(A'A) = \sqrt{\sum (a_{k\ell})^2}. \quad (3.20) $$

Let $1 \le i \le r$. Then, for all $S \in \mathbb{T}_0^i$,

$$ \sum_{i=1}^{r} \|S\|_{F}^{2} = 4 \sum_{1 \le i < j \le r} \|S^{(i,j)}\|_{F}^{2}. \quad (3.21) $$

Proof. This identity is an easy consequence of Lemma 6 and (3.20). $\square$

## 3.3 Preliminaries for estimating the flag of PS's

The flag $F^I(\Sigma)$ of PS's is estimated by the flag of eigenprojections $F^I(\hat{\Sigma}_n)$, of type I, introduced in paragraph 3.3.1 below. Thus, we need a distance on $\mathcal{F}^I$ to measure the deviation between flags of type I. However, there is no known metric on $\mathcal{F}^I$ for which the geodesic distance on $\mathcal{F}^I$ is available in closed form. To overcome this lack, we introduce in (3.24) below, an extrinsic distance $\mathcal{D}_{\text{ex}}$ on $\mathcal{F}^I$ i.e. that is the restriction to $\mathcal{F}^I$ of a distance on $G^I$, which has an explicit expression.

### 3.3.1 Flag of eigenprojections

Let $\Sigma \in \text{Sym}_d$, whose flag of eigenspaces, denoted by $F^I(\Sigma)$, is of type I := $(q_i)_{1\le i\le r}$. For $S \in \text{Sym}_d^\ne$, let $\Gamma$ and $C$ be respective orthogonal matrices of eigenvectors of $\Sigma$ and $S$, associated to eigenvalues in decreasing order.

For all $1 \le i \le r$, $\pi^i(C^{(i)})$ is the projector onto a $q_i$-linear subspace spanned by eigenvectors of $C$. This projector is independent of the choice of $C$ i.e. depends only on $S$. Thus it is called the *i-th eigenprojection* of $S$ wrt I, denoted by $P_i(S)$. Then, the collection $F^I(S) := (P_i(S))_i$ form a flag, called the *flag of eigenprojections* of $S$ wrt I. By (3.7),

$$ F^I(S) = (P_i(S))_i = (\pi^i(C^{(i)}))_i = \pi^I(C). \quad (3.22) $$

Assume that $C \in W_{F^I(\Sigma)}$. Set $E := \Gamma'C$. Since $\Gamma$ spans $F^I(\Sigma)$, by Lemma 5, the geodesic paramatrization of $C$ wrt $(F^I(\Sigma), \Gamma)$ is fully determined by that of $E$ wrt $(\mathcal{P}_0^I, I_d)$, for which

$$ \pi^i(E^{(i)}) = \Gamma' \cdot \pi^i(C^{(i)}) = \Gamma' \cdot P_i(S), \quad 1 \le i \le r. \quad (3.23) $$

### 3.3.2 Extrinsic distance

Let $\mathcal{R} = (R_i)_i \in \mathcal{F}^I$ and $W_{\mathcal{R}}$ defined in (3.18). Then,

$$ \pi^I(W_{\mathcal{R}}) = \{\mathcal{P} = (P_i)_i \in \mathcal{F}^I : \forall 1 \le i \le r, P_i \notin \text{Cut}(R_i)\}. $$

Now, for $\mathcal{P} \in \pi^I(W_{\mathcal{R}})$, set

$$ \mathcal{D}_{\text{ex}}(\mathcal{P}, \mathcal{R}) := \left( \sum_{i=1}^{r} \delta^i(P_i, R_i)^2 \right)^{1/2} \quad \text{where } \delta^i(P_i, R_i) := \left\| \text{Log}_{P_i}^{G^i}(R_i) \right\|_F. \quad (3.24) $$
**Corollary 3.** For any $\alpha \in (0, 1)$, a confidence region for $F^1(\Sigma)$ of asymptotic level $(1 - \alpha)$ is given by

$$R_{n,\alpha} := \{\pi^1(Q) \in \mathcal{F}^1 : \hat{T}_n(Q) \le \chi_{\mathrm{D}^1}^2(1-\alpha)\},$$

where $\chi_{\mathrm{D}^1}^2(1 - \alpha)$ is the quantile of order $(1 - \alpha)$ of the $\chi_{\mathrm{D}^1}^2$ distribution.

**Corollary 4.** Consider the following null hypothesis assumption.

$$H_0: \pi^1(Q_0) = F^1(\Sigma), \text{ for } Q_0 \in O(d).$$

For any $\alpha \in (0, 1)$, consider the test which accepts $H_0$ when $\hat{T}_n(Q_0) \le \chi_{\mathrm{D}^1}^2(1 - \alpha)$ and reject $H_0$ else. Then, this test is of asymptotic level $\alpha$.

## 3.6 Conclusion

Given a normally distributed random vector $X$ valued in $\mathbb{R}^d$, $d \ge 1$, a geometric framework allows us to develop an asymptotic procedure to infer the set of all its PS's. In addition, we provide easily implementable tests concerning the complete collection of principal subspaces and confidence regions for them. These results opens many questions which could lead to useful extensions. Among them, two of the most striking questions are:

* How to estimate or learn the type I?

* Is it possible to relax the Gaussian assumption?

## Acknowledgements

The authors have received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation program (grant agreement G-Statistics No 786854). It was also supported by the French government through the 3IA Côte d'Azur Investments ANR-19-P3IA-0002 managed by the National Research Agency.

## References

Anderson, T.W. (1963). *Asymptotic theory for principal components analysis*. Ann. Math. Stat. (1963), Vol. 34, 122-148.

Batzies, E., Huper, K., Machado, L. and Silva Leite, F. (2015). *Geometric mean and geodesic regression on Grassmannians*. Linea Algebra Appl., 466:83-101, 2015.

Bhattacharya, R. and Patrangenaru, V. (2005). *Large sample theory of intrinsic and extrinsic sample means on manifolds-II*. The Annals of Statistics (2005), Vol. 33, No. 3, 1225-1259.

Bendokat, T., Zimmermann, R., Absil, P.A. (2020). *A Grassmann Manifold Handbook: Basic Geometry and Computational Aspects*. arXiv preprint arXiv:2011.13699, 2020.

Lee, J.M. (2018). *Introduction to Riemannian manifolds*. volume 176 of Graduate Texts in Mathematics. Springer, Cham, 2018.

Tyler, D.E. (1981). *Asymptotic inference for eigenvectors*. The Annals of Statistics 9 (4), 725-736, 1981.

Van der Vaart, A.W. (2000). *Asymptotic Statistics*. Cambridge Series in Statistical and Probabilistic Mathematics, Volume 3.

Ye, K., Wong, K.S.W., Lim, L.H. (2022). *Optimization on flag manifolds*. Mathematical Programming (Series A), 2022, 194(1-2), pp. 621-660.
Thus, $\delta^i$ is the geodesic distance on $G^i$. The distance $\mathcal{D}_{\text{ex}}$ is extrinsic, i.e. it is the restriction to $\mathcal{F}^I$ of a distance on $G^I$, and has an explicit expression, provided by that of $\text{Log}_{P_i}^{G^i}$. Clearly, the action of $O(d)$ on $\mathcal{F}^I$ preserves $\mathcal{D}_{\text{ex}}$, i.e. for all $Q \in O(d)$,

$$ \mathcal{D}_{\text{ex}}(Q * \mathcal{P}, Q * \mathcal{R}) = \mathcal{D}_{\text{ex}}(\mathcal{P}, \mathcal{R}). \qquad (3.25) $$

In particular, with the notations of paragraph 3.2.4, we deduce from (3.25) that for $E := \Gamma'C$,

$$ \mathcal{D}_{\text{ex}}(\pi^{\text{I}}(E), \mathcal{P}_{0}^{\text{I}}) = \mathcal{D}_{\text{ex}}(\Gamma' * \pi^{\text{I}}(C), \Gamma' * F^{\text{I}}(\Sigma)) = \mathcal{D}_{\text{ex}}(F^{\text{I}}(S), F^{\text{I}}(\Sigma)). $$

### 3.3.3 A discrepancy between flags

In Section 4, to obtain a pivotal statistic for the flag $F^I(\Sigma)$ of PS's, we need to standardize some Gaussian rv's. Due to such renormalizations, the resulting confidence regions for $F^I(\Sigma)$ are expressed in terms of a discrepancy between $F^I(\Sigma)$ and $F^I(\hat{\Sigma}_n)$, instead of the extrinsic distance $\mathcal{D}_{\text{ex}}$ between them. In fact, this discrepancy, defined in (3.27) below, is a deformation of $\mathcal{D}_{\text{ex}}$. This means that, without these renormalizations, this discrepancy is equal to $\mathcal{D}_{\text{ex}}$. Let $\mathcal{D}(I)$ be the set of all block-diagonal matrices whose structure agrees with the type I, i.e.

$$ M \in \mathcal{D}(I) \iff M = \text{Diag} \left( M^{(1,1)}, \dots, M^{(i,i)}, \dots, M^{(r,r)} \right). $$

In this paragraph, $\Gamma$ denotes a variable matrix in $O(d)$. Let $\mathbb{K} = (K^i)_{1 \le i \le r}$ be a collection of matrices such that for all $1 \le i \le r$, $K^i \in \mathcal{D}(I)$. Let $\mathcal{P} = (P_i)_i \in \mathcal{F}^I$. By Lemma 4,

$$ \Gamma \in \mathbb{W}_{\mathcal{P}} \iff \forall 1 \le i \le r, \Gamma' P_i \Gamma \notin \text{Cut}(P_0^i). \qquad (3.26) $$

Then, for $\Gamma \in \mathbb{W}_{\mathcal{P}}$, set

$$ \mathcal{D}_{\mathbb{K}}(\Gamma, \mathcal{P}) := \sqrt{\sum_{i=1}^{r} (\delta_{K^i}^i (\Gamma, P_i))^2} \quad \text{where} \quad \delta_{K^i}^i (\Gamma, P_i) := \left\| K^i \text{Log}_{P_0^i}^{G^i} (\Gamma' \cdot P_i) K^i \right\|_{\mathcal{F}}. \quad (3.27) $$

**Lemma 7.** Let $\mathcal{P} = (P_i)_i \in \mathcal{F}^I$ and $\Gamma \in \mathbb{W}_{\mathcal{P}}$. Then, $\mathcal{D}_{\mathbb{K}}(\Gamma, \mathcal{P})$ is independent of the class of $\Gamma$ modulo $O(I)$, i.e. for all $H \in O(I)$,

$$ \Gamma H \in \mathbb{W}_{\mathcal{P}} \quad \text{and} \quad \mathcal{D}_{\mathbb{K}}(\Gamma H, \mathcal{P}) = \mathcal{D}_{\mathbb{K}}(\Gamma, \mathcal{P}). \qquad (3.28) $$

*Proof.* Recall that $O(I)$ is the isotropy group of $\mathcal{P}_0^I$ under the action of $O(d)$ on $\mathcal{F}^I$, i.e. for all $H \in O(I)$ and $1 \le i \le r$, $H \cdot P_0^i = P_0^i$. So, by (3.26), for all $H \in O(I)$,

$$ \Gamma \in \mathbb{W}_{\mathcal{P}} \iff \Gamma H \in \mathbb{W}_{\mathcal{P}}. \qquad (3.29) $$

Then, by Corollary 1, for any $\mathcal{P} = (P_i)_i \in \mathcal{F}^I$, $\Gamma \in \mathbb{W}_{\mathcal{P}}$ and $H \in O(I)$,

$$ \text{Log}_{P_0^i}^{G^i} ((\Gamma H)' \cdot R_i) = H' (\text{Log}_{P_0^i}^{G^i} (\Gamma' \cdot R_i)) H. $$

This, combined to Remark 7 below, implies that for all $1 \le i \le r$,

$$ \delta_{K^i}^i (\Gamma H, P_i) = \left\| K^i (H' \text{Log}_{P_0^i}^{G^i} (\Gamma' \cdot P_i) H) K^i \right\|_{\mathcal{F}} = \left\| H' (K^i \text{Log}_{P_0^i}^{G^i} (\Gamma' \cdot P_i) K^i) H \right\|_{\mathcal{F}}. $$

Finally, by the properties of $\|\cdot\|_F$, the rhs hereabove is equal to $\delta_{K^i}^i(\Gamma, P_i)$. $\square$

**Remark 7.** For all $H \in O(I)$ and $K \in \mathcal{D}(I)$, $H$ and $K$ commute, i.e. $HK = KH$.

Now, we may introduce a discrepancy between flags of type I as follows.
**Definition 10.** Let $K = (K^i)_i \in (\mathcal{D}(I))^r$. Let $\mathcal{P}, \mathcal{R} \in \mathcal{F}^1$ such that $\mathcal{R} \in \pi^1(\mathbb{W}_\mathcal{P})$. Then, by Lemma 7, we may define the $K$-discrepancy between $\mathcal{R}$ and $\mathcal{P}$ as follows.

$$ \tilde{\mathcal{D}}_{\mathcal{K}}(\mathcal{R}, \mathcal{P}) := \mathcal{D}_{\mathcal{K}}(\Gamma, \mathcal{P}), \qquad (3.30) $$

for any $\Gamma \in O(d)$ such that $\pi^1(\Gamma) = \mathcal{R}$, which is defined modulo $O(1)$. In particular, it is easy to check that, for $\mathbb{I} := (I_d, \dots, I_d) \in (\mathcal{D}(I))^r$, we have that $\tilde{\mathcal{D}}_{\mathbb{I}}(\mathcal{R}, \mathcal{P}) = \mathcal{D}_{\text{ex}}(\mathcal{R}, \mathcal{P})$. Thus, the $K$-discrepancy $\tilde{\mathcal{D}}_{\mathcal{K}}$ is interpreted as a deformation of the extrinsic distance $\mathcal{D}_{\text{ex}}$.

**Remark 8.** In general, $\tilde{\mathcal{D}}_{\mathcal{K}}$ is not a distance on $\mathcal{F}^1$. However, we readily prove that $\tilde{\mathcal{D}}_{\mathcal{K}}$ fulfills the following properties. (i) The separation property holds, i.e. for all $\mathcal{R}, \mathcal{P} \in \mathcal{F}^1$

$$ \tilde{\mathcal{D}}_{\mathcal{K}}(\mathcal{R}, \mathcal{P}) = 0 \iff \mathcal{R} = \mathcal{P}. \qquad (3.31) $$

(ii) The $K$-discrepancy preserves the action of $O(d)$, i.e. for all $Q \in O(d)$,

$$ \tilde{\mathcal{D}}_{\mathcal{K}}(Q * \mathcal{R}, Q * \mathcal{P}) = \tilde{\mathcal{D}}_{\mathcal{K}}(\mathcal{R}, \mathcal{P}), \quad \mathcal{R}, \mathcal{P} \in \mathcal{F}^1. $$

# 4 Limit theorem and estimation of the flag of PS's

If $E_n \in \mathbb{W}_{\mathcal{P}_0^I}$, then we establish in Theorem 5 the limiting distribution of its geodesic parametrization wrt $(\mathcal{P}_0^I, I_d)$. In fact, $\text{Pr}(E_n \notin \mathbb{W}_{\mathcal{P}_0^I}) \to 0$ as $n \to \infty$: see (i) of Proposition 7 in Appendix. Thus, we obtain Gaussian and Haar limiting distributions, as in Theorem 1, where $E_n$ is parametrized by its usual coordinates. Then, we derive from (i) of Theorem 5, confidence regions for the flag $F^1(\Sigma)$ of PS's, in terms of a discrepancy (in the sense of Definition 10) between $F^1(\Sigma)$ and the flag $F^1(\hat{\Sigma}_n)$ of eigenprojections of $\hat{\Sigma}_n$, both of type I.

## 4.1 Statement of main theorem

First, for $1 \le i \le r$, we define the limit function $g_0^i$ valued in $\mathbb{T}_0^i$ involved in (4.1) below. For $u \in \text{Sym}_d$, set

$$ \mathcal{G}^{i,j}(u) := \frac{1}{\lambda_i - \lambda_j} u^{(i,j)}, \quad 1 \le i, j \le r \text{ with } i \neq j. $$

Then, consider the matrices $\mathcal{G}_<^i(u)$ and $\mathcal{G}_>^i(u)$ of respective sizes $q_i \times (q_1 + \dots + q_{i-1})$ and $q_i \times (q_{i+1} + \dots + q_r)$ defined by

$$ \mathcal{G}_<^i(u) := (\mathcal{G}^{i,1}(u) \quad \dots \quad \mathcal{G}^{i,i-1}(u)) \quad \text{and} \quad \mathcal{G}_>^i(u) := (\mathcal{G}^{i,i+1}(u) \quad \dots \quad \mathcal{G}^{i,r}(u)) $$

Finally, define the map $g_0^i : \text{Sym}_d \to \mathbb{T}_0^i$ by

$$ g_0^1(u) = \begin{pmatrix} 0_{q_1} & \mathcal{G}_>^1(u) \\ (\mathcal{G}_>^1(u))' & 0 \end{pmatrix}, \quad g_0^r(u) = \begin{pmatrix} 0 & (\mathcal{G}_<^r(u))' \\ (\mathcal{G}_<^r(u))' & 0_{q_r} \end{pmatrix} $$

$$ g_0^i(u) = \begin{pmatrix} 0 & (\mathcal{G}_<^i(u))' & 0 \\ \mathcal{G}_<^i(u) & 0_{q_i} & \mathcal{G}_>^i(u) \\ 0 & (\mathcal{G}_>^i(u))' & 0 \end{pmatrix}, \quad i \neq 1, r. $$

Now, we can state our main result hereafter, where the rv $U$ is defined in Theorem 2.

**Theorem 5.** (i) For any $1 \le i \le r$, the following CLT holds in $\mathbb{T}_0^i$.

$$ G_n^i := \sqrt{n} \text{Log}_{\mathcal{P}_0^i} (\pi^i (E_n^{(i)})) \xrightarrow{n \to \infty} g_0^i(U). \qquad (4.1) $$
(ii) For any $1 \le i \le r$, the following convergence holds in $O(q_i)$.

$$H_n^i := (I_d^{(i)})' \left( \mathcal{H}_{(\pi^i(E_n^{(i)}), P_0^i)}^{q_i} (E_n^{(i)}) \right) \xrightarrow{n \to \infty} \mathcal{E}^{i,i}(U), \quad (4.2)$$

where the distribution of $\mathcal{E}^{i,i}(U)$ is the conditional Haar invariant distribution.

**Remark 9.** In conclusion, we conjecture a possible application of (4.2).

## 4.2 Convergence to a $\chi^2$ distribution

In Corollary 4, we derive a pivotal statistic for the matrix $\Gamma$, i.e. which depends only on the sample, and whose limit distribution does not depend on any unknown parameter. Thus, Lemma 8 allows to concatenate the CLT's of (4.1) for $1 \le i \le r$, providing the convergence of a pivotal statistic to a $\chi^2$ distribution. Namely, by Theorem 2, for $i \ne j$, the entries of $(g_0^i(U))^{(i,j)} = \frac{1}{\lambda_i - \lambda_j} U^{(i,j)}$ are real i.i.d rv's $\mathcal{N}(0, \sigma_{i,j}^2)$, of standard deviation $\sigma_{i,j} = \frac{\sqrt{\lambda_i \lambda_j}}{|\lambda_i - \lambda_j|}$. We normalize these entries by setting

$$\bar{g}_0^i(U) := K^i g^i(U) K^i \quad \text{where} \quad K^i := \operatorname{Diag} \left( \frac{1}{\sigma_{i,1}} I_{q_1}, \dots, I_{q_i}, \dots, \frac{1}{\sigma_{i,r}} I_{q_r} \right). \quad (4.3)$$

**Lemma 8.** (i) The non-null entries of the matrix $\bar{g}_0^i(U)$ are real iid rv's $\mathcal{N}(0, 1)$.

(ii) The blocks $\{ (\bar{g}_0^i(U))^{(i,j)} : 1 \le i < j \le r \}$ are mutually independent.

(iii) The real rv $\frac{1}{4} \sum_{i=1}^r \|\bar{g}_0^i(U)\|_F^2$ is distributed as a $\chi_{\text{D}^1}^2$, where

$$\text{D}^I := \frac{1}{2} \left( d^2 - \sum_{i=1}^{r} q_i^2 \right).$$

*Proof.* Clearly, $\bar{g}_0^i(U)$ is a random matrix valued in $\mathbb{T}_0^i$ and for $i \ne j$,

$$\left(\bar{g}_0^i(U)\right)^{(i,j)} = \frac{1}{\sigma_{i,j}} \left( \frac{1}{\lambda_i - \lambda_j} U^{(i,j)} \right). \quad (4.4)$$

This proves (i) and (ii). Then, (iii) holds, by Corollary 3 applied with $S = \bar{g}_0^i(U)$. $\square$

**Proposition 4.** For $(K^i)_{1 \le i \le r}$ defined in (4.3),

$$\frac{n}{4} \sum_{i=1}^{r} \left\| K^i \operatorname{Log}_{P_0^i}^{G^i} \left( \Gamma' \cdot P_i(\hat{\Sigma}_n) \right) K^i \right\|_F^2 \xrightarrow{n \to \infty} \chi_{\text{D}^1}^2, \quad (4.5)$$

*Proof.* By (i) of Theorem 5 and (iii) of Lemma 8,

$$\frac{n}{4} \sum_{i=1}^{r} \left\| K^i \operatorname{Log}_{P_0^i}^{G^i} \left( \pi^i(E_n^{(i)}) \right) K^i \right\|_F^2 \xrightarrow{n \to \infty} \frac{1}{4} \sum_{i=1}^{r} \left\| \bar{g}_0^i(U) \right\|_F^2 \sim \chi_{\text{D}^1}^2. \quad (4.6)$$

Now, by (3.23) applied to $S = \hat{\Sigma}_n \in \operatorname{Sym}_d^\ne$ a.s., we have that $\pi^i(E_n^{(i)}) = \Gamma' \cdot P_i(\hat{\Sigma}_n)$ a.s. So, by injecting this relation in the rhs of (4.6), we deduce that (4.5) holds. $\square$

However, in the rhs of (4.5), the matrix $K^i$ still depends on the unknown parameters $(\lambda_i)_i$. Thus, in $K_i$, we replace $\lambda_i$ by $\tilde{\lambda}_n^i$, where (with the notations of paragraph 2.1), denoting by $(\mu_k(\hat{\Sigma}_n))_{1 \le k \le d}$ the eigenvalues of $\hat{\Sigma}_n$ and by $(\beta_i)_i$ the partition of $\{1, ..., d\}$ wrt I,

$$\tilde{\lambda}_n^i := \frac{1}{q_i} \sum_{k \in \beta_i} \mu_k(\hat{\Sigma}_n) \xrightarrow{n \to \infty} \lambda_i.$$
**Corollary 4.** The statistic $\hat{S}_n$ defined hereafter is a pivotal statistic for $\Gamma$.

$$ \hat{S}_n := \frac{n}{4} \sum_{i=1}^{r} \left\| \hat{K}_n^i \operatorname{Log}_{P_0^i} (\Gamma' \cdot P_i (\hat{\Sigma}_n)) \hat{K}_n^i \right\|_F^2 \xrightarrow[n \to \infty]{d} \chi_{D^1}^2, \quad (4.7) $$

where, for $1 \le i \le r$, $\hat{K}_n^i$ is derived from $K^i$ by replacing $\lambda_i$ by $\hat{\lambda}_n^i$.

*Proof.* Since $\hat{K}_n^i \to K^i$ a.s. as $n \to \infty$, this follows from Proposition 4. □

### 4.3 Confidence regions

For $n \ge 1$, set $\hat{K}_n := (\hat{K}_n^i)_{1 \le i \le r}$. Thus, by definition,

$$ \hat{S}_n = \frac{n}{4} \left[ \mathfrak{D}_{\hat{K}_n} (\Gamma, F^1(\hat{\Sigma}_n)) \right]^2, \quad (4.8) $$

where $\mathfrak{D}_{\hat{K}_n}$ is defined in (3.27). Thus, $\tilde{\mathfrak{D}}_{\hat{K}_n}$ is the $\hat{K}_n$-discrepancy, in the sense of Definition 10. Since $\pi^1(\Gamma) = F^1(\Sigma)$, Corollary 4 and (4.8) imply that

$$ \hat{S}_n = \frac{n}{4} \left[ \tilde{\mathfrak{D}}_{\hat{K}_n} (F^1(\Sigma), F^1(\hat{\Sigma}_n)) \right]^2 \xrightarrow[n \to \infty]{d} \chi_{D^1}^2. \quad (4.9) $$

We deduce hereafter confidence regions for $F^1(\Sigma)$ of the desired form.

**Proposition 5.** For any $\alpha \in (0, 1)$, let $\chi_{D^1}^2(1 - \alpha)$ be the quantile of order $(1 - \alpha)$ of the $\chi_{D^1}^2$ distribution. For $n \ge 1$, set

$$ R_{n,\alpha} := \left\{ \mathcal{P} \in \mathcal{F}^1 : \frac{n}{4} \left[ \tilde{\mathfrak{D}}_{\hat{K}_n} (\mathcal{P}, F^1(\hat{\Sigma}_n)) \right]^2 \le \chi_{D^1}^2(1-\alpha) \right\}. $$

Then, by (4.9), $R_{n,\alpha}$ is a confidence region for $F^1(\Sigma)$ of asymptotic level $(1 - \alpha)$. By (3.31),

$$ \tilde{\mathfrak{D}}_{\hat{K}_n} (\mathcal{P}, F^1(\hat{\Sigma}_n)) = 0 \iff Q = F^1(\hat{\Sigma}_n). $$

So, $R_{n,\alpha}$ is interpreted as a full deformed ellipsoid in $\mathcal{F}^1$, whose center is $F^1(\hat{\Sigma}_n)$.

**Corollary 5.** Fix $Q_0 \in O(d)$. Consider the following null hypothesis assumption. $H_0: \pi^1(Q_0) = F^1(\Sigma)$. For any $\alpha \in (0, 1)$, consider the test which accepts $H_0$ when $\pi^1(Q_0) \in R_{n,\alpha}$ and rejects $H_0$ else. Then, this test is of asymptotic level $\alpha$.

### 4.4 Simulation

Figure 1 illustrates the convergence in distribution of $\hat{T}_n$ to a $\chi^2$ distribution. The parameters for this simulation are $d = 4$ and $I = (1, 1, 1, 1)$. Then, $D^1 = 6$. For the sample size, we take $n = 10000$. The histogram in blue represents the distribution of $\hat{T}_n$ and the curve in red is that of the probability distribution function of the $\chi_6^2$ distribution. We see that the distribution of $\hat{T}_n$ is indeed very close to that of the $\chi_6^2$ one.

## 5 Conclusion

Given a Gaussian random vector $X$ whose covariance matrix $\Sigma$ has possibly repeated eigenvalues, we develop a geometric method to estimate the flag $F^1(\Sigma)$ of its PS's, for which we provide confidence regions and implementable tests. These results open many questions, among which the following ones, which we present with their motivations.
<figure><img src="image_5.png" /><figcaption>Figure 1: Illustration of the convergence of $\hat{T}_n$ to the $\chi_6^2$ distribution</figcaption></figure>

(i) Are our geometric CLT's of Theorem 5 valid when the distribution of $X$ is elliptic? If so, then, for such a distribution, the estimation of $F^I(\Sigma)$ would be derived as in paragraphs 4.2 and 4.3, which are only based on such CLT's. This question is motivated by [14]. Therein, when $X$ is elliptic, a CLT for each $P_i(\hat{\Sigma}_n)$ is obtained, but not of the form of Eq.(4.1)

(ii) Does $H_n^i$ converge in $O(q_i)$ to $E^{i,i}$ wrt the Kullback-Leibler (KL) divergence?

By Remark 5, among all $q_i$-frames generating $rg(P_0^i)$, we have that $\mathcal{H}_{(\Gamma' P_i(\hat{\Sigma}_n)\Gamma, P_0^i)}^{q_i}(E_n^{(i)})$ minimizes the geodesic distance in $\text{St}_{q_i}$ to $E_n^{(i)}$. We conjecture that this geometric optimality should imply a stronger mode of convergence of $H_n^i$ to $E^{i,i}$, i.e. wrt the KL divergence in the compact group $O(q_i)$. Our conjecture is motivated by the case of iid convolutions on compact groups, for which the convergence in distribution to the Haar measure has been first obtained: See [6] and references therein. Later, it was proved in [5] or [7], by information theoretic methods, that they converge wrt the KL divergence.

# 6 Appendix: Method for the proof of Theorem 5

This Appendix is devoted to the proof of Theorem 5 and to some tools developed for it. The order of appearance of the proofs does not follow that of the preceding sections. Instead, the most technical parts are presented at the end. Thus, the proof of the generalized $\delta$-method of Theorem 3 is deferred to the end of this Appendix, i.e. to subsection 6.5.

Now, we describe the proof of Theorem 5. By Theorem 2, $U_n$ converges in distribution to $U$. Thus, we prove Theorem 5 by expressing, in (6.7), the $\text{lhs}$ of (4.1) and (4.2), i.e. $G_n^i$ and $H_n^i$, as functions of $U_n$ and then by deriving, from the generalized $\delta$-method of Theorem 3, that these functions of $U_n$ converge to functions of $U$ which are the $\text{rhs}$ of (4.1) and (4.2). However, $G_n^i$ and $H_n^i$ are defined only when $\pi^i(E_n^{(i)}) \notin \text{Cut}(P_0^i)$. So, instead of applying Theorem 3 with $U_n$, we apply it with truncations of $U_n$, in the sense of paragraph 6.2.2.

## 6.1 Preliminaries from Anderson's proof

First, we present results from the proof of Theorem 1 which will be used later for ours. Consider the open subset of $\text{Sym}_d$ defined by

$$ \mathbb{S}_0 := \{ S \in \text{Sym}_d : S^{(i,i)} \in \text{Sym}_{q_i}^{\neq}, 1 \le i \le r \}. \qquad (6.1) $$
Let $(u_n)_{n \ge 1}$ be a sequence valued in $((\mathbb{S}_n^A)^*)_{n \ge 1}$ which converges to $u \in \mathbb{S}_0$. Set

$$ \mathcal{E}_n := e_n(u_n). \qquad (6.2) $$

**Proposition 6.** For $1 \le i, j \le r$ with $i \ne j$, we have that $\mathcal{E}_n^{(i,i)}$ and $\sqrt{n}\mathcal{E}_n^{(i,j)}$ converge to finite limits denoted respectively by $\mathcal{E}^{i,i}(u)$ and $\mathcal{F}^{i,j}(u)$ which satisfy

$$ \lambda_i \mathcal{G}^{i,j}(u) + \lambda_j \mathcal{L}^{i,j}(u) = u^{(i,j)} \quad \text{and} \quad \mathcal{G}^{i,j}(u) + \mathcal{L}^{i,j}(u) = 0, \qquad (6.3) $$

with $\mathcal{G}^{i,j}(u) := \mathcal{E}^{i,i}(u) (\mathcal{F}^{j,i}(u))'$ and $\mathcal{L}^{i,j}(u) := \mathcal{F}^{i,j}(u) (\mathcal{E}^{j,j}(u))'$. Thus, by (6.3),

$$ \mathcal{G}^{i,j}(u) = \frac{1}{\lambda_i - \lambda_j} u^{(i,j)}. \qquad (6.4) $$

Recall that $U$ is the limit rv of Theorem 2. Then, the following Lemma is obtained in [1].

**Lemma 9.** $\mathcal{E}^{i,i}(U)$ is distributed as $E^{i,i}$ in (i) of Theorem 1 and $\mathcal{F}^{i,j}(U)$ is distributed as $F^{i,j}$ in (ii) of Theorem 1.

## 6.2 Expression of $G_n^i$ and $H_n^i$ as functions of a truncation of $U_n$

### 6.2.1 Expression of $G_n^i$ and $H_n^i$ as functions of $U_n$

First, in (6.5), we express $E_n$ in function of $U_n$. Recall that $E_n = \psi(T_n)$, where $\psi: \text{Sym}_d^\mp \to O(d)$ is the eigenvector map of Definition 1. In view of truncations, we extend the map $\psi$ as follows: see Remark 10.

**Definition 11.** Let $\bar{\psi}: \text{Sym}_d^\mp \cup \{\Delta\} \to O(d)$ be the map such that for $S \in \text{Sym}_d^\mp$, $\bar{\psi}(S) = \psi(S)$ and $\bar{\psi}(\Delta) = I_d$. Then, $\bar{\psi}$ is called the extended eigenvector map.

(i) By definition, $U_n := \sqrt{n}(T_n - \Delta)$. Then, $T_n = \phi_n(U_n)$, where

$$ \phi_n(u) = \Delta + n^{-1/2}u, \quad u \in \text{Sym}_d. $$

So, $E_n = \psi(\phi_n(U_n))$ provided that $U_n \in \phi_n^{-1}(\text{Sym}_d^\mp)$. Thus, consider the sets

$$ (\mathbb{S}_n^A)^* := \phi_n^{-1}(\text{Sym}_d^\mp) \quad \text{and} \quad \mathbb{S}_n^A := (\mathbb{S}_n^A)^* \cup \{0\}. $$

Since $\phi_n(0) = \Delta$, we may define the map $e_n: \mathbb{S}_n^A \to O(d)$ by $e_n(u) = \bar{\psi}(\phi_n(u))$. Then,

$$ U_n \in (\mathbb{S}_n^A)^* \implies E_n = e_n(U_n). \qquad (6.5) $$

(ii) Now, for $1 \le i \le r$, we express $G_n^i$ and $H_n^i$ as functions of $U_n$. Define the map $p_n^i$ by

$$ p_n^i: \mathbb{S}_n^A \to \mathcal{G}^i \quad \text{and} \quad p_n^i(u) = \pi^i(e_n(u)^{(i)}). $$

If $p_n^i(u) \notin \text{Cut}(P_0^i)$, then $g_n^i(u)$ and $h_n^i(u)$ introduced hereafter are well-defined.

$$ g_n^i(u) := \sqrt{n} \text{Log}_{P_0^i} G_n^i(p_n^i(u)) \in \mathbb{T}_0^i \quad \text{and} \quad h_n^i(u) := \mathcal{H}_{(p_n^i(u), P_0^i)}^i (e_n(u)^{(i)}) \in \text{St}^i \qquad (6.6) $$

Thus, we consider the sets

$$ \mathbb{W}_n^i := \{u \in \text{Sym}_d : p_n^i(u) \notin \text{Cut}(P_0^i)\} \quad \text{and} \quad (\mathbb{S}_n^i)^* := (\mathbb{S}_n^A)^* \cap \mathbb{W}_n^i. $$

Now, Equation (6.5) implies that $\pi^i(E_n^{(i)}) = p_n^i(U_n)$ provided that $U_n \in (\mathbb{S}_n^A)^*$. So,

$$ U_n \in (\mathbb{S}_n^i)^* \implies G_n^i = g_n^i(U_n) \quad \text{and} \quad H_n^i = (I_d^{(i)})' h_n^i(U_n). \qquad (6.7) $$
### 6.2.2 Truncations of rv's

Let $M$ be a metric space with Borel $\sigma$-algebra $\mathcal{B}$.

**Definition 12.** Let $Z: (\Omega, \mathcal{A}) \to (M, \mathcal{B})$ be a rv and $A \in \mathcal{A}$. If $M$ is a vector space, let $Z1_A^+$ be the rv such that $(Z1_A^+)(\omega) = Z(\omega)$ if $\omega \in A$ and $(Z1_A^+)(\omega) = 0 \in M$ else. If $M$ is a subgroup of $GL_q(\mathbb{R})$, denote by $Z1_A^\times$ the rv such that, $(Z1_A^\times)(\omega) = Z(\omega)$ if $\omega \in A$ and $(Z1_A^\times)(\omega) = I_q$ else. Then, we say that the rv's $Z1_A^+$ and $Z1_A^\times$ are truncations of $Z$ wrt $A$.

**Lemma 10.** Let $(X_n)_{n \ge 1}$ be a sequence of rv's valued in the metric space $M$ and $(A_n)_{n \ge 1}$ a sequence of events in $\mathcal{A}$. Assume that $X_n 1_{A_n}^+ \xrightarrow{d} X$ or $X_n 1_{A_n}^\times \xrightarrow{d} X$ as $n \to \infty$, where $X$ is a rv valued in $M$. If $P(A_n) \xrightarrow{n \to \infty} 1$, then $X_n \xrightarrow{n \to \infty} X$ also.

*Proof.* The proof follows readily from the definitions. $\square$

### 6.2.3 Expression of $G_n^i$ and $H_n^i$ as functions of $V_n^i$

The rv $U_n$ is valued in the vector space $\text{Sym}_d$. So, for $1 \le i \le r$, we define below a truncation of $U_n$. Thus, set

$$V_n^i := U_n 1_{A_n^i} \quad \text{where} \quad A_n^i := \{U_n \in (\mathbb{S}_n^i)^*\}.$$

So, $V_n^i$ is valued in the set $\mathbb{S}_n^i$ defined hereafter, which holds for all $\omega \in \Omega$ and not only a.s.

$$\mathbb{S}_n^i := (\mathbb{S}_n^i)^* \cup \{0\}.$$

**Remark 10.** By definition of the map $\bar{\psi}$, $e_n(0) = \bar{\psi}(\phi_n(0)) = \bar{\psi}(\Delta) = I_d$, so that $p_n^i(0) = \pi^i(I_d^{(i)}) = P_0^i$. So, $g_n^i(u)$ and $h_n^i(u)$ are also well-defined when $u = 0$, for which

$$g_n^i(0) = 0 \in \mathbb{T}_0^i \quad \text{and} \quad (I_d^{(i)})' h_n^i(0) = I_{q_i}. \qquad (6.8)$$

So, we may consider the maps $g_n^i: \mathbb{S}_n^i \to \mathbb{T}_0^i$ and $h_n^i: \mathbb{S}_n^i \to \text{St}^i$ defined by (6.6). Then, by (6.7) and (6.8), we express hereafter $G_n^i 1_{A_n^i}$ and $H_n^i 1_{A_n^i}$ as functions of $V_n^i$:

$$G_n^i 1_{A_n^i} = g_n^i(V_n^i) \quad \text{and} \quad H_n^i 1_{A_n^i} = (I_d^{(i)})' h_n^i(V_n^i). \qquad (6.9)$$

## 6.3 Proof of Theorem 5

We prove Theorem 5 by applying Theorem 3. Thus, we define in Table 1 the elements to which we apply Theorem 3. In Table 1, $\mathbb{S}_0$ is the subset defined in (6.1) and the map $\mathcal{E}^{i,i}: \mathbb{S}_0 \to O(q_i)$ is defined in Proposition 6.

<table><caption>Table 1: Elements to which Theorem 3 is applied</caption><thead><tr><th></th><th>proof of (4.1)</th><th>proof of (4.2)</th></tr></thead><tbody><tr><td>Metric spaces S and T</td><td>S = Sym<sub>d</sub> and T = T<sub>0</sub><sup>i</sup></td><td>S = Sym<sub>d</sub> and T = O(q<sub>i</sub>)</td></tr><tr><td>Function g<sub>n</sub></td><td>g<sub>n</sub><sup>i</sup> : S<sub>n</sub><sup>i</sup> → T<sub>0</sub><sup>i</sup></td><td>(I<sub>d</sub><sup>(i)</sup>)' h<sub>n</sub><sup>i</sup> : S<sub>n</sub><sup>i</sup> → O(q<sub>i</sub>)</td></tr><tr><td>Domain S<sub>n</sub> and subset S<sub>n</sub><sup>*</sup></td><td>S<sub>n</sub><sup>i</sup> and (S<sub>n</sub><sup>i</sup>)<sup>*</sup></td><td></td></tr><tr><td>Random variable V<sub>n</sub></td><td>V<sub>n</sub><sup>i</sup> = U<sub>n</sub> 1<sub>{U<sub>n</sub>∈(S<sub>n</sub><sup>i</sup>)<sup>*</sup>}</sub></td><td></td></tr><tr><td>Limit V of V<sub>n</sub></td><td>V<sub>n</sub><sup>i</sup> d n→∞ U</td><td></td></tr><tr><td>Limit function g<sub>0</sub></td><td>g<sub>0</sub><sup>i</sup> : S<sub>0</sub> → T<sub>0</sub><sup>i</sup></td><td>ε<sup>i,i</sup> : S<sub>0</sub> → O(q<sub>i</sub>)</td></tr></tbody></table>
### 6.3.1 Conditions to apply Theorem 3

We prove hereafter that all the conditions of Theorem 3 hold with the elements of Table 1. Firstly, by (ii) of Proposition 7 below, for all $1 \le i \le r$, $V_n^i \xrightarrow{d} U$ as $n \to \infty$. Then, by (iii) of Proposition 7, Assumption (2.2) holds. Finally, by Propositions 8 and 9 below, Assumption (2.3) holds.

**Proposition 7.** (i) For all $1 \le i \le r$, $\Pr(A_n^i) \to 1$ as $n \to \infty$, where we recall that $A_n^i := \{U_n \in (\mathbb{S}_n^i)^*\}$. In particular, we deduce that

$$ \Pr(E_n \in \mathbb{W}_{P_0^1}) = \Pr(\forall 1 \le i \le r, U_n \in (\mathbb{S}_n^i)^*) \xrightarrow{n \to \infty} 0. $$

(ii) For all $1 \le i \le r$, the sequence $(V_n^i)_n$ converges in distribution to $U$ as $n \to \infty$.

(iii) For all $1 \le i \le r$, $\Pr(V_n^i \in (\mathbb{S}_n^i)^*) \to 1$ as $n \to \infty$.

*Proof.* See paragraph 6.4.1. □

**Proposition 8.** Let $(u_n)_{n \ge 1}$ be a sequence valued in $((\mathbb{S}_n^i)^*)_{n \ge 1}$ which converges to $u \in \mathbb{S}_0$. Then, $g_n^i(u_n) \to g_0^i(u)$ as $n \to \infty$.

*Proof.* See paragraph 6.4.5. □

**Proposition 9.** Let $(u_n)_{n \ge 1}$ be as in Proposition 8. Then, $(I_d^{(i)})' h_n^i(u_n) \to \mathcal{E}^{i,i}(u)$ as $n \to \infty$, where $\mathcal{E}^{i,i}(u)$ is defined in Proposition 6.

*Proof.* See paragraph 6.4.3. □

### 6.3.2 End of proof of Theorem 5

We deduce that we may apply Theorem 3, from which (6.9) implies that

$$ G_n^i \mathbf{1}_{A_n^i}^+ = g_n^i(V_n^i) \xrightarrow{n \to \infty} g_0^i(U) \quad \text{and} \quad H_n^i \mathbf{1}_{A_n^i}^* = (I_d^{(i)})' h_n^i(V_n^i) \xrightarrow{n \to \infty} \mathcal{E}^{i,i}(U). \quad (6.10) $$

Finally, by Proposition 7, $\Pr(A_n^i) \to 1$ as $n \to \infty$. So, by Lemma 10 and (6.10), we conclude that (4.1) and (4.2) hold, which proves Theorem 5.

## 6.4 Proof of Proposition 7, 8 and 9

### 6.4.1 Proof of Proposition 7

*Proof.* (i) By definition, $U_n \in (\mathbb{S}_n^i)^*$ iff $U_n \in (\mathbb{S}_n^A)^*$ and $p_n^i(U_n) \notin \text{Cut}(P_0^i)$. Recall that $U_n \in (\mathbb{S}_n^A)^* \implies p_n^i(U_n) = \pi^i(E_n^{(i)})$ and that $U_n \in (\mathbb{S}_n^A)^*$ a.s. So, since $P_0^i = I_d^{(i)}(I_d^{(i)})'$,

$$ \Pr(U_n \in (\mathbb{S}_n^i)^*) = \Pr(\pi^i(E_n^{(i)}) \notin \text{Cut}(P_0^i)) = \Pr(E_n^{(i)}(E_n^{(i)})' \notin \text{Cut}(I_d^{(i)}(I_d^{(i)})').) $$

Now, by the description of the cut locus in Grassmannians given in Lemma 4,

$$ E_n^{(i)}(E_n^{(i)})' \notin \text{Cut}(P_0^i) \iff \text{rk}((I_d^{(i)})' E_n^{(i)}) = q_i \iff \text{rk}(E_n^{(i,i)}) = q_i. $$

So, (i) holds provided that $\Pr(\text{rk}(E_n^{(i,i)}) = q_i) \to 1$ as $n \to \infty$, which we derive by applying Lemma 11 below. We check hereafter that its assumptions hold. First, by Theorem
1, $E_n^{(i,i)}$ converges in distribution to $E^{i,i}$, where $rk(E^{i,i}) = q_i$ a.s. Furthermore, for all $n \ge 1$, $rk(E_n^{(i,i)}) \le q_i$. So, the assumptions of Lemma 11 hold, from which (i) is deduced.

(ii) and (iii): By definition of $(\mathbb{S}_n^i)^*$ and $V_n^i$,

$$ \mathrm{Pr} \left( U_n \in (\mathbb{S}_n^i)^* \right) \le \mathrm{Pr} \left( V_n^i = U_n \text{ and } U_n \in (\mathbb{S}_n^i)^* \right). \quad (6.11) $$

On the one hand, the rhs of (6.11) is bounded by $\mathrm{Pr}(V_n^i = U_n)$. So, by (i), $\mathrm{Pr}(V_n^i = U_n) \to 1$ as $n \to \infty$. Since $(U_n)$ converges in distribution to $U$, we deduce that (ii) holds. On the other hand, the rhs of (6.11) is bounded by $\mathrm{Pr}(V_n^i \in (\mathbb{S}_n^i)^*)$. So (i) implies that (iii) holds. $\square$

The following Lemma is used in the proof of Proposition 7. It is proved in [14].

**Lemma 11.** Let $B_n$ and $B$ be random matrices of same size. Assume that $B_n \xrightarrow{n \to \infty} B$, that, a.s., $rk(B) = b$ and that $\mathrm{Pr}(\mathrm{rk}(B_n) \le b) \xrightarrow{n \to \infty} 1$. Then, $\mathrm{Pr}(\mathrm{rk}(B_n) = b) \xrightarrow{n \to \infty} 1$.

### 6.4.2 Preliminaries for the proofs of Propositions 8 and 9

Let $(u_n)_{n \ge 1}$ be a sequence valued in $(\mathbb{S}_n^i)^*$ which converges to $u \in \mathbb{S}_0$, where $\mathbb{S}_0$ is defined in (6.1). Set

$$ \mathcal{E}_n := e_n(u_n) \quad \text{and} \quad \mathcal{P}_n^i := p_n^i(u_n), \quad n \ge 1, 1 \le i \le r. \quad (6.12) $$

Since $(\mathbb{S}_n^i)^* \subset (\mathbb{S}_n^A)^*$, Proposition 6 still holds with this sequence $(u_n)_{n \ge 1}$.

**Lemma 12.** For $n \ge 1$, there exist matrices $\mathcal{A}_n$ and $\mathcal{B}_n$ such that

$$ \mathcal{E}_n = \mathcal{A}_n + n^{-1/2} \mathcal{B}_n, \quad (6.13) $$

where $\mathcal{A}_n \in \mathcal{D}(I)$ i.e. $\mathcal{A}_n = \mathrm{Diag}(\mathcal{A}_n^{(1,1)}, \dots, \mathcal{A}_n^{(i,i)}, \dots, \mathcal{A}_n^{(r,r)})$ and the sequences $(\mathcal{A}_n)_{n \ge 1}$ and $(\mathcal{B}_n)_{n \ge 1}$ converge respectively to finite limits $\mathcal{A}(u)$ and $\mathcal{B}(u)$ such that

$$ \mathcal{A}(u) \in \mathcal{D}(I), \quad \mathcal{A}(u)^{(i,i)} = \mathcal{E}^{i,i}(u) \quad \text{and} \quad \mathcal{B}(u)^{(i,j)} = \mathcal{F}^{i,j}(u). \quad (6.14) $$

*Proof.* Set $\mathcal{A}_n := \mathrm{Diag}(e_n^1(u_n), \dots, e_n^i(u_n), \dots, e_n^r(u_n)) \in \mathcal{D}(I)$ and $\mathcal{B}_n$ defined by: $\mathcal{B}_n^{(i,j)} = f_n^{i,j}(u_n)$ if $i \ne j$ and $\mathcal{B}_n^{(i,i)} = 0$. So, $e_n^i(u_n) = \mathcal{E}_n^{(i,i)}$ and for $j \ne i$, we have that $f_n^{i,j}(u_n) = \sqrt{n}\mathcal{E}_n^{(i,j)}$. By Proposition 6, $\mathcal{A}_n$ and $\mathcal{B}_n$ converge to limits satisfying (6.14). $\square$

**Corollary 6.** With the notations of Lemma (12), we have that

$$ \mathcal{E}_n^{(i)} \xrightarrow{n \to \infty} (0 \quad \dots \quad \mathcal{E}^{i,i}(u) \quad \dots \quad 0)' \quad (6.15) $$

Recall that we have set $\mathcal{P}_n^i := p_n^i(u_n)$. Then,

$$ \mathcal{P}_n^i = (\mathcal{E}_n^{(i)})' (\mathcal{E}_n^{(i)})' = \mathcal{A}_n^{(i)} (\mathcal{A}_n^{(i)})' + n^{-1/2} \Gamma_n^i + n^{-1} \Phi_n^i, \quad (6.16) $$

where $\Gamma_n^i := \mathcal{A}_n^{(i)} (\mathcal{B}_n^{(i)})' + \mathcal{B}_n^{(i)} (\mathcal{A}_n^{(i)})'$ and $\Phi_n^i := \mathcal{B}_n^{(i)} (\mathcal{B}_n^{(i)})'$. Furthermore,

$$ \mathcal{P}_n^i \xrightarrow{n \to \infty} P_0^i \quad \text{as } n \to \infty. \quad (6.17) $$

*Proof.* (6.15) holds by Lemma 12 and (6.16) by (6.13). Finally, by (6.16), since $\mathcal{A}(u) \in \mathcal{D}(I)$,

$$ \mathcal{P}_n^i \xrightarrow{n \to \infty} (\mathcal{A}(u)^{(i)})' (\mathcal{A}(u)^{(i)})' = \mathrm{Diag} \left( 0_{q_1}, \dots, (\mathcal{A}(u)^{(i,i)})' (\mathcal{A}(u)^{(i,i)})', \dots, 0_{q_r} \right). \quad (6.18) $$

By (6.14) and (6.3), $\mathcal{A}(\mathcal{U})^{(i,i)} = \mathcal{E}^{i,i}(u) \in \mathcal{O}(q_i)$, so that the rhs of (6.18) is equal to $P_0^i$. $\square$
### 6.4.3 Proof of Proposition 9

We compute the limit of $h_n^i(u_n)$. By definition,

$$h_n^i(u_n) = \mathcal{H}_{(\mathcal{P}_n^i, P_0^i)}^{i} (\mathcal{E}_n^{(i)}) = (\mathcal{E}_n^{(i)} \quad \Delta_n^i \mathcal{E}_n^{(i)}) (\exp_m \begin{pmatrix} 0 & -\mathcal{C}_n^i \\ I_{q_i} & 0 \end{pmatrix}) \begin{pmatrix} I_{q_i} \\ 0_{q_i} \end{pmatrix}, \quad (6.19)$$

where $\Delta_n^i := \text{Log}_{P_n^i}^{G_n^i}(P_0^i)$ and $\mathcal{C}_n^i := (\mathcal{E}_n^{(i)})' (\Delta_n^i)^2 (\mathcal{E}_n^{(i)})$. By (6.17), $\Delta_n^i$ converges to $\text{Log}_{P_0^i}^{G_0^i}(P_0^i) = 0 \in \mathbb{T}_0^i$, i.e. $\mathcal{C}_n^i$ converges to 0. So, setting $h^i(u) := (0 \quad \dots \quad \mathcal{E}^{i,i}(u) \quad \dots \quad 0)' \in \text{St}^i$,

$$h_n^i(u_n) \underset{n \to \infty}{\longrightarrow} (h^i(u) \quad 0) \begin{pmatrix} \exp_m \begin{pmatrix} 0 & 0 \\ I_{q_i} & 0 \end{pmatrix} \end{pmatrix} \begin{pmatrix} I_{q_i} \\ 0_{q_i} \end{pmatrix}. \quad (6.20)$$

We notice that the rhs of (6.20) is equal to $\mathcal{H}_{(P_0^i, P_0^i)}^i(h^i(u)) = h^i(u)$. Therefore, by (6.20), $(I_d^{(i)})' h_n^i(u_n)$ converges to $(I_d^{(i)})' h^i(u) = \mathcal{E}^{i,i}(u)$ as $n \to \infty$. This concludes the proof.

### 6.4.4 Technical preliminaries for the proof of Proposition 8

**Lemma 13.** Let $1 \le i \le r$. Then, for all $A, B, C \in \mathcal{D}_0^i(\mathbb{R})$ and $M \in M_d(\mathbb{R})$,

$$[A, P_0^i] = 0 \quad \text{and} \quad [BMC, P_0^i] = 0.$$

*Proof.* This follows from elementary calculations. $\square$

**Lemma 14.** $M_d(\mathbb{R})$ is endowed with a norm $\|\cdot\|$ such that, for $A, B \in M_d(\mathbb{R})$, $\|AB\| \le \|A\| \cdot \|B\|$. Let $(u_n)_{n \ge 1}$ be a sequence valued in $M_d(\mathbb{R})$ and $m \ge 1$.

(i) If $u_n = o(1)$, then $\sum_{k=m}^{\infty} (u_n)^k = o(1)$ as $n \to \infty$.

(ii) Assume that $u_n = O(a_n)$, where $(a_n)_{n \ge 1}$ is a sequence valued in $(0, \infty)$ such that $a_n \to \infty$ as $n \to \infty$. Then, $\sum_{k=m}^{\infty} (u_n)^k = O((a_n)^m)$ as $n \to \infty$.

*Proof.* This follows readily from the properties of sums of geometric series. $\square$

### 6.4.5 Proof of Proposition 8

*Proof.* Set $\mathcal{P}_n^i := p_n^i(u_n)$. Then, $g_n^i(u_n) = \sqrt{n}\text{Log}_{P_0^i}^{G_n^i}(\mathcal{P}_n^i)$. By the explicit formula for the Riemannian Logarithm in Grassmannians given by (3.15),

$$\text{Log}_{P_0^i}^{G_n^i}(\mathcal{P}_n^i) = \frac{1}{2} [\log_m((I_d - 2\mathcal{P}_n^i)\mathcal{I}_i), P_0^i], \quad \text{where } \mathcal{I}_i := (I_d - 2P_0^i). \quad (6.21)$$

By Corollary 6, $\mathcal{P}_n^i \to P_0^i$ as $n \to \infty$. So, by (6.21), $\text{Log}_{P_0^i}^{G_n^i}(\mathcal{P}_n^i)$ converges to 0. Thus, we need to prove that the rate of this convergence is in $O(n^{-1/2})$. Since $\mathcal{I}_i^2 = I_d$, we have that

$$(I_d - 2\mathcal{P}_n^i) \mathcal{I}_i = I_d + \mathcal{K}_n^i, \quad \text{where } \mathcal{K}_n^i = o(1). \quad (6.22)$$

Then, we expand the $\log_m(\cdot)$ in (6.21). Thus, by (6.22), for all $n$ large enough,

$$\text{Log}_{P_0^i}^{G_n^i}(\mathcal{P}_n^i) = \frac{1}{2} [\log_m(I_d + \mathcal{K}_n^i), P_0^i] = \frac{1}{2} \sum_{k=1}^{\infty} \frac{(-1)^{k+1}}{k} [(\mathcal{K}_n^i)^k, P_0^i]. \quad (6.23)$$

In the rhs of (6.23), for $k \ge 1$, the rate of convergence to 0 of $(\mathcal{K}_n^i)^k = o(1)$ is not controlled. This explains why we needed a power series for the $\log_m(\cdot)$ and not only a
Taylor expansion. To prove Proposition 8, we need to compute the limit of $g_n^i(u_n)$. By (6.23),

$$g_n^i(u_n) = \sqrt{n}\mathrm{Log}_{P_0^i} (P_n^i) = \frac{1}{2} \sum_{k=1}^{\infty} \frac{(-1)^{k+1}}{k} \left[ (\mathcal{K}_n^i)^k, P_0^i \right] \quad (6.24)$$

where $P_n^i := p_n^i(u_n)$ and by (6.16), after calculations,

$$\mathcal{K}_n^i := (I_d - 2\mathcal{P}_n^i)(I_d - 2P_0^i) - I_d = \mathfrak{A}_n^i - 2n^{-1/2}\Gamma_n^i \mathcal{I}_i - 2n^{-1}\Phi_n^i \mathcal{I}_i, \quad (6.25)$$

with $\mathfrak{A}_n^i := -2\mathrm{Diag}(0_{q_1}, \dots, I_{q_i} - (\mathcal{A}_n^{(i,i)}) (\mathcal{A}_n^{(i,i)})'$, \dots, 0_{q_r}) = o(1)$. We split the sum in (6.24) by isolating the term for $k=1$. Thus, by Lemma 15 below,

$$\mathrm{Log}_{P_0^i}^{G_n^i}(\mathcal{P}_n^i) = \frac{1}{2}[\mathcal{K}_n^i, P_0^i] + \sum_{k=2}^{\infty} \frac{(-1)^{k+1}}{k} [(\mathcal{K}_n^i)^k, P_0^i] = \frac{1}{2}[\mathcal{K}_n^i, P_0^i] + o(n^{-1/2}). \quad (6.26)$$

Now, we deal with the first term in the rhs of (6.26). By Lemma 16 below,

$$\frac{\sqrt{n}}{2} [\mathcal{K}_n^i, P_0^i] = g_0^i(u) + O(n^{-1/2}) = g_0^i(u) + o(1). \quad (6.27)$$

$$\text{Finally, by (6.26) and (6.27), } g_n^i(u_n) = \sqrt{n}\mathrm{Log}_{P_0^i}^{G_n^i}(\mathcal{P}_n^i) = g_0^i(u) + o(1). \quad \square$$

<figure><img src="image_6.png" /></figure>

### 6.4.6 Auxiliary results for the proof of Proposition 8

**Lemma 15.** Recall that the sequence $(\mathcal{K}_n^i)_n$ is defined in (6.25). Then,

$$\sum_{k=2}^{\infty} \frac{(-1)^{k+1}}{k} [(\mathcal{K}_n^i)^k, P_0^i] = o(n^{-1/2}). \quad (6.28)$$

*Proof. First part:* By the expression of $\mathcal{K}_n^i$ in (6.25),

$$\mathcal{K}_n^i = \mathfrak{A}_n^i + \mathcal{M}_n^i, \quad \text{where } \mathfrak{A}_n^i \in \mathcal{D}_0^i(\mathcal{I}), \quad \text{with } \mathfrak{A}_n^i = o(1) \text{ and } \mathcal{M}_n^i = O(n^{-1/2}). \quad (6.29)$$

In the rhs of (6.24), $(\mathcal{K}_n^i)^k$ still contains some $(\mathfrak{A}_n^i)^\ell = o(1)$, whose rate of convergence to 0 is not controlled. Lemma 13 provides simplifications. Thus, by (6.29), for all $k \ge 2$,

$$(\mathcal{K}_n^i)^k = (\mathfrak{A}_n^i)^k + \left( \sum_{\ell=0}^{k-1} (\mathfrak{A}_n^i)^\ell (\mathcal{M}_n^i)^{k-\ell} \right) + \left( \sum_{\ell=1}^k (\mathcal{M}_n^i)^\ell (\mathfrak{A}_n^i)^{k-\ell} \right) + \mathcal{R}_n^i, \quad (6.30)$$

where $\mathcal{R}_n^i$ is a sum of terms of the form $(\mathfrak{A}_n^i)^\ell (\mathcal{M}_n^i)^{\ell'} (\mathfrak{A}_n^i)^{\ell''}$, with $\ell + \ell' + \ell'' = k$ and $\ell > 1$, $\ell'' > 1$. Since $\mathfrak{A}_n^i \in \mathcal{D}_0^i(\mathcal{I})$, Lemma 13 implies that

$$[(\mathfrak{A}_n^i)^k, P_0^i] = 0 \quad \text{and} \quad [\mathcal{R}_n^i, P_0^i] = 0. \quad (6.31)$$

By combining (6.30) and (6.31),

$$\sum_{k=2}^{\infty} \frac{(-1)^{k+1}}{k} [(\mathcal{K}_n^i)^k, P_0^i] = \frac{1}{2}[\alpha_n^i, P_0^i] + \frac{1}{2}[\beta_n^i, P_0^i], \quad (6.32)$$

where, setting $\epsilon_k := \frac{(-1)^{k+1}}{k}$,

$$\alpha_n^i = \sum_{k=2}^{\infty} \sum_{\ell=0}^{k-1} \epsilon_k (\mathfrak{A}_n^i)^\ell (\mathcal{M}_n^i)^{k-\ell} \quad \text{and} \quad \beta_n^i = \sum_{k=2}^{\infty} \sum_{\ell=1}^{k} \epsilon_k (\mathcal{M}_n^i)^\ell (\mathfrak{A}_n^i)^{k-\ell}.$$
**Second part:** We estimate $\alpha_n^i$ and $\beta_n^i$ as $n \to \infty$. In the double sum defining $\alpha_n^i$, we wish to factorize $(\mathfrak{A}_n^i)^\ell$ in sums indexed by $k$, by inverting the summation indices. Namely, we invert on the following triangular domain $\tau := \{(k, \ell) \in \mathbb{N}^2 : k \ge 2, 1 \le \ell \le k-1\}$. Thus, for fixed $k$, we split the sum in $\ell$ into two parts: the terms corresponding to $\ell = 0$ and those to $\ell \ge 1$:

$$ \alpha_n^i = \left( \sum_{k=2}^{\infty} \epsilon_k (\mathcal{M}_n^i)^k \right) + \left( \sum_{k=2}^{\infty} \sum_{\ell=1}^{k-1} \epsilon_k (\mathfrak{A}_n^i)^\ell (\mathcal{M}_n^i)^{k-\ell} \right). \quad (6.33) $$

Then, in the second part of the rhs in (6.33), we invert the indices, which lie in $\tau$.

$$ \sum_{k=2}^{\infty} \sum_{\ell=1}^{k-1} \epsilon_k (\mathfrak{A}_n^i)^\ell (\mathcal{M}_n^i)^{k-\ell} = \sum_{\ell=1}^{\infty} \sum_{k=\ell+1}^{\infty} \epsilon_k (\mathfrak{A}_n^i)^\ell (\mathcal{M}_n^i)^{k-\ell} = \sum_{\ell=1}^{\infty} (\mathfrak{A}_n^i)^\ell \left( \sum_{k'=1}^{\infty} \epsilon_{\ell+k'} (\mathcal{M}_n^i)^{k'} \right). \quad (6.34) $$

Since $\mathcal{M}_n^i = O(n^{-1/2})$ and $\mathfrak{A}_n^i = o(1)$, Lemma 14 implies that

$$ \sum_{k'=1}^{\infty} \epsilon_{\ell+k'} (\mathcal{M}_n^i)^{k'} = O(n^{-1/2}), \quad \text{uniformly in } \ell, \quad \text{and} \quad \sum_{\ell=1}^{\infty} (\mathfrak{A}_n^i)^\ell = o(1). \quad (6.35) $$

Therefore, by combining (6.34) and (6.35),

$$ \sum_{k=2}^{\infty} \sum_{\ell=1}^{k-1} \epsilon_k (\mathfrak{A}_n^i)^\ell (\mathcal{M}_n^i)^{k-\ell} = o(n^{-1/2}). \quad (6.36) $$

By Lemma 14, since $\mathcal{M}_n^i = O(n^{-1/2})$, we obtain for the first part of the rhs in (6.33):

$$ \sum_{k=2}^{\infty} \epsilon_k (\mathcal{M}_n^i)^k = O\left(\left(n^{-1/2}\right)^2\right) = O(n^{-1}). \quad (6.37) $$

By (6.33), (6.36) and (6.37), we deduce that $\alpha_n^i = o(n^{-1/2})$. Similarly, $\beta_n^i = o(n^{-1/2})$. $\square$

**Lemma 16.** We have that

$$ \frac{\sqrt{n}}{2} [\mathcal{K}_n^i, P_0^i] = g_0^i(u) + O(n^{-1/2}). \quad (6.38) $$

*Proof.* Since $[\mathfrak{A}_n^i, P_0^i] = 0$, the expression of $\mathcal{K}_n^i$ in (6.25) yields that

$$ \frac{1}{2} [\mathcal{K}_n^i, P_0^i] = \frac{1}{2} [\mathfrak{A}_n^i - 2n^{-1/2}\Gamma_n^i \mathcal{I}_i - 2n^{-1}\Phi_n^i \mathcal{I}_i, P_0^i] = -n^{-1/2} [\Gamma_n^i \mathcal{I}_i, P_0^i] + O(n^{-1}). \quad (6.39) $$

Thus, it is enough to prove that $-[\Gamma_n^i \mathcal{I}_i, P_0^i] \to g_0^i(u)$ as $n \to \infty$. Consider the matrices

$N^i := (\mathcal{A}_n^{(i,i)} (\mathcal{B}_n^{(1,i)})' \quad \dots \quad \mathcal{A}_n^{(i,i)} (\mathcal{B}_n^{(i-1,i)})')$ and $N^i := (\mathcal{A}_n^{(i,i)} (\mathcal{B}_n^{(i+1,i)})' \quad \dots \quad \mathcal{A}_n^{(i,i)} (\mathcal{B}_n^{(r,i)})')$ of respective sizes $q_i \times (q_1 + \dots + q_{i-1})$ and $q_i \times (q_{i+1} + \dots + q_r)$. After calculations:

$$ \begin{aligned} -[\Gamma_n^1 \mathcal{I}_1, P_0^1] &= \begin{pmatrix} 0_{q_1} & N_>^1 \\ (N_>^1)' & 0 \end{pmatrix}, & -[\Gamma_n^r \mathcal{I}_r, P_0^r] &= \begin{pmatrix} 0 & (N_<^r)' \\ N_<^r & 0_{q_r} \end{pmatrix} \\ -[\Gamma_n^i \mathcal{I}_i, P_0^i] &= \begin{pmatrix} 0 & (N_<^i)' & 0 \\ N_<^i & 0_{q_i} & N_>^i \\ 0 & (N_>^i)' & 0 \end{pmatrix}, & i \neq 1, r. \end{aligned} $$

By (6.14) in Lemma 12, for $j \neq i$,

$$ \mathcal{A}_n^{(i,i)} (\mathcal{B}_n^{(j,i)})' \underset{n \to \infty}{\longrightarrow} \mathcal{E}^{i,i}(u) (\mathcal{F}^{j,i}(u))' := \mathcal{G}^{i,j}(u) = \frac{1}{\lambda_i - \lambda_j} u^{(i,j)}. \quad (6.40) $$

By definition of the functions $g_0^i$, this concludes the proof. $\square$
## 6.5 Proof of Theorem 3

*Proof.* It is enough to prove that for any closed set $T \subset \mathbb{T}$,

$$ \overline{\lim} \Pr(V_n \in g_n^{-1}(T)) \le \Pr(V \in g_0^{-1}(T)). \quad (6.41) $$

Indeed, if (6.41) holds, then the Portmanteau Theorem allows to conclude the proof. First,

$$ \overline{\lim} \Pr(V_n \in g_n^{-1}(T)) \le \overline{\lim} \Pr(V_n \in g_n^{-1}(T) \cap \mathbb{S}_n^*) + \overline{\lim} \Pr(V_n \in g_n^{-1}(T) \cap (\mathbb{S}_n^*)^c). \quad (6.42) $$

By (2.2), the second term of the rhs in (6.42) vanishes. So,

$$ \overline{\lim} \Pr(V_n \in g_n^{-1}(T)) \le \overline{\lim} \Pr(V_n \in g_n^{-1}(T) \cap \mathbb{S}_n^*). \quad (6.43) $$

Consider the following sequence of sets $(R_m)_{m \ge 1}$ and the set $R^*$ defined by

$$ R_m := \bigcup_{n \ge m} \{g_n^{-1}(T) \cap \mathbb{S}_n^*\}, \quad m \ge 1 \quad \text{and} \quad R^* := \bigcap_{m \ge 1} \overline{R_m}. $$

Fix $m \ge 1$. Then, for all $n \ge m$, $\{g_n^{-1}(T) \cap \mathbb{S}_n^*\} \subset \overline{R_m}$. Since $V_n \xrightarrow{d} V$,

$$ \overline{\lim} \Pr(V_n \in g_n^{-1}(T) \cap \mathbb{S}_n^*) \le \overline{\lim} \Pr(V_n \in \overline{R_m}) \le \Pr(V \in \overline{R_m}). \quad (6.44) $$

Since the rhs of (6.44) is independent of $m$, taking the limit as $m \to \infty$ in (6.44) yields that

$$ \overline{\lim} \Pr(V_n \in g_n^{-1}(T) \cap \mathbb{S}_n^*) \le \Pr(V \in R^*). \quad (6.45) $$

Indeed, the sequence $(R_m)_{m \ge 1}$ is decreasing. Now, we claim that

$$ R^* \subset g_0^{-1}(T) \cup (\mathbb{S}_0)^c. \quad (6.46) $$

Indeed, let $v \in R^*$. Then, by definition of $R^*$, we can construct a sequence $(v_m)_{m \ge 1}$ valued in $(R_m)_{m \ge 1}$ which converges to $v$. Now, by definition of $R_m$, we may build a sub-sequence $(v_{\phi(m)})_{m \ge 1}$ of $(v_m)_m$ such that for all $m \ge 1$, $v_{\phi(m)} \in g_{\phi(m)}^{-1}(T) \cap \mathbb{S}_{\phi(m)}^*$. Thus, $g_{\phi(m)}(v_{\phi(m)}) \in T$ and, since $v_{\phi(m)} \in \mathbb{S}_{\phi(m)}^*$, assumption (2.3) implies that $g_{\phi(m)}(v_{\phi(m)}) \to g_0(v)$ as $m \to \infty$. Since $T$ is closed, we deduce that $v \in g_0^{-1}(T)$. This proves the claim (6.46). Thus, by (6.46) and the assumption $\Pr(V \in \mathbb{S}_0) = 1$,

$$ \Pr(V \in R^*) \le \Pr(V \in g_0^{-1}(T) \cup (\mathbb{S}_0)^c) \le \Pr(V \in g_0^{-1}(T)). \quad (6.47) $$

We combine (6.43), (6.45) and (6.47) to conclude that (6.41) holds. $\square$

## Acknowledgements

The authors have received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation program (grant agreement G-Statistics No 786854). It was also supported by the French government through the 3IA Côte d'Azur Investments ANR-19-P3IA-0002 managed by the National Research Agency.

## References

[1] Anderson, T.W. (1963). *Asymptotic theory for principal components analysis*. Ann. Math. Stat. (1963), Vol. 34, 122-148.

[2] Batzies, E., Huper, K., Machado, L. and Silva Leite, F. (2015). *Geometric mean and geodesic regression on Grassmannians*. Linea Algebra Appl., 466:83-101, 2015.
[3] Bhattacharya, R. and Patrangenaru, V. (2005). *Large sample theory of intrinsic and extrinsic sample means on manifolds-II*. The Annals of Statistics (2005), Vol. 33, No. 3, 1225-1259.

[4] Bendokat, T., Zimmermann, R., Absil, P.A. (2024). *A Grassmann manifold handbook: basic geometry and computational aspects*. Adv Comput Math 50, 6 (2024).

[5] Harremoes, P. (2009). *Maximum Entropy on Compact Groups*. Entropy 2009, 11, 222-237.

[6] Johnson, O. (2004). *Information Theory and Central Limit Theorem*. Imperial College Press: London, 2004.

[7] Johnson, O.T., Suhov, Y.M. (2000). *Entropy and convergence on compact groups*. J. Theoret. Probab. 2000, 13, 843-857.

[8] Kato, T. (1995). *Perturbation theory for linear operators*. Springer-Verlag Berlin Heidelberg, 1995.

[9] Mankovich, N., Camps-Valls, G., Birdal, T. (2024). *Fun with Flags: Robust Principal Directions via Flag Manifolds*. Proc. of the IEEE/CVF Conference on Computer Vision and Pattern Recognition (CVPR), 2024, pp. 330-340.

[10] Michor, P.W. (2008). *Topics in Differential Geometry*. Volume 93 of Graduate Studies in Mathematics. American Mathematical Soc., 2008.

[11] Pennec, X. (2018). *Barycentric subspace analysis on manifolds*. Ann. Statist. 46 (6A) 2711 - 2746, December 2018.

[12] Rabenoro, D., Pennec, X. (2024). *Effective formulas for the geometry of normal homogeneous spaces. Application to flag manifolds*. arXiv preprint arXiv: 2302.14810, 2024.

[13] Szwagier, T., Pennec, X. (2024). *The curse of isotropy: from principal components to principal subspaces*. arXiv preprint arXiv: 2307.15348, 2024.

[14] Tyler, D.E. (1981). *Asymptotic inference for eigenvectors*. The Annals of Statistics 9 (4), 725-736, 1981.

[15] Van der Vaart, A.W. (2000). *Asymptotic Statistics*. Cambridge Series in Statistical and Probabilistic Mathematics, Volume 3.

[16] Ye, K., Wong, K.S.W., Lim, L.H. (2022). *Optimization on flag manifolds*. Mathematical Programming (Series A), 2022, 194(1-2), pp. 621-660.