Add new spatially dependent lorenz-forcing calibration comparison between methods

[odunbar]

Purpose

Close #424

To-do

• Note: results might be different than @rgjini python results because EKI and ETKI don't have the prior included as regularization within the cost function yet

Content

• Included spin-up period in the statistics function: (i) Look at some trajectories of x(t) from prior, how long is the time-to-steady-state (make the spin up for each ensemble > 2*max of this) (ii) optional create a heatmap of these
• Initialized approach (for initial state condition) - (may want to change the method later).
• Included RMSE calculation within the iteration loop: get_error(ekiobj) , apply sqrt and /sqrt(length(y))
• introduced break/continue in loop when RMSE < 1.
• Add outer for loop to run different ensemble sizes for the problem
• Add testing different random seeds for each of the algorithms

---

• I have read and checked the items on the review checklist.

[rgjini]

Work in progress (WIP)

[rgjini]

The spin-up period in the statistics function is working properly:

• I checked the trajectory of x(t) and it seems like the initial spinup is a little less than 2 total time steps. Multiplied by 2, gives 4 time units for spinup
• I looked at heat maps in python, attached is an example of a trajectory
  

[rgjini]

I finished up the example and everything seems to be working! I fixed the random seed and ensemble size for loops, so those are both working. I kept the amount of different ensemble sizes and random seeds low for now, the code does take a while to run (maybe 10 minutes). I also have a few plots I made at the end. Two plots are of the final optimized solution compared to the true parameters for different ensemble sizes and the last one is supposed to be how many iterations each algorithm took for each ensemble size. For some reason, the last plot isn't looking how I would expect, but it might be related to an issue I ran into that I will describe in another comment. Here are some of the plots I made!:

[rgjini]

I finished including all the changes we talked about into the example script and I believe everything is ready to go! Attached are some new plots with the parameter estimates, the model output from the ensemble means compared with the data (and associated errors), and a plot showing how many forward runs each of the algorithms take to converge:

[odunbar]

Hi @rgjini These look awesome, I have rerun the simulations, and got the results below! I have also rebased to main so merging is possible again!

If you are pleased with this then I am happy to merge

PS: I think it is still worth discussing the exact observations (i.e. the form of "G") at a later date, as I think this experiment does very long averages to eliminate an observational noise that EKP methods don't need to avoid and then we add artificial noise back on to the observations. In particular, the current setup way we have G and the noise I think only works with T~500, but really we could do maybe have good learning with something as small as T~10.

[odunbar]

[rgjini]

@odunbar Thank you! The plots look great to me, so I'll go ahead and merge :)

I agree, I think it's worth discussing more. If you are still free at 10:30 before we meet with Matti and Tapio, maybe we can discuss the observations more then? I'm also happy to cancel that meeting, and we can either discuss it during our meeting at 11 or in a couple of weeks?