Intrinsic Riemannian Functional Data Analysis
Zhenhua Lin, Fang Yao
TL;DR
This paper introduces a new geometric framework for analyzing Riemannian functional data, enabling principal component analysis and linear regression directly on manifolds using intrinsic methods.
Contribution
It develops tensor Hilbert spaces along curves on manifolds and derives Karhunen-Loeve expansion for Riemannian processes, advancing intrinsic geometric analysis tools.
Findings
Framework applies to Euclidean and non-Euclidean manifolds.
Develops intrinsic PCA and regression methods with asymptotic analysis.
Demonstrates effectiveness through simulations and real data examples.
Abstract
In this work we develop a novel and foundational framework for analyzing general Riemannian functional data, in particular a new development of tensor Hilbert spaces along curves on a manifold. Such spaces enable us to derive Karhunen-Loeve expansion for Riemannian random processes. This framework also features an approach to compare objects from different tensor Hilbert spaces, which paves the way for asymptotic analysis in Riemannian functional data analysis. Built upon intrinsic geometric concepts such as vector field, Levi-Civita connection and parallel transport on Riemannian manifolds, the developed framework applies to not only Euclidean submanifolds but also manifolds without a natural ambient space. As applications of this framework, we develop intrinsic Riemannian functional principal component analysis (iRFPCA) and intrinsic Riemannian functional linear regression (iRFLR)…
Click any figure to enlarge with its caption.
Figure 1
Figure 2
Figure 3
Figure 4
Figure 5
Figure 6| 0.244 (0.056) | 0.135 (0.029) | 0.085 (0.019) | ||||
|---|---|---|---|---|---|---|
| iRMISE | aRMISE | iRMISE | aRMISE | iRMISE | aRMISE | |
| 0.279 (0.073) | 0.331 (0.078) | 0.147 (0.037) | 0.180 (0.042) | 0.086 (0.022) | 0.106 (0.027) | |
| 0.478 (0.133) | 0.514 (0128) | 0.264 (0.064) | 0.287 (0.061) | 0.147 (0.044) | 0.167 (0.042) | |
| Manifold | FPC | iRFPCA | RFPCA | ||||
|---|---|---|---|---|---|---|---|
| iRFLR | RFLR | |||||||
|---|---|---|---|---|---|---|---|---|
| Estimation | normal | |||||||
| SPD | ||||||||
| SPD | ||||||||
| Prediction | normal | |||||||
| SPD | ||||||||
| SPD | ||||||||
Peer Reviews
No public reviews on file for this paper yet. If you reviewed it on a platform where reviews are public (OpenReview, ICLR, NeurIPS, ICML), you can paste yours below so the community can read it here.
Videos
No videos yet. Explain this paper in a talk, walkthrough, or lecture? Add one.
Intrinsic Riemannian Functional Data Analysis
Zhenhua Linlabel=e1][email protected] [
Fang Yaolabel=e2][email protected] [ National University of Singapore\thanksmarkm1 and Peking University\thanksmarkm2
Z. Lin
Department of Statistics and Applied Probability
National University of Singapore
Singapore 117546
F. Yao
Corresponding author
Department of Probability and Statistics
Center for Statistical Science
Peking University, Beijing, China
(2019; 10 2017; 10 2018)
Abstract
In this work we develop a novel and foundational framework for analyzing general Riemannian functional data, in particular a new development of tensor Hilbert spaces along curves on a manifold. Such spaces enable us to derive Karhunen–Loève expansion for Riemannian random processes. This framework also features an approach to compare objects from different tensor Hilbert spaces, which paves the way for asymptotic analysis in Riemannian functional data analysis. Built upon intrinsic geometric concepts such as vector field, Levi–Civita connection and parallel transport on Riemannian manifolds, the developed framework applies to not only Euclidean submanifolds but also manifolds without a natural ambient space. As applications of this framework, we develop intrinsic Riemannian functional principal component analysis (iRFPCA) and intrinsic Riemannian functional linear regression (iRFLR) that are distinct from their traditional and ambient counterparts. We also provide estimation procedures for iRFPCA and iRFLR, and investigate their asymptotic properties within the intrinsic geometry. Numerical performance is illustrated by simulated and real examples.
62G05,
62J05,
Functional principal component,
functional linear regression,
intrinsic Riemannian Karhunen–Loève expansion,
parallel transport,
tensor Hilbert space.,
doi:
10.1214/18-AOS1787
keywords:
[class=AMS]
keywords:
††volume: 47††issue: 6
\arxiv
arxiv:1812.01831
MnLargeSymbols’164 MnLargeSymbols’171
and
t2Fang Yao’s research is partially supported by National Natural Science Foundation of China Grant 11871080, a Discipline Construction Fund at Peking University and Key Laboratory of Mathematical Economics and Quantitative Finance (Peking University), Ministry of Education. Data were provided by the Human Connectome Project, WU-Minn Consortium (Principal Investigators: David Van Essen and Kamil Ugurbil; 1U54MH091657) funded by the 16 NIH Institutes and Centers that support the NIH Blueprint for Neuroscience Research, and by the McDonnell Center for Systems Neuroscience at Washington University.
1 Introduction
Functional data analysis (FDA) advances substantially in the past two decades, as the rapid development of modern technology enables collecting more and more data continuously over time. There is rich literature spanning more than seventy years on this topic, including development on functional principal component analysis such as Rao (1958); Kleffe (1973); Dauxois, Pousse and Romain (1982); Silverman (1996); Yao, Müller and Wang (2005a); Hall and Hosseini-Nasab (2006); Zhang and Wang (2016), and advances on functional linear regression such as Yao, Müller and Wang (2005b); Hall and Horowitz (2007); Yuan and Cai (2010); Kong et al. (2016), among many others. For a thorough review of the topic, we refer readers to the review article Wang, Chiou and Müller (2016) and monographs Ramsay and Silverman (2005); Ferraty and Vieu (2006); Hsing and Eubank (2015); Kokoszka and Reimherr (2017) for comprehensive treatments on classic functional data analysis. Although traditionally functional data take values in a vector space, more data of nonlinear structure arise and should be properly handled in a nonlinear space. For instance, trajectories of bird migration are naturally regarded as curves on a sphere which is a nonlinear Riemannian manifold, rather than the three-dimensional vector space . Another example is the dynamics of brain functional connectivity. The functional connectivity at a time point is represented by a symmetric positive-definite matrix (SPD). Then the dynamics shall be modelled as a curve in the space of SPDs that is endowed with either the affine-invariant metric (Moakher, 2005) or the Log-Euclidean metric (Arsigny et al., 2007) to avoid the “swelling” effect (Arsigny et al., 2007). Both metrics turn SPD into a nonlinear Riemannian manifold. In this paper, we refer this type of functional data as Riemannian functional data, which are functions taking values on a Riemannian manifold and modelled by Riemannian random processes, that is, we treat Riemannian trajectories as realizations of a Riemannian random process.
Analysis of Riemannian functional data is not only challenged by the infinite dimensionality and compactness of covariance operator from functional data, but also obstructed by the nonlinearity of the range of functions, since manifolds are generally not vector spaces and render many techniques relying on linear structure ineffective or inapplicable. For instance, if the sample mean curve is computed for bird migration trajectories as if they were sampled from the ambient space , this naïve sample mean in general does not fall on the sphere of earth. For manifolds of tree-structured data studied in Wang and Marron (2007), as they are naturally not Euclidean submanifolds which refer to Riemannian submanifolds of a Euclidean space in this paper, the naïve sample mean can not even be defined from ambient spaces, and thus a proper treatment of manifold structure is necessary. While the literature for Euclidean functional data is abundant, works involving nonlinear manifold structure are scarce. Chen and Müller (2012) and Lin and Yao (2017) respectively investigate representation and regression for functional data living in a low-dimensional nonlinear manifold that is embedded in an infinite-dimensional space, while Lila and Aston (2017) focuses principal component analysis on functional data whose domain is a two-dimensional manifold. None of these deal with functional data that take values on a nonlinear manifold, while Dai and Müller (2017) is the only endeavor in this direction for Euclidean submanifolds.
As functional principal component analysis (FPCA) is an essential tool for FDA, it is of im-portance and interest to develop this notion for Riemannian functional data. Since manifolds are in general not vector spaces, classic covariance functions/operators do not exist naturally for a Riemannian random process. A strategy that is often adopted, for example, Shi et al. (2009) and Cornea et al. (2017), to overcome the lack of vectorial structure is to map data on the manifold into tangent spaces via Riemannian logarithm map defined in Section 2.2. As tangent spaces at different points are different vector spaces, in order to handle observations from different tangent spaces, some existing works assume a Euclidean ambient space for the manifold and identify tangent vectors as Euclidean vectors. This strategy is adopted by Dai and Müller (2017) on Riemannian functional data such as compositional data modelled on the unit sphere for the first time. Specifically, they assume that functional data are sampled from a time-varying geodesic submanifold, where at a given time point, the functions take values on a geodesic submanifold of a common manifold. Such a common manifold is further assumed to be a Euclidean submanifold that allows to identify all tangent spaces as hyperplanes in a common Euclidean space (endowed with the usual Euclidean inner product). Then, with the aid of Riemannian logarithm map, Dai and Müller (2017) are able to transform Riemannian functional data into Euclidean ones while accounting for the intrinsic curvature of the underlying manifold.
To avoid confusion, we distinguish two different perspectives to deal with Riemannian manifolds. One is to work with the manifold under consideration without assuming an ambient space surrounding it or an isometric embedding into a Euclidean space. This perspective is regarded as completely intrinsic, or simply intrinsic. Although generally difficult to work with, it can fully respect all geometric structure of the manifold. The other one, referred to as ambient here, assumes that the manifold under consideration is isometrically embedded in a Euclidean ambient space, so that geometric objects such as tangent vectors can be processed within the ambient space. For example, from this point of view, the local polynomial regression for SPD proposed by Yuan et al. (2012) is intrinsic, while the aforementioned work by Dai and Müller (2017) takes the ambient perspective.
Although it is possible to account for some of geometric structure in the ambient perspective, for example, the curved nature of manifold via Riemannian logarithm map, several issues arise due to manipulation of geometric objects such as tangent vectors in the ambient space. First, the essential dependence on an ambient space restricts potential applications. It is not immediately applicable to manifolds that are not a Euclidean submanifold or do not have a natural isometric embedding into a Euclidean space, for example, the Riemannian manifold of () SPD matrices endowed with the affine-invariant metric (Moakher, 2005) which is not compatible with the -dimensional Euclidean metric. Second, although an ambient space provides a common stage for tangent vectors at different points, operation on tangent vectors from this ambient perspective can potentially violate the intrinsic geometry of the manifold. To illustrate this, consider comparison of two tangent vectors at different points (this comparison is needed in the asymptotic analysis of Section 3.2; see also Section 2.4). From the ambient perspective, taking the difference of tangent vectors requires moving a tangent vector parallelly within the ambient space to the base point of the other tangent vector. However, the resultant tangent vector after movement in the ambient space is generally not a tangent vector for the base point of the other tangent vector; see the left panel of Figure 1 for a geometric illustration. In another word, the ambient difference of two tangent vectors at different points is not an intrinsic geometric object on the manifold, and the departure from intrinsic geometry can potentially affect the statistical efficacy and/or efficiency. Lastly, since manifolds might be embedded into more than one ambient space, the interpretation of statistical results crucially depends on the ambient space and could be misleading if one does not choose the ambient space appropriately.
In the paper, we develop a completely intrinsic framework that provides a foundational theory for general Riemannian functional data that paves the way for the development of intrinsic Riemannian functional principal component analysis and intrinsic Riemannian functional linear regression, among other potential applications. The key building block is a new concept of tensor Hilbert space along a curve on the manifold, which is described in Section 2. On one hand, our approach experiences dramatically elevated technical challenges relative to the ambient counterparts. For example, without an ambient space, it is nontrivial to perceive and handle tangent vectors. On the other hand, the advantages of the intrinsic perspective are at least three-fold, in contrast to ambient approaches. First, our results immediately apply to many important Riemannian manifolds that are not naturally a Euclidean submanifold but commonly seen in statistical analysis and machine learning, such as the aforementioned SPD manifolds and Grassmannian manifolds. Second, our framework features a novel intrinsic proposal for coherent comparison of objects from different tensor Hilbert spaces on the manifold, and hence makes the asymptotic analysis sensible. Third, results produced by our approach are invariant to embeddings and ambient spaces, and can be interpreted independently, which avoids potential misleading interpretation in practice.
As important applications of the proposed framework, we develop intrinsic Riemannian functional principal component analysis (iRFPCA) and intrinsic Riemannian functional linear regression (iRFLR). Specifically, estimation procedures of intrinsic eigenstructure are provided and their asymptotic properties are investigated within the intrinsic geometry. For iRFLR, we study a Riemannian functional linear regression model, where a scalar response intrinsically and linearly depends on a Riemannian functional predictor through a Riemannian slope function, a concept that is formulated in Section 4, along with the concept of linearity in the context of Riemannian functional data. We present an FPCA-based estimator and a Tikhonov estimator for the Riemannian slope function and explore their asymptotic properties, where the proposed framework of tensor Hilbert space again plays an essential role.
The rest of the paper is structured as follows. The foundational framework for intrinsic Riemannian functional data analysis is laid in Section 2. Intrinsic Riemannian functional principal component analysis is presented in Section 3, while intrinsic Riemannian functional regression is studied in Section 4. In Section 5, numerical performance is illustrated through simulations, and an application to Human Connectome Project analyzing functional connectivity and behavioral data is provided.
2 Tensor Hilbert space and Riemannian random process
In this section, we first define the concept of tensor Hilbert space and discuss its properties, including a mechanism to deal with two different tensor Hilbert spaces at the same time. Then, random elements on tensor Hilbert space are investigated, with the proposed intrinsic Karhunen–Loève expansion for the random elements. Finally, practical computation with respect to an orthonormal frame is given. Throughout this section, we assume a -dimensional, connected and geodesically complete Riemannian manifold equipped with a Riemannian metric , which defines a scalar product for the tangent space at each point . This metric also induces a distance function on . A preliminary for Riemannian manifolds can be found in the appendix. For a comprehensive treatment on Riemannian manifolds, we recommend the introductory text by Lee (1997) and also Lang (1995).
2.1 Tensor Hilbert spaces along curves
Let be a measurable curve on a manifold and parameterized by a compact domain equipped with a finite measure . A vector field along is a map from to the tangent bundle such that for all . It is seen that the collection of vector fields along is a vector space, where the vector addition between two vector fields and is a vector field such that for all , and the scalar multiplication between a real number and a vector field is a vector field such that for all . Let be the collection of (equivalence classes of) measurable vector fields along such that with identification between and in (or equivalently, and are in the same equivalence class) when . Then is turned into an inner product space by the inner product , with the induced norm given by . Moreover, we have that
Theorem 1**.**
For a measurable curve on , is a separable Hilbert space.
We call the space the tensor Hilbert space along , as tangent vectors are a special type of tensor and the above Hilbertian structure can be defined for tensor fields along . The above theorem is of paramount importance, in the sense that it suggests to serve as a cornerstone for Riemannian functional data analysis for two reasons. First, as shown in Section 2.2, via Riemannian logarithm maps, a Riemannian random process may be transformed into a tangent-vector-valued random process (called log-process in Section 2.2) that can be regarded as a random element in a tensor Hilbert space. Second, the rigorous theory of functional data analysis formulated in Hsing and Eubank (2015) by random elements in separable Hilbert spaces fully applies to the log-process.
One distinct feature of the tensor Hilbert space is that, different curves that are even parameterized by the same domain give rise to different tensor Hilbert spaces. In practice, one often needs to deal with two different tensor Hilbert spaces at the same time. For example, in the next subsection we will see that under some conditions, a Riemannian random process can be conceived as a random element on the tensor Hilbert space along the intrinsic mean curve . However, the mean curve is often unknown and estimated from a random sample of . Since the sample mean curve generally does not agree with the population one, quantities of interest such as covariance operator and their sample versions are defined on two different tensor Hilbert spaces and , respectively. For statistical analysis, one needs to compare the sample quantities with their population counterparts and hence involves objects such as covariance operators from two different tensor Hilbert spaces.
In order to intrinsically quantify the discrepancy between objects of the same kind from different tensor Hilbert spaces, we utilize the Levi–Civita connection (p. 18, Lee 1997) associated with the Riemannian manifold . The Levi–Civita connection is uniquely determined by the Riemannian structure. It is the only torsion-free connection compatible with the Riemannian metric. Associated with this connection is a unique parallel transport operator that smoothly transports tangent vectors at along a curve to and preserves the inner product. We shall emphasize that the parallel transportation is performed intrinsically. For instance, tangent vectors being transported always stay tangent to the manifold during transportation, which is illustrated by the right panel of Figure 1. Although transport operator depends on the curve connecting and , there exists a canonical choice of the curve connecting two points, which is the minimizing geodesic between and (under some conditions, almost surely the minimizing geodesic is unique between two points randomly sampled from the manifold). The smoothness of parallel transport also implies that if and are not far apart, then the initial tangent vector and the transported one stays close (in the space of tangent bundle endowed with the Sasaki metric (Sasaki, 1958)). This feature is desirable for our purpose, as when sample mean approaches to , one expects a tangent vector at converges to its transported version at . Owing to these nice properties of parallel transport, it becomes an ideal tool to construct a mechanism of comparing objects from different tensor Hilbert spaces as follows.
Suppose and are two measurable curves on defined on . Let be a family of smooth curves that is parameterized by the interval (the way of parameterization here does not matter) and connects to , that is, and , such that is measurable for all . Suppose and let be a smooth vector field along such that and , where denotes the Levi–Civita connection of the manifold . The theory of Riemannian manifolds shows that such a vector field uniquely exists. This gives rise to the parallel transporter along , defined by . In other words, parallelly transports to along the curve . As the parallel transporter determined by the Levi–Civita connection, preserves the inner product of tangent vectors along transportation, that is, for . Then we can define the parallel transport of vector fields from to , denoted by , for all and . One immediately sees that is a linear operator on tensor Hilbert space. Its adjoint, denoted by , is a map from to and is given by for and . Let denote the set of all Hilbert–Schmidt operators on , which is a Hilbert space with the Hilbert–Schmidt norm . We observe that the operator also gives rise to a mapping from to , defined by for and . The operator is called the parallel transporter of Hilbert–Schmidt operators. Below are some important properties of and , where for .
Proposition 2**.**
Suppose , , and for two measurable curves and on . Let and be parallel transporters along a family of smooth curves defined above, such that is smooth and is measurable. Then the following statements regarding and hold.
The operator is a unitary transformation from to . 2. 2.
. Also, . 3. 3.
. 4. 4.
If is invertible, then . 5. 5.
, where are scalar constants, and . 6. 6.
.
We define for and , and for operators and . To quantify the discrepancy between an element in and another one in , we can use the quantity . Similarly, we adopt as discrepancy measure for two covariance operators and . These quantities are intrinsic as they are built on intrinsic geometric concepts. In light of Proposition 2, they are symmetric under the parallel transport, that is, transporting to yields the same discrepancy measure as transporting to . We also note that, when , the difference operators and reduce to the regular vector and operator difference, that is, becomes , while becomes . Therefore, and can be viewed as generalization of the regular vector and operator difference to a Riemannian setting. One shall note that and depend on the choice of the family of curves , a canonical choice of which is discussed in Section 3.2.
2.2 Random elements on tensor Hilbert spaces
Let be a Riemannian random process. In order to introduce the concept of intrinsic mean function for , we define a family of functions indexed by :
[TABLE]
For a fixed , if there exists a unique that minimizes over all , then is called the intrinsic mean (also called Fréchet mean) at , denoted by , that is,
[TABLE]
As required for intrinsic analysis, we assume the following condition.
A.1
The intrinsic mean function exists.
We refer readers to Bhattacharya and Patrangenaru (2003) and Afsari (2011) for conditions under which the intrinsic mean of a random variable on a general manifold exists and is unique. For example, according to Cartan–Hadamard theorem, if the manifold is simply connected and complete with nonpositive sectional curvature, then intrinsic mean function always exists and is unique as long as for all , for some .
Since is geodesically complete, by Hopf–Rinow theorem (p. 108, Lee 1997), its exponential map at each is defined on the entire . As might not be injective, in order to define its inverse, we restrict to a subset of the tangent space . Let denote the set of all tangent vectors such that the geodesic fails to be minimizing for for each . Now, we define only on . The image of , denoted by , consists of points in , such that for some . In this case, the inverse of exists and is called Riemannian logarithm map, which is denoted by and maps to . We shall make the following assumption.
A.2
.
Then, is almost surely defined for all . The condition is superfluous if is injective for all , like the manifold of SPDs endowed with the affine-invariant metric.
In the sequel we shall assume satisfies conditions A.1 and A.2. An important observation is that the log-process (denoted by for short) is a random vector field along . If we assume continuity for the sample paths of , then the process is measurable with respect to the product -algebra and the Borel algebra , where is the -algebra of the probability space. Furthermore, if , then according to Theorem 7.4.2 of Hsing and Eubank (2015), can be viewed as a tensor Hilbert space valued random element. Observing that according to Theorem 2.1 of Bhattacharya and Patrangenaru (2003), the intrinsic covariance operator for is given by . This operator is nuclear and self-adjoint. It then admits the following eigendecomposition (Theorem 7.2.6, Hsing and Eubank, 2015)
[TABLE]
with eigenvalues and orthonormal eigenelements that form a complete orthonormal system for . Also, with probability one, the log-process of has the following Karhunen–Loève expansion
[TABLE]
with being uncorrelated and centered random variables. Therefore, we obtain the intrinsic Riemannian Karhunen–Loève (iRKL) expansion for given by
[TABLE]
The elements are called intrinsic Riemannian functional principal component (iRFPC), while the variables are called intrinsic iRFPC scores. This result is summarized in the following theorem whose proof is already contained in the above derivation and hence omitted. We shall note that the continuity assumption on sample paths can be weakened to piece-wise continuity.
Theorem 3** (Intrinsic Karhunen–Loève Representation).**
Assume that satisfies assumptions A.1 and A.2. If sample paths of are continuous and , then the intrinsic covariance operator of admits the decomposition (2), and the random process admits the representation (4).
In practice, the series at (4) is truncated at some positive integer , resulting in a truncated intrinsic Riemannian Karhunen–Loève expansion of , given by with . The quality of the approximation of for is quantified by , and can be shown by a method similar to Dai and Müller (2017) that if the manifold has non-negative sectional curvature everywhere, then . For manifolds with negative sectional curvatures, such inequality in general does not hold. However, for Riemannian random process that almost surely lies in a compact subset of , the residual can be still bounded by up to a scaling constant.
Proposition 4**.**
Assume that conditions A.1 and A.2 hold, and the sectional curvature of is bounded from below by . Let be a subset of . If , we let , and if , we assume that is compact. Then, for some constant , for all . Consequently, if almost surely, then .
2.3 Computation in orthonormal frames
In practical computation, one might want to work with specific orthonormal bases for tangent spaces. A choice of orthonormal basis for each tangent space constitutes an orthonormal frame on the manifold. In this section, we study the representation of the intrinsic Riemannian Karhunen–Loève expansion under an orthonormal frame and formulas for change of orthonormal frames.
Let be a continuous orthonormal frame, that is, each is a vector field of such that and for and all . At each point , form an orthonormal basis for . The coordinate of with respect to is denoted by , with the subscript indicating its dependence on the frame. The resulting process is called the -coordinate process of . Note that is a regular valued random process defined on , and classic theory in Hsing and Eubank (2015) applies to . For example, its norm is defined by , where denotes the canonical norm on . One can show that . Therefore, if , then the covariance function exists and is matrix-valued, quantified by (Kelly and Root, 1960; Balakrishnan, 1960), noting that as for all . Also, the vector-valued Mercer’s theorem implies the eigendecomposition
[TABLE]
with eigenvalues and corresponding eigenfunctions . Here, the subscript in is to emphasize the dependence on the chosen frame. One can see that is a coordinate representation of , that is, .
The coordinate process admits the vector-valued Karhunen–Loève expansion
[TABLE]
under the assumption of mean square continuity of , according to Theorem 7.3.5 of Hsing and Eubank (2015), where . While the covariance function and eigenfunctions of depend on frames, and in (4) and (6) are not, which justifies the absence of in their subscripts and the use of the same notation for eigenvalues and iRFPC scores in (2), (4), (5) and (6). This follows from the formulas for change of frames that we shall develop below.
Suppose is another orthonormal frame. Change from to can be characterized by a unitary matrix . For example, and hence for all . Then the covariance function of is given by
[TABLE]
and consequently,
[TABLE]
From the above calculation, we immediately see that are also eigenvalues of . Moreover, the eigenfunction associated with for is given by
[TABLE]
Also, the variable in (4) and (6) is the functional principal component score for associated with , as seen by . The following proposition summarizes the above results.
Proposition 5** (Invariance Principle).**
Let be a -valued random process satisfying conditions A.1 and A.2. Suppose and are measurable orthonormal frames that are continuous on a neighborhood of the image of , and denotes the -coordinate log-process of . Assume is the unitary matrix continuously varying with such that for .
The -norm of for , defined by , is independent of the choice of frames. In particular, for all orthonormal frames . 2. 2.
If , then the covariance function of exists for all and admits decomposition of (5). Also, (2) and (5) are related by for all , and the eigenvalues coincide. Furthermore, the eigenvalues of and the principal component scores of Karhunen–Loève expansion of do not depend on . 3. 3.
The covariance functions and of respectively and , if they exist, are related by (7). Furthermore, their eigendecomposions are related by (8) and for all . 4. 4.
If and sample paths of are continuous, then the scores (6) coincide with the iRFPC scores in (4).
We conclude this subsection by emphasizing that the concept of covariance function of the log-process depends on the frame , while the covariance operator, eigenvalues, eigenelements and iRFPC scores do not. In particular, the scores , which are often the input for further statistical analysis such as regression and classification, are invariant to the choice of coordinate frames. An important consequence of the invariance principle is that, these scores can be safely computed in any convenient coordinate frame without altering the subsequent analysis.
2.4 Connection to the special case of Euclidean submanifolds
Our framework applies to general manifolds that include Euclidean submanifolds as special examples to which the methodology of Dai and Müller (2017) also applies. When the underlying manifold is a -dimensional submanifold of the Euclidean space with , we recall that the tangent space at each point is identified as a -dimensional linear subspace of . For such Euclidean manifolds, Dai and Müller (2017) treat the log-process of as a -valued random process, and derive the representation for the log-process (equation (5) in their paper) within the ambient Euclidean space. This is distinctly different from our intrinsic representation (3) based on the theory of tensor Hilbert spaces, despite their similar appearance. For instance, equation. (5) in their work can only be defined for Euclidean submanifolds, while ours is applicable to general Riemannian manifolds. Similarly, the covariance function defined in Dai and Müller (2017), denoted by , is associated with the ambient log-process , that is, . Such an ambient covariance function can only be defined for Euclidean submanifolds but not general manifolds.
Nevertheless, there are connections between the ambient method of Dai and Müller (2017) and our framework when is a Euclidean submanifold. For instance, the mean curve is intrinsically defined in the same way in both works. For the covariance structure, although our covariance function is a matrix-valued function while is a matrix-valued function, they both represent the intrinsic covariance operator when is a Euclidean submanifold. To see so, first, we observe that the ambient log-process as defined in Dai and Müller (2017) at the time , although is ambiently -dimensional, lives in a -dimensional linear subspace of . Second, the orthonormal basis for the tangent space at can be realized by a full-rank matrix by concatenating vectors . Then is the -coordinate process of . This implies that . On the other hand, since , one has . Thus, and determine each other and represent the same object. In light of this observation and the invariance principle stated in Proposition 5, when is a Euclidean submanifold, can be viewed as the ambient representation of the intrinsic covariance operator , while is the coordinate representation of with respect to the frame . Similarly, the eigenfunctions of are the ambient representation of the eigenelements of . The above reasoning also applies to sample mean functions and sample covariance structure. Specifically, when is a Euclidean submanifold, our estimator for the mean function discussed in Section 3 is identical to the one in Dai and Müller (2017), while the estimators for the covariance function and eigenfunctions proposed in Dai and Müller (2017) are the ambient representation of our estimators stated in Section 3.
However, when quantifying the discrepancy between the population covariance structure and its estimator, Dai and Müller (2017) adopt the Euclidean difference as a measure. For instance, they use to represent the discrepancy between the sample eigenfunctions and the population eigenfunctions, where is the sample version of . When , the sample version of , is not equal to , and belong to different tangent spaces. In such case, the Euclidean difference is a Euclidean vector that does not belong to the tangent space at either or , as illustrated in the left panel of Figure 1. In other words, the Euclidean difference of ambient eigenfunctions does not obey the geometry of the manifold, hence might not properly measure the intrinsic discrepancy. In particular, the measure might completely result from the departure of the ambient Euclidean geometry from the manifold, rather than the intrinsic discrepancy between the sample and population eigenfunctions, as demonstrated in the left panel of Figure 1. Similar reasoning applies to . In contrast, we base on Proposition 2 to propose an intrinsic measure to characterize the intrinsic discrepancy between a population quantity and its estimator in Section 3.2.
3 Intrinsic Riemannian functional principal component analysis
3.1 Model and estimation
Suppose admits the intrinsic Riemannian Karhunen–Loève expansion (4), and are a random sample of . In the sequel, we assume that trajectories are fully observed. In the case that data are densely observed, each trajectory can be individually interpolated by using regression techniques for manifold valued data, such as Steinke, Hein and Schölkopf (2010), Cornea et al. (2017) and Petersen and Müller (2017). This way the densely observed data could be represented by their interpolated surrogates, and thus treated as if they were fully observed curves. When data are sparse, delicate information pooling of observations across different subjects is required. The development of such methods is substantial and beyond the scope of this paper.
In order to estimate the mean function , we define the finite-sample version of in (1) by
[TABLE]
Then, an estimator for is given by
[TABLE]
The computation of depends on the Riemannian structure of the manifold. Readers are referred to Cheng et al. (2016) and Salehian et al. (2015) for practical algorithms. For a subset of , denotes the set , where is the ball with center and radius in . We use to denote the set . In order to define , at least with a dominant probability for a large sample, we shall assume a slightly stronger condition than A.2:
A.2′**
There is some constant such that .
Then, combining the fact that we will show later, we conclude that for a large sample, almost surely, for all . Therefore, under this condition, is well-defined almost surely for a large sample.
The intrinsic Riemannian covariance operator is estimated by its finite-sample version
[TABLE]
This sample intrinsic Riemannian covariance operator also admits an intrinsic eigendecomposion for . Therefore, the estimates for the eigenvalues are given by , while the estimates for are given by . These estimates can also be conveniently obtained under a frame, due to the invariance principle stated in Proposition 5. Let be a chosen orthonormal frame and be the sample covariance function based on , where is the coordinate process of under the frame with respect to . We can then obtain the eigendecomposition , which yields for . Finally, the truncated process for is estimated by
[TABLE]
where are estimated iRFPC scores. The above truncated iRKL expansion can be regarded as generalization of the representation (10) in Dai and Müller (2017) from Euclidean submanifolds to general Riemannian manifolds.
3.2 Asymptotic properties
To quantify the difference between and , it is natural to use the square geodesic distance as a measure of discrepancy. For the asymptotic properties of , we need the following regularity conditions.
B.1
The manifold is connected and complete. In addition, the exponential map is surjective at every point .
B.2
The sample paths of are continuous.
B.3
is finite. Also, for all compact subsets , .
B.4
The image of the mean function is bounded, that is, the diameter is finite, .
B.5
For all , .
To state the next condition, let . The calculus of manifolds suggests that , where denotes the gradient operator at . For each , let denote the Hessian of the real function , that is, for vector fields and on ,
[TABLE]
B.6
, where denotes the smallest eigenvalue of an operator or matrix.
B.7
and , where for a real function on .
The assumption B.1 regarding the property of manifolds is met in general, for example, the -dimensional unit sphere , SPD manifolds, etc. By the Hopf–Rinow theorem, the condition also implies that is geodesically complete. Conditions similar to B.2, B.5, B.6 and B.7 are made in Dai and Müller (2017). The condition B.4 is a weak requirement for the mean function and is automatically satisfied if the manifold is compact, while B.3 is analogous to standard moment conditions in the literature of Euclidean functional data analysis and becomes superfluous when is compact. If , then B.3 is equivalent to , a condition commonly made in the literature of functional data analysis, for example, Hall and Hosseini-Nasab (2006). The following theorem summarizes asymptotic properties of for general Riemannian manifolds. Its proof can be obtained by mimicking the techniques in Dai and Müller (2017), with additional care of the case that is noncompact. For completeness, we provide its proof in the Appendix. As noted by Dai and Müller (2017), the root- rate can not be improved in general.
Theorem 6**.**
Under conditions A.1, A.2′** and B.1-B.4, the following holds.
* is uniformly continuous, and is uniformly continuous with probability one.* 2. 2.
* is a uniformly strong consistent estimator for , that is, .* 3. 3.
If additional conditions B.5-B.7 are assumed, then converges in distribution to a Gaussian measure on the tensor Hilbert space . 4. 4.
If additional conditions B.5-B.7 are assumed, then and .
For asymptotic analysis of the estimated eigenstructure, as discussed in Section 2.1 and 2.4, and are not defined for intrinsic Riemannian functional data analysis, since they are objects originated from different tensor Hilbert spaces and induced by different curves and . Therefore, we shall adopt the intrinsic measure of discrepancy developed in Section 2.1, namely, and , where the dependence of and on and is suppressed for simplicity. As mentioned at the end of Section 2.1, both and depend on the choice of the family of curves . Here, a canonical choice for each is the minimizing geodesic between and . The existence and uniqueness of such geodesics can be deduced from assumptions A.1 and A.2′**. Also, the continuity of and implies the continuity of and hence the measurability of for each . By Proposition 2, one sees that , recalling that is a vector field along . It can also be seen that are eigenpairs of . These identities match our intuition that the transported sample covariance operator ought to be an operator derived from transported sample vector fields, and that the eigenfunctions of the transported operator are identical to the transported eigenfunctions.
To state the asymptotic properties for the eigenstructure, we define
[TABLE]
[TABLE]
Theorem 7**.**
Assume that every eigenvalue has multiplicity one, and conditions A.1, A.2′** and B.1-B.7 hold. Suppose tangent vectors are parallel transported along minimizing geodesics for defining the parallel transporters and . If , then . Furthermore, and for all ,
[TABLE]
If is replaced by , then (10) holds with probability 1.
In this theorem, (10) generalizes Lemma 4.3 of Bosq (2000) to the Riemannian setting. Note that the intrinsic rate for is optimal. Also, from (10) one can deduce the optimal rate for a fixed . We stress that these results apply to not only Euclidean submanifolds, but also general Riemannian manifolds.
4 Intrinsic Riemannian functional linear regression
4.1 Regression model and estimation
Classical functional linear regression for Euclidean functional data is well studied in the literature, that is, the model relating a scalar response and a functional predictor by , where is the intercept, is the slope function and represents measurement errors, for example, Cardot, Ferraty and Sarda (2003), Cardot, Mas and Sarda (2007), Hall and Horowitz (2007) and Yuan and Cai (2010), among others. However, for Riemannian functional data, both and take values in a manifold and hence the product is not well defined. Rewriting the model as , where is the canonical inner product of the square integrable functions, we propose to replace by the inner product on the tensor Hilbert space , and define the following Riemannian functional linear regression model:
[TABLE]
where we require conditions A.1 and A.2. Note that is a manifold valued function defined on , namely the Riemannian slope function of the model (11), and this model is linear in terms of . We stress that the model (11) is intrinsic to the Riemannian structures of the manifold.
According to Theorem 2.1 of Bhattacharya and Patrangenaru (2003), the process is centered at its mean function, that is, , which also implies that . Since the focus is the Riemannian slope function , in the sequel, without loss of generality, we assume that is centered and hence . In practice, we may use sample mean of for centering. With the same reasoning as Cardot, Ferraty and Sarda (2003), one can see that and hence are identifiable if A.1, A.2 and the following condition hold.
C.1
The random pair satisfies , where are eigenpairs of the covariance operator of .
To estimate the Riemannian slope function , by analogy to the FPCA approach in traditional functional linear regression (Hall and Horowitz, 2007), we represent by . Then with and . The empirical iRFPCA estimate of is then given by
[TABLE]
where , , and . We can also obtain a Tikhonov estimator as follows. For , define , where is the identity operator on . The Tikhonov estimator is given by
[TABLE]
In practice, it is convenient to conduct computation with respect to an orthonormal frame . For each , the -coordinate along the curve is denoted by . Note that by Theorem 2.1 of Bhattacharya and Patrangenaru (2003), implies that . Thus, the empirical covariance function is computed by . Then, and , recalling that are eigenfunctions of . Finally, the -coordinate of is given by . For the Tikhonov estimator, we first observe that shares the same eigenfunctions with , while the eigenvalues of is related to by . Therefore, and its -coordinate is given by . We emphasize that both estimators and are intrinsic and hence independent of frames, while their coordinate representations depend on the choice of frames. In addition, as a special case, when is a Euclidean submanifold, an argument similar to that in Section 2.4 can show that, if one treats as an ambient random process and adopts the FPCA and Tikhonov regularization approaches (Hall and Horowitz, 2007) to estimate the slope function in the ambient space, then the estimates are the ambient representation of our estimates and in (12) and (13), respectively.
4.2 Asymptotic properties
In order to derive convergence of the iRFPCA estimator and the Tikhonov estimator, we shall assume the sectional curvature of the manifold is bounded from below by to exclude pathological cases. The compact support condition on in the case might be relaxed to weaker assumptions on the tail decay of the distribution of . Such weaker conditions do not provide more insight for our derivation, but complicate the proofs significantly, which is not pursued further.
C.2
If , is assumed to lie in a compact subset almost surely. Moreover, errors are identically distributed with zero mean and variance not exceeding a constant .
The follow conditions concern the spacing and the decay rate of eigenvalues of the covariance operator, as well as the strength of the signal . They are standard in the literature of functional linear regression, for example, Hall and Horowitz (2007).
C.3
For , .
C.4
, and .
Let be the collection of distributions of satisfying conditions C.2-C.4. The following theorem establish the convergence rate of the iRFPCA estimator for the class of models in .
Theorem 8**.**
Assume that conditions A.1, A.2′**, B.1-B.7 and C.1-C.4 hold. If , then
[TABLE]
For the Tikhonov estimator , we have a similar result. Instead of conditions C.3-C.4, we make the following assumptions, which again are standard in the functional data literature.
C.5
.
C.6
, and .
Let be the class of distributions of that satisfy C.2 and C.5-C.6. The following result shows that the convergence rate of is at least in terms of integrated square geodesic distance.
Theorem 9**.**
Assume that conditions A.1, A.2′**, B.1-B.7, C.1-C.2 and C.5-C.6 hold. If , then
[TABLE]
It is important to point out that the theory in Hall and Horowitz (2007) is formulated for Euclidean functional data and hence does not apply to Riemannian functional data. In particular, their proof machinery depends on the linear structure of the sample mean function for Euclidean functional data. However, the intrinsic empirical mean generally does not admit an analytic expression, which hinges derivation of the optimal convergence rate. We leave the refinement on minimax rates of iRFPCA and Tikhonov estimators to future research. Note that model (11) can be extended to include a finite and fixed number of scalar predictors with slight modification, and the asymptotic properties of and remain unchanged.
5 Numerical examples
5.1 Simulation studies
We consider two manifolds that are frequently encountered in practice111Implementation of iRFPCA on these manifolds can be found on Github: https://github.com/linulysses/iRFDA.. The first one is the unit sphere which is a compact nonlinear Riemannian submanifold of for a positive integer . The sphere can be used to model compositional data, as exhibited in Dai and Müller (2017) which also provides details of the geometry of . Here we consider the case of . The sphere consists of points satisfying . Since the intrinsic Riemannian geometry of is the same as the one inherited from its ambient space (referred to as ambient geometry hereafter), according to the discussion in Section 2.4, the ambient approach to FPCA and functional linear regression yields the same results as our intrinsic approach.
The other manifold considered is the space of symmetric positive definite matrices, denoted by . The space includes nonsingular covariance matrices which naturally arise from the study of DTI data (Dryden, Koloydenko and Zhou, 2009) and functional connectivity (Friston, 2011). Transitionally, it is treated as a convex cone of the linear space of symmetric matrices. However, as discussed in Arsigny et al. (2007), this geometry suffers from the undesired swelling effect. To rectify the issue, several intrinsic geometric structures have been proposed, including the affine-invariant metric (Moakher, 2005) which is also known as the trace metric (Chapter XII, Lang 1999). This metric is invariant to affine transformation (i.e., change of coordinates) on and thus suitable for covariance matrices; see Arsigny et al. (2007) as well as Fletcher and Joshib (2007) for more details. Moreover, the affine-invariant metric has a negative sectional curvature, and thus the Fréchet mean is unique if it exists. In our simulation, we consider . We emphasize that the affine-invariant geometry of is different from the geometry inherited from the linear space . Thus, the ambient RFPCA of Dai and Müller (2017) might yield inferior performance on this manifold.
We simulate data as follows. First, the time domain is set to . The mean curves for and are, respectively, with and , and , that is a matrix. The Riemannian random processes are produced in accordance to , where for and for . We set iRFPCs , where is an orthonormal frame over the path , with being orthonormal Fourier basis functions on , and is an orthonormal matrix that is randomly generated but fixed throughout all simulation replicates. We take for all manifolds. Each curve is observed at regular design points . The slope function is with . Two different types of distribution for in (11) are considered, namely, normal and Student’s distribution with degree of freedom . Note that the latter is a heavy-tailed distribution, with a smaller suggesting a heavier tail and ensuring the existence of variance. In addition, the noise is scaled to make the signal-to-noise ratio equal to . Three different training sample sizes are considered, namely, , and , while the sample size for test data is . Each simulation setup is repeated independently times.
First, we illustrate the difference between the intrinsic measure and the ambient counterpart for the discrepancy of two random objects residing on different tangent spaces, through the examples of the sphere manifold and the first two iRFPCs. Recall that the metric of agrees with its ambient Euclidean geometry, so that both iRFPCA and RFPCA essentially yield the same estimates for iRFPCs. We propose to use the intrinsic root mean integrated squared error (iRMISE) to characterize the difference between and its estimator , while Dai and Müller (2017) adopt the ambient RMISE (aRMISE) , as discussed in Section 2.4. The numerical results of iRMISE and aRMISE for and , as well as the RMISE for , are showed in Table 1. We see that, when is small and hence is not sufficiently close to , the difference between iRMISE and aRMISE is visible, while such difference decreases as the sample size grows and converges to . In particular, aRMISE is always larger than iRMISE since aRMISE contains an additional ambient component that is not intrinsic to the manifold, as illustrated on the left panel of Figure 1.
We now use iRMISE to assess the performance of iRFPCA by comparing to the ambient counterpart RFPCA proposed by Dai and Müller (2017). Table 2 presents the results for the top eigenelements. The first observation is that iRFPCA and RFPCA yield the same results on the manifold , which numerically verifies our discussion in Section 2.4. We notice that in Dai and Müller (2017) the quality of estimation of principal components is not evaluated, likely due to the lack of a proper tool to do so. In contrast, our framework of tensor Hilbert space provides an intrinsic gauge (e.g., iRMISE) to naturally compare two vector fields along different curves. For the case of which is not a Euclidean submanifold, the iRFPCA produces more accurate estimation than RFPCA. In particular, as sample size grows, the estimation error for iRFPCA decreases quickly, while the error of RFPCA persists. This coincides with our intuition that when the geometry induced from the ambient space is not the same as the intrinsic geometry, the ambient RFPCA incurs loss of statistical efficiency, or even worse, inconsistent estimation. In summary, the results numerically demonstrate that the RFPCA proposed by Dai and Müller (2017) does not apply to manifolds that do not have an ambient space or whose intrinsic geometry differs from its ambient geometry, while our iRFPCA are applicable to such Riemannian manifolds.
For functional linear regression, we adopt iRMISE to quantify the quality of the estimator for slope function , and assess the prediction performance by prediction RMSE on independent test dataset. For comparison, we also fit the functional linear model using the principal components produced by RFPCA (Dai and Müller, 2017), and hence we refer to this competing method as RFLR. For both methods, the tuning parameter which is the number of principal components included for , is selected by using an independent validation data of the same size of the training data to ensure fair comparison between two methods. The simulation results are presented in Table 3. As expected, we observe that on both methods produce the same results. For the SPD manifold, in terms of estimation, we see that iRFLR yields far better estimators than RFLR does. Particularly, we again observe that, the quality of RFLR estimators does not improve significantly when sample size increases, in contrast to estimators based on the proposed iRFLR. For prediction, iRFLR outperforms RFLR by a significant margin. Interestingly, comparing to estimation of slope function where the RFLR estimator is much inferior to the iRFLR one, the prediction performance by RFLR is relatively closer to that by iRFLR. We attribute this to the smoothness effect brought by the integration in model (11). Nevertheless, although the integration cancels out certain discrepancy between the intrinsic and the ambient geometry, the loss of efficiency is inevitable for the RFLR method that is bound to the ambient spaces. In addition, we observe that the performance of both methods for Gaussian noise is slightly better than that in the case of heavy-tailed noise.
5.2 Data application
We apply the proposed iRFPCA and iRFLR to analyze the relationship between functional connectivity and behavioral data from the HCP 900 subjects release (Essen et al., 2013). Although neural effects on language (Binder et al., 1997), emotion (Phana et al., 2002) and fine motor skills (Dayan and Cohen, 2011) have been extensively studied in the literature, scarce is the exploration on human behaviors that do not seem related to neural activities, such as endurance. Nevertheless, a recent research by Raichlen et al. (2016) suggests that endurance can be related to functional connectivity. Our goal is to study the endurance performance of subjects based on their functional connectivity.
The data consists of subjects who are healthy young adults, in which each subject is asked to walk for two minutes and the distance in feet is recorded. Also, each subject participates in a motor task, where participants are asked to act according to presented visual cues, such as tapping their fingers, squeezing their toes or moving their tongue. During the task, the brain of each subject is scanned and the neural activities are recorded at 284 equispaced time points. After preprocessing, the average BOLD (blood-oxygen-level dependent) signals at 68 different brain regions are obtained. The details of experiment and data acquisition can be found in the reference manual of WU-Minn HCP 900 Subjects Data Release that is available on the website of human connectome project.
Our study focuses on regions that are related to the primary motor cortex, including precentral gyrus, Broca’s area, etc. At each design time point , the functional connectivity of the th subject is represented by the covariance matrix of BOLD signals from regions of interest (ROI). To practically compute , let be an -dimensional column vector that represents the BOLD signals at time from the ROIs of the th subject. We then adopt a local sliding window approach (Park et al., 2017) to compute by
[TABLE]
where is a positive integer that represents the length of the sliding window to compute for . Without loss of generality, we reparameterize each from the domain to . In practice, the window length is required to be sufficiently large to make nonsingular, but not too large in order to avoid significant bias. We found that the analysis below is not sensitive to the choice of as long as it is in a reasonable range, for example, . Thus, we simply set throughout our study.
The functional predictors are then the functions whose values are covariance matrices, or SPDs, and hence is regarded as a curve on . The heat maps of the first two iRFPCs along 8 time points are depicted in the top two rows of Figure 2, which shows the modes of variation in the functional connectivity with contrasts between early and middle/later times in both iRFPCs. The scalar response is the distance that subjects walked in two minutes. A scatterplot of the first two principal component scores labeled by the response variable for 10 randomly selected subjects is given in Figure 3 respectively for iRFPCA and RFPCA, to visualize the difference of scores produced by these methods. Gender, age, gait speed and strength are included as baseline covariates, selected by the forward-stepwise selection method (Section 3.3, Hastie, Tibshirani and Friedman 2009). Among these covariates, gender and age are in accordance with the common sense about endurance, while gait speed and muscle strength could be influential since endurance is measured by the distance walked in two minutes. Our primary interest is to assess the significance of the functional predictor when effect of the baseline covariates is controlled.
To fit the intrinsic functional linear model, we adopt the cross-validation procedure to select the number of components to be included in representing the Riemannian functional predictor and the Riemannian slope function . For assessment, we conduct 100 runs of 10-fold cross-validation, where in each run we permute the data independently. In each run, the model is fitted on 90% data and the MSE for predicting the walking distance is computed on the other 10% data for both iRFLR and RFLR methods. The fitted intrinsic Riemannian slope function displayed in the bottom panel of Figure 2 shows the pattern of weight changes. The MSE for iRFLR is reduced by around 9.7%, compared to that for RFLR. Moreover, the for iRFLR is , with a p-value based on a permutation test of 1000 permutations, which is significant at level 5%. In contrast, the for RFLR drops to and the p-value is that does not spell significance at all.
Appendix A Background on Riemannian Manifold
We introduce geometric concepts related to Riemannian manifolds from an intrinsic perspective without referring to any ambient space.
A smooth manifold is a differentiable manifold with all transition maps being differentiable. For each point on the manifold , there is a linear space of tangent vectors which are derivations. A derivation is a linear map that sends a differentiable function on into and satisfies the Leibniz property. For example, if is the derivation associate with the tangent vector at , then for any , where is a collection of real-valued differentiable functions on . For submanifolds of a Euclidean space for some , tangent vectors are often perceived as vectors in that are tangent to the submanifold surface. If one interprets a Euclidean tangent vector as a directional derivative along the vector direction, then Euclidean tangent vectors coincide with our definition of tangent vectors on a general manifold. The linear space is called the tangent space at . The disjoint union of tangent spaces at each point constitutes the tangent bundle, which is also equipped with a smooth manifold structure induced by . The tangent bundle of is conventionally denoted by . A (smooth) vector field is a map from to such that for each . It is also called a smooth section of . Noting that a tangent vector is a tensor of type , a vector field can be viewed as a kind of tensor field, which assigns a tensor to each point on . A vector field along a curve on is a map from an interval to such that . For a smooth function from a manifold and to another manifold , the differential of at is a linear map from to , such that for all and .
An affine connection on is a bilinear mapping that sends a pair of smooth vector fields to another smooth vector field and satisfies and the Leibniz rule in the second argument, that is, for smooth real-valued functions on . An interesting fact is that is determined by and a neighborhood of . A vector field is parallel if for all vector fields . When a connection is interpreted as the covariant derivative of a vector field with respect to another one, a parallel vector field is intuitively a constant vector field in some sense, as vanishing derivative indicates constant. Similarly, a vector field along a curve is parallel if for all , where denotes the differential of at . Given a tangent vector at , there exists a unique vector field along such that and . In this case, is the parallel transport of along from to . In particular, if , then is called a geodesic of the connection . Also, for any , there exists a unique geodesic such that and . Then the exponential map at , denoted by , is defined by .
A Riemannian manifold is a smooth manifold with a metric tensor , such that for each , the tensor defines an inner product on . For Riemannian manifold, the fundamental theorem of Riemannian geometry says that there is a unique connection that is 1) compatible with the metric tensor, that is, for every pair of vector fields and , and every vector , ; 2) torsion-free, that is, for any vector fields and , , where denotes the Lie bracket of a pair of vector fields , and is defined by for all (by definition, is also a vector field). This connection is called Levi–Civita connection, which is uniquely determined by the Riemannian metric tensor. The Riemannian curvature tensor is a map that sends a triple of smooth vector fields to another vector field . Then the sectional curvature is defined for a 2-plane spanned by two linearly independent tangent vector at the same point , given by
[TABLE]
The metric tensor also induces a distance function on that turns a connected Riemannian manifold into a metric space. The function is defined in the following way. For a continuously differentiable curve , the length of is given by . Then is the infimum of over all continuously differentiable curves joining and . For a connected and complete Riemannian, given two points on the manifold, there is a minimizing geodesic connecting these two points.
Appendix B Implementation for
Given -valued functional data , below we briefly outline the numerical steps to perform iRFPCA. The computation details for can be found in Dai and Müller (2017).
Step 1.
Compute the sample Fréchet mean . As there is no analytic solution, the recursive algorithm developed by Cheng et al. (2016) can be used.
Step 2.
Select an orthonormal frame along . For , at each , the tangent space is isomorphic to . This space has a canonical linearly independent basis with , defined in the following way. For an integer , let be the largest integer such that . Let . Then is defined as the matrix that has 1 at , 1 at and 0 elsewhere. Because the inner product on the space is given by
[TABLE]
for , in general this basis is not orthonormal in . To obtain an orthonormal basis of for any given , we can apply the Gram-Schmidt procedure on the basis . The orthonormal bases obtained in this way smoothly vary with and hence form an orthonormal frame of along .
Step 3.
Compute the -coordinate representation of each . For , the logarithm map at a generic is given by for , where denotes the matrix logarithm function. Therefore,
[TABLE]
Using the orthonormal basis obtained in the previous step, one can compute the coefficient representation of for any given .
Step 4.
Compute the first eigenvalues and eigenfunctions of the empirical covariance function . This step is generic and does not involve the manifold structure. For , the classic univariate FPCA method such as Hsing and Eubank (2015) can be employed to derive the eigenvalues and eigenfunctions of . When , each observed coefficient function is vector-valued. FPCA for vector-valued functional data can be performed by the methods developed in Happ and Greven (2018) or Wang (2008).
Step 5.
Compute the scores . Finally, compute the approximations of by the first principal components using
[TABLE]
where for , the exponential map at a generic is given by
[TABLE]
for , where denotes the matrix exponential function.
Appendix C Proofs of Main Theorems
Proof of Theorem 1.
We first show that is a Hilbert space. It is sufficient to prove that the inner product space is complete. Suppose is a Cauchy sequence in . We will later show that there exists a subsequence such that
[TABLE]
Since is complete, the limit is -a.s. well defined and in . Fix any and choose such that implies . Fatou’s lemma applying to the function implies that if , then . This shows that . Since , we see that . The arbitrariness of implies that . Because , we conclude that converges to in .
It remains to show (14). To do so, we choose so that . This is possible since is a Cauchy sequence. Let . By Cauchy–Schwarz inequality, . Thus, Then (14) follows, because otherwise, if the series diverges on a set with , then a choice of such that for contradicts the above inequality.
Now let be a measurable orthonormal frame. For every element , the coordinate representation of with respect to is denoted by . One can see that is an element in the Hilbert space of square integrable -valued measurable functions with norm for . If we define the map by , we can immediately see that is a linear map. It is also surjective, because for any , the vector field along given by for is an element in , since . It can be also verified that preserves the inner product. Therefore, it is a Hilbertian isomorphism. Since is separable, the isomorphism between and implies that is also separable. ∎
Proof of Proposition 2.
The regularity conditions on , and ensure that and are measurable. Part 1, 2 and 6 can be deduced from the fact that is a unitary operator between two finite-dimensional real Hilbert spaces and its inverse is . To reduce notational burden, we shall suppress the subscripts from and below. For Part 3,
[TABLE]
To prove Part 4, assume . Then, noting that and , we have
[TABLE]
and
[TABLE]
Part 5 is seen by the following calculation: for ,
[TABLE]
∎
Proof of Proposition 4.
The case is already given by Dai and Müller (2017) with . Suppose . The second statement follows from the first one if we let , and for any fixed and note that is independent of .
For the first statement, the inequality is clearly true if , or . Now suppose , and are distinct points on . The minimizing geodesic curves between these points form a geodesic triangle on . By Toponogov’s theorem (the hinge version), , where is the model space with constant sectional curvature . For , it is taken as the hyperboloid with curvature . Let , and . The interior angle of geodesics connecting to and to is denoted by . Denote , the law of cosine on gives
[TABLE]
By definition that and , . By Taylor expansion , we have
[TABLE]
where with denotes the norm of tangent vectors. The last inequality is due to the fact that , .
For the second term , we first have . Also, the Euclidean law of cosine implies that and . Therefore,
[TABLE]
where is the diameter of . Combining with (15),
[TABLE]
where and constants uniform for . Therefore
[TABLE]
It can be shown that . Thus,
[TABLE]
where , or in other words, . ∎
Proof of Proposition 5.
Part 1 follows from a simple calculation. To lighten notations, let denote for a valued function defined on . Suppose is the coordinate of under . Because
[TABLE]
one concludes that and hence is an eigenfunction of corresponding to the eigenvalue . Other results in Part 2 and 3 have been derived in Section 3. The continuity of and , in conjunction with , implies that is a mean square continuous random process and the joint measurability of passes to . Then can be regarded a random element of the Hilbert space that is isomorphic to . Also, the isomorphism maps to for each in the sample space. Then, Part 4 follows from Theorem 7.4.3 of Hsing and Eubank (2015). ∎
Proof of Theorem 6.
The strong consistency stated in Part 2 is an immediate consequence of Lemma 12. For Part 1, to prove continuity of , fix . Let be compact. By B.3, . Thus,
[TABLE]
The continuity assumption of sample paths implies as . Then by condition B.5, as , and the the continuity of follows. The uniform continuity follows from the compactness of . Given Lemma 11 and 12, the a.s. continuity of can be derived in a similar way. The first statement of Part 4 is a corollary of Part 3, while the second statement follows from the first one and the compactness of . It remains to show Part 3 in order to conclude the proof, as follows.
Let and be the minimizing geodesic from to at unit time. The first-order Taylor series expansion at yields
[TABLE]
where an expression for is provided in the proof of Lemma 10.
Since , we deduce from (16) that
[TABLE]
By LLN, in probability, while is invertible for all by condition B.6. In light of Lemma 10, this result suggests that with probability tending to one, for all , is invertible, and also
[TABLE]
and according to (17),
[TABLE]
where the terms do not depend on . Given this, we can now conclude the proof of Part 3 by applying a central limit theorem in Hilbert spaces (Aldous, 1976) to establish that the process converges to a Gaussian measure on tensor Hilbert space with for a random element in the tensor Hilbert space . ∎
Proof of Theorem 7.
Note that
[TABLE]
For , it is seen that
[TABLE]
With smoothness of , continuity of and compactness of , one can show that . By the uniform consistency of , with the same Taylor series expansion in (16) and the technique in the proof of Lemma 10, it can be established that . Also note that by LLN, . Then, with Part 4 of Theorem 6,
[TABLE]
Similar calculation shows that and . Now, by Dauxois, Pousse and Romain (1982), . Thus, . According to Part 1 & 5 of Proposition 2, are also eigenvalues of . The results for and follow from Bosq (2000). Those for are due to . ∎
Proof of Theorem 8.
In this proof, both and are understood to be uniform for the class . Let . Then
[TABLE]
The first term is of order uniform for the class , according to a technique similar to the one in the proof of Lemma 10, as well as Theorem 6 (note that the results in Theorem 6 are uniform for the class ). Then the convergence rate is established if one can show that
[TABLE]
which follows from
[TABLE]
and Proposition 4. It remains to show (18).
We first observe that because ,
[TABLE]
Define
[TABLE]
Then
[TABLE]
It is clear that the term is asymptotically dominated by and . Note that the compactness of in condition C.2 implies . Then, by Theorem 7, for , we have the bound
[TABLE]
It is easy to verify that, because , is asymptotically dominated by . Thus,
[TABLE]
Now we focus on the term . Note that . Then
[TABLE]
With Theorem 6 and CLT, it can be shown that Note that which implies that uniformly for all . In conjunction with the choice of , one can conclude that
[TABLE]
Thus,
[TABLE]
Finally, (18) follows from the above equation, the equation (19) and . ∎
Proof of Theorem 9.
In this proof, both and are understood to be uniform for the class . First, . The first term is of order uniform for the class , according to a technique similar to the one in the proof of Lemma 10, as well as Theorem 6 (note that the results in Theorem 6 are uniform for the class ). For the second term, using Proposition 2, one can show that , where and denote the identity operators on and , respectively. It remains to show .
Define , , and
[TABLE]
It is easy to see that the dominant term is and for and , respectively. With , it has been shown in Hall and Horowitz (2007) that . Denote . Then
[TABLE]
By Theorem 7, . Also, one can see that , with the assumption that . Also, . Using the similar technique in Hall and Horowitz (2005), we can show that , and hence conclude that . For ,
[TABLE]
Combining all results above, we deduce that and thus
[TABLE]
according to condition C.2 and Proposition 4. ∎
Appendix D Ancillary Lemmas
Lemma 10**.**
, where is as in (16).
Proof.
With the continuity of and compactness of , the existence of local smooth orthnormal frames (e.g., Proposition 11.17 of Lee (2002)) suggests that we can find a finite open cover for such that there exists a smooth orthonormal frame for the th piece of , where denotes topological closure of a set . For fixed , by mean value theorem, it can be shown that
[TABLE]
for and , where is the orthonormal frame extended by parallel transport of along minimizing geodesic.
Take as . For each , by the same argument of Lemma 3 of Kendall and Le (2011), together with continuity of and the continuity of the frame , we can choose a continuous positive such that, and for where denotes the ball on centered at with radius ,
[TABLE]
In the above, plays a role of in (20). Let and . We then have
[TABLE]
For (21), the Lipschitz condition of B.7 and smoothness of imply that
[TABLE]
As , could be chosen so that as . Thus, with probability tending to one, by Markov inequality, we deduce that
[TABLE]
For the first term in (22), LLN shows that
[TABLE]
or
[TABLE]
For the second term in (22), the compactness of the Lipschitz condition of B.7 and smoothness of also imply that . Consequently, by LLN,
[TABLE]
Combining (23), (24) and (25), with , one concludes that . ∎
Lemma 11**.**
Suppose conditions A.1 and B.1-B.3 hold. For any compact subset , one has
[TABLE]
Proof.
By applying the uniform SLLN to , for a given ,
[TABLE]
Therefore, there exists a set such that , and for all ,
[TABLE]
since by condition B.3. Fix . By the inequality , for all and ,
[TABLE]
with . Now, let be chosen such that if and . Suppose is a -net in with . Applying uniform SLLN again, there exists a set such that , for all , and
[TABLE]
for all with . Then, for all , for all , we have
[TABLE]
and this concludes the proof. ∎
Lemma 12**.**
Assume conditions A.1 and B.1-B.5 hold. Given any , there exists such that and for all , and for all , .
Proof.
Let and . It is sufficient to show that there exists and for all , such that for all ,
[TABLE]
This is because the above two inequalities suggest that for all , is not attained at with , and hence .
Let . We first show that there exists a compact set and for some such that , and both and are greater than for all , and . This is trivially true when is compact, by taking . Now assume is noncompact. By the inequality , one has
[TABLE]
and by Cauchy–Schwarz inequality,
[TABLE]
Similarly,
[TABLE]
Now, we take at a sufficiently large distance from such that on for all , where (Heine–Borel property yields compactness of , since it is bounded and closed). Since converges to uniformly on a.s. by Lemma 11, we can find a set such that and for , and on for all and .
Finally, let and Then is compact and by condition B.5, for some constant . By Lemma 11, one can find a set with and for , such that for all , (i) and (ii) . Since on for all with , we conclude that for all if for . The proof is completed by noting that can serve the we are looking for. ∎
The reference list from the paper itself. Each links out to its DOI / PubMed record.
- 1Afsari (2011) {barticle} [author] \bauthor \bsnm Afsari, \bfnm Bijan \binits B. ( \byear 2011). \btitle Riemannian L p superscript 𝐿 𝑝 L^{p} center of mass: existence, uniqueness, and convexity. \bjournal Proceedings of the American Mathematical Society \bvolume 139 \bpages 655–673. \endbibitem
- 2Aldous (1976) {barticle} [author] \bauthor \bsnm Aldous, \bfnm David J. \binits D. J. ( \byear 1976). \btitle A Characterisation of Hilbert Space Using the Central Limit Theorem. \bjournal Journal London Mathematical Society \bvolume 14 \bpages 376–380. \endbibitem
- 3Arsigny et al. (2007) {barticle} [author] \bauthor \bsnm Arsigny, \bfnm Vincent \binits V., \bauthor \bsnm Fillard, \bfnm Pierre \binits P., \bauthor \bsnm Pennec, \bfnm Xavier \binits X. and \bauthor \bsnm Ayache, \bfnm Nicholas \binits N. ( \byear 2007). \btitle Geometric Means in a Novel Vector Space Structure on Symmetric Positive-Definite Matrices. \bjournal SIAM Journal of Matrix Analysis and Applications \bvolume 29 \bpages 328–347. \endbibitem
- 4Balakrishnan (1960) {barticle} [author] \bauthor \bsnm Balakrishnan, \bfnm A. V. \binits A. V. ( \byear 1960). \btitle Estimation and Detection Theory for Multiple Stochastic Processes. \bjournal Journal of Mathematical Analysis and Applications \bvolume 1 \bpages 386–410. \endbibitem
- 5Bhattacharya and Patrangenaru (2003) {barticle} [author] \bauthor \bsnm Bhattacharya, \bfnm Rabi \binits R. and \bauthor \bsnm Patrangenaru, \bfnm Vic \binits V. ( \byear 2003). \btitle Large sample theory of intrinsic and extrinsic sample means on manifolds. I. \bjournal The Annals of Statistics \bvolume 31 \bpages 1–29. \bdoi 10.1214/aos/1046294456 \endbibitem
- 6Binder et al. (1997) {barticle} [author] \bauthor \bsnm Binder, \bfnm Jeffrey R. \binits J. R., \bauthor \bsnm Frost, \bfnm Julie A. \binits J. A., \bauthor \bsnm Hammeke, \bfnm Thomas A. \binits T. A., \bauthor \bsnm Cox, \bfnm Robert W. \binits R. W., \bauthor \bsnm Rao, \bfnm Stephen M. \binits S. M. and \bauthor \bsnm Prieto, \bfnm Thomas \binits T. ( \byear 1997). \btitle Human Brain Language Areas Identified by Functional Magnetic Resonance Imaging. \bjournal Journal of
- 7Bosq (2000) {bbook} [author] \bauthor \bsnm Bosq, \bfnm Denis \binits D. ( \byear 2000). \btitle Linear Proceses in Function Spaces. \bseries Lecture Notes in Statistics. \bpublisher Springer. \endbibitem
- 8Cardot, Ferraty and Sarda (2003) {barticle} [author] \bauthor \bsnm Cardot, \bfnm Hervé \binits H., \bauthor \bsnm Ferraty, \bfnm Frédéric \binits F. and \bauthor \bsnm Sarda, \bfnm Pascal \binits P. ( \byear 2003). \btitle Spline estimators for the functional linear model. \bjournal Statistica Sinica \bvolume 13 \bpages 571–591. \endbibitem