Rare Event Sampling for the Earth and Planetary Sciences

By Robert J. Webber

In the Earth and planetary sciences, high-resolution predictive models can reveal the probability of a major hurricane, a collision between two planets, and other rare and extreme events. However, unlocking the full potential of these high-resolution models requires long runtimes on high-performance computing clusters. Due to these intense computational requirements, ensembles of 10 to 100 simulations stretch the limits of the best available supercomputers.

Unfortunately, running small ensembles with up to 100 members only leads to a rough estimate of the likelihood of a rare, extreme event. For example, the return time for a Category 3 hurricane to make landfall in New England might be 10, 20, or 50 years; a simulation on that timescale will generate only a few events at best. Moreover, many models do not give a single, definitive answer because they depend sensitively on the initial conditions or they are stochastic and evolve in different ways based on the random seed. To identify return times, we therefore need to run a large ensemble of costly simulations.

Extreme weather and climate predictions are further complicated by the uncertainty of global climate change, since the warming atmosphere may lead to more extreme heat waves [6] and more intense tropical cyclones [5]. Our climate’s uncertain future heightens the importance of accurate extreme event predictions, which can help policymakers understand the risk levels that are associated with different climate change scenarios.

Rather than simply accept the current limitations of a small ensemble size, computational scientists are inventing creative ways to obtain more accurate predictions. Our research group is investigating a computational approach that is inspired by evolutionary biology; we use “survival of the fittest” to steer our simulations toward the exhibition of rare, extreme events. We created and tested a new survival-of-the-fittest algorithm that increases the accuracy of rare event estimation in the Earth and planetary sciences [1, 8].

The survival-of-the-fittest approach for rare event sampling is not new. Stanislaw Ulam and John Von Neumann suggested this approach when experimenting with the first computers in the 1950s [4], and scientists have since applied the technique to areas such as quantum chemistry [3] and nuclear physics [2]. Building on this past research, our current efforts—as well as other ongoing studies—are exploring the hypothesis that these algorithms can be adapted to calculate extreme event probabilities for the Earth and planetary sciences [1, 7, 8].

The survival-of-the-fittest algorithm is based on a fitness function \(\phi(\mathbf{X}_t)\) for each simulated trajectory \(\mathbf{X}_t\); this function reaches a high value when the rare event is likely to occur. We then “split” or “kill” our trajectories at regular intervals, with the number of progeny dependent on the fitness level \(\phi(\mathbf{X}_t)\). We create identical copies of the split trajectories and perturb their states slightly or assign different random seeds so that the “children” evolve along unique paths in the future. Killing the other trajectories frees up valuable computing resources that allow for more splitting at a later time.

This approach can inspire extreme instances of a rare event, like a raging Category 5 hurricane in Boston or a collision between Mercury and Venus. However, is it possible to learn the probability of a rare event based on these simulations?

To calculate rare event probabilities, we assign probability weights \(w_t^{(1)}, w_t^{(2)}, \ldots, w_t^{(N)}\) to each of our simulations \(\mathbf{X}_t^{(1)}, \mathbf{X}_t^{(2)}, \ldots, \mathbf{X}_t^{(N)}\). Initially, all of the weights are equal and sum to \(1\). We then evolve the weights whenever we split or kill a trajectory. For example, if we split the trajectory \(\mathbf{X}_t^{(1)}\) into two children, we assign each child the probability weight \(w_t^{(1)}/2\). Conversely, if \(\mathbf{X}_t^{(1)}\) and \(\mathbf{X}_t^{(2)}\) are two equally weighted trajectories and we randomly kill one of them, we assign the survivor the sum of the probability weights \(w_t^{(1)} + w_t^{(2)}\). After running all of the simulations up to the target time \(t = T\), we can simply add up the weights for all of the surviving simulations that exhibited the rare event. This approach leads to unbiased rare event probability statistics, and mathematical analysis demonstrates that our algorithm can be far more efficient than the brute-force approach of running trajectories without any splitting or killing.

We tested our survival-of-the-fittest algorithm by simulating intense hurricanes that approach Boston [8] and rare solar system trajectories wherein Mercury is ejected from its orbit [1]. Figures 1 and 2 depict the results of these extreme event simulations. Figure 1 shows simulations of Hurricane Earl at high levels of intensity (measured by minimum sea surface pressure); the simulations with our quantile diffusion Monte Carlo (QDMC) algorithm in Figure 1b are much more intense than the simulations that use direct sampling without any splitting or killing (see Figure 1a). Figure 2 presents QDMC simulations of Mercury’s instability, which leads to close encounters between Mercury and Venus. In both cases, our algorithm yielded high-accuracy predictions that would ordinarily have required 10 times as much computational effort, if not more.

However, we acknowledge the challenge of extreme event probability calculations and the hurdles that complicate widespread adoption of our approach. One such obstacle is that our rare event sampling algorithm requires a fitness function \(\phi\) that approximates the likelihood of the rare event’s occurrence. We would ideally approximate the optimal fitness function \(\phi\) via machine learning, but we have always defined \(\phi\) in our case studies through physically motivated suggestions and testing. As another hurdle, our approach is only viable with many simulations; our algorithm has limited potential with just 10. Although we have run it with 100 simulations, it performs even better with 1,000.

Despite these hurdles, we envision a future in which researchers routinely utilize survival-of-the-fittest algorithms to study extreme events in the Earth and planetary sciences. This outlook seems likely for several reasons. First, extreme events are fascinating, fearsome, and potentially quite destructive; as such, there will always be a need to better understand and prepare for them. Second, Earth and planetary scientists acknowledge the importance of creative new approaches to extreme event simulations that ensure reliable statistics while also accounting for uncertainties in the initial conditions and parameters. And third, the adaptability of survival-of-the-fittest algorithms to different models and computational architectures make them versatile computational tools.

Looking forward, we believe that survival-of-the-fittest algorithms will emerge as a mainstream approach for the calculation of extreme event probabilities. We are hence working to develop the most practical, elegant, and theoretically well-founded algorithm for this task.

---

Robert J. Webber delivered a minisymposium presentation on this research at the 2022 SIAM Conference on Mathematics of Data Science, which took place in San Diego, Ca., last year.

References
[1] Abbot, D.S., Webber, R.J., Hadden, S., Seligman, D., & Weare, J. (2021). Rare event sampling improves Mercury instability statistics. Astrophys. J., 923(2), 236. 
[2] Carter, L.L., & Cashwell, E.D. (1975). Particle transport simulation with the Monte Carlo method (Report No. TID-26607). Los Alamos, NM: Los Alamos Scientific Laboratory.
[3] Grimm, R.C., & Storer, R.G. (1971). Monte-Carlo solution of Schrödinger’s equation. J. Comput. Phys., 7(1), 134-156.
[4] Kahn, H. (1955). Use of different Monte Carlo sampling techniques (Tech. Report P-766). Santa Monica, CA: RAND Corporation.
[5] Knutson, T.R., Chung, M.V., Vecchi, G., Sun, J., Hsieh, T.-L., & Smith, A.J.P. (2021). Climate change is probably increasing the intensity of tropical cyclones. In C. Le Quéré, P. Liss, & P. Forster (Eds.), Critical issues in climate change science. ScienceBrief Reviews. 
[6] Mora, C., Dousset, B., Caldwell, I.R., Powell, F.E., Geronimo, R.C., Bielecki, C.R., … Trauernicht, C. (2017). Global risk of deadly heat. Nat. Clim. Change, 7, 501-506.
[7] Ragone, F., Wouters, J., & Bouchet, F. (2018). Computation of extreme heat waves in climate models using a large deviation algorithm. Proc. Natl. Acad. Sci., 115(1), 24-29.
[8] Webber, R.J., Plotkin, D.A., O’Neill, M.E., Abbot, D.S., & Weare, J. (2019). Practical rare event sampling for extreme mesoscale weather. Chaos, 29(5), 053109. 

Robert J. Webber is a postdoctoral scholar in the Computing & Mathematical Sciences Department at the California Institute of Technology, where he is supervised by Joel A. Tropp. In his research, Webber analyzes and applies Monte Carlo algorithms and algorithms for randomized numerical linear algebra.