Time Series Analysis Henrik Madsen hm@imm.dtu.dk Informatics and - PowerPoint PPT Presentation

H. Madsen, Time Series Analysis, Chapmann Hall Time Series Analysis Henrik Madsen hm@imm.dtu.dk Informatics and Mathematical Modelling Technical University of Denmark DK-2800 Kgs. Lyngby Henrik Madsen 1

H. Madsen, Time Series Analysis, Chapmann Hall Outline of the lecture State space models, 1st part: Model: Sec. 10.1 The Kalman filter: Sec. 10.3 Cursory material: Sec. 10.3.2 (Empirical-Bayesian description) Henrik Madsen 2

H. Madsen, Time Series Analysis, Chapmann Hall State space models System model; A full description of the dynamical system (i.e. including the parameters) Observations; Noisy measurements on some parts (states) of the system Goal; reconstruct and predict the state of the system System Output Input State: X t u Y t t Henrik Madsen 3

H. Madsen, Time Series Analysis, Chapmann Hall State space models; examples Estimate the temperature inside a solid block of material when we measure the temperature on the surface (with noise) Noisy measurements of the position of a ship; give a better estimate of the current position A model of a car engine: Input; fuel. State; Fuel and temperature in various parts. Observarions: Sensor output PK/PD-modeling: State: Amount of drug in blood, liver, muscules, . . . Observations: Amount in blood (with noise), Input: Drug. Henrik Madsen 4

H. Madsen, Time Series Analysis, Chapmann Hall Determining the model structure The system model is often based on physical considerations; this often leads to dynamical models consisting of differential equations An m ’th order differential equation can be formulated as m 1st order differential equations Sampling such a system leads to a linear state space model and there exist a way of coming from the coefficients in continuous time to the coefficients in discrete time Henrik Madsen 5

H. Madsen, Time Series Analysis, Chapmann Hall The linear stochastic state space model System equation: = AX t − 1 + Bu t − 1 + e 1 ,t X t Observation equation: = CX t + e 2 ,t Y t dim( X t ) = m is called X : State vector the order of the system Y : Observation vector { e 1 ,t } and { e 2 ,t } mutually u : Input vector independent white noise e 1 : System noise V [ e 1 ] = Σ 1 , V [ e 2 ] = Σ 2 e 2 : Observation noise A , B , C , Σ 1 , and Σ 2 are known matrices The state vector contains all information available for future evaluation; the state vector is a Markov process Henrik Madsen 6

H. Madsen, Time Series Analysis, Chapmann Hall Example – a falling body Height above ground: z ( t ) Initial conditions: Position z ( t 0 ) and velocity z ′ ( t 0 ) Physical considerations: d 2 z dt 2 = − g States: Position x 1 ( t ) = z ( t ) and velocity x 2 ( t ) = z ′ ( t ) Only the position is measured y ( t ) = x 1 ( t ) Continuous time description x ( t ) = [ x 1 ( t ) x 2 ( t )] T : � 0 � � � 1 0 x ′ ( t ) = x ( t ) + g 0 0 − 1 � � y ( t ) = 1 0 x ( t ) Henrik Madsen 7

H. Madsen, Time Series Analysis, Chapmann Hall Example – a falling body (cont’nd) Solving the equations: x 2 ( t ) = − g ( t − t 0 ) + x 2 ( t 0 ) − g 2( t − t 0 ) 2 + ( t − t 0 ) x 2 ( t 0 ) + x 1 ( t 0 ) x 1 ( t ) = Sampling: t = sT , t 0 = ( s − 1) T , and T = 1 � 1 � � − 1 / 2 � 1 = x s − 1 + g x s 0 1 − 1 � � y s = 1 0 x s Adding disturbances and measurement noise: � 1 � � − 1 / 2 � 1 = x s − 1 + g + e 1 ,s x s 0 1 − 1 � � = 1 0 x s + e 2 ,s y s Henrik Madsen 8

H. Madsen, Time Series Analysis, Chapmann Hall Example – a falling body (cont’nd) Given measurements of the position at time points 1 , 2 , . . . , s we could: Predict the future position and velocity x s + k | s ( k > 0 ) Reconstruct the current position and velocity from noisy measurements x s | s Interpolate to find the best estimate of the position and velocity at a previous time point x s + k | s ( k < 0 ) (estimate the path in the state space; vary k so that s + k varied from 1 to s ) we will focus on reconstruction and prediction Henrik Madsen 9

H. Madsen, Time Series Analysis, Chapmann Hall Requirement In order to predict, reconstruct or interpolate the m -dimensional state in the system = AX t − 1 + Bu t − 1 + e 1 ,t X t = CX t + e 2 ,t Y t the system must be observable, i.e. � � CA m − 1 � T � C T . . ( CA ) T . . · · · . . . . rank . = m. For the falling body (S-PLUS): > qr( cbind(t(C), t(C %*% A)) )$rank [1] 2 Where A and C is taken from the discrete time description of the system. Henrik Madsen 10

H. Madsen, Time Series Analysis, Chapmann Hall The Kalman filter Initialization: � X 1 | 0 = E [ X 1 ] = µ 0 , Σ xx 1 | 0 = V [ X 1 ] = V 0 , and 1 | 0 C T + Σ 2 thereby Σ yy 1 | 0 = C Σ xx For: t = 1 , 2 , 3 , . . . t | t − 1 C T � � − 1 Σ yy Σ xx = K t t | t − 1 � � Reconstruction: � � Y t − C � = X t | t − 1 + K t X t | t X t | t − 1 t | t − 1 − K t Σ yy Σ xx Σ xx t | t − 1 K T = t | t t � A � = X t | t + Bu t X t +1 | t t | t A T + Σ 1 Σ xx A Σ xx Prediction: = t +1 | t t +1 | t C T + Σ 2 Σ yy C Σ xx = t +1 | t Henrik Madsen 11

H. Madsen, Time Series Analysis, Chapmann Hall Multi step predictions Not part of the Kalman filter as stated above Can be calculated recursively for a given t starting with k = 1 for which � X t + k | t and Σ t + k | t are calculated as part of the Kalman filter � A � = X t + k | t + Bu t + k X t + k +1 | t t + k | t A T + Σ 1 Σ xx A Σ xx = t + k +1 | t The future input must be decided Henrik Madsen 12

H. Madsen, Time Series Analysis, Chapmann Hall Naming and history The filter is named after Rudolf E. Kalman, though Thorvald Nicolai Thiele and Peter Swerling actually developed a similar algorithm earlier. It was during a visit of Kalman to the NASA Ames Research Center that he saw the applicability of his ideas to the problem of trajectory estimation for the Apollo program, leading to its incorporation in the Apollo navigation computer. From http://en.wikipedia.org/wiki/Kalman_filter Henrik Madsen 13

H. Madsen, Time Series Analysis, Chapmann Hall The Foundation of the Kalman filter Theorem 2.6 (Linear projection) The theorem is concerned with the random vectors X and Y for which the means, variances and covariances are used The state is called X t and the observation is called Y t and we could write down the theorem for these t − 1 = ( Y T 1 , . . . , Y T We have additional information; Y T t − 1 ) We include this information by considering the random vectors X t |Y t − 1 and Y t |Y t − 1 instead E [( X t |Y t − 1 ) | ( Y t |Y t − 1 )] = E [ X t | Y t , Y t − 1 ] = E [ X t |Y t − 1 ] + C [ X t , Y t |Y t − 1 ] V − 1 [ Y t |Y t − 1 ]( Y t − E [ Y t |Y t − 1 ]) V [( X t |Y t − 1 ) | ( Y t |Y t − 1 )] = V [ X t | Y t , Y t − 1 ] = V [ X t |Y t − 1 ] − C [ X t , Y t |Y t − 1 ] V − 1 [ Y t |Y t − 1 ] C T [ X t , Y t |Y t − 1 ] Henrik Madsen 14

H. Madsen, Time Series Analysis, Chapmann Hall The Foundation of the Kalman filter (cont’nd) E [ X t | Y t , Y t − 1 ] = E [ X t |Y t − 1 ] + C [ X t , Y t |Y t − 1 ] V − 1 [ Y t |Y t − 1 ]( Y t − E [ Y t |Y t − 1 ]) V [ X t | Y t , Y t − 1 ] = V [ X t |Y t − 1 ] − C [ X t , Y t |Y t − 1 ] V − 1 [ Y t |Y t − 1 ] C T [ X t , Y t |Y t − 1 ] � � − 1 � � � X t | t − 1 + Σ xy � Σ yy Y t − � X t | t = Y t | t − 1 t | t − 1 t | t − 1 � � − 1 � � T t | t − 1 − Σ xy Σ yy Σ xy Σ xx Σ xx = t | t t | t − 1 t | t − 1 t | t − 1 � � − 1 K t = Σ xy Σ yy t | t − 1 t | t − 1 K t is called the Kalman gain , because it determine how much the 1-step prediction error influence the update of the state estimate Henrik Madsen 15

H. Madsen, Time Series Analysis, Chapmann Hall The Foundation of the Kalman filter (cont’nd) The 1-step predictions are obtained directly from the state space model: � A � = X t | t + Bu t X t +1 | t � C � = Y t +1 | t X t +1 | t Which results in the prediction errors: � X t +1 − � X t +1 | t = A � X t +1 | t = X t | t + e 1 ,t +1 � Y t +1 − � Y t +1 | t = C � = X t +1 | t + e 2 ,t +1 Y t +1 | t t | t A T + Σ 1 And in therefore: Σ xx A Σ xx = t +1 | t t +1 | t C T + Σ 2 Σ yy C Σ xx = t +1 | t Σ xy t +1 | t can also be calculated Henrik Madsen 16

H. Madsen, Time Series Analysis, Chapmann Hall Kalman filter applied to a falling body Description of the system: � 1 � � − 1 / 2 � � � 1 A = B = C = 1 0 0 1 − 1 � 2 . 0 � � � 0 . 8 Σ 1 = Σ 2 = 10000 0 . 8 1 . 0 Initialization: Released 10000 m above ground at 0 m/s � 10000 � � 0 � � � 0 � Σ yy Σ xx X 1 | 0 = 1 | 0 = 1 | 0 = 10000 0 0 0 Henrik Madsen 17

H. Madsen, Time Series Analysis, Chapmann Hall Kalman filter applied to a falling body (cont’nd) 1st observation ( t = 1 ): y 1 = 10171 � � T Reconstruction: K 1 = 0 0 � 10000 � � 0 � 0 � Σ xx X 1 | 1 = 1 | 1 = 0 0 0 Prediction: � 9995 . 09 � � � � � 2 0 . 8 � Σ yy Σ xx X 2 | 1 = 2 | 1 = 2 | 1 = 10002 − 9 . 82 0 . 8 1 Henrik Madsen 18