Numerical Linear Algebra issues in Singular Spectrum Analysis of - PowerPoint PPT Presentation

Numerical Linear Algebra issues in Singular Spectrum Analysis of time series Dario Fasino — with E. Bozzo, R. Carniel University of Udine (Italy) Genova, 2ggALN’12 D. Fasino (Udine) NLA issues in SSA 17/02/2012 1 / 18

Singular Spectrum Analysis (SSA): Introduction SSA is a quite recent technique for the analysis of experimental time series, based on the SVD of certain Hankel matrices. Let x = ( x 1 , x 2 , . . . , x ℓ ) T a finite time series, ℓ = m + n − 1. The m × n Hankel matrix   x 1 x 2 · · · x n   x 2 x 3 · · · x n + 1   X m , n =   . . . . . . . .   . . . . · · · x m x m + 1 x ℓ is known as trajectory matrix, X m , n = T m , n ( x ) . D. Fasino (Udine) NLA issues in SSA 17/02/2012 2 / 18

Singular Spectrum Analysis (SSA): Introduction SSA is a quite recent technique for the analysis of experimental time series, based on the SVD of certain Hankel matrices. Let x i = p k ( i ) , Let z ∈ C . Then   a k -degree algebraic poly. Then z n − 1 1 · · ·     z n z · · · x 1 · · · x n    = 1 . rank   . . .   . . . . . . . . .   = k . rank  . . . . . . z m − 1 z ℓ − 2 · · · x m · · · x ℓ Time series made up by trigonometric, algebraic, exponential terms have small rank trajectory matrices. D. Fasino (Udine) NLA issues in SSA 17/02/2012 2 / 18

Singular Spectrum Analysis (SSA): Introduction SSA idea Use SVD of trajectory matrices to decompose time series into constant terms, trends, oscillatory components, noise. . . Given x = ( x 1 , . . . , x ℓ ) and I 1 , . . . , I k a partition of { 1 , . . . , n } do: Build up trajectory matrix: X = T m , n ( x ) . 1 Compute SVD X = U Σ V T ; singular triples: ( u i , σ i , v i ) . 2 Group triples: X ( k ) = � i ∈I k σ i u i v T i . 3 Note: � k X ( k ) = X . Hankelization (diagonal averaging): H ( k ) = H ( X ( k ) ) . 4 Note: � k H ( k ) = X . Extract components: x ( k ) = T − 1 m , n ( H ( k ) ) . 5 Note: � k x ( k ) = x . D. Fasino (Udine) NLA issues in SSA 17/02/2012 3 / 18

Example Figure: SSA of a mixture of time series. ℓ = 200, n = 10 , . . . , 30. Above: individual time series and respective SVs. Below: composite time series and respective SVs. D. Fasino (Udine) NLA issues in SSA 17/02/2012 4 / 18

Example Figure: SSA of a mixture of time series. ℓ = 200, n = 30. Left: first original component (blue) and its reconstruction (red). Right: second original component (blue) and its reconstruction (red). D. Fasino (Udine) NLA issues in SSA 17/02/2012 5 / 18

References N. Golyandina, V. Nekrutkin, A. Zhigljavsky. Analysis of time series structure. SSA and related techniques. Chapman & Hall/CRC, 2001. R. Carniel et al . On the singular values decoupling in the Singular Spectrum Analysis of volcanic tremor at Stromboli. Nat. Hazards Earth Syst. Sci. , 6 (2006), 903–909. E. Bozzo, R. Carniel, D. F . Relationship between SSA and Fourier analysis: Theory and application to the monitoring of volcanic activity. Comp. Math. Appl. 60 (2010), 812–820. V. Busoni. Risultati di tipo perturbativo nell’analisi dello spettro singolare di serie temporali. Tesi di Laurea in Matematica, Universit` a di Udine, 2011. D. Fasino (Udine) NLA issues in SSA 17/02/2012 6 / 18

A motivating problem (and our answer) April 5, 2003: a major paroxism destroyed a seismic station on the Stromboli volcano. SSA analysis of the sismogram suggests the presence of a consistent preparatory phase. What information is conveyed in the SVs of (not so large) Figure: SSA of a volcanic tremor trajectory matrices coming from sismogram. chaotic time series? Singular values of 7200 trajectory matrices X 3000 , 10 (smoothed plot; x-axis in hours) D. Fasino (Udine) NLA issues in SSA 17/02/2012 7 / 18

A motivating problem (and our answer) Our result: The behaviour of SVs mirrors qualitative modifications in the power spectrum of the time series. Figure: SSA of a volcanic tremor sismogram. Singular values of 7200 trajectory matrices X 3000 , 10 (smoothed plot; x-axis in hours) D. Fasino (Udine) NLA issues in SSA 17/02/2012 7 / 18

Stationary time series Definition The infinite time series � x = ( x 1 , x 2 , . . . ) is called stationary if m � 1 ∀ i , j � 0 lim x i + k x j + k = R ( i − j ) , m m →∞ k = 1 where R : Z �→ R is the covariance function of x . Equivalently, the covariance matrix 1 mX T T n = lim m , n X m , n m →∞ exists for all n (and is a Toeplitz matrix). D. Fasino (Udine) NLA issues in SSA 17/02/2012 8 / 18

Stationary time series Definition The infinite time series � x = ( x 1 , x 2 , . . . ) is called stationary if m � 1 ∀ i , j � 0 lim x i + k x j + k = R ( i − j ) , m m →∞ k = 1 where R : Z �→ R is the covariance function of x . By Herglotz theorem, � x is a stationary time series iff there exists a (unique, bounded) nondecreasing function µ ( t ) on I = [ 0 , 2 π ] such that � R ( k ) = 1 e i 2 π kt d µ ( t ) , k ∈ Z . 2 π I D. Fasino (Udine) NLA issues in SSA 17/02/2012 8 / 18

Stationary time series Definition The infinite time series � x = ( x 1 , x 2 , . . . ) is called stationary if m � 1 ∀ i , j � 0 lim x i + k x j + k = R ( i − j ) , m m →∞ k = 1 where R : Z �→ R is the covariance function of x . A special case: � x is called aperiodic or chaotic whenever � R ( k ) = 1 e i 2 π kt f ( t ) d t , k ∈ Z . 2 π I The function f ∈ L 1 ( I ) , f � 0, is the spectral density of � x and the symbol of the Toeplitz matrix sequence { T n } . D. Fasino (Udine) NLA issues in SSA 17/02/2012 8 / 18

Asymptotic distributions Definition Two triangular sequences { ξ ( n ) } i = 1 ... n and { ζ ( n ) } i = 1 ... n , with n ∈ N , i i are equally distributed (or asymptotically equidistributed ), ξ ( n ) ∼ ζ ( n ) i i if, for all continuous functions F having bounded support, � �� n � � � � 1 ξ ( n ) ζ ( n ) lim F − F = 0 . i i n n →∞ i = 1 Theorem Let { T n } be a sequence of Toeplitz matrices whose symbol f ∈ L 1 ( I ) is Riemann-integrable. Then λ i ( T n ) ∼ f ( 2 π i / n ) . D. Fasino (Udine) NLA issues in SSA 17/02/2012 9 / 18

Asymptotic distributions Definition Two triangular sequences { ξ ( n ) } i = 1 ... n and { ζ ( n ) } i = 1 ... n , with n ∈ N , i i are equally distributed (or asymptotically equidistributed ), ξ ( n ) ∼ ζ ( n ) i i if, for all continuous functions F having bounded support, � �� n � � � � 1 ξ ( n ) ζ ( n ) lim F − F = 0 . i i n n →∞ i = 1 Theorem [Tyrtyshnikov ’96] Let A n , B n be m ( n ) × n matrices, with m ( n ) ≥ n . If 1 n � A n − B n � 2 lim F = 0 = ⇒ σ i ( A n ) ∼ σ i ( B n ) . n →∞ D. Fasino (Udine) NLA issues in SSA 17/02/2012 9 / 18

SSA and Fourier analysis Let � x = ( x 1 , x 2 . . . ) be a stationary time series,   x 1 x 2 · · · x n   . . . . . .   x 2 . . .   1 1   . . . � � . . . X m , n = √ mX m , n X m , n = √ m .   . . . x 1    . .  . .   . x m · · · . · · · x m x 1 x n − 1 Lemma For any integer sequence m ( n ) such that n / m ( n ) → 0 as n → ∞ , the singular values of � X m ( n ) , n and � X m ( n ) , n are equally distributed. D. Fasino (Udine) NLA issues in SSA 17/02/2012 10 / 18

SSA and Fourier analysis Indeed: n � 1 1 n − k n � � X m ( n ) , n − � X m ( n ) , n � 2 ( x k − x m ( n )+ k ) 2 = F m ( n ) n k = 1 n � 2 x 2 k + x 2 < m ( n )+ k m ( n ) k = 1 � 1 � n n � � 2 n + 1 x 2 x 2 = → 0 . k m ( n )+ k m ( n ) n n � �� � k = 1 k = 1 � �� � � �� � → 0 → R ( 0 ) → R ( 0 ) Note: � x stationary ⇒ any trailing subsequence is stationary. D. Fasino (Udine) NLA issues in SSA 17/02/2012 11 / 18

SSA and Fourier analysis For any fixed window length m , let x = ( x 1 , . . . , x m ) T , x m − 1 ) T = F m x be its Fourier transform. and let ˆ x = (ˆ x 0 , . . . , ˆ Let Λ m = Diag ( 1 , e i 2 π/ m , . . . , e i 2 ( m − 1 ) π/ m ) and P m = F H m Λ m F m . The matrix T m , n = � m , n � X T X m , n is Toeplitz, with symbol n − 1 n − 1 � � f m , n ( z ) = 1 m x e i kz = 1 x T P − k x H Λ − k x e i kz . ˆ m ˆ m m k = − n + 1 k = − n + 1 Due to the stationarity assumption, lim m →∞ T m , n = T n , hence for “fastly growing” m ( n ) we have lim n →∞ 1 n � T m ( n ) , n − T n � 2 F = 0. X m ( n ) , n ) 2 ∼ σ i ( � X m ( n ) , n ) 2 = λ i ( T m ( n ) , n ) ∼ λ i ( T n ) ∼ f ( 2 π i / n ) . σ i ( � D. Fasino (Udine) NLA issues in SSA 17/02/2012 12 / 18