Representation and Characterization of Nonstationary Processes by Dilation Operators and Induced Shape Space Manifolds

We proposed in this work the introduction of a new vision of stochastic processes through geometry induced by dilation. The dilation matrices of a given process are obtained by a composition of rotation matrices built in with respect to partial correlation coefficients. Particularly interesting is the fact that the obtention of dilation matrices is regardless of the stationarity of the underlying process. When the process is stationary, only one dilation matrix is obtained and it corresponds therefore to Naimark dilation. When the process is nonstationary, a set of dilation matrices is obtained. They correspond to Kolmogorov decomposition. In this work, the nonstationary class of periodically correlated processes was of interest. The underlying periodicity of correlation coefficients is then transmitted to the set of dilation matrices. Because this set lives on the Lie group of rotation matrices, we can see them as points of a closed curve on the Lie group. Geometrical aspects can then be investigated through the shape of the obtained curves, and to give a complete insight into the space of curves, a metric and the derived geodesic equations are provided. The general results are adapted to the more specific case where the base manifold is the Lie group of rotation matrices, and because the metric in the space of curve naturally extends to the space of shapes; this enables a comparison between curves’ shapes and allows then the classification of random processes’ measures.


Introduction
The analysis and/or the representation of nonstationary processes has been tackled for four or five decades now by time-scale/time-frequency analysis [1,2], by Fourier-like representation when the processes belong to the periodically correlated (PC) subclass [3,4], or by partial correlation coefficients (parcors) series [5,6], to cite a few. One of the advantages of dealing with parcors resides in their strong relation to the measure of the process by the one-to-one relation with correlation coefficients [7,8]. They consequently appear explicitly in the Orthogonal Polynomial on the Real Line/Unit Circle decomposition of the measure [9,10], on the Matrices Orthogonal Polynomials on the Unit Circle [11] and its applications [12], are the elements for the construction of dilation matrices that appear in the Cantero Moral, and Velazquez (CMV)/Geronimus, Gragg, and Teplyaev (GGT) matrices [13], for the Schur flows problem with upper Hessenberg matrices [14] that are also seen in the literature as evolution operators [10] or shift operator [15], and finally appear in the state-space representation [7,16]. The dilation theory takes its roots from the operator theory [17], which bridges the process's measure and unitary operators. In its simplest version, the dilation theory corresponds to Naimark dilation [17,18], and states that given a sequence of correlation coefficients, there exists a unitary matrix W such that R n (1 0 0 · · · )W n (1 0 0 · · · ) T where · T denotes the transposition. When the process is not stationary, its associated correlation matrix is no more Toeplitz structured, a set of matrices is required [16] and the previous expression becomes R i,j (1 0 0 · · · )W i+1 W i+2 · · · W j (1 0 0 · · · ) T . The matrices W i are theoretically understood as infinite rotation matrices, which become finite when the correlation coefficients sequence is itself finite. In that particular case, the matrices W i belong to SO(n) or SU(n), the special orthogonal or unitary group, respectively, and the process's measure is totally described by the set of W i . As a consequence, the measure of the process is beautifully characterized for the nonstationary case, by a sampled trajectory induced by the dilation matrices on the appropriate Lie group. When the process is periodically correlated, the sequence of parcors is periodic, as well as the sequence of dilation matrices, which yields a closed path as illustrated in Figure 1. This raises the question of comparison of processes by means of their dilation matrices. Many efforts have been made in the last decade to exploit the hyperbolic geometric structure not of the correlation matrices directly but of the related parcors when obtained in stationary conditions [8,[19][20][21][22], as well as the information geometry structure that is closely related to it [23][24][25]. As the Kullback-Leibler divergence let do, the comparison of stationary processes is then made by comparing curves, whose sampled points are parcors sequences, defined on several copies of the Poincaré disk through geodesics deformation. Treatment of the nonstationary case has not been tackled to our knowledge with the previously mentioned approaches. In this paper, we hope to initiate interest in filling this gap by extending the representation and the characterization of the processes's measure in a nonstationary context, inspired on the one hand by information geometric application and interpretation of parcors [19,20] or correlation matrices [26][27][28], and on the other hand by theories and applications dealing with curves on manifolds [29,30], closely related to some aspects of Euler's equations [31]. Therefore, we combine Constantinescu's approach to dilation and shape analysis to the propose of seeing stochastic processes as elements of a Lie group. Characterizing the time-varying measure of the process is now tackled by studying curves (or sampled curves) on special groups.
To support the reader, some insights on dilation theory are given in Section 2. Practical implementations of dilation matrices according to the operator theory approach [16,18] or the lattice filter structure approach [32,33] are also discussed and the strong connection between parcors and the dilation matrices is emphasized. Section 3 focuses on the geometry of the curves induced by the dilation on particular manifolds. The general framework is first introduced by recalling concepts of distances and shape of curves when the ambient space is not flat. Next, the square root velocity (SRV) functions are developed and adapted to the Lie group, and a procedure to compare nonstationary processes through their time evolution trajectory is presented. Finally, a conclusion is drawn in Section 4 and the reader will find some technical tools in Appendices A and B. Figure 1. Illustration of a sampled closed trajectory drawn in SO(n) or SU(n) that materializes the time varying of the Periodically Correlated (PC) measure for a stochastic process. Each W i is a dilation matrix built through the parcors. Recall that a PC process is a process such that R s,t = R s+T,t+T for a certain T, where R ·,· stands for the correlation function of the process.

Outline of the Dilation Theory
Let us give some insights into the dilation theory. In its fundamental definition, the dilation theory consists of a Hilbert space H and an operator-valued function f , i.e., an L(H)-valued function, to find a larger Hilbert space H and another application F such that f is the orthogonal projection of F : where P H denotes the orthogonal projection onto the Hilbert space H. The ideas of the dilation theory are: • there exists a larger space from which the original function (or matrix) is deduced; • we can choose the "dilated" function to be simpler. For instance, when dealing with matrices, each of its coefficients can be expressed as the projection of a larger unitary matrix. In this case, we obtain a unitary dilation. This approach has been for example developed in [34][35][36] for the stationary dilation of periodically-correlated processes.

Dilation and Rotation of Contractions
For an operator T on a Hilbert space H, we denote by T * the adjoint operator, i.e., the operator on H such that Tx, y = x, T * y for all x, y ∈ H. An operator T ∈ L(H) is said to be a contraction if || T ||≤ 1 where || · || is the operator norm. We deduce the expression for the defect operator D T = (I − T * T) 1/2 and its adjoint D T * = (I − TT * ) 1/2 .
One of the easiest results is that, given a contraction Γ, the following unitary operator satisfies, for all n ∈ N Γ n = 1 0 J(Γ) n 1 0 .
In other words, the elementary rotation of a contraction also corresponds to the unitary dilation operator of the contraction. This operator is called the Julia operator, and corresponds to the Halmos extension [15] of a contraction.

Dilation and Isometries
Following the idea and the formulation of Naimark, the dilation theory can be restated in terms of dilation of the sequence of operators or sequence of numbers when the dimension of the underlying Hilbert space is 1. A sequence of operators {R n } ∞ n=1 acting on H is said to be positive if Assuming now that R * n = R −n and R 0 = I, leads to the following Toeplitz matrix: which is positive-definite. Remark that this matrix can be seen as the correlation matrix of a stationary process, as it is positive and Toeplitz [37]. Owing to this property, we obtain the existence of an operator U such that [18,38] and Theorem 1.1 in [39]: R n = P H U n | H , for all n ≥ 0 and U an isometry on K (6) as a result of the Naimark dilation theorem [17]. Furthermore, if K = n≥0 U n H, where denotes the linear span, then U is unique up to an isomorphism.

Dilation and Measure
From Bochner's theorem [37], we know a matrix of type (5) can be seen as the Fourier coefficient of a given positive Borelian measure. This is also known as the moment or trigonometric problem [16]. Therefore, we can restate the dilation problem in terms of measure. If we denote by E λ an operator-valued distribution function on [0, 2π], then the function is positive-definite. This shows the strong correspondence between the spectral measure and the dilation theory. There hence exists a unitary operator on a Hilbert space K such that R n = P H U(n) where P H stands for the orthogonal projection. With the spectral representation of unitary operators, U = 2π 0 e iλ dE λ and we have or, in an equivalent form: Note that the operator-valued measure F λ is in fact an orthogonal projection-valued measure because all its increments are orthogonal. With dilation matrices having been introduced, we give in the next section a methodology to understand how they are obtained.

Construction of Dilation Matrices
As mentioned previously, given an SPD matrix R = R i,j i,j∈N , it is possible to build a sequence of matrices {W i } i∈N such that R i,j = 1 0 0 · · · 0 W i W i+1 · · · W j−1 1 0 0 · · · 0 T . Let the general framework where R i,j is a complex operator satisfying R i,j ∈ L(H j , H i ) with {H n } n a sequence of Hilbert spaces and L the set of linear applications. For example, consider the stochastic process {X n } n , where X n ∈ L 2 (P) is a squared integrable random variable with respect to the probability space (Ω, F , P). Then the stochastic process can be viewed as an operator: ∼ X n : C → L 2 (P), ∼ X n λ = λX n , and the correlation kernel becomes To start the construction, let first the following theorem Theorem 1 (Structure of a positive-definite block matrix). Let X and Z be positive operators in L(H X ) and L(H Z ) respectively. Then the following are equivalent: There exists a unique contraction Γ in L(R(Z), R(X)) such that

Proof. Appendix A
Let us now apply this relation repeatedly on an SPD matrix. To fix ideas, let the R be the 3 × 3 (block-)matrix : and apply Theorem 1 to and finally to R 1,2 R 1,3 . Note that when a square root of a (block-)matrix has to be chosen, it is done according to the Schur decomposition given in Appendix A. For example, we have R 1,2 = R 1/2 1,1 Γ 1,2 R 1/2 1,1 . At each step, a contraction Γ i,j is generated with respect to the indices of the upper and lower (block-)matrices of the main diagonal, e.g., Γ 1,2 for the first We thus obtain a one-to-one correspondence between the SPD matrix R and the set of contractions Γ i,j i=1,2 j=3 . Regarding the huge work of Constantinescu [16], we will call these contractions the Schur-Constantinescu parameters. We consider now unit variance and arbitrary size n × n for the SPD matrix, which allows us to write the correspondence as: Once (12) is established, each dilation matrix W i is built-up as a product of Givens rotations of a sequence of Schur-Constantinescu parameters in the following way: where G Γ i,i+l denotes the Givens rotation of Γ i,i+l as follows: where the "non-identity" part, consisting of a Julia operator is located at the entry (i, i). When the SPD matrix is Toeplitz, which corresponds to a stationary underlying process, then all dilation matrices W i are identical and take the form which is nothing less than the Naimark dilation introduced in the first part, i.e., For the sake of completeness, we give the correspondence between the coefficients of the SPD matrix (the left-hand side of (12)) and the Schur-Constantinescu parameters: there exists a family {Γ k,j | k, j = 1, · · · n, k j} of contraction such that where B k,k is the Cholesky's square-root of R k,k . and a row contraction associated to the set of parameters {Γ k,m | k < m ≤ j}, a column contraction associated to the set of parameters {Γ m,j | m = j − 1, · · · k}, and finally Proof. This theorem is proved in ( [16], Theorem 5.3).
A different approach leading to the same results can be found in [40], using directly the Kolmogorov decomposition. In [32] the Naimark dilation is constructed using the lattice filter and finally applications of this decomposition in quantum mechanics are to be found in [41,42] for example.
We now give some remarks to conclude this part: is a semi-positive definite complex-valued Toeplitz kernel, then all the {Γ n } are complex-valued and respect |Γ i | < 1. The structure and the construction procedure for obtaining such a complex-valued parameter is identical whether the kernel is real or complex.

•
The framework proposed by Constantinescu and recalled previously is quite general and can be extended to Matrix Orthogonal Polynomial on the Unit Circle (MOPUC) development. By referring to ( [43], Section 3.1), matrix polynomials stem from a Szëgo recursion akin to the scalar case and thus provide a sequence of Verblunsky coefficients, or parcors that are matrices. Again in ( [43], Section 3.11), a correspondence between MOPUC and CMV matrices [13] is provided, which are equivalent to dilation matrices. The construction procedure remains the same for the dilation matrices, but the parcors become in that case matrices (they are matrix-valued Verblunsky coefficients), which can be obtained by a matrix version of the Schur/Geronimus algorithm [44].
Dilation matrices being now fully introduced, we focus the attention of the reader on the hidden information contained in their timely geometrical dissemination.

Analysis of Curves on a Manifold Induced by the Dilation
Parcors, composing dilation matrices, have already been given a geometrical point of view, as, for example, in [8] where the sequence of parcors associated with a stationary process is seen as a point onto the Poincaré polydisk P n , that is, the product of the Poincaré disk. To give geometrical settings, a distance to characterize individual parcors is then proposed and discussed. In [45], a stochastic process is studied under the local stationarity assumption. To each stationary slice of the process corresponds a sequence of parcors, represented as a point in the Poincaré polydisk P n as well. A trajectory is then generated on that space which materializes a curve on the manifold P n . The underlying computations are quite intricate because of the product manifolds, and the question of nonstationarity arises. Based on the works of Le Brigant [45,46], Celledoni et al. [47] and Zhang et al. [48], we propose then to give particular attention to this question. We first make use of the dilation theory introduced in Section 2. When the process under study is nonstationary, a set of matrices W i is obtained. The basic idea for having geometric information on the nonstationary process is therefore to characterize the trajectory formed by the set of dilation matrices. These matrices are theoretically operators of infinite dimension, but as we dispose of only a finite set of parcors, the theoretical matrices of (15) are truncated. Matrices respecting (15) are general rotation matrices that become perfect rotation operators belonging to SO(n) for real processes and SU(n) when dealing with complex processes, when their dimensions are reduced to n × n. Our aim is finally to analyse those curves living on the Lie group of rotation matrices and emphasize the geometry or, more precisely, the intrinsic geometry formulation of these objects. For example, we aim at comparing different curves coming from different processes or at resuming many realizations of a stochastic process (multiple measurements) through the computation of the mean of the associated several curves. The question as to computation complexity still exists, but many results have been proposed recently to overcome this difficulty and to propose closed-form formulations [30,49]. In particular, it is predicated to extract the shape of the trajectory for it contains the essential information, in a topologic sense.
To allow the curves comparison, we have based our development on the works of Le Brigant [45] and Celledoni et al. [47]. First, we define the manifold M given by the set of all curves in the base manifold. This leads to another space, the shape space, for which the manifold M will be a principal bundle. We dispose then of a metric in M from which a metric on the shape space is deduced. These steps are now explained in the following.

Preliminaries on Lie Groups
As we are going to deal with curves on a Lie group, we start with some preliminaries. A metric ·, · on a Lie group is said to be left invariant if: where (dL a ) b is the derivative in the manifold field sense (so the tangent map) of the left translation L a at b. A left-invariant metric gives the same number whenever the vectors are translated on the left. It is straightforward to adapt this definition to a right-invariant metric. A metric that is both left and right invariant is called a bi-invariant metric. A Lie group endowed with a bi-invariant metric has plenty of import properties that can be exploited for our study of curves on shape spaces. We list some of them in the following [50].
• The geodesics through e (the identity element) are the integral curves t → exp(tu), u ∈ g, that is, the one-parameter groups. In addition, because left and right are isometries and isometries maps geodesics to geodesics, the geodesics through any point a ∈ G are the left (right) translates of the geodesics through e Of course, we have • The Levi-Civita connection is given by where [·, ·] denotes the Lie bracket. We can now link these formulas to our base manifold SO(n). The Killing form, B, of a Lie algebra is the symmetric bilinear form B : where tr denotes the trace operator and ad denotes the adjoint representation of the group, namely, the map ad : G −→ GL(g) such that, for all a ∈ G ad a : g −→ g is the linear isomorphism defined by ad a = d(R −1 a • L a ) e . If we now assume B to be negative-definite, then -B is an inner product and is adjoint invariant. Thus, it is a classical result of the compact semi-simple Lie theory that -B induces a bi-invariant metric on G.
The Lie algebra of SO(n) is the set of skew-symmetric matrices which verifies M T = −M. The Killing form on SO(n) is given by B so(n) = (n − 2)tr(XY), and as a result of the skew symmetry, we have −B so(n) = (n − 2)tr(XY T ). Therefore, it induces a bi-invariant metric and the previous formula can be plugged into the expression of the metric on the space of curves. In the sequel, the manifold that supports the curves is SO(n) endowed with its bi-invariant metric.

Basic Outline of Geometry
Curves of interest are those living in the Lie group of real rotation matrices; this yields c : [0, 1] → SO(n). For the sake of clarity, assume that c ∈ C ∞ ([0, 1], SO(n)) , we will come back to the case of discrete curves later. To study the geometrical features of such curves, we interest ourselves with the set of all curves lying in SO(n) (where SO(n) is seen as a manifold) with nonvanishing velocity, i.e., M = {c ∈ C ∞ ([0, 1], SO(n)) : c (t) = 0 ∀t}, this is in fact a sub-manifold of C ∞ ([0, 1], SO(n)). A curve c is thus a particular point in M. The tangent space at a curve c is given by where TSO(n) denotes the tangent bundle of the base manifold SO(n). Note that a tangent vector is a curve in the tangent space of SO(n) (see Theorem 5.6 in [51]). When comparing two curves, it is natural that the distance between these two curves should remain the same if the curves are only reparametrized, that is, if we define other curves that pass through the same points than the original curves but at different speeds. When the curve is discretized as we will see in the sequel, doing a reparametrization is equivalent to changing the chosen points (see Figure 2). A reparametrization is represented by increasing diffeomorphism φ ∈ D : [0, 1] → [0, 1] acting on the right of the curve by composition. In other words, we required that the Riemannian metric g on M satisfies the following property: for all c ∈ M, u, v ∈ T c M and φ ∈ D. This property is called reparametrization invariance. We insist on the fact that g is the metric on M, the space of all curves on SO(n) and not on SO(n) itself. In terms of distances, this gives where d M denote the distance on M corresponding to the metric g. The reparametrization introduced above induces an equivalence relation between points in M such that with this equivalence relation, a quotient space can be constructed as the collection of equivalence classes; it is named the shape space and has the following writing: A distance function on the shape space is obtained from the distance on M as follows: where [c 0 ] and [c 1 ] are representatives of the equivalence classes of c 0 and c 1 respectively. It can be shown that this distance is independent of the choice of the representatives. Some precautions has to be taken here: whereas M is a submanifold of the Fréchet manifold C ∞ ([0, 1], SO(n)), has proven by [52], Theorem 10.4, the shape space S is not a manifold and the principal bundle structure π = M → S is not formally defined. However, a manifold structure can be obtained if we only consider free immersion [53]. As the metric defined on the shape space is reparametrization-invariant, it is constant along the "fibers" (the origin point is fixed). Further explanations on the Riemannian submersion can be found in [54]. Closed curves being of main interest in this work, we can also define the set Basically, the closure of a curve just imposes the equality of the first and the last point of it, and not of their first derivative. Consequently, M C turns into We need now to introduce the Square Root Velocity function (SRV function) [55], in which a curve is represented by its starting point and its normalized velocity at each time t. There are several possibilities to define the SRV of a curve. The more general definition is the following However, we can go further and benefit from the specific case of the Lie group. In this section, we will denote the base manifold G = SO(n) to emphasize its group structure, and g an element of the group. As in [47], we consider only curves that start at the identity; this is because other curves can be reduced to this case by right or left translation. In these settings, it is interesting to turn the SRV function into the Transported SRV function (TSRV). This is basically the SRV that has been parallel transported to a reference point. Different versions have been given in [47,48,56] which differ in the choice of their reference point. For our case of study, the identity is our natural curve starting point and is thus a particularly good choice for being the reference point. In a Lie group, a parallel transport operation can be defined, here again, by the right (or left) translation. This justifies that we can take, as suggested in [47], a TSRV function of the following form: where g is the Lie algebra, R is the right translation on the group, R g 1 (g 2 ) = g 2 g 1 , R g * = T e R g is the tangent map at the identity, || · || is a norm induced by a right-invariant metric on G, and T c(t)→I c denotes the parallel transport from c(t) to the identity according to the curve c. A curve is now represented pointwise as an element of the tangent bundle (c(0), q(t)) ∈ M × g (recall that q draws a curve in the tangent bundle), and c(0) is the identity element of the Lie group. The inverse of the SRV function is then straightforward: for every q ∈ C ∞ ([0, 1], TM), there exists a unique curve c such that F(c i ) = q i and c(t) = t 0 q(r) || q(r) || dr where || · || is the norm in SO(n).

Metric and Distance over M and S
We now give insights on a relevant metric that should be used on M to compare different closed trajectories. The following development and expression of metrics and distances can be found in [45]. The distance on the shape space is used to compare how the curves are intrinsically different. It has been seen in [57] that the simple L 2 metric on M given by where ·, · is the Riemannian metric on SO(n), induced a vanishing metric on the shape space, that is, we cannot differentiate shape with this metric. To overcome this difficulty, the family of elastic metric, derived from the Sobolev metric [58,59], has been investigated for it is non-vanishing on the shape space. In the case of a Euclidean space R n , it admits the expression: where D l u = h / || c ||, D l u T = D l u, w w, with w = c / || c || and D l u N = D l u − D l u T . Here, we are only interested in the special metric that has been proposed in [45], and which is an adaptation of the elastic metric for the Riemannian manifold. In our case this gives: where ∇ is the Levi-Civita connection that corresponds to ·, · ; ∇ l u = 1 || c || ∇ c h, ∇ l u T = ∇ l u, w w, w = c / || c ||. For the computations being done now in a manifold space, the Levi-Civita connection replaced the ordinary derivative of R n . Once geometry has been settled in M, the geometry of the shape space can be derived from its quotient structure. Let the tangent bundle be decomposed into a vertical and a horizontal subspace: TM = H M ⊕ V M , with V M = ker (T c π) and T c the tangent map, π : M → S the principal bundle, and H M = (V M ) ⊥ , see Figure 3. This metric is reparametrisation invariant, that is, constant along the fibers, hence we have where [g] denotes the metric on the shape space. A similar result in a different (but still close) context is used in Lemma 1 of [60]. In terms of distances, this can be understood in the following sense. The geodesic s → [c](s) between [c 0 ] and [c 1 ] in the shape space is the projection of the horizontal geodesic linking c 0 to the fiber containing c 1 . In fact, the horizontal geodesic between c 0 of c 1 intersects the fiber at c 1 at the reparametrized version of c 1 , c 1 • φ which gives the distance in the shape space: where [d] denotes the distance in S, and d g denotes the distance on the space of curves induced by the aforementioned Riemannian metric. In the TSRV formulation, the distance problem of Equation (37) yields an optimisation problem: which is solved by a traditional gradient descent algorithm or a dynamic linear programming [47]. Finally, we have to mention that in a practical situation, the above formula has to be discretized. This is the object of [46]. Formulae are essentially similar, but in this setting, a curve is now represented by a set of points c disc (x 0 , x 1 , · · · , x n ) and the tangent space turns into Concerning the metric on the space of curves, it becomes where, as before, for a u ∈ T c disc M, we define a path of piecewise geodesic curves (s, t) → c u (s, t) such that the following traditional initial conditions are fulfilled c u 0, k n = x k , and (∂c u /∂t) 0, k n = n log x k (x k+1 ). This is the discrete analogue of the tangent vector of a continuous curve at time t. The log function is the inverse of the exponential map on the base manifold, SO(n) for us, and here c u (s, ·) must be a geodesic on SO(n) between x k/n and x (k+1)/n . The SRV function that appears in the formula refer to the SRV function of the piecewise geodesics c u (s, ·). Then, the discretized version of the SRV function,

The Geodesic Equation
Let us now give the geodesic equation, relative to our chosen measure. As a result of the TSRV, the geodesic equation takes a much simpler form than what can be found in [45,46]. The formula can be found in [47]. For the sake of completeness, we give a reformulated proof in Appendix B. Recall that a geodesic is a particular path of curves. A path of curves is a continuous set of curve s → c(s, ·) such that for each s, c(s, ·) is a point in M, or, equivalently, a curve in M, (see Figure A1). Thus, for each curve of the path of curves, we can defined its TSRV function. Then for all s ∈ [0, 1], we have (we omit the letter 's' for clarity): q = T c(t)→I c (∂c/∂t) || ∂c/∂t || Theorem 3. A path of curves [0, 1] s → (c(s, 0), q(s, t)) (t is the parameter of the curve c(s, ·))is a geodesic on M if and only if ∇ ∂c/∂s (∇ ∂c/∂s q(s, t)) (s, t) = 0 ∀s, t Proof. Appendix B Thus, we have a quite familiar expression for the geodesic interpolation between two curves c 0 and c 1 , expressed in their TSRV domain: for s ∈ [0, 1]. This expression is nothing but a linear interpolation on the transported tangent space. We now have all the ingredients to give the procedure for nonstationary processes characterization and comparison: 1. Input: a set of rotation matrices {W i } i , seen as a partial observation of a closed trajectory on SO(n). 2. Map the set of matrices {W i } into the a set of matrices in the Lie algebra {V i } using the inverse exponential map. 3. Interpolate with splines between matrices V i [61,62]. 4. Go back in the base manifold SO(n) with the exponential map. 5. Shift the interpolated curve in order to fulfill the condition c(0) = e and compute the SRV transformation given by (41) 6. Compute the distance defined by (38). The optimization is carried out by dynamic programming. 7. Output: distance between two curves in the manifold defined by the set of curves in SO(n), and geodesic path between the curves.
The interpolation computations are carried out in the Lie algebra, which is a vector space, and thus it does not demand great computational resources. The discretization step, which amounts to choosing certain values among the continuous curve, is also done in the Lie algebra. In this way, we avoid the calculation of matrices exponential that would have been discarded at the end. We finally note that geodesic shooting [45,63] or other path straightening methods could be applied to obtain a geodesic path between two curves, and between the shapes of the two curves.

Results
In order to expose how the approach of this work gives interesting results for PC processes understanding, we propose to compare four PC processes, displayed along with Figure 4. We also bring their corresponding SO(3) representation on Figures 5 and 6 with 200 interpolated points and 50 interpolated points respectively. For this scenario we have generated four PC processes with 1000 samples each. A classical amplitude modulated model a(t) cos(2π f / f e t) where a(t) is a zero mean and unit variance stationary random process with a period of 20 points, a periodic AR(2) with a period of 20 points, a periodic AR(2) (AutoRegressive) with a period of 54 points, and a periodic ARMA(2,1) (AutoRegressive Moving Average) with a period of 20 points have been generated. We have used the R package PerARMA to generate the periodic ARMA and AR signals and we finally used the PerPACF function of this package to estimate the 20 (or 54) sequences of three parcors each. The analysis of Figure 4 with Figure 5 shows that the spectral measure of the amplitude modulated signal of Figure 4a has dilation matrices which do not spread a lot; we could think that this process is almost stationary due to the weak distance between each matrices. A contrario, whereas the temporal form of the PARMA(2,1) signal of Figure 4d is quite identical to the amplitude modulated signal of Figure 4a, their representation on SO (3) is very different. The spectral measure of the PARMA(2,1) signal spread much more. Lastly, when we observe the Figure 4b,c which are generated with the same model but with a different period, we can see that the more the number of points per period is important, the more the curve wraps. We also note one of the advantages of using spline interpolation. Along with Figure 5, we can remark that the curvature is well approximated owing to more and closer interpolated points. For this figure, we had approximatively four times more interpolated points than original ones whereas for Figure 6 we computed roughly twice more interpolated points than the original ones. Actually, for such a number, a change in the grid could be thought but it seems again that the spline interpolation gives good points repartition. This is mainly illustrated by Figure 7 for which the difference between geodesic interpolation of the curve associated with the periodic AR(2) signal of Figure 4c for 200 interpolated points and 50 interpolated points is very weak.  To end this analysis by the example, we have computed the distance defined by (38) between the PC process of Figure 8 and all the PC processes studied and displayed on Figures 4 and 5. The distances are reported inside Table 1. Clearly, the distances between the shapes of the curves characterizing the spectral measure of each PC process, reveal some spectral proximity between the PC processes benchmarked. We have gray colored the row of the PAR(2) signal model indexed by letter (c). Whatever the number of interpolated points and the dimension of the base manifold, the spectral representation through dilation operators of this signal is the nearest on SO(3) and the second nearest for SO(4) and SO(5) to the PAR(2) signal of reference. The second interesting signal is the one indexed by letter (b) and stands for a PAR(2) signal with exactly the same model parameters as that of the signal of reference but with 54 points of periodicity instead of 20. We have lightly gray colored its associated row on Table 1 when it had the shortest distance. The curve associated with this signal has many wraps on its representation, and it has consequently the greatest distance on SO(3), but increasing the dimension improves the comparison. Indeed, it finally has the shortest distance on SO(4) and SO (5). This is particularly interesting to see that there is a competition between the curves of a PAR(2) with different model parameters but the same periodicity and a PAR(2) with the same model parameters but different periodicity.  Figure 4 to the gold standard PC process of Figure 8 through the distance of their curves' shapes on SO(3), SO(4) and SO (5). We have interpolated with roughly twice and four times the number of original points. We also applied here a DP to solve the optimization assignment problem.    We end by noticing that the PARMA(2,1) has the second shortest distance to the PAR(2) signal model of reference on SO(3). As Figure 5 shows, their spectral measure evolves in a similar way with one major loop and a second less important loop. However, once the dimension of the base manifold increases, the assumption that the two processes may be close is strongly rejected by the fact that the PARMA(2,1) signal model has the longest distance. Finally, these observations leave open besides the question of the topology of these curves and how it could be used for the classification.

Conclusions
We have introduced a new vision of stochastic processes through geometry induced by dilation. The dilation matrices of given processes were obtained by a composition of rotations whose angles correspond to the well-known parcors, reflexion coefficients or Verblunski coefficients. The advantage of working with these particular matrices is that they are strongly related to the stochastic measure of the process, and thus, to its spectra. Furthermore, the dilation theory is independent of the stationarity of the underlying process; when the signal is stationary, its dilation operator is related to the Naimark dilation whereas when the signal is nonstationary, a set of dilation matrices is obtained and it is related to the Kolmogorov decomposition. Rigorously, dilation matrices are infinite dimensional, although we turn them into rotation matrices by truncation. Each of them belongs to the Special Orthogonal Group SO(n) or the Special Unitary Group SU(n) depending on the real-or complex-valued process under study. We focused our attention on the Periodically Correlated (PC) class of nonstationary processes for which a timely ordered set of dilation matrices describes the process measure. This set draws a closed curve on the Lie group of rotation matrices, and describing or classifying the different PC processes is made by curves comparison. We use for that the Square Root Velocity (SRV) function which represents a curve by its starting point and by its normed velocity vector on the space or curves. The metric in the space of curve naturally extends to the space of shapes. It is then possible to compare the shape of curves when the metric is translated into the Lie algebra, achieving therefore a closed-form expression and easy computation. Nonstationary processes are then characterized via their embedded curves.

Conflicts of Interest:
The authors declare no conflict of interest.

Appendix A. Defect Operator, Elementary Rotation
Introducing the defect operator of a contraction T as being D T = (I − T * T) 1/2 , we have the following factorisation: where X and Y are positive matrices. Note that this is a Cholesky factorisation-type result. This type of decomposition is used as the square root of matrices in the construction of the dilation. A corollary is that the operator I T T * I is positive if and only if T is a contraction.
Theorem A1. Let X and Y be operators in z. The following statements are equivalent: • There exists a contraction Γ in z such that X = ΓY, Proof. This result can be proved by taking the contraction Γ with respect to ΓXh = Yh. [41].
As a corollary If, X * X = Y * Y, then there exists a partial isometry V such that VX = Y. It is easy to see that we can choose V to be the contraction Γ defined above. Isometry V can also be assumed unitary. For a positive operator A ∈ L(H), if we denote by A 1/2 its unique positive square root, then every L such that L * L = A is related to A 1/2 by A 1/2 = VL (or A 1/2 = L * V * ).
Let us state another theorem that intervenes much in Constantinescu's factorisation of positive-definite kernel. Note that in the following, R(Γ) will denote the close range of the operator Γ. We start with a basic case: Theorem A2 (row contraction). Let T = [T 1 T 2 ] ∈ L(H 1 ⊕ H 2 , H), then || T || 0 if and only if there exists contractions Γ 1 ∈ L(H 1 , H) and Γ 2 ∈ L(H 2 , H) such that Proof. The proof is a simple application of Theorem A1. For the if part, it is obvious that we can take Γ 1 to be T 1 . Then || T || 1 implies Hence, there exists ∆ such that ∆D Γ * 1 = T * 2 . Choosing Γ 2 = ∆ * finishes the argument.
In the same way as that of the Cholesky factorisation, we can write down the defect operator for the whole contraction T = [T 1 T 2 ] [41] to be Therefore, with Theorem A2, we have an operator α such that Similarly, and the general case is: Theorem A3 (Structure of row contraction). The following are equivalent: • The operator T n = [T 1 T 2 · · · T n ] in L(⊕ n k=1 H k , H ) is a contraction • T 1 = Γ 1 is a contraction and, for k > 2, there exists uniquely determined contractions Γ k ∈ L(H k , R(γ k )) such that T k = D Γ * 1 D Γ * 2 · · · D Γ * k−1 Γ k .
Furthermore, the defect operators of the whole contraction T are of the form and D 2 T * = D Γ * 1 · · · D Γ * n D Γ * n · · · D Γ *  + R (q(s, t, τ), ∇ ∂ s q(s, t, τ)) (∂ s c(s, 0, τ)), ∂ τ c(s, 0, τ + ∇ ∂c/∂s ∇ ∂c/∂τ q(s, t, 0), ∇ ∂c/∂s q(s, t, 0) ds. (A15) Geodesic corresponds to minimal energy. It means that every other path that starts and ends at the same points should require more energy to travel than the geodesic. We then have to solve E (0) = 0 for every ∂ τ c(s, 0, τ) and every ∇ ∂ τ (q(s, t, τ)). This gives the result. Now when the framework is given by the TSRV and not by the SRV, only the second part of the geodesic equation remains as a result of the fixed starting point which corresponds to the identity element. This very much simplifies the equation, even though the derivation is the same.