Joint state and parameter estimation with an iterative ensemble - - PowerPoint PPT Presentation

joint state and parameter estimation with an iterative
SMART_READER_LITE
LIVE PREVIEW

Joint state and parameter estimation with an iterative ensemble - - PowerPoint PPT Presentation

Joint state and parameter estimation with an iterative ensemble Kalman smoother Marc Bocquet 1 , 2 , Pavel Sakov 3 1 Universit e Paris-Est, CEREA, joint lab Ecole des Ponts ParisTech and EdF R&D, France 2 INRIA, Paris-Rocquencourt


slide-1
SLIDE 1

Joint state and parameter estimation with an iterative ensemble Kalman smoother

Marc Bocquet 1,2, Pavel Sakov 3

1Universit´

e Paris-Est, CEREA, joint lab ´ Ecole des Ponts ParisTech and EdF R&D, France

2INRIA, Paris-Rocquencourt Research center, France 3Bureau of Meteorology, Australia

(bocquet@cerea.enpc.fr)

  • M. Bocquet

9th EnKF workshop, Bergen, Norway, 22-24 June 2014 1 / 22

slide-2
SLIDE 2

Context

◮ New methods called ensemble variational methods that mix variational and ensemble approaches (see Lorenc, 2013 for an almost perfect definition): Hybrid methods, 4D-Var-Ben, 4D-En-Var, Ensemble of data assimilation (EDA) and IEnKF/IEnKS.

Lorenc A. 2013. Recommended nomenclature for EnVar data assimilation methods. In Research Activities in Atmospheric and Oceanic Modelling , WGNE.

◮ The IEnKF/IEnKS differ from the other ones in that they are more natural (simple?), regardless of the numerical cost. ◮ The IEnKS has a great potential for parameter estimation, as it is variational but avoids the derivation of the adjoint.

  • M. Bocquet

9th EnKF workshop, Bergen, Norway, 22-24 June 2014 2 / 22

slide-3
SLIDE 3

Context

The IEnKS: at the crossroad between the EnKF and 4D-Var

◮ The IEnKS follows the scheme of the EnKF: Analysis in ensemble space → Posterior ensemble generation → Ensemble forecast ◮ Except that The analysis in ensemble space is variational [e.g. Zupanski, 2005] over a finite time

  • windows. It may require several iterations in strongly nonlinear conditions [Gu &

Oliver, 2007; Sakov et al., 2012; Bocquet and Sakov, 2012-2014].

The gradient of the 4D cost function is computed with the ensemble [Gu & Oliver,

2007;Liu et al., 2008]: no need for the tangent linear/adjoint.

◮ It generalises the iterative extended Kalman filter/smoother [Wishner et al., 1969;

Jazwinski, 1970; Bell, 1994] to ensemble methods.

◮ It is a unified/straightforward scheme (no hybridization so to speak).

  • M. Bocquet

9th EnKF workshop, Bergen, Norway, 22-24 June 2014 3 / 22

slide-4
SLIDE 4

The iterative ensemble Kalman smoother

The IEnKS: the cycling

◮ L: length of the data assimilation window, ◮ S: shift of the data assimilation window in between two updates. S∆ S∆ t t t

1

t t t t t

L L−1 L+1

yL ∆ L t t t t

L+1 L+2 L−2 L−3

t tL−1

L

y y y

L+2 L+1 L−2 L−3

y yL−1

  • M. Bocquet

9th EnKF workshop, Bergen, Norway, 22-24 June 2014 4 / 22

slide-5
SLIDE 5

The iterative ensemble Kalman smoother

The IEnKS: a variational standpoint

◮ Analysis IEnKS cost function in state space p(x0|yL) ∝ exp(−J (x0)): J (x0) =

L

k=1

1 2 (yk −Hk ◦Mk←0(x0))T βkR−1

k (yk −Hk ◦Mk←0(x0))

+ 1 2 (x0 −x0)P−1

0 (x0 −x0) .

(1) {β0,β1,...,βL} weight the observations impact within the window. ◮ Reduced scheme in ensemble space, x0 = x0 +A0w, where A0 is the ensemble anomaly matrix:

  • J (w) = J (x0 +A0w).

(2) ◮ IEnKS cost function in ensemble space [Hunt et al., 2007; Bocquet and Sakov, 2012]:

  • J (w) =1

2

L

k=1

(yk −Hk ◦Mk←0 (x0 +A0w))T βkR−1

k (yk −Hk ◦Mk←0 (x0 +A0w))

+ 1 2(N −1)wTw. (3)

  • M. Bocquet

9th EnKF workshop, Bergen, Norway, 22-24 June 2014 5 / 22

slide-6
SLIDE 6

The iterative ensemble Kalman smoother

The IEnKS: minimisation scheme

◮ As a variational reduced method, one can use Gauss-Newton [Sakov et al., 2012], Levenberg-Marquardt [Bocquet and Sakov, 2012; Chen and Oliver, 2013], quasi-Newton, etc., minimisation schemes. ◮ Gauss-Newton scheme (the Hessian is approximate): w(p+1) = w(p) − H −1

(p) ∇

J(p)(w(p)), x(p) = x(0) +A0w(p) , ∇ J(p) = −

L

k=1

YT

k,(p)βkR−1 k

  • yk −Hk ◦Mk←0(x(p)

0 )

  • +(N −1)w(p) ,
  • H(p) = (N −1)IN +

L

k=1

YT

k,(p)βkR−1 L Y(p) ,

Yk,(p) = [Hk ◦Mk←0]′

|x(p)

0 A0 .

(4) ◮ One solution to compute the 4D sensitivities: the bundle scheme. It simply mimics the action of the tangent linear by finite difference: Yk,(p) ≈ 1 ε Hk ◦Mk←0

  • x(p)1T +εA0
  • IN − 11T

N

  • .

(5)

  • M. Bocquet

9th EnKF workshop, Bergen, Norway, 22-24 June 2014 6 / 22

slide-7
SLIDE 7

The iterative ensemble Kalman smoother

The IEnKS: ensemble update and the forecast step

◮ Generate an updated ensemble from the previous analysis: E⋆

0 = x⋆ 01T +

√ N −1A0 H −1/2

U where U1 = 1. (6) ◮ Just propagate the updated ensemble from t0 to tS: ES = MS←0(E0). (7) In the quasi-static case: S = 1.

  • M. Bocquet

9th EnKF workshop, Bergen, Norway, 22-24 June 2014 7 / 22

slide-8
SLIDE 8

The iterative ensemble Kalman smoother

IEnKS: single vs multiple data assimilation

S∆ S∆ t t t

1

t t t t t

L L−1 L+1

yL ∆ L t t t t

L+1 L+2 L−2 L−3

t tL−1

L

y y y

L+2 L+1 L−2 L−3

y yL−1 y y y

−2 2 β β β β β β β β

L L

β0

L−1 L−1 L−1 L

◮ SDA IEnKS: The observation vector are assimilated once and for all. Exact scheme. ◮ MDA IEnKS: The observation vector are assimilated several times and poundered with weights βk within the window. Exact scheme in the linear/Gaussian limit. More stable for long windows than the SDA scheme.

  • M. Bocquet

9th EnKF workshop, Bergen, Norway, 22-24 June 2014 8 / 22

slide-9
SLIDE 9

Numerical experiments

Application to the Lorenz-95 model

◮ Weakly nonlinear case: Lorenz-95, M = 40, N = 20, ∆t = 0.05, R = I. ◮ Comparison of 4D-Var S = 1, EnKS S = 1, SDA IEnKS S = 1, SDA IEnKS S = L, and MDA IEnKS S = 1.

1 5 10 15 20 25 30 35 40 45 50

DAW length L

0.160 0.165 0.170 0.175 0.180 0.185 0.190 0.155 0.200 0.195 0.210 0.220

Filtering analysis RMSE

4D-Var S=1 EnKS-N S=1 SDA IEnKS-N S=1 SDA IEnKS-N S=L MDA IEnKS-N S=1 1 5 10 15 20 25 30 35 40 45 50

DAW length L

0.04 0.05 0.06 0.07 0.08 0.09 0.10 0.12 0.14 0.16 0.18 0.20 0.22

Smoothing analysis RMSE

4D-Var S=1 EnKS-N S=1 SDA IEnKS-N S=1 SDA IEnKS-N S=L MDA IEnKS-N S=1

  • M. Bocquet

9th EnKF workshop, Bergen, Norway, 22-24 June 2014 9 / 22

slide-10
SLIDE 10

Numerical experiments

Application to the Lorenz-95 model

◮ Strongly nonlinear case: Lorenz-95, M = 40, N = 20, ∆t = 0.20, R = I. ◮ Comparison of 4D-Var S = 1, EnKS S = 1, SDA IEnKS S = 1, SDA IEnKS S = L, and MDA IEnKS S = 1.

1 2 3 4 5 6 7 8 9 10

DAW length L

0.30 0.32 0.34 0.28 0.36 0.29 0.31 0.33 0.35 0.37 0.38 0.39 0.40 0.42 0.44 0.46 0.27

Filtering analysis RMSE

4D-Var S=1 EnKS-N S=1 SDA IEnKS-N S=1 SDA IEnKS-N S=L MDA IEnKS-N S=1 1 2 3 4 5 6 7 8 9 10

DAW length L

0.08 0.10 0.12 0.14 0.16 0.18 0.20 0.22 0.24 0.26 0.28 0.30 0.32 0.34 0.36 0.38 0.40

Smoothing analysis RMSE

4D-Var S=1 EnKS-N S=1 SDA IEnKS-N S=1 SDA IEnKS-N S=L MDA IEnKS-N S=1

  • M. Bocquet

9th EnKF workshop, Bergen, Norway, 22-24 June 2014 10 / 22

slide-11
SLIDE 11

Localisation

IEnKF/IEnKS: Localisation

◮ Localisation in an EnVar context is non-trivial because localisation and the evolution model do not commute: Mk←0 (C◦B0)MT

k←0 = C◦

  • Mk←0B0MT

k←0

  • .

(8) ◮ Local analysis of IEnKF, and comparison with a non-scalable optimal approach.

0.1 0.2 0.3 0.4 0.5 0.6

Time interval between updates

0.15 0.2 0.3 0.4 0.5 0.6 0.7

Analysis rmse

CL IEnKF (bundle, opt. infl., c=10, non-scalable) N=10 IEnKF-N (bundle) N=20 LA IEnKF-N (bundle, c=10, εN=1) N=10

  • M. Bocquet

9th EnKF workshop, Bergen, Norway, 22-24 June 2014 11 / 22

slide-12
SLIDE 12

Localisation

IEnKF/IEnKS: Localisation

◮ Local analysis of IEnKS, and comparison with a non-scalable optimal approach (filtering performance).

5 10 15 20 30 40 50 1 2 3 4 6 7 8 9

DAW length L

0.16 0.18 0.20

Analysis RMSE

SDA IEnKS-N filtering N=20 MDA IEnKS-N filtering N=20 LA EnKS-N filtering N=10 l=10 LA MDA IEnKS-N filtering N=10 l=10 NSCL SDA IEnKS opt.infl. filtering N=10 l=10 LA SDA IEnKS-N filtering N=10 l=10 EnKS-N filtering N=20

  • M. Bocquet

9th EnKF workshop, Bergen, Norway, 22-24 June 2014 12 / 22

slide-13
SLIDE 13

Augmented state formalism

IEnKF/IEnKS: Augmented state formalism

◮ IEnKS treats parameters the way both 4D-Var and EnKF treat them. ◮ The state space is augmented from x ∈ RM to a vector z = x θ

  • ∈ RM+P ,

(9) Technically, there is nothing more to the joint state and parameter IEnKS than in the state IEnKS. ◮ A forward model needs to be introduced for the parameters: For instance, the persistence model (θk+1 = θk),

  • r some jittering such as a Brownian motion (θk+1 = θk +εk).
  • M. Bocquet

9th EnKF workshop, Bergen, Norway, 22-24 June 2014 13 / 22

slide-14
SLIDE 14

Augmented state formalism

Estimation of the Lorenz-95 forcing parameter F

◮ F is static but unknown.

1000 2000 3000 4000 5000

Time

8.10 8.05 8 7.95 7.90

Analysis of parameter F

EnKF-N EnKS-N L=50 IEnKF-N IEnKS-N L=5 IEnKS-N L=10 IEnKS-N L=30

◮ Augmented state vector ∈ R41, N = 20. The forcing of the true model is F = 8.

  • M. Bocquet

9th EnKF workshop, Bergen, Norway, 22-24 June 2014 14 / 22

slide-15
SLIDE 15

Augmented state formalism

Estimation of the Lorenz-95 forcing parameter F

◮ Setup: Lorenz-95, M = 40, N = 20, ∆t = 0.05, R = I. ◮ Comparison of 4D-Var S = 1, EnKS S = 1, SDA IEnKS S = 1, SDA IEnKS S = L, and MDA IEnKS S = 1.

0 1 5 10 20 30 50

Data assimilation window length (in ∆t)

0.25 0.05 0.10 0.15 0.07 0.125 0.20 0.30

Analysis RMSE (state variables)

4D-Var filtering 4D-Var smoothing EnKF-N/EnKS-N filtering EnKS-N smoothing MDA IEnKS-N filtering MDA IEnKS-N smoothing

1 5 10 30 50 100 20

Data assimilation window length (in ∆t)

0.001 0.002 0.004 0.008 0.016 0.032 0.064 0.128

Analysis RMSE (parameter F)

4D-Var filtering/smoothing EnKF-N/EnKS-N filtering EnKS-N smoothing MDA IEnKS-N filtering/smoothing

  • M. Bocquet

9th EnKF workshop, Bergen, Norway, 22-24 June 2014 15 / 22

slide-16
SLIDE 16

Augmented state formalism

Estimation of the Lorenz-95 forcing parameter F

◮ The forcing parameter F is time-varying. ◮ Internalised model error (F is in the augmented state) + unaccounted external model error (the true F is time-varying = persistence assumption). Method / F profile Sinusoidal Step-wise EnKF-N 0.063 0.079 EnKS-N L=50 0.040 0.063 4D-Var L=50 0.030 0.045 MDA IEnKS-N L=50 0.020 0.031

1000 2000 3000

Time

7,4 7,6 7,8 8 8,2 8,4 8,6 8,8

Retrospective analysis of parameter F

EnKF-N EnKS-N L=50 4D-Var L=50 MDA IEnKS-N L=50 Truth

1000 2000 3000

Time

7,5 8 8,5 9

Retrospective analysis of parameter F

EnKF-N L=50 EnKS-N L=50 4D-Var L=50 MDA IEnKS-N L=50 Truth

  • M. Bocquet

9th EnKF workshop, Bergen, Norway, 22-24 June 2014 16 / 22

slide-17
SLIDE 17

Towards using the IEnKF in a chemistry and transport model

Extending the Lorenz-95 model

◮ An online tracer model: Lorenz-95 (wind field) + tracer xm−1 cm− 1

2

xm cm+ 1

2

xm+1

  • Φm−1

Em− 1

2

Φm Em+ 1

2

Φm+1 ◮ The tracer is advected by the wind field of the Lorenz-95 model. We have chosen to use the simple Godunov/upwind scheme which is positive and conservative. dxm dt = (xm+1 −xm−2)xm−1 −xm +F , (10) dcm+ 1

2

dt = Φm −Φm+1 −λcm+ 1

2 +Em+ 1 2 ,

(11) where Φm = xmcm− 1

2

if xm ≥ 0, (12) = xmcm+ 1

2

if xm < 0. (13) ◮ Uniform emission of the tracer with the flux Em+ 1

2 . Deposited on the whole domain,

using a simple scavenging scheme parametrised by a scavenging ratio λ.

  • M. Bocquet

9th EnKF workshop, Bergen, Norway, 22-24 June 2014 17 / 22

slide-18
SLIDE 18

Towards using the IEnKF in a chemistry and transport model

Extending the Lorenz-95 model

100 200 300 400 500 600 5 10 15 20 25 30 35

  • 7.5
  • 5.0
  • 2.5

0.0 2.5 5.0 7.5 10.0 12.5

100 200 300 400 500 600 5 10 15 20 25 30 35

5.0 10.0 15.0 20.0 25.0 30.0 35.0 40.0 45.0

◮ Time evolution of the wind (top) and concentration (bottom) fields of the coupled Lorenz-95 - tracer model.

  • M. Bocquet

9th EnKF workshop, Bergen, Norway, 22-24 June 2014 18 / 22

slide-19
SLIDE 19

Towards using the IEnKF in a chemistry and transport model

Extending the Lorenz-95 model

1 5 10 20 40

Data assimilation window length (in ∆t)

0.25 0.05 0.10 0.15 0.07 0.12 0.20 0.30 0.04 0.03 0.02

Analysis RMSE (wind)

4D-Var filtering 4D-Var smoothing EnKS-N filtering EnKS-N smoothing SDA IEnKS-N filtering SDA IEnKS-N smoothing MDA IEnKS-N filtering MDA IEnKS-N smoothing

1 5 10 20 40

Data assimilation window length (in ∆t)

0.25 0.05 0.10 0.15 0.07 0.12 0.30 0.20 0.04 0.40

Analsyis RMSE (concentration)

4D-Var filtering 4D-Var smoothing EnKS-N filtering EnKS-N smoothing SDA IEnKS-N filtering SDA IEnKS-N smoothing MDA IEnKS-N filtering MDA IEnKS-N smoothing

◮ Mean filtering and smoothing analysis rmses of the wind variables (left) and concentration variables (right) of the online tracer model, as a function of the daw length for the IEnKS (finite-size variant), the EnKF/EnKS, and 4D-Var (with optimal inflation of the prior).

  • M. Bocquet

9th EnKF workshop, Bergen, Norway, 22-24 June 2014 19 / 22

slide-20
SLIDE 20

Towards using the IEnKF in a chemistry and transport model

Extending the Lorenz-95 model

  • 20
  • 15
  • 10
  • 5

5 10 15 20

Distance

  • 0.4
  • 0.2

0.2 0.4 0.6 0.8 1

Correlation coefficient

wind-wind concentration-concentration wind-concentration

10 20 30 40 50 60 70 80 10 20 30 40 50 60 70 80

  • 0.30
  • 0.15

0.00 0.15 0.30 0.45 0.60 0.75 0.90

◮ Structure functions of the mean correlation of the errors of the initial condition from the IEnKS applied to the online tracer model. ◮ The full error covariance matrix obtained from the IEnKS can be used to help 4D-Var by building better background statistics. It does help 4D-Var in the estimation of parameters, but does little to the estimation of the state variables whose error covariance matrix is quite dynamical.

  • M. Bocquet

9th EnKF workshop, Bergen, Norway, 22-24 June 2014 20 / 22

slide-21
SLIDE 21

Conclusions

Conclusions

The iterative ensemble Kalman smoother (IEnKS) is a method to seamlessly combine the advantages of variational and ensemble Kalman filtering. The IEnKS is a generalisation of the iterative ensemble Kalman filter (IEnKF). It is an EnVar method. It is flow-dependent, tangent linear and adjoint free. The IEnKF/IEnKS have the potential to (significantly) outperform both the EnKF and the 4D-Var in all regimes. IEnKS already does so with toy-models. IEnKS is very well suited for parameter (or joint state/parameter) estimation, and does so in a very simple way via the augmented state formalism. More complex reactive air quality toy-model under development in order to test the IEnKS on challenging atmospheric chemistry problems.

  • M. Bocquet

9th EnKF workshop, Bergen, Norway, 22-24 June 2014 21 / 22

slide-22
SLIDE 22

References

References

◮ Bocquet, M., Sakov, P., 2012. Combining inflation-free and iterative ensemble Kalman filters for strongly nonlinear systems. Nonlin. Processes Geophys. 19, 383–399. ◮ Bocquet, M., Sakov, P., 2014. An iterative ensemble Kalman smoother. Q. J. R.

  • Meteorol. Soc., in press, doi: 10.1002/qj.2236.

◮ Bocquet, M., Sakov, P., 2013. Joint state and parameter estimation with an iterative ensemble Kalman smoother. Nonlin. Processes Geophys., 20, 803–818. ◮ Gu, Y., Oliver, D. S., 2007. An iterative ensemble Kalman filter for multiphase fluid flow data assimilation. SPE Journal 12, 438–446. ◮ Hunt, B. R., Kostelich, E. J., Szunyogh, I., 2007. Efficient data assimilation for spatiotemporal chaos: A local ensemble transform Kalman filter. Physica D 230, 112–126. ◮ Sakov, P., Oliver, D., Bertino, L., 2012. An iterative EnKF for strongly nonlinear

  • systems. Mon. Wea. Rev. 140, 1988–2004.

◮ Yang, S.-C., E. Kalnay, and B. Hunt, 2012. Handling nonlinearity and non-Gaussianity in the Ensemble Kalman Filter: Experiments with the three-variable Lorenz model. Mon. Wea. Rev., 140, 2628–2646.

  • M. Bocquet

9th EnKF workshop, Bergen, Norway, 22-24 June 2014 22 / 22