{'fraction_non_alphanumeric': 0.05153498663224855, 'fraction_numerical': 0.01373651700931133, 'mean_word_length': 4.138904680689786, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 6, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Video-Based Action Recognition Using Rate-Invariant Analysis of Covariance Trajectories

Statistical classification of actions in videos is mostly performed by extracting relevant features, particularly covariance features, from image frames and studying time series associated with temporal evolutions of these features. A natural mathematical representation of activity videos is in form of parameterized trajectories on the covariance manifold, i.e. the set of symmetric, positive-definite matrices (SPDMs). The variable execution-rates of actions implies variable parameterizations of the resulting trajectories, and complicates their classification. Since action classes are invariant to execution rates, one requires rate-invariant metrics for comparing trajectories. A recent paper represented trajectories using their transported square-root vector fields (TSRVFs), defined by parallel translating scaled-velocity vectors of trajectories to a reference tangent space on the manifold. To avoid arbitrariness of selecting the reference and to reduce distortion introduced during this mapping, we develop a purely intrinsic approach where SPDM trajectories are represented by redefining their TSRVFs at the starting points of the trajectories, and analyzed as elements of a vector bundle on the manifold. Using a natural Riemannain metric on vector bundles of SPDMs, we compute geodesic paths and geodesic distances between trajectories in the quotient space of this vector bundle, with respect to the re-parameterization group. This makes the resulting comparison of trajectories invariant to their re-parameterization. We demonstrate this framework on two applications involving video classification: visual speech recognition or lip-reading and hand-gesture recognition. In both cases we achieve results either comparable to or better than the current literature.


Introduction
The problem of classification of human actions or activities in video sequences is both important and challenging. It has applications in video surveillance, lip reading, pedestrian tracking, hand-gesture recognition, manufacturing quality control, human-machine interfaces, and so on. Since the size of video data is generally very high, the task of video classification is often performed by extracting certain low-dimensional features of interest -geometric, motion, colorimetric features, etc -from each frame and then forming temporal sequences of these features for full videos. Consequently, analysis of videos get replaced by modeling and classification of longitudinal observations in a certain feature space. (Some papers discard temporal structure by pooling all feature together but that represents a severe loos of information.) Since many features are naturally constrained to lie on nonlinear manifolds, the corresponding video representations form parameterized trajectories on these manifolds. Examples of these manifolds include unit spheres, Grassmann manifolds, tensor manifolds, and the space of probability distributions.
One of the most commonly used and effective feature in image analysis is a covariance matrix, as shown via applications in medical imaging [1,2] and computer vision [3,4,5,6,7,8]. These matrices are naturally constrained to be symmetric positive-definite matrices (SPDMs) and have also played a prominent role as region descriptors in texture classification, object detection, object tracking, action recognition and face recognition. Tuzel et al. [3] introduced the concept of covariance tracking where they extracted a covariance matrix for each video frame and studied the temporal evolution of this matrix in the context of pedestrian tracking in videos. Since the set of SPDMs is a nonlinear manifold, denoted byP, a whole video segment can be represented as a (parameterized) trajectory onP. In this paper we focus on the problem of classification of actions or activities by treating them as parameterized trajectories onP. The two specific applications we will study are: visual-speech recognition and hand-gesture classification. Fig. 1 illustrates examples of video frames for these two applications.
One challenge in characterizing activities as trajectories comes from the variability in execution rates. The execution rate of an activity dictates the parameterization of the corresponding trajectory. Even for the same activity performed by the same person, the execution rates can potentially differ a lot. Different execution rate implies that the corresponding trajectories go through the same sequences of points inP but have different parameterizations. Directly analyzing such trajectories without alignment, e.g. comparing the difference, calculating point-wise mean and covariance, can be misleading (the mean is not representative of individual trajectories, and the variance is artificially inflated).
To make these issues precise, we develop some notation first. Let α : [0, 1] →P be a trajectory and let γ : [0, 1] → [0, 1] be a positive diffeomorphism such that γ(0) = 0 and γ(1) = 1. This γ plays the role of a time-warping function, or a re-parameterization function, so that the composition α • γ is now a time-warped or re-parameterized version of α. In other words, the trajectory α•γ goes through the same set of points as α but at a different rate (speed).
1. Pairwise Registration: Now, let α 1 , α 2 : [0, 1] →P be two trajectories onP. The process of registration of α 1 and α 2 is to find a warping γ such that α 1 (t) is optimally registered to α 2 (γ(t)) for all t ∈ [0, 1]. In order to ascribe a meaning to optimality, we need to develop a precise criterion. Figure 1: Examples of video frames in visual-speech recognition (first two rows) and handgesture classification (last two).
2. Groupwise or Multiple Registration: This problem can be extended to more than two trajectories: let α 1 , α 2 , . . . , α n be n trajectories onP, and we want to find out time warpings γ 1 , γ 2 , . . . , γ n such that for all t, the variables {α i (γ i (t))} n i=1 are optimally registered. A solution for pairwise registration can be extended to the multiple alignment problem as follows -for the given trajectories, first define a template trajectory and then align each given trajectory to this template in a pairwise fashion. One way of defining this template is to use the mean of given trajectories under an appropriately chosen metric.
Notice that the problem of comparisons of trajectories is different from the problem of curve fitting or trajectory estimation from noisy data. Many papers have studied spline-type solutions for fitting curves to discrete, noisy data points on manifolds [9,10,11,12,13] but in this paper we assume that the trajectories are already available through some other means.

Past Work & Their Limitations
There are very few papers in the literature for analyzing, in the sense of comparing, averaging or clustering, trajectories on nonlinear manifolds. Let dP denote the geodesic distance resulting from the chosen Riemannian metric onP. It can be shown that the quantity 1 0 dP(α 1 (t), α 2 (t))dt forms a proper distance on the setP [0,1] , the space of all trajectories onP. For example, [14] uses this metric, combined with the arc-length distance on S 2 , to cluster hurricane tracks. However, this metric is not immune to different temporal evolutions of hurricane tracks. Handling this variability requires performing some kind of temporal alignment. It is tempting to use the following modification of this distance to align two trajectories: but this can lead to degenerate solutions (also known as the pinching problem, described for real-valued functions in [15]). Pinching implies that a severely distorted γ is used to eliminate (or minimize) those parts of α 2 that do not match with α 1 , which can be done even when α 2 is mostly different from α 1 . While this degeneracy can be avoided using a regularization penalty on γ, some of the other problems remain, including the fact that the solution is not symmetric. A recent solution, presented in [16,4], develops the concept of elastic trajectories to deal with the parameterization variability. It represents each trajectory by its transported square-root vector field (TSRVF) defined as: where c is pre-determined reference point onP and → denotes a parallel transport of the vectorα(t) from the point α(t) to c along a geodesic path. This way a trajectory can be mapped into the tangent space T c (P) and one can compare/align them using the L 2 norm on that vector space. More precisely, the quantity inf γ h α 1 − h α 2 •γ provides not only a criterion for optimality of γ but also approximates a proper metric for averaging and other statistical analyses. (The exact metric is based on the use of semigroupΓ, the set of all absolutely continuous, weakly increasing functions, rather than Γ. We refer the reader to [16] for details.) This TSRVF representation is an extension of the SRVF used for elastic shape analysis of curves in Euclidean spaces [17]. One limitation of this framework is that the choice of reference point, c, is left arbitrary. The results can potentially change with c and make it difficult to interpret the results. A related, and bigger issue, is that the transport of tangent vectorsα(t) to c, along geodesics, can introduce large distortion, especially when the trajectories are far from c on the manifold.

Our Approach
We present a different approach that does not require a choice of c. Here the trajectories are represented by their transported vector fields but without a need to transport them to a global reference point. For each trajectory α i , the reference point is chosen to be its starting point α i (0), and the transport is performed along the trajectory itself. In other words, for each t, the velocity vectorα i (t) is transported along α to the tangent space of the starting point α i (0). This idea has been used previously in [11] and others for mapping trajectories into vector spaces, and results in a relatively stable curve, with smaller distortion than the TSRVFs of [16]. We then develop a metric-based framework for comparing, averaging, and modeling such curves in a manner that is invariant to their re-parameterizations. Consequently, this framework provides a natural solution for removal of rate, or phase, variability from trajectory data. Another issue is the choice of Riemannian metric onP. Although a larger number of Riemannian structures and metrics have been used forP in past papers and books [18,19], they do not provide all the mathematical tools we will need for ensuing statistical analysis. Consequently, we use a different Riemannian structure than those previous papers, a structure that allows relevant mathematical tools for applying the proposed framework.
The rest of this paper paper is organized as follows. In Section 3, we introduce a framework of aligning, averaging and comparing of trajectories on a general manifold M . Since we mainly focus onP, in Section 2, we introduce the details of the Riemannian structure onP used to implement the comparison framework described in Section 3. In Section 4, we demonstrate the proposed work with real-world action recognition data involving two applications: lip-reading and hand-gesture recognition.

Riemannian Structure onP
Here we will discuss the geometry ofP and impose a Riemannian structure that facilitates our analysis of trajectories onP. Most of the background material is derived in Appendix A with only the final expressions noted here. For a beginning reader in differential geometry, we strongly recommend reading Appendix A first.
We start by choosing an appropriate Riemannian metric onP. Then, from the resulting structure, we derive expressions for the following: (1) geodesic paths between arbitrary points onP; (2) parallel transport of tangent vectors along geodesic paths; (3) exponential map; (4) inverse exponential map; and (5) Riemannian curvature tensor. Several past papers have studied the space of SPDMs as a nonlinear manifold and imposed a metric structure on that manifold [2,1,18,20]. While they mostly focus on defining distances, a few of them originate from a Riemannian structure with expressions for geodesics and exponential maps. However, they do not provide expressions for all desired items (parallel transport and Riemannian curvature tensor). In this section, we utilize a particular Riemannian structure onP, which is similar but not identical to the one in [2], and is particularly convenient for our purposes. This Riemannian structure has been used previously for other applications such as spline-fitting in [13].
LetP be the space of n×n SPDMs, and let P be its subset of matrices with determinant one. The tangent spaces ofP and P at I, where I is the identity matrix, are T I (P) = {A|A t = A} and T I (P) ≡ p(n) = {A|A t = A and tr(A) = 0}. The exponential map at I can be shown to be the standard matrix exponential: for any P ∈ P (P ∈P), there is A ∈ T I (P) (Ã ∈ T I (P)) such that P = e A (P = eÃ), where e A denotes the matrix (Lie) exponential; this relationship is oneto-one. Our approach is first to identify the space P with the quotient space SL(n)/SO(n) and borrow the Riemannian structure from the latter directly. Then, we straightforwardly extend the Riemannian structure on P toP. The Riemannian geometries of SL(n) and its quotient space SL(n)/SO(n) are discussed in Appendix A. As described there, the Riemannian metric at any point G is defined by pulling back the tangent vectors under G −1 to I, and then using the trace metric (see Eqn. 8). This definition leads to expressions for exponential map, its inverse, parallel transport of tangent vectors, and the Riemannian curvature tensor on SL(n). It also induces a Riemannian structure on the quotient space SL(n)/SO(n) in a natural way because it is invariant to the action of SO(n) on SL(n).

Riemannian Structure on P
To make the connection with P, we state a useful result termed the polar decomposition of square matrices. Recall that for any square matrix G ∈ SL(n), one can decompose it uniquely as G = P S where P is a SPDM with determinant one and S ∈ SO(n). We note in passing that this fact makes P a section of SL(n) under the action of SO(n). (If a group acts on a manifold, then a section of that action is defined to be a subset of the manifold that intersects each orbit of the action in at most one point. It is easy to see that P satisfies that condition for the action of SO(n) on SL(n). For any SO(n) orbit [G], P intersects that orbit in only one point given by P such that G = P S.) We also note that this section is not an orthogonal section since it is not perpendicular to the orbit, i.e., the tangent space T P (P) is not perpendicular to the tangent space T P ([G]). We will identify P with the quotient space SL(n)/SO(n) via a map π defined as: π : SL(n)/SO(n) → P, π([G]) = GGt , for anyG ∈ [G]. One can check that this map is well defined and is a diffeomorphism. This square-root is the symmetric, positive-definite square-root of a symmetric matrix. One can verify that π([G]) lies in P by lettingG = P S (polar decomposition), and then π([G]) = GGt = √ P SS t P = P . The inverse map of π is given by: π −1 (P ) = [P ] ≡ {P S|S ∈ SO(n)} ∈ SL(n)/SO(n). This establishes a one-to-one correspondence between the quotient space SL(n)/SO(n) and P. We will use the map π to push forward the chosen Riemannian metric from the quotient space SL(n)/SO(n) to P.
Geodesic between two points on P: With that induced Riemannian metric, we can derive the geodesic path and the geodesic distance between any P 1 and P 2 in P. The idea is to pullback these points into the quotient space using π −1 , compute the geodesic there and then map the result back to P using π. As mentioned earlier π −1 (P 1 ) = [P 1 ] and π −1 (P 2 ) = [P 2 ]. The expression for geodesic between these points in the quotient space is given in Appendix A. We compute A 12 ∈ p(n) such that e A 12 = P −1 1 P 2 S 12 for some S 12 ∈ SO(n). (Let P 12 = P −1 1 P 2 S 12 . Note that P 12 ∈ P, and the rotation matrix S 12 brings P −1 1 P 2 to P.) The corresponding geodesic in the quotient space (SL(n)/SO(n)) is given by t → [P 1 e tA 12 ] and, therefore, the desired geodesic in P is t → π([P 1 e tA 12 ]) = P 1 e 2tA 12 P 1 .
Parallel transport of tangent vectors along geodesic paths on P: To determine the parallel transport along the geodesic path t → π([P 1 e tA 12 ]), we recall from Appendix A that the two orthogonal subspaces of sl(n): so(n) and p(n), which we call the vertical and the horizontal tangent subspaces, respectively. We identify the tangent space T P (P) with the horizontal subspace in T P (SL(n)), so that T P (P) = {P A|A ∈ p(n)}. Now let X ∈ T P 1 (P) be a tangent vector that needs to be translated along a geodesic path from P 1 to P 2 , given by t → π([P 1 e tA 12 ]). Similar to the case of quotient space in Appendix A, let B ∈ p(n) such that X is identified with P 1 B. In the quotient space, the parallel transport of X along the geodesic is a vector field t → P 1 e tA 12 B. In P, the geodesic is obtained by the forward map π, t → π([P 1 e tA 12 ]). Therefore, the parallel transport vector field in P along the geodesic is the image of P 1 e tA 12 B under dπ: where T (t) ∈ SO(n) and π([P 1 e tA 12 ])T (t) = P 1 e tA 12 .
Riemannian curvature tensor on P: Let X, Y and Z be three vectors on T P (P), then the Riemannian curvature tensor on P with these arguments can be calculated in the follows. Let A, B and C be elements in p(n) such that X = P A, Y = P B and Z = P C, the tensor is given by: We summarize these mathematical tools for trajectories analysis on the manifold P: 1. Exponential map: Give a point P ∈ P and a tangent vector V ∈ T P (P), the exponential map is given as: exp P (V ) = √ P e 2P −1 V P .

2.
Geodesic distance: For any two points P 1 , P 2 ∈ P, the geodesic distance between them in is given by: d P (P 1 , P 2 ) = d P (I, P 12 ) = A 12 where e A 12 = P 12 ∈ P and P 12 = 3. Inverse exponential map: For any P 1 , P 2 ∈ P, the inverse exponential map can be calculated using the formula: 4. Parallel transport: Given P 1 , P 2 ∈ P and a tangent vector V ∈ T P 1 (P), the tangent vector at P 2 which is the parallel transport of V along the shortest geodesic from P 1 to P 2 is: 5. Riemannian curvature tensor: For any point P ∈ P and tangent vectors X, Y and Z ∈ T P (P), the Riemannian curvature tensor R(X, Y )(Z) is given as:

Extension of Riemannian Structure toP
We now extend the Riemannian structure on P toP. Since for anyP ∈P we have det(P ) > 0, we can expressP = (P, 1 n log(det(P ))) with P =P det(P ) 1/n ∈ P. Thus,P is identified with the product space of P × R. Moreover, for any smooth function ψ on P, the Riemannian metric g (P,x) = ψ 2 dx 2 + g P , whereg (P,x) and g P denote the Riemannian metrics onP and P respectively, givesP the structure of a "warped" Riemannian product. In this paper we set ψ ≡ √ n.
Parallel transport of tangent vectors along geodesic paths onP: IfṼ is a tangent vector tõ P atP , whereP is identified withP = (P, x) and x = 1 n log(det(P )), we can expressṼ as V = (V, v) with V being a tangent vector of P at P . The parallel transport ofṼ is the parallel transport of each of its two components in the corresponding space, and the transport of v is itself. The other part, i.e. the parallel translation of V in P has been dealt earlier.
Riemannian curvature tensor onP: By identifying the tangent vectorṼ = (V, v), the nonzero curvature tensors are just the part in P. Therefore, the curvature tensor forX = (X, x), Y = (Y, y) andZ = (Z, z) on TP (P) is the same as R(X, Y )Z. Again, we summarize these mathematical tools onP: 1. Exponential map: GiveP ∈P and a tangent vectorṼ ∈ TP (P).

Differences From Past Riemannian Structures
Pennec et al. [2] and others [21,22,23] and even others before them have utilized a Riemannian structure on P that is commonly known as the Riemanian metric in the literature. Here we clarify the difference between our Riemannian structure and the structure used in these past papers.
We have used a natural bijection π : SL(n)/SO(n) → P given by π([G]) = √ GG T to inherit the Riemannian structure from the quotient space to P. Its inverse is given by π −1 (P ) = [P ]. There is another bijection between these same two sets, defined by π B ([G]) = GG T , and with the inverse π −1 B (P ) = [ (Note that an element of P can be raised to any positive power in a well-defined way by diagonalizing.) If we use π B to transfer the Riemannian structure from the quotient space SL(n)/SO(n) onto P, we will reach the one used in earlier papers. In other words, the structure on the quotient space is the same for both cases, only the inherited structure on P is different due to the different bijections. Note that for all P ∈ P, π B (π −1 (P )) = P 2 . It follows that the map from P → P defined by P → P 2 forms an isometry with our metric on one side and the older metric on the other.
The Riemannian metric on SL(n) (trace metric on pulled-back matrices at identity) is invariant under right translation by elements of SO(n) and under left translation by elements of SL(n). It follows that the induced metric on the quotient space is invariant under the natural left action of SL(n) on SL(n)/SO(n). (This action is given by (g, [G]) → [gG]).) It is a wellknown fact that, up to a fixed scalar multiple, this is the only metric on SL(n)/SO(n) that is invariant under the action by SL(n). This aspect is same for both the cases. If we push forward the action of SL(n) on the quotient space via π, the resulting action of SL(n) on P is given by (G, P ) → √ GP 2 G T , but if we push forward the action via π B , the resulting action of SL(n) on P is given by Given a simple relationship between the two Riemannian structures, one may conjecture that the relevant formulae (parallel transport, curvature tensor, etc) can be mapped from one setup to another. While this is true in principle, the practice is much harder since it involves the expression for the differential of this map which is somewhat complicated.

Analysis of Trajectories onP
Now that we have expressions for calculating relevant geometric quantities onP, we return to the problem of analyzing trajectories onP. In the following we derive the framework for a general Riemannian manifold M , keeping in mind that M =P in our applications. If α is a trajectory on M , then α • γ is a trajectory that follows the same sequence of points as α but at the evolution rate governed by γ. More technically, the group Γ acts on F, F × Γ → F, according to

Representation of Trajectories
We introduce a new representation of trajectories that will be used to compare and register them. We assume that for any two points α(τ 1 ), α(τ 2 ) ∈ M, τ 1 = τ 2 , we have a mechanism for parallel transporting any vector v ∈ T α(τ 1 ) (M ) along α from α(τ 1 ) to α(τ 2 ), denoted by Definition 1 Let α : [0, 1] → M denote a smooth trajectory starting with p = α(0). Given a trajectory α, and the velocity vector fieldα, define its transported square-root vector field (TSRVF) to be a scaled parallel-transport of the vector field along α to the starting point p according to: where | · | denotes the norm that is defined through the Riemannian metric on M .
Note the difference in this definition from the one in [16] where the parallel transport was along geodesics to a reference point c. Here the parallel transport is along α and to the starting point p. The concept of parallel transportation of velocity vector along trajectories has also been used previously in [11] and others. This reduces distortion in representation relative to the parallel transport of [16] to a faraway reference point.
This TSRVF representation maps a trajectory α on M to a curve q in T p (M ). What is the range space of all such mapping? For any point p ∈ M , denote the set of all square-integrable curves in T p (M ) as C p ≡ L 2 ([0, 1], T p (M )). The space of interest, then, becomes an infinitedimensional vector bundle C = p∈M C p , which is the indexed union of C p for every p ∈ M . Note that the TSRVF representation is bijective: any trajectory α is uniquely represented by a pair (p, q(·)) ∈ C, where p ∈ M is the starting point, q ∈ C p is its TSRVF. We can reconstruct the trajectory from (p, q) using the covariant integral defined later (in Algorithm 3).

Riemannian Structure on C
In order to compare trajectories, we will compare their corresponding representations in C and that requires a Riemannian structure on the latter space. Let α 1 , α 2 be two trajectories on M , with starting points p 1 and p 2 , respectively, and let the corresponding TSRVFs be q 1 and q 2 . Now α 1 , α 2 are represented as two points in the vector bundle (p 1 , q 1 ), (p 2 , q 2 ) ∈ C over M . This representation space is an infinite-dimensional bundle, whose fiber over each point p in M is C p .
We impose the following Riemannian structure on C. For an element ( To see this, suppose we have a curve in C given by (x(s), v(s, τ )), s, τ ∈ [0, 1]. The velocity vector to this curve at s = 0 is given by (x s (0), ∇ xs v(0, ·)) ∈ T x (M ) ⊕ C x , where x s denotes dx/ds, and ∇ xs denotes covariant differentiation of tangent vectors. The Riemannian inner product on C is defined in an obvious way: If (u 1 , w 1 (·)) and (u 2 , w 2 (·)) are both elements of where the inner products on the right denote the original Riemannian metric in T x (M ). Next, for given two points (p 1 , q 1 ) and (p 2 , q 2 ) on C, we are interested in finding the geodesic path connecting them. Let (x(s), v(s)) , s ∈ [0, 1] be a path with (x(0), s(0)) = (p 1 , q 1 ) and (x(1), s(1)) = (p 2 , q 2 ). We have the following characterization of geodesics on C.
Theorem 1 A parameterized path [0, 1] → C given by s → (x(s), v(s, τ )) on C (where the variable τ corresponds to the parametrization in C x ), is a geodesic on C if and only if: Here R(·, ·)(·) denotes the Riemannian curvature tensor, x s denotes dx/ds, and ∇ xs denotes the covariant differentiation of tangent vectors on tangent space T x(s) (M ).
Proof: We will prove this theorem in two steps.
and, again, the inner products on the right denote the original Riemannian metric on T x (M ). Suppose we have a path in I → T M given by s → (x(s), v(s)). We define the energy of this path by The integrand is the inner product of the velocity vector of the path with itself. It is a standard result that a geodesic on T M can be characterized as a path that is a critical point of this energy function on the set of all paths between two fixed points in T M . To derive local equations for this geodesic, we now assume we have a parameterized family of paths denoted by (x (s, t), v(s, t)), where s is the parameter of each individual path in the family (as above) and the variable t tells us which path in the family we are in. Assume 0 ≤ s ≤ 1 and t takes values on (−δ, δ) for some small δ. We want all the paths in this family to start and end at the same points of T M , so assume that (x(0, t), v(0, t)) and (x(1, t), v(1, t)) are constant functions of t.
The energy of the path with index t is given by: To simplify notation in what follows, we will write ∇ s for ∇ xs and ∇ t for ∇ xt . To establish conditions for (x, v) to be critical, we take the derivative of E(t) with respect to t at t = 0: We will use two elementary facts: , without presenting their proofs. Plugging these facts into the above calculation then makes E (0) equal to: The second equality comes from using integration by parts on the first and third term, taking into account the fact that x t and ∇ t v vanish at s = 0, 1, (since all the paths begin and end at the same point). Now, using the standard identities Thus we have derived the geodesic equations for T M .
(2) Now we consider the case of the infinite dimensional vector bundle C → M whose fiber over x ∈ M is L 2 (I, T x (M )), I = [0, 1]. A point in C is denoted by (x, v(τ )), where the variable τ corresponds to the I-parameter in L 2 (I, T x (M )). The tangent space to C at (x, v(τ )) is T x (M ) ⊕ L 2 (I, T x (M )). Suppose (u 1 , w 1 (τ )) and (u 2 , w 2 (τ )) are elements of this tangent space and we use the Riemannian metric: Now we want to work out the local equations for geodesics in C. A path in C is denoted by (x(s), v(s, τ )). The energy calculation is basically the same as above but surround everything with integration with respect to τ . So, it starts out with (Of course x s · x s does not involve the parameter τ , but surrounding it with 1 0 . . . dτ does not change its value!) In order to do the variational calculation, we now consider a parametrized family of such paths, denoted by (x(s, t), v(s, t, τ )) where we assume that x(0, t) and x(1, t) are constant functions of t, and for each τ , v(0, t, τ ) and v(1, t, τ ) are constant functions of t, since we want every path in our family to start and end at the same points of C.
Then, following through the computation exactly as in earlier case, we obtain In order for our path (x(s), v(s, τ )) to be critical for E, E (0) must vanish for every variation x t (s) of x(s) and ∇ t (v(s, τ )) of v(s, τ ), which is clearly true if and only if  2) The second equation leads to a fact that v is covariant linear, i.e. v(s, τ ) = a(s, τ ) + sb(s, τ ) and ∇ xs a = ∇ xs b = 0 for every s and τ . For a geodesic path connecting (p 1 , q 1 ) and (p 2 , q 2 ), it is natural to let a(s, τ ) = q 1 (τ ) x(0)→x(s) and b(s, τ ) = w(τ ) x(0)→x(s) , where q 1 (τ ) x(0)→x(s) and w(τ ) x(0)→x(s) represent the parallel transport of q 1 (τ ) and w(τ ) along x(0) to x(s), and w is the difference between the TSRVFs q 2 and q 1 in T x(0) (M ), defined as (q 2 ) x(1)→x(0) − q 1 . In Fig. 2, we illustrate the geodeisc path between two trajectories on a simpler manifold M = S 2 . In each case, the yellow solid line denotes the baseline x(s) and the intermediate lines are the covariant integrals (in Algorithm 3) of v(s) with starting point x(s). As comparison, the dash yellow line shows the geodesic between the starting points p 1 and p 2 on S 2 .
Theorem 1 is only a characterization of geodesics but does not provide explicit expressions for them. In the following section, we develop a numerical solution for constructing geodesics on C.

Numerical Computation of Geodesic in C
There are two main approaches in numerical construction of geodesic paths on manifolds. The first approach, called path-straightening, initializes with an arbitrary path between the given two points on the manifold and then iteratively "straightens" it until a geodesic is reached. The second approach, called the shooting method, tries to "shoot" a geodesic from the first point, iteratively adjusting the shooting direction, so that the resulting geodesic passes through the second point. In this paper, we use the shooting method to obtain the geodesic paths on C.
In order to implement the shooting method, we need the exponential map on C. Given a point (p, q) ∈ C and a tangent vector (u, w) ∈ T (p,q) (C), the exponential map exp (p,q) (s(u, w)) for s ∈ [0, 1] gives a geodesic path (x(s), v(s)) on C. Equation 3 helps us with this construction as follows. The two equations essentially provide expressions for second-order covariant derivatives of x and v components of the path. Therefore, using numerical techniques, we can perform covariant integration of these quantities to recover the path itself.

Obtain x(
Once we have a numerical procedure for the exponential map, we can establish the shooting method for finding geodesics. Let (p 1 , q 1 ) be the starting point and (p 2 , q 2 ) be the target point, the shooting method iteratively updates the tangent or shooting vector (u, w) on T (p 1 ,q 1 ) (C) such that exp (p 1 ,q 1 ) ((u, w)) = (p 2 , q 2 ). Then, the geodesic between (p 1 , q 1 ) and (p 2 , q 2 ) is given by The key step here is to use the current discrepancy between the point reached, exp (p 1 ,q 1 ) ((u, w)), and the target, (p 2 , q 2 ), to update the shooting vector (u, w), at each iteration. There are several possibilities for performing the updates and we discuss one here. Since we have two components to update, u and w, we will update them separately: (1) Fix w and update u. For the u component, the increment can come from parallel translation of the vector exp −1 p (p 2 ) (the difference between the reached pointp and the target point p 2 ) fromp to p 1 , wherep is the first component of reached point exp (p 1 ,q 1 ) ((u, w)). (2) Fix u and update w. For the w component, we can take the difference between q 2 and the second component of the point reached (denoted asq) as the increment. This is done by parallel translatingq to T p 2 (M ) (the same space as q 2 ) and calculate the difference, and then parallel translate the difference to T p 1 (M ) to update w.

Algorithm 2 Shooting algorithm for calculating geodesic on C
Given (p 1 , q 2 ), (p 2 , q 2 ) ∈ C, select one point, say (p 1 , q 1 ), as the starting point and the other, (p 2 , q 2 ), as the target point. the shooting algorithm for calculating the geodesic from (p 1 , q 1 ) to (p 2 , q 2 ) is: 1. Initialize the shooting direction: find the tangent vector u at p 1 such that the exponential map exp p 1 (u) = p 2 on the manifold M . Parallel transport q 2 to the tangent space of p 1 along the shortest geodesic between p 1 and p 2 , denoted as q 2 . Initialize w = q 2 − p 1 . Now we have a pair (u, w) ∈ T (p 1 ,q 1 ) (C).
2. Construct a geodesic starting from (p 1 , q 1 ) in the direction (u, w) using the numerical exponential map in Algorithm 1. Let us denote this geodesic path as (x(s), v(s)), where s is the time parameter for the geodesic flow.
4. Iteratively, update the shooting direction (u, w) to reduce the discrepancy to zero. This update can be done using a two-stage approach: (1) fix u and update w until converge; (2) fix w and update u until converge.
Recall that trajectories on M and their representations in C are bijective. For each pair (p, q) ∈ C, one can reconstruct the corresponding trajectory α using covariant integration. A numerical implementation of this procedure is summarized in Algorithm 3.
Algorithm 2 allows us to calculate the geodesic between two points in C. So, for each point along the geodesic (x(s), v(s)) in C, one can easily reconstruct the trajectory on M using Algorithm 3. Here, one sets x(s) as the starting point and v(s) as the TSRVF of the trajectory. Fig. 3 shows one example of calculating geodesic using the numerical method in Algorithm 2, where M = P, the set of 3 × 3 SPDMs. In this plot, each SPDM matrix is visualized by an ellipsoid and a trajectory onP by a sequence of ellipsoids. The left panel shows two original trajectories α 1 and α 2 (their representations are (p 1 , q 1 ) and (p 2 , q 2 )), the end point of the trajectory shot (exp (p 1 ,q 1 ) (u, w)), and baseline path x(s). In this case, we selected (p 1 , q 1 ) as the starting point and computed the shooting direction (u, w) such that exp (p 1 ,q 1 ) (u, w) = (p 2 , q 2 ). The bottom panel shows the evolution of L 2 norm between the shot trajectory and the target (p 2 , q 2 ) during the shooting algorithm.

Geodesic Distance on C
Using the natural Riemannian metric on C (defined in Eqn. 2), the geodesic distance between the two points is defined as the following.
Definition 2 Given two trajectories α 1 , α 2 and their representations (p 1 , q 1 ), (p 2 , q 2 ) ∈ C, and let (x(s), v(s)) ∈ C, s ∈ [0, 1] be the geodesic between (p 1 , q 1 ) and (p 2 , q 2 ) on C, the geodesic distance is given as: This distance has two components: (1) the length between the starting points on M , l x = 1 0 |ẋ(s)|ds; and (2) the standard L 2 norm on C p 2 between the TSRVFs of the two trajectories, where q 1 represents the parallel transport of q 1 ∈ C p 1 along x to C p 2 . Since we have a numerical approach for approximating the geodesic, this same algorithm can also provide an estimate for the geodesic distance.

Geodesic Distance on Quotient Space C/Γ
The main motivation of using TSRVF representation for trajectories on M and constructing the distance d c to compare two trajectories comes from the following. If a trajectory α is warped by γ, resulting in α • γ, what is the TSRVF of the time-warped trajectory? The new TSRVF is given by: Theorem 2 For any two trajectories α 1 , α 2 ∈ F and their representations (p 1 , q 1 ), (p 2 , q 2 ) ∈ C, the metric d c satisfies d c ((p 1 , q α 1 •γ ), (p 2 , q α 2 •γ )) = d c ((p 1 , q 1 ), (p 2 , q 2 )), for any γ ∈ Γ Proof: First, if a trajectory is warped by γ ∈ Γ, the resulting trajectory is α • γ, and the starting point does not change. If the original representation of α is (p, q) ∈ C, the resulting representation is (p, q α•γ ). Second, the baseline x(s) connecting the starting points of given two trajectories will not change with respect to arbitrary γ's since the starting points are fixed, so for l x . Starting with the left side, we simplify d 2 c ((p 1 , q α 1 •γ ), (p 2 , q α 2 •γ )) as follows: where we have used s = γ(t). Theorem 2 reveals the advantage of using TSRVF representation: the action of Γ on C under the metric d c is by isometries. The isometry property of time-warping action under the metric d c allows us to compare trajectories in a manner that the comparison is invariant to the time warping. This is achieved through defining a distance in the quotient space of reparameterization group.
To form the quotient space of C modulo the re-parameterization group, we first introduceΓ as the set of all non-decreasing, absolutely continuous functions γ on [0, 1] such that γ(0) = 0 and γ(1) = 1. This set is a semigroup with the composition operation (it does not have a welldefined inverse). It can be shown that Γ is a dense subset ofΓ. Consequently, the orbit of a TSRVF q under the action ofΓ is exactly the same as the closure of the orbit of q under the action of Γ. Since the orbits underΓ form closed sets [16], while those under Γ do not, we choose to work with the former, at least for the formal development. But in practice, we will approximate the theoretical solutions using the elements of Γ. We define the quotient space C/Γ as the set of all orbits under the action ofΓ, with each orbit being: To understand the orbit, one can view it as an equivalence class. For any two trajectories α 1 , α 2 and their representations in C, (p 1 , q 1 ), (p 2 , q 2 ), we define them to be equivalent when: (1) p 1 = p 2 ; and (2) there exists a sequence γ i ∈Γ such that q α 2 •γ i converges to q 1 . In other words, if two trajectories have the same starting point, and the TSRVF of one can be time-warped into the TSRVF of the other, using a sequence of time-warpings, then these two trajectories are deemed equivalent to each other. Theorem 2 indicates that if two trajectories are warped by the same γ function, the distance d c between them remains the same. In other word, the orbits in C are "parallel" to each other.
Our goal is to define a distance such that it is invariant to the time-warping of trajectories. This can be achieved by comparing trajectories through their equivalence classes (or the orbits). That is, define a metric on the quotient space C/Γ using the inherent Riemannian metric from C. The geodesic distance on C/Γ is defined as follows. Definition 3 The geodesic distance d q on C/Γ is the shortest distance between two orbits in C, given as The geodesic d q between (p 1 , [q 1 ]) and (p 2 , [q 2 ]) is obtained by forming geodesics between all possible cross pairs in sets (p 1 , [q 1 ]) and (p 2 , [q 2 ]). Since the group action is by isometries, we can fixed one point, say (p 1 , q 1 ), and search over all (p 2 , [q 2 ]) that minimizes Eqn. 6.
Next, we focus on the problem of pairwise temporal registration between two trajectories. Eqn. 6 not only defines a metric on the quotient space C/Γ but also provides an objective function for registering two trajectories: the optimal γ * for Eqn. 6 aligns α 2 to α 1 . That is the point α 1 (t) is optimally matched with the point α 2 (γ * (t)). To solve Eqn. 6, it is equivalent to optimize the following equation: where (x, v) is the geodesic between (p 1 , q 1 ) and (p 2 , q α 2 •γ ), and q 1,x means parallel transport q 1 along x to C p 2 . Note that the time-warping γ acting on α 2 changes the underlying geodesic (x, v) between two trajectories. Algorithm 4 describes a numerical solution for optimizing Eqn. 7 on a general manifold M .
3. Parallel transportq 2 to the tangent space T p 1 (M ) along x(s), denoted asq 2 . Alignq 2 to q 1 using Dynamic Programming Algorithm and obtain the optimal warping function γ.

Note that
Step 2 corresponds to the first argument (x, v) and Step 3 corresponds to the second argument γ in Eqn. 7, respectively. The optimization over the warping function in Step 3 is achieved using the Dynamic Programming Algorithm [24]. Here one samples the interval [0, 1] using N discrete points and then restricts to only piecewise linear γ's that pass through that N × N grid. In Fig. 4, we present one example of aligning two trajectories α 1 and α 2 on space of 3 × 3 SPDMs.
Fast Approximation: Since pairwise temporal registration algorithm involves multiple evaluations of the exponential map and dynamic programming alignment, it is not computationally efficient. One way to speed up this optimization is to use an approximate method: approximate the baseline x(s) connecting two trajectories first (using geodesic between the starting points of trajectories) and then align their TSRVFs. Another way is to find an explicit expression for the exponential map. This seems possible only for a simple manifold, such as M = S 2 , but for a Before: α 1 and α 2 After: α 1 and α 2 • γ * complicated manifold, such as M =P, the analytical expressions are not known. In the experiment section, we use the approximate method to speed up the registration and comparison.
If we compare the proposed framework with that in [16], we see several advantages. The proposed work perseveres the invariance properties achieved in [16], but does not require choosing a reference point. Also, the proposed framework naturally includes the difference between the starting points of two trajectories that are ignored in [16]. Since the velocity vectors here are transported to the starting point of a trajectory, along that trajectory, as opposed to a transport to an arbitrary reference point in [16], this representation is more stable.

Statistical Summarization of Multiple Trajectories
Since d q defines a metric in the quotient space C/Γ, this framework allows us to perform statistical analysis of multiple trajectories in C/Γ. Given a set of trajectories {α i , i = 1 . . . k}, we are interested in computing the average of these trajectories and using it as a template for registering these trajectories. This sample average is calculated using the notion of Karcher mean [25]. In the space of C/Γ, the Karcher mean is defined to be: Note that (µ p , [µ q ]) is an orbit (equivalence class of trajectories) and one can select any element of this orbit as a template to help to align multiple trajectories.

Video-Based Activities Recognition
Now we turn to evaluation of the framework developed so far on test datasets involving action recognition using videos. While there are several application secures involving videobased action recognition, we try two applications here: (1) hand-gesture recognition, where one classifies the hand motion into pre-determined classes using video data, and (2) visual speech recognition, where one classifies phrases uttered using videos of lip movements.

Representations of Videos as Trajectories onP
The video data is typically extremely high-dimensional, a short video with 50-frames could have more than a million pixel values. To represent a video or a frame in this video, researchers extract few relevant features from the data, forming low-dimensional representations and devise statistical analysis of these features directly. One common technique is to extract a covariance matrix of relevant features to represent each frame in the video [3,26,8]. The definition of a covariance feature for a frame is as follows: At each pixel (or patch) of the frame image, extract an d-dimension local feature (e.g. pixel location, intensity, spatial derivatives, HOG features, etc), and then calculate a covariance matrix of those local features over the whole frame, i.e. sum over the spatial coordinates. More precisely , let f x ∈ R d denote a d-dimensional feature at location x in the image I. The empirical estimate of the covariance matrix of S is given by: x∈I f x is the empirical mean feature vector. The covariance matrix provides a natural way to fuse multiple local features, across the whole image. The dimension of the representation space, (d 2 + d)/2, is much smaller than the image size |I|.
A video of human activities is now replaced by a sequence of covariance matrices. Therefore, each video is represented by a parameterized trajectory on the space ofP, and we can utilize the previous framework to analyze these trajectories. Sometime, as in the hand-gesture recognition, we can localize features better by dividing the image frame into four quadrants and computing a covariance matrix for each quadrant. Then, each frame is represented as an element ofP 4 , and the whole video as t → α(t) ∈P 4 .

Hand Gesture Recognition
Hand gesture recognition using videos is an important research area since gesture is a natural way of communication and expressing intentions. People use gestures to depict sign language for deaf, convey messages in loud environment and to interface with computers. Therefore, an accurate and automated gesture recognition system would broadly enhance human-computer interaction and enrich our daily lives. In this section, we are interested in applying our framework in video-base (dynamic) hand gesture recognition. We use the Cambridge hand-gesture dataset [7] which has 900 video sequences with nine different hand gestures: 100 video sequences for each gesture. The nine gestures result from 3 primitive hand shapes and 3 primitive motions, and as collected under different illumination conditions. Some example gestures are shown in Fig. 6. The gestures are imaged under five different illuminations, labeled as Set1, Set2, . . . , Set5.  In addition to the illumination variability, the main challenge here comes from the fact that hands in this database are not well aligned, e.g. the proportion of a hand in an image and the location of the hand are different in different video sequences. To reduce these effects we evenly split one image into four quadrants (upper-left, upper-right, bottom-left, bottom-right) with some overlaps. Each of the four quadrants is represented by a sequence of covariance matrices (on the SPDMs manifoldP). In this experiment, we use HOG features [27] to form a covariance matrix per image quadrant as follows. We use 2 × 2 blocks of 8 × 8 pixel cells with 7 histogram channels to form HOG features. Those HOG features are then used to generate 7 × 7 covariance matrix for each quadrant of each frame. Thus, our representation of a video is now given by t → α(t) ∈P 4 .
Since we have split each hand gesture into four dynamic parts, the total distance between any two hand gestures is a composite of four corresponding distances. For each corresponding dynamic part, e.g. the upper-left part, we first the parts across videos (using Algorithm 4) individually and then compare them using the metric d q , denoted by d upl . The final distance is obtained using an weighted average of the four parts: d = λ 1 d upl +λ 2 d upr +λ 3 d downl +λ 4 d downr and 4 i=1 λ i = 1. For each illumination set, we use different weights, and these weights are trained using randomly selected half of the data (90 video sequences) in that set, and the other half of the data are used for the testing. Table 1 shows our results using the nearest neighbor classifier on all five sets. One can see that after the alignment, the recognition rate has significant improvement on every set. Also, we have reported the state-of-art results on this database [28,29]. One can see that our method outperforms these methods.  [6] 81% 81% 78% 86 % -RLPP [5] 86% 86% 85% 88 % -PM 1-NN [29] 89% 86% 89% 87 % -PMLSR [28] 93% 89% 91% 94 % -Our before alignment 94% 91% 90% 88% 77% after alignment 99% 97% 97% 96% 98% Improvement 5% 6% 7% 8% 21%

Visual Speech Recognition
This application is concerned with visual speech recognition (VSR), or lip-reading, using closeup videos of human facial movements. Speech recognition is important because it allows computers to interpret human speech and take appropriate actions. It also has applications in biometric security, human-machine interaction, manufacturing and so on. Speech recognition is performed through multiple modalities -the common speech data consists of both audio and visual components. In the case the audio information is either not available or it is corrupted by noise, it becomes important to understand the speech using the visual data only. This motivates the need for VSR. The process of visual speech recognition is to understand the words uttered by speakers, derived from the visual cues. Movements of the tongue, lips, jaw and other speech related articulators are involved to make sound. To represent a video of such dynamic process, we extract a covariance matrix for each frame. Now each video becomes a parameterized trajectory on the space ofP and we can utilize the proposed framework to analysis these trajectories.
In the experiments reported here, we utilize the commonly used OuluVS dataset [30] which includes 20 speakers, each uttering 10 everyday greetings five times: Hello, Excuse me, I am sorry, Thank you, Good bye, See you, Nice to meet you, You are welcome, How are you, Have a good time. Thus, totally the database has 1000 videos; all the image sequences are segmented, having the mouth region determined by manually labeled eye positions in each frame [31]. Some examples of the segmented mouth images are shown in Fig. 7.
For each video frame I, we extract a 7 × 7 covariance matrix representing seven features: {x, y, I(x, y), | ∂I ∂x |, | ∂I ∂y |, | ∂ 2 I ∂ 2 x |, | ∂ 2 I ∂ 2 y |}. The resulting trajectories inP are aligned using Algorithm 4 and compared using distance d q defined in Eqn. 6. In Fig. 8 (a), we show some optimal γ's obtained to align one video of phrase ("excuse me") to other videos of the same phrase spoken by the same person. One can see that there exist temporal differences in the original videos and they need to be aligned before further analysis. In (b), we show the histogram of (d c − d q )'s (differences between distances before and after alignment). In this case, each person has 50 videos, and we can calculate (50 × 49)/2 pairwise distances before and after alignment, and their differences. For all 20 persons in this dataset, we have 20 × (50 × 49)/2 = 24500 such differences. From the histogram of these differences, one can see that after the alignment, the distances (d q 's) consistently become smaller. Note that one can choose other features to obtain a better representation perhaps, but the main point here is the improvement in temporal alignment and reduced distances between trajectories.
To compare classification performance with previous methods [30,4], we perform the experiment on a subset of the whole dataset, which contains 800 video sequences by removing some short videos due to the restriction of the method in [30]. (Although our method do not have this restriction, to compare fairly, we perform the experiment in the same subset). Then, we perform the Speaker-Dependent Test (SDT, see [30] for details) on this subset. The recognition rate is calculated based on Nearest Neighbor (NN) classifier. Table 2 shows the average first nearest neighbor (1NN) classification rate of our method and previous methods. One can see that even with a simple classifier, our method has the classification rate of 78.6%, which is 8.1% better than [4]'s. This indicates the advantages of proposed intrinsic method compared with the one in [4]. One can also see that, there are 37.6% percent of improvement after the alignment (registration), which demonstrates the importance of removing temporal difference in comparing of the dynamic systems in computer vision. Several other papers have reported higher classification rates on this dataset by they generally use advanced classier (e.g. SVM, leave-one-out cross valuation), additional information (e.g. audio, image depth) and machine learning techniques [32,33,34]. Thus, their results are not directly comparable with our results.

Conclusion
In summary, we have proposed metric-based approach for simultaneous alignment and comparisons of trajectories onP, the Riemannian manifold of covariance matrices (SPDMs). In order to facilitate our analysis, we impose a Riemannian structure on this manifold, induced by the action of SO(n) and SL(n), resulting in explicit expressions for geometric quantities such as parallel transport and Riemannian curvature tensor. Returning to the trajectories, the basic idea is to represent each trajectory by a starting pointP ∈P and a TSRVF which is a curve in the tangent space TP (P). The metric for comparing these elements is a composed of: (a) the length of the path between the starting points and (b) the distortion introduced in parallel translation TSRVFs along that path. The search for optimal path, or a geodesic, is based on a shooting method, that in itself uses geodesic equations for computing the exponential map. Using a numerical implementation of the exponential map, we derive numerical solutions for pairwise alignment of covariance trajectories and to quantify their differences using a rate-invariant distance. We have applied this framework to covariance tracking in video data, with applications to hand-gesture recognition and visual-speech recognition, and obtain state of the art results in each case.