Basic 2

Multiple synthetic benchmark experiments¶

Imports¶

NB: If you're on Gooble Colab, then replace %matplotlib notebook below by
!python -m pip install git+https://github.com/nansencenter/DAPPER.git
Also note that liveplotting does not work on Colab.

In [1]:

%matplotlib notebook
import dapper as dpr  # noqa: I001
import dapper.da_methods as da

%matplotlib notebook import dapper as dpr # noqa: I001 import dapper.da_methods as da

DA method configurations¶

In [2]:

from dapper.mods.Lorenz63.sakov2012 import HMM  # Expected rmse.a:

xps = dpr.xpList()
xps += da.Climatology()                                      # 7.6
xps += da.OptInterp()                                        # 1.25
#xps += da.Persistence()                                      # 10.7
#xps += da.PreProg(lambda k, xx, yy: xx[k])                   # 0
xps += da.Var3D(xB=0.1)                                      # 1.03
xps += da.ExtKF(infl=90)                                     # 0.87
xps += da.EnKF('Sqrt'   , N=3   , infl=1.30)                 # 0.82
xps += da.EnKF('Sqrt'   , N=10  , infl=1.02, rot=True)       # 0.63
xps += da.EnKF('PertObs', N=500 , infl=0.95, rot=False)      # 0.56
xps += da.EnKF_N(         N=10             , rot=True)       # 0.54
xps += da.iEnKS('Sqrt'  , N=10  , infl=1.02, rot=True)       # 0.31
xps += da.PartFilt(       N=100 , reg=2.4  , NER=0.3)        # 0.38
xps += da.PartFilt(       N=800 , reg=0.9  , NER=0.2)        # 0.28
# xps += da.PartFilt(       N=4000, reg=0.7  , NER=0.05)       # 0.27
# xps += da.PFxN(xN=1000,   N=30  , Qs=2     , NER=0.2)        # 0.56

from dapper.mods.Lorenz63.sakov2012 import HMM # Expected rmse.a: xps = dpr.xpList() xps += da.Climatology() # 7.6 xps += da.OptInterp() # 1.25 #xps += da.Persistence() # 10.7 #xps += da.PreProg(lambda k, xx, yy: xx[k]) # 0 xps += da.Var3D(xB=0.1) # 1.03 xps += da.ExtKF(infl=90) # 0.87 xps += da.EnKF('Sqrt' , N=3 , infl=1.30) # 0.82 xps += da.EnKF('Sqrt' , N=10 , infl=1.02, rot=True) # 0.63 xps += da.EnKF('PertObs', N=500 , infl=0.95, rot=False) # 0.56 xps += da.EnKF_N( N=10 , rot=True) # 0.54 xps += da.iEnKS('Sqrt' , N=10 , infl=1.02, rot=True) # 0.31 xps += da.PartFilt( N=100 , reg=2.4 , NER=0.3) # 0.38 xps += da.PartFilt( N=800 , reg=0.9 , NER=0.2) # 0.28 # xps += da.PartFilt( N=4000, reg=0.7 , NER=0.05) # 0.27 # xps += da.PFxN(xN=1000, N=30 , Qs=2 , NER=0.2) # 0.56

With Lorenz-96 instead¶

In [3]:

# from dapper.mods.Lorenz96.sakov2008 import HMM   # Expected rmse.a:
# xps = dpr.xpList()
# xps += da.Climatology()                                     # 3.6
# xps += da.OptInterp()                                       # 0.95
# xps += da.Var3D(xB=0.02)                                    # 0.41
# xps += da.ExtKF(infl=6)                                     # 0.24
# xps += da.EnKF('PertObs', N=40, infl=1.06)                  # 0.22
# xps += da.EnKF('Sqrt'   , N=28, infl=1.02, rot=True)        # 0.18
# # More sophisticated:
# xps += da.EnKF_N(         N=24, rot=True)                   # 0.21
# xps += da.EnKF_N(         N=24, rot=True, xN=2)             # 0.18
# xps += da.iEnKS('Sqrt'  , N=40, infl=1.01, rot=True)        # 0.17
# # With localisation:
# xps += da.LETKF(          N=7 , infl=1.04, rot=True, loc_rad=4)  # 0.22
# xps += da.SL_EAKF(        N=7 , infl=1.07, rot=True, loc_rad=6)  # 0.23

# from dapper.mods.Lorenz96.sakov2008 import HMM # Expected rmse.a: # xps = dpr.xpList() # xps += da.Climatology() # 3.6 # xps += da.OptInterp() # 0.95 # xps += da.Var3D(xB=0.02) # 0.41 # xps += da.ExtKF(infl=6) # 0.24 # xps += da.EnKF('PertObs', N=40, infl=1.06) # 0.22 # xps += da.EnKF('Sqrt' , N=28, infl=1.02, rot=True) # 0.18 # # More sophisticated: # xps += da.EnKF_N( N=24, rot=True) # 0.21 # xps += da.EnKF_N( N=24, rot=True, xN=2) # 0.18 # xps += da.iEnKS('Sqrt' , N=40, infl=1.01, rot=True) # 0.17 # # With localisation: # xps += da.LETKF( N=7 , infl=1.04, rot=True, loc_rad=4) # 0.22 # xps += da.SL_EAKF( N=7 , infl=1.07, rot=True, loc_rad=6) # 0.23

Other models (suitable xp's listed in HMM files):¶

In [4]:

# from dapper.mods.LA           .evensen2009 import HMM
# from dapper.mods.KS           .bocquet2019 import HMM
# from dapper.mods.LotkaVolterra.settings101 import HMM

# from dapper.mods.LA .evensen2009 import HMM # from dapper.mods.KS .bocquet2019 import HMM # from dapper.mods.LotkaVolterra.settings101 import HMM

Launch¶

Write some more non-arg parameters to the xps. In this case we set the seed, so that repeat experiments produce exactly the same result.

In [5]:

for xp in xps:
    xp.seed = 3000

for xp in xps: xp.seed = 3000

Adjust experiment duration

In [6]:

HMM.tseq.T = 50

HMM.tseq.T = 50

Run/assimilate (for each xp in xps)

In [7]:

save_as = xps.launch(HMM, liveplots=False)

save_as = xps.launch(HMM, liveplots=False)

Experiment gets stored at /Users/para/data/dpr_data/noname/run_2024-10-31__17-12-44

Truth & Obs: 100%|█████████████| 5000/5000 [00:00<00:00, 41765.45it/s]
Climatology: 100%|█████████████| 5000/5000 [00:00<00:00, 58936.74it/s]
Truth & Obs: 100%|█████████████| 5000/5000 [00:00<00:00, 42079.80it/s]
OptInterp: 100%|███████████████| 5000/5000 [00:00<00:00, 25958.39it/s]
Truth & Obs: 100%|█████████████| 5000/5000 [00:00<00:00, 42096.02it/s]
Var3D: 100%|███████████████████| 5000/5000 [00:00<00:00, 26433.50it/s]
Truth & Obs: 100%|█████████████| 5000/5000 [00:00<00:00, 38184.81it/s]
ExtKF infl:90: 100%|███████████| 5000/5000 [00:00<00:00, 19740.78it/s]
Truth & Obs: 100%|█████████████| 5000/5000 [00:00<00:00, 41851.22it/s]
EnKF infl:1.3 upd_a:Sqrt N:3 rot:False: 100%|█| 5000/5000 [00:00<00:00
Truth & Obs: 100%|█████████████| 5000/5000 [00:00<00:00, 40812.37it/s]
EnKF infl:1.02 upd_a:Sqrt N:10 rot:True: 100%|█| 5000/5000 [00:00<00:0
Truth & Obs: 100%|█████████████| 5000/5000 [00:00<00:00, 42888.91it/s]
EnKF infl:0.95 upd_a:PertObs N:500 rot:False: 100%|█| 5000/5000 [00:00
Truth & Obs: 100%|█████████████| 5000/5000 [00:00<00:00, 40611.08it/s]
EnKF_N infl:1 N:10 rot:True xN:1: 100%|█| 5000/5000 [00:00<00:00, 1325
Truth & Obs: 100%|█████████████| 5000/5000 [00:00<00:00, 40498.23it/s]
iEnKS infl:1.02 upd_a:Sqrt N:10 rot:True: 100%|█| 201/201 [00:03<00:00
Truth & Obs: 100%|█████████████| 5000/5000 [00:00<00:00, 42317.38it/s]
PartFilt N:100 reg:2.4 NER:0.3: 100%|█| 5000/5000 [00:00<00:00, 13981.
Truth & Obs: 100%|█████████████| 5000/5000 [00:00<00:00, 42299.37it/s]
PartFilt N:800 reg:0.9 NER:0.2: 100%|█| 5000/5000 [00:00<00:00, 8458.1

Print results¶

In [8]:

print(xps.tabulate_avrgs(statkeys=["rmse.a"]))

print(xps.tabulate_avrgs(statkeys=["rmse.a"]))

    da_method     infl  upd_a      N  rot    xN  reg  NER  |  rmse.a  1σ  
--  -----------  -----  -------  ---  -----  --  ---  ---  -  ------------
 0  Climatology                                            |    7.4  ±0.3 
 1  OptInterp                                              |    1.23 ±0.04
 2  Var3D                                                  |    1.07 ±0.06
 3  ExtKF        90                                        |    0.88 ±0.09
 4  EnKF          1.3   Sqrt       3  False                |    0.8  ±0.1 
 5  EnKF          1.02  Sqrt      10  True                 |    0.59 ±0.06
 6  EnKF          0.95  PertObs  500  False                |    1.1  ±0.3 
 7  EnKF_N        1               10  True    1            |    0.71 ±0.07
 8  iEnKS         1.02  Sqrt      10  True                 |    0.23 ±0.02
 9  PartFilt                     100             2.4  0.3  |    0.32 ±0.03
10  PartFilt                     800             0.9  0.2  |    0.20 ±0.03