Nonlinear Data Assimilation using an extremely efficient Particle Filter

1
Nonlinear Data Assimilation using an extremely
efficient Particle Filter

• Peter Jan van Leeuwen
• Data-Assimilation Research Centre
• University of Reading

2
The Agulhas System
3
In-situ observations
4
In-situ observations
5
In situ observationsTransport through Mozambique
Channel
6
Data assimilation

• Uncertainty points to use of probability density
  functions.

P(u)
0.5
0.0
1.0
u (m/s)
7
Data assimilation general formulation
Bayes theorem
Solution is pdf! NO INVERSION !!!
8
How is this used today?

• Present-day data-assimilation systems are based
  on linearizations and search for one optimal
  state
• (Ensemble) Kalman filter assumes Gaussian pdfs
• 4DVar smoother assumes Gaussian pdf for initial
  state and observations (no model errors)
• Representer method as 4DVar but with Gaussian
  model errors
• Combinations of these

9
Prediction smoothers vs. filters

• The smoother solves for the mode of the
  conditional joint pdf p( x0T d0T) (modal
  trajectory).
• The filter solves for the mode of the conditional
  marginal pdf p( xT d0T).
• For linear dynamics these give the same
  prediction.

10

• Filters maximize the marginal pdf

• Smoothers maximize the joint pdf

These are not the same for nonlinear problems !!!
11
Example
Nonlinear model
2 xn
xn1 0.5 xn _________ nn
1 e (xn - 7)
Initial pdf
x0 N(-0.1, 10)
Model noise
nn N(0, 10)
12
Example marginal pdfs
Note mode is at x - 0.1
Note mode is at x8.5
x0
xn
13
Example joint pdf
Mode joint pdf
x0
xn
Modes marginal pdfs
14
And what about the linearizations?

• Kalman-like filters solve for the wrong state
  gives rise to bias.
• Variational methods use gradient methods, which
  can end up in local minima.
• 4DVar assumes perfect model gives rise to
  bias.

15
Where do we want to go?

• Represent pdf by an ensemble of model states
• Fully nonlinear

Time
16
How do we get there? Particle filter?
Use ensemble
with
the weights.
17
What are these weights?

• The weight w_i is the pdf of the observations
  given the model state i.
• For M independent Gaussian distributed
  observation errors

18
Standard Particle filter
19
Particle Filter degeneracy resampling

• With each new set of observations the old weights
  are multiplied with the new weights.
• Very soon only one particle has all the weight
• Solution
• Resampling duplicate high-weight particles are
  abandon low-weight particles

20
Problems

• Probability space in large-dimensional systems is
  empty the curse of dimensionality

u(x1)
u(x2)
T(x3)
21
Standard Particle filter
Not very efficient !
22
Specifics of Bayes Theorem I
We know from Bayes Theorem
Now use
in which we introduced the transition density
23
Specifics of Bayes Theorem II
This can be rewritten as
q is the proposal transition density, which might
be conditioned on the new observations!
This leads finally to
24
Specifics of Bayes Theorem III
How do we use this? A particle representation of

Giving
Now we choose from the proposal transition
density
for each particle i.
25
Particle filter with proposal density
Stochastic model
Proposed stochastic model
Leads to particle filter with weights
26
Meaning of the transition densities
the probability of this specific value for the
random model error. For Gaussian model errors we
find
A similar expression is found for the proposal
transition
27
Particle filter with proposal transitiondensity
28
Experiment Lorentz 1963 model (3 variables
x,y,z, highly nonlinear)
x-value
Measure only X-variable
y-value
29
Standard Particle filter with resampling 20
particles
Typically 500 particles needed !
X-value
Time
30
Particle filter with proposal transition density
3 particles
X-value
Time
31
Particle filter with proposal transition density
3 particles
Y-value (not observed)
Time
32
However degeneracy

• For large-scale problems with lots of
  observations this method is still degenerate
• Only a few particles get high weights the other
  weights are negligibly small.
• However, we can enforce almost equal weight for
  all particles

33
Equal weights

1. Write down expression for each weight with q
   deterministic

Prior transition density
Likelihood
2. When H is linear this is a quadratic function
in for each particle. 3. Determine
the target weight
34
Almost Equal weights I
1
4
3
Target weight
2
5
4. Determine corresponding model states, e.g.
solving alpha in
35
Almost equal weights II

• But proposal density cannot be deterministic
• Add small random term to model equations from a
  pdf with broad wings e.g. Gauchy
• Calculate the new weights, and resample if
  necessary

36
Application Lorenz 1995
N40 F8 dt 0.005 T 1000 dt Observe
every other grid point Typically 10,000
particles needed
37
Ensemble mean after 500 time steps20 particles
Position
38
Ensemble evolution at x2020 particles
Time step
39
Ensemble evolution at x35(unobserved) 20
particles
40
Isnt nudging enough?
Only nudged
Nudged and weighted
41
Isnt nudging enough?
Unobserved variable
Only nudged
Nudged and weighted
42
Conclusions

• The nonlinearity of our problem is growing

• Particle filters with proposal transition
  density
• solve for fully nonlinear solution
• very flexible, much freedom
• application to large-scale problems
  straightforward

43
Future

• Fully nonlinear filtering (smoothing) forces us
  to concentrate on the transition densities, so on
  the errors in the model equations.
• What is the sensitivity to our choice of the
  proposal?
• What can we learn from studying the statistics of
  the nudging terms?
• How do we use the pdf???