Particle filters for infinite-dimensional systems: combining localization and optimal transportation Sebastian Reich University of Potsdam and University of Reading MCMSki 2014 Chamonix 7 January 2014 Sebastian Reich (UP and UoR) Particle filters for infinite-dimensional systems 7 January 2014 1 / 22
Talk Synopsis Introduction from a “weather prediction perspective” Ensemble prediction McKean data analysis cycle Linear ensemble transform filters (LETFs) Sequential importance resampling (SIR) particle filter Ensemble Kalman filter (EnKF) Ensemble transform particle filter (ETPF) Optimal transportation Example: Lorenz-63 Spatially extended systems Localization Example: Lorenz-96 Future work and references Sebastian Reich (UP and UoR) Particle filters for infinite-dimensional systems 7 January 2014 2 / 22
Introduction Ensemble prediction We consider a state space model in form of an iteration z n + 1 = Ψ( z n ) , n ≥ 0 , with z n ∈ R N z . The initial states z 0 are realizations of a random variable (RV) Z 0 : Ω → R N z with PDF π Z 0 . Hence all z n ’s become realizations of associated RVs Z n with marginal PDFs π Z n . Ensemble prediction relies on M independent realizations z 0 i = Z 0 ( ω i ) (MC or quasi-MC) from the initial RV and associated trajectories z n + 1 = Ψ( z n i ) , n ≥ 0 , i = 1 , . . . , M . i Sebastian Reich (UP and UoR) Particle filters for infinite-dimensional systems 7 January 2014 3 / 22
Introduction Ensemble prediction We consider a state space model in form of an iteration z n + 1 = Ψ( z n ) , n ≥ 0 , with z n ∈ R N z . The initial states z 0 are realizations of a random variable (RV) Z 0 : Ω → R N z with PDF π Z 0 . Hence all z n ’s become realizations of associated RVs Z n with marginal PDFs π Z n . Ensemble prediction relies on M independent realizations z 0 i = Z 0 ( ω i ) (MC or quasi-MC) from the initial RV and associated trajectories z n + 1 = Ψ( z n i ) , n ≥ 0 , i = 1 , . . . , M . i Sebastian Reich (UP and UoR) Particle filters for infinite-dimensional systems 7 January 2014 3 / 22
Introduction Ensemble prediction We consider a state space model in form of an iteration z n + 1 = Ψ( z n ) , n ≥ 0 , with z n ∈ R N z . The initial states z 0 are realizations of a random variable (RV) Z 0 : Ω → R N z with PDF π Z 0 . Hence all z n ’s become realizations of associated RVs Z n with marginal PDFs π Z n . Ensemble prediction relies on M independent realizations z 0 i = Z 0 ( ω i ) (MC or quasi-MC) from the initial RV and associated trajectories z n + 1 = Ψ( z n i ) , n ≥ 0 , i = 1 , . . . , M . i Sebastian Reich (UP and UoR) Particle filters for infinite-dimensional systems 7 January 2014 3 / 22
Introduction Ensemble prediction Ensemble prediction for Lorenz-63 model with initial distribution N ( 0 , 0 . 01 I ) and M = 1000 ensemble members. initial conditions solutions at t = 10 solutions at t = 100 50 50 50 40 40 40 30 30 30 20 20 20 10 10 10 0 0 0 � 10 � 10 � 10 40 40 40 20 20 20 0 0 0 � 20 � 20 � 20 20 20 20 0 0 0 � 40 � 40 � 40 � 20 � 20 � 20 Sebastian Reich (UP and UoR) Particle filters for infinite-dimensional systems 7 January 2014 4 / 22
Introduction McKean data analysis cycle We now aim to shadow or track an unknown reference solution z n + 1 = Ψ( z n ref ) , ref which is accessible to us only through partial and noisy observations y n obs = h ( z n ref ) + ξ n , n ≥ 1 . Here h : R N z → R N y is the forward operator and the ξ n ’s are realizations of independent Gaussian RVs with mean zero and covariance matrix R . Bayesian inference is used to merge ensemble predictions and partial observations in order to obtain an optimal shadowing of the reference trajectory. This is called the analysis step in weather forecasting. Sebastian Reich (UP and UoR) Particle filters for infinite-dimensional systems 7 January 2014 5 / 22
Introduction McKean data analysis cycle We now aim to shadow or track an unknown reference solution z n + 1 = Ψ( z n ref ) , ref which is accessible to us only through partial and noisy observations y n obs = h ( z n ref ) + ξ n , n ≥ 1 . Here h : R N z → R N y is the forward operator and the ξ n ’s are realizations of independent Gaussian RVs with mean zero and covariance matrix R . Bayesian inference is used to merge ensemble predictions and partial observations in order to obtain an optimal shadowing of the reference trajectory. This is called the analysis step in weather forecasting. Sebastian Reich (UP and UoR) Particle filters for infinite-dimensional systems 7 January 2014 5 / 22
Introduction McKean data analysis cycle Summary of the McKean approach to the analysis step : PDFs RVs MC Ref.: Del Moral (2004), CJC & SR (2013), YC & SR (2014). Sebastian Reich (UP and UoR) Particle filters for infinite-dimensional systems 7 January 2014 6 / 22
Introduction McKean data analysis cycle (A) Fit a Gaussian N (¯ z f , P f ) to the forecast ensemble { z f i } and assume that the forward operator can be linearized. Then the analysis is also Gaussian N (¯ z a , P a ) with z a = ¯ z f − K ( H ¯ z f − y obs ) , P a = P f − KHP f . ¯ Here K denotes the Kalman gain matrix . (B) Use the empirical measure M π f ( z ) = 1 � δ ( z − z f i ) M i = 1 to define the analysis measure M π a ( z ) = � w i δ ( z − z f i ) i = 1 with importance weights − 1 � 2 ( h ( z f i ) − y obs ) T R − 1 ( h ( z f � exp i ) − y obs ) w i = � � � M − 1 2 ( h ( z f j ) − y obs ) T R − 1 ( h ( z f j = 1 exp j ) − y obs ) Sebastian Reich (UP and UoR) Particle filters for infinite-dimensional systems 7 January 2014 7 / 22
Introduction McKean data analysis cycle (A) Fit a Gaussian N (¯ z f , P f ) to the forecast ensemble { z f i } and assume that the forward operator can be linearized. Then the analysis is also Gaussian N (¯ z a , P a ) with z a = ¯ z f − K ( H ¯ z f − y obs ) , P a = P f − KHP f . ¯ Here K denotes the Kalman gain matrix . (B) Use the empirical measure M π f ( z ) = 1 � δ ( z − z f i ) M i = 1 to define the analysis measure M π a ( z ) = � w i δ ( z − z f i ) i = 1 with importance weights − 1 � 2 ( h ( z f i ) − y obs ) T R − 1 ( h ( z f � exp i ) − y obs ) w i = � � � M − 1 2 ( h ( z f j ) − y obs ) T R − 1 ( h ( z f j = 1 exp j ) − y obs ) Sebastian Reich (UP and UoR) Particle filters for infinite-dimensional systems 7 January 2014 7 / 22
Introduction McKean data analysis cycle Implementation of the McKean approach then either requires coupling two Gaussians (approach A) or two empirical measures (approach B). Approach A leads to the family of ensemble Kalman filters (Evensen, 2006), while Approach B leads to particle filters (Doucet et al, 2001). Optimal couplings in the sense of minimizing some cost function are known in both cases (CJC & SR, 2013). We now provide a unifying mathematical framework in form of linear ensemble transform filters (LETFs) (YC & SR, 2014) instead of focusing on the McKean Markov process aspect. Sebastian Reich (UP and UoR) Particle filters for infinite-dimensional systems 7 January 2014 8 / 22
Linear ensemble transform filters The analysis steps of a standard particle filter as well as an ensemble Kalman filter (EnKF) are of the form M � z a z f j = i s ij , i = 1 where { z f i } M i = 1 is the forecast ensemble and { z a i } M i = 1 is the analysis ensemble. The matrix S ∈ R M × M with entries s ij depends on y obs and the forecast ensemble. In some cases, S is the realization of a matrix-valued RV S : Ω → R M × M , i.e. S = S ( ω ) . Sebastian Reich (UP and UoR) Particle filters for infinite-dimensional systems 7 January 2014 9 / 22
Linear ensemble transform filters SIR particle filter Sequential importance resampling (SIR) particle filters rely on S : Ω → R M × M satisfying the following properties s ij ∈ { 0 , 1 } , � M i = 1 s ij = 1, � � � M 1 = w i . E j = 1 s ij M One often adds a rejuvenation step : M � z a z f i s ij + ξ a j = j . i = 1 The ξ a j ’s are realizations of M independent Gaussian RVs with mean zero and appropriate covariance matrix. Sebastian Reich (UP and UoR) Particle filters for infinite-dimensional systems 7 January 2014 10 / 22
Linear ensemble transform filters SIR particle filter Sequential importance resampling (SIR) particle filters rely on S : Ω → R M × M satisfying the following properties s ij ∈ { 0 , 1 } , � M i = 1 s ij = 1, � � � M 1 = w i . E j = 1 s ij M One often adds a rejuvenation step : M � z a z f i s ij + ξ a j = j . i = 1 The ξ a j ’s are realizations of M independent Gaussian RVs with mean zero and appropriate covariance matrix. Sebastian Reich (UP and UoR) Particle filters for infinite-dimensional systems 7 January 2014 10 / 22
Linear ensemble transform filters Ensemble Kalman filter (EnKF) The EnKF with perturbed observations leads to 1 M − 1 ( z f z f ) T H T ( H T P f H + R ) − 1 ( y obs − Hz f i − ¯ s ij = δ ij + j + η j ) with η j a realization of N ( 0 , R ) . The coefficients satisfy � M i = 1 s ij = 1, z a = � M w i z f ¯ j = 1 ˜ i with “importance weights” M w i = 1 � ˜ s ij , M j = 1 The coefficients s ij ’s can take negative values ! Ensemble inflation or rejuvenation is often applied. Sebastian Reich (UP and UoR) Particle filters for infinite-dimensional systems 7 January 2014 11 / 22
Recommend
More recommend
