# Chaos in a simple model of a delta network

See allHide authors and affiliations

Edited by Andrea Rinaldo, École Polytechnique Fédérale de Lausanne, Lausanne, Switzerland, and approved September 17, 2020 (received for review May 22, 2020)

## Significance

River deltas are threatened with rapid change due to human activity, highlighting the importance of predicting their dynamics. Based on a model of a simple delta network, we find that delta dynamics can be chaotic, implying a fundamental constraint on detailed long-term predictions of deltas. However, our estimates of the fundamental predictability horizon suggest substantial room for improvement in delta modeling. We test our model by using it to generate stratigraphy, translating our measure of predictability from time to vertical distance, and identifying the time scale required for long-term averages to become predictable. Although chaos precludes long-term prediction of deltas, the statistical behavior of deltas and their response to external changes such as human activity are predictable.

## Abstract

The flux partitioning in delta networks controls how deltas build land and generate stratigraphy. Here, we study flux-partitioning dynamics in a delta network using a simple numerical model consisting of two orders of bifurcations. Previous work on single bifurcations has shown periodic behavior arising due to the interplay between channel deepening and downstream deposition. We find that coupling between upstream and downstream bifurcations can lead to chaos; despite its simplicity, our model generates surprisingly complex aperiodic yet bounded dynamics. Our model exhibits sensitive dependence on initial conditions, the hallmark signature of chaos, implying long-term unpredictability of delta networks. However, estimates of the predictability horizon suggest substantial room for improvement in delta-network modeling before fundamental limits on predictability are encountered. We also observe periodic windows, implying that a change in forcing (e.g., due to climate change) could cause a delta to switch from predictable to unpredictable or vice versa. We test our model by using it to generate stratigraphy; converting the temporal Lyapunov exponent to vertical distance using the mean sedimentation rate, we observe qualitatively realistic patterns such as upwards fining and scale-dependent compensation statistics, consistent with ancient and experimental systems. We suggest that chaotic behavior may be common in geomorphic systems and that it implies fundamental bounds on their predictability. We conclude that while delta “weather” (precise configuration) is unpredictable in the long-term, delta “climate” (statistical behavior) is predictable.

Deltas are landforms arising from the deposition of sediment by a river entering a standing body of water. Deltas worldwide are highly populated and widely relied upon for agriculture and navigation but are under threat of collapse due to relative sea-level rise and sediment-supply reduction (1⇓–3). Additionally, delta deposits are important reservoirs for groundwater and hydrocarbons, and better understanding their stratigraphic architecture would improve prediction of subsurface fluid flow (4).

Delta networks deliver sediment to different parts of the delta, controlling where land is built. Deltas have been proposed to self-organize their flux partitioning for a given network topology (5). The flux partitioning in the network can change over time through the process of avulsion, which may result in the abandonment of old channels and creation of new ones or merely a shift in the flux partitioning without modification of the number of channels (termed “soft avulsion”) (6). A wide body of research has aimed to better understand avulsions in order to better predict their timing and location and to explore how deltas record themselves in stratigraphy (7⇓⇓⇓⇓⇓⇓–14).

Our ability to accurately predict the future evolution of river deltas is hampered by our incomplete representation of the physical processes in our models, as well as uncertainty in boundary conditions. A similar problem arises in the prediction of subsurface stratigraphy based on incomplete data. Is there a fundamental limit to our ability to accurately predict delta evolution, which would apply even if we could perfectly simulate the relevant physical processes and boundary conditions?

Complex and irregular behavior, typically assumed to be stochastic, can be found everywhere in geomorphology. Examples include the motions of bedload particles (15), the evolution of dune fields (16), the organization of drainage networks (17, 18), the dynamics of braided (19) and meandering (20) rivers, and the kinematics of delta surfaces (21). Reproducing in detail the irregular dynamics of geomorphic systems is beyond the reach of existing models, but it is unknown to what extent detailed prediction of these systems is impossible rather than “merely” difficult.

In chaotic systems, long-term prediction is fundamentally impossible. Strogatz (22) proposed the following definition of chaos: “Chaos is aperiodic long-term behavior in a deterministic system that exhibits sensitive dependence on initial conditions.” Sensitive dependence on initial conditions refers to the idea that two systems started with nearly identical initial conditions diverge exponentially. This means that any error in our measurement of a system’s initial condition will grow exponentially through time. Because we can never measure a system with infinite precision, this implies a fundamental limit to long-term prediction. In contrast to stochastic systems, unpredictability arises in chaotic systems despite their lack of randomness; chaotic systems are by definition deterministic.

While many geomorphic systems display highly complex dynamics, demonstrating chaotic behavior is difficult. Major challenges include the lack of very long time series, which are typically required to detect chaos (23⇓–25), and the fact that some methods to detect chaos can be fooled by certain types of noise (26). In the study of geomorphic systems, some have attempted to use spatial information to supplement the lack of temporal information needed to detect chaos (27), for example, from a time series of ripple migration (28).

Despite the complex behavior of natural geomorphic systems, a growing body of research on so-called “reduced-complexity models” has shown that complex dynamics that are at least qualitatively similar to those of natural systems can be generated even by models obeying relatively simple rules (16, 19, 29⇓⇓⇓⇓–34). Demonstration of chaotic dynamics in an idealized and simplified model of a natural system can be a strong suggestion that the system is chaotic; for example, the Lorenz equations are an extreme simplification of atmospheric convection but are believed to indicate that the atmosphere is itself chaotic (35). Chaotic behavior has been uncovered in models of geological phenomena such as earthquakes (36) and the geodynamo (37). Dynamic stratigraphic models are typically nonlinear dissipative systems and, therefore, are potentially susceptible to exhibiting chaos (38). More generally, the same could be said of many geomorphological models. However, to date, demonstrations of chaos in geomorphological models are rare: one example is from Pelletier (39), who identified deterministic chaos in a model of landform evolution, albeit with a definition of chaos that differs from the one we use in this paper. Another example is Phillips (40), who found chaos in a difference equation modeling regolith cover on a hillslope. Chaos has also been proposed for numerical models of river meandering (41, 42).

In this paper, we consider a numerical model of a simple delta network. The network consists of three coupled bifurcations: a channel splits into two branches at a bifurcation, each of which split into two additional branches. Our model is a further extension of the quasi-one dimensional (1D) bifurcation model developed by Bolla Pittaluga et al. (43) (hereafter, BRT) and extended by Salter et al. (44, 45) (hereafter, SPV) to include the effect of deposition in the downstream branches. We find that the coupling between the upstream and downstream bifurcations leads to chaotic dynamics in the flux partitioning of the network. This suggests that delta networks that have coupled bifurcations behave chaotically and that it is possible to quantitatively estimate the time scale over which long-term detailed prediction of their avulsion dynamics becomes impossible.

## Model Overview

Our starting point is the BRT (43) quasi-1D bifurcation model, which has been widely used and built upon (46⇓⇓⇓⇓⇓–52). This model allows for a transverse bed slope immediately upstream of the bifurcation, controlling the steering of water and sediment to the downstream branches. This is implemented by dividing the reach immediately upstream of the bifurcation into two laterally interacting cells, which have a length of

Whereas prior studies chose downstream boundary conditions resulting in sediment bypass (i.e., sediment output equal to sediment input), river deltas are landforms resulting from sediment deposition. SPV (44) imposed deposition via a specified bypass fraction F, which is the fraction of the supplied sediment that exits the downstream end of the network. A bypass fraction *F* = 1) results in a fixed point, i.e., unchanging discharge partitioning, deposition creates the possibility of ongoing avulsion dynamics. The basic mechanism is as follows: because more sediment is routed to the dominant branch than the subordinate branch, its slope tends to become lower relative to the subordinate branch over time. Eventually, this slope difference causes a switch in the discharge partitioning. As with bypass, a symmetric fixed point always exists but may be stable or unstable. However, when it is unstable, rather than evolving toward a stable asymmetric fixed point, the bifurcation undergoes repeating avulsion dynamics, either soft avulsion (6) or full avulsion. These avulsion cycles are stable limit cycles, meaning that for a range of initial conditions, the bifurcation tends toward a periodic solution that is independent of the initial conditions.

To summarize, when the evolution of the bifurcation is decoupled from the downstream slopes, the system evolves toward a stable fixed point (43). In contrast, the coupling between the positive feedback in the bifurcation and deposition in the downstream branches results in stable limit cycles (44).

Here, we consider the scenario in which bifurcations are coupled in a simple network, with an upstream bifurcation coupled to bifurcations in each of its two downstream branches, as shown in Fig. 1. The model consists of a total of seven branches coupled via three bifurcations. The widths (

For simplicity, we assume that the geometry of the network (as described by **13**. This allows us to calculate the characteristic avulsion time scale, which we define as **12**), which, in this paper, we will confine to **4**); and the model parameters α and r, which specify the length of the divided reach and strength of bed load steering from a lateral bed slope, respectively. The parameter values used in this paper are listed in Table 1.

## Dynamics

We begin by describing the types of dynamics that can be observed in our network model. At any given time, the model consists of N bed elevations, where N is set by the discretization of the branches. Because the system remains unchanged following an arbitrary vertical translation, in fact only

For the delta network, a fixed point corresponding to a symmetric partitioning at all three bifurcations and uniform aggradation across the network always exists. However, as with the single-bifurcation models discussed above, the fixed point may be stable or unstable, depending on model parameters. When the fixed point is unstable, then the flux partitioning of the network becomes dynamic. At any given moment in time, the flux partitioning tends to be uneven, causing faster sedimentation in some areas than others. However, persistent unequal sedimentation would cause the relative elevation of some areas to become increasingly low, a situation favoring increased sediment supply to those areas. Therefore, flux redistribution over time leads to uniform aggradation everywhere on average, with the rate

For some parameter values, we observe periodic dynamics associated with stable limit cycles (Fig. 2 *A* and *C*). We note that these limit cycles are often asymmetric, implying the existence of another limit cycle with the asymmetry reversed (altogether, there are eight equivalent labelings of the bifurcation network; for each of these symmetries, either a limit cycle must be itself symmetric, or otherwise a reciprocal limit cycle exists). Note that even when a limit cycle is asymmetric, the average sedimentation rate over the course of a limit cycle is equal everywhere; this can be achieved through different combinations of avulsion magnitude and duration.

For other choices of parameters, we observe aperiodic dynamics (i.e., the time series of flux partitioning never repeats itself exactly). In some cases, the dynamics correspond to full avulsion, which we define as avulsion behavior where, at times, at least one branch is completely abandoned (44). Alternatively, the network may remain in the soft-avulsion regime (*SI Appendix*, Fig. S1), in which all branches maintain a positive water discharge indefinitely. Despite the aperiodicity of the dynamics, the system remains within a bounded region of phase space, and the statistical behavior is stationary when averaged over sufficient duration. Aperiodicity does not by itself imply chaos; in fact, Fig. 2*B* shows an example of a quasiperiodic attractor (trajectories are confined to the surface of a 2-torus but never repeat exactly). However, in our system, quasiperiodicity appears to be rare relative to true chaos (Fig. 2 *D* and *E*). Next, we will show how the existence of chaos in our model can be established.

## Sensitive Dependence on Initial Conditions

The hallmark signature of chaos is the sensitive dependence on initial conditions. This is defined by the exponential growth of arbitrarily small perturbations: the trajectories of two simulations with ever-so-slightly different initial conditions diverge exponentially until the difference between simulations is on the scale of the size of the attractor. In other words, given an initial separation *A* and *B* and *SI Appendix*, Fig. S2). In contrast, for parameter values resulting in a stable limit cycle, the MLE is zero. We observe that the MLE appears to vary smoothly across transitions from a stable limit cycle to a strange attractor (*SI Appendix*, Fig. S3) and find that at least some of these transitions follow a period-doubling route to chaos, similar to the logistic map (*SI Appendix*, Figs. S4 and S5).

## Characterizing the Attractor

When viewed in phase space, strange attractors are typically fractal. We can quantify the fractal dimension of the attractor using the correlation dimension (23). This technique involves selecting random points on the attractor and measuring the distance ℓ between all possible pairs of points. The correlation integral *C*). We have observed fractal dimensions in excess of four, although the fractal dimension varies as a function of parameter values. Typically, at parameter values near transitions to a limit cycle, the fractal dimension tends to be lower (slightly larger than two). The structure of these lower-dimensional chaotic attractors can be visualized using a 3D plot of the phase space, as shown in Fig. 2, whereas the higher-dimensional attractors cannot be adequately represented through three variables.

The existence of a strange attractor does not imply that avulsion is random. Although the avulsion process is not perfectly periodic, there may still be a dominant frequency. The simplest way to check for dominant frequencies is through a Fourier transform of the time series of water-discharge partitioning through one of the four branches. When we take the Fourier transform of the time series of a limit cycle, we obtain a sharp peak at the frequency associated with the periodicity of the signal, with decaying peaks at the harmonics. Compared to the limit cycle, we find that power is distributed across a broad range of frequencies in the case of a strange attractor (Fig. 3*D*). Nevertheless, a broad peak can typically be observed in the Fourier spectrum. We observe that this peak frequency roughly corresponds to the average time between discharge peaks in the time series, but for a chaotic time series, there is enough variability in the timing between discharge peaks to make long-term prediction impossible. As a point of reference, we can find the frequency at which the upstream bifurcation avulses independently of the downstream bifurcations and the frequency at which the downstream bifurcations avulse independently of the upstream bifurcation. This is done by leaving either the upstream or downstream bifurcations unperturbed; the bifurcations are unstable in that any arbitrarily small perturbation to the symmetry grows, but without the initial perturbation, the bifurcations remain symmetric indefinitely. We do not find a consistent relationship between the dominant frequency of the strange attractor and the frequencies of the unstable limit cycles; it may line up closely with one of the unstable limit-cycle frequencies, or it may be higher or lower than either unstable limit-cycle frequency (*SI Appendix*, Fig. S6).

## Application to Stratigraphy

The lack of long observational time series presents a challenge to testing our model against field data. However, long records of delta dynamics are potentially accessible from stratigraphy. Our model does not generate stratigraphy directly, because we assume a single grain size, and we assume a fixed channel network (i.e., all avulsions are reoccupational). As a workaround, we can plot the Shields stress at the time of deposition, which should manifest itself via grain size and/or sedimentary structures. Additionally, statistics based on time variation of the bed-surface elevation are straightforward to compute using the model and have already been measured for experimental deltas and ancient systems recorded in stratigraphy (58).

We synthesize an image of a longitudinal slice of stratigraphy generated by one of the four branches, colored by Shields stress at the time of deposition (Fig. 4*A*). Assuming that Shields stress can be considered a proxy for grain size, we observe features that are qualitatively consistent with observed stratigraphy: most beds fine upwards (59), we observe both stratigraphically transitional and abrupt avulsions (9), and some “mud plugs” (deposits formed at low Shields stress) are preserved (60).

For a quantitative comparison, we can calculate a dimensionless compensation statistic *B*, the deposition rate is not constant through time and incorporates periods of erosion and/or stasis (58). However, averaged over sufficiently long time, the mean deposition rate at any given point approaches the average,

## Discussion

Our results based on a model of a simple delta network indicate that two-way coupling between upstream and downstream bifurcations in a network can lead to chaotic dynamics in the flux partitioning of the network. This is in contrast to the behavior of a single bifurcation, which may produce limit cycles but never chaos. Although through a different mechanism, the possibility for interaction between bifurcations in a network was previously noted for braided rivers (64). In our system, the physical mechanisms behind coupling are as follows: firstly, the more intuitive mechanism is the upstream-to-downstream coupling. A change in the water partitioning in the upstream bifurcation instantaneously changes the upstream boundary conditions for the downstream bifurcations by affecting both the width-to-depth ratio and the Shields stress. As shown by previous studies (43, 51, 53), bifurcation asymmetry is sensitive to both of these parameters. Secondly, downstream-to-upstream coupling occurs because avulsions in the downstream bifurcations propagate a signal upstream through the bed-elevation profiles of their feeder channels. This alters the branch slopes at the upstream bifurcation, thereby affecting the flow partitioning. We can conceptually understand how avulsion affects the bed elevation of the feeder channel as follows: typically, following an avulsion, the new channel is steeper than the old one. This locally steeper region of the bed-elevation profile produces a higher sediment-transport rate. Via the Exner equation, this causes erosion (or at least reduced deposition) that propagates upstream. In an extreme example, following the 1855 avulsion of the Yellow River, 5 to 10 m of erosion propagated up to 100 km upstream (8). We note that backwater effects are not included in our model; therefore, the downstream-to-upstream signal propagation is purely morphodynamic, a variant on the morphodynamic backwater concept (65, 66). In contrast, the upstream-to-downstream transfer of information is effectively instantaneous, because it takes place through the water discharge. Similarly, if we were to incorporate backwater effects in our model, this would provide a mechanism for instantaneous coupling from downstream to upstream through the hydrodynamics (67).

The two-way coupling between bifurcations in a delta network is a core concept of this work. Through appropriate analysis of experiments, high-fidelity models, or satellite/historical/stratigraphic records, we believe that it should be possible to obtain evidence of this coupling in physical delta-network systems. For example, a validating experiment could be based on a fixed-width/location distributary network, ideally under flow conditions that promote bifurcation stability but limit the formation of alternate bars (48), which have been shown to complicate discharge-partitioning dynamics in single-bifurcation experiments (45). We also expect two-way coupling between bifurcations to occur under more general conditions, e.g., within delta networks where channel widths and network geometry can change through time.

Our study is not the first to use the BRT model for a system involving multiple coupled bifurcations; Kleinhans et al. (49) modeled a network consisting of 15 coupled bifurcations and confluences. In their study, the modeled river network evolved toward a frozen state with highly asymmetric flux partitionings at bifurcations. Similarly, in a Delft3D model, delta networks became frozen when they reached the edge of the computational domain, which caused deposition to cease (68). These studies highlight the importance of long-term deposition in sustaining avulsion dynamics (44). Indeed, when we run our simple delta-network model with a bypass fraction of 1, the system inevitably evolves toward a frozen state, i.e., a fixed point, with a potentially asymmetric flux partitioning.

We find that differing yet physically reasonable parameters can result in qualitatively different dynamics. Whereas the stability of the symmetric fixed point depends on model parameters such as the width-to-depth ratio in a straightforward way, the occurrence of periodic vs. chaotic dynamics is less intuitive; for instance, we find that stable periodic behavior can occur at parameter values that are book-ended by chaotic behavior. So-called “periodic windows” are common in chaotic systems. Our results suggest that a change to the forcing conditions of a delta (e.g., due to climate or tectonics) could cause a change from predictable to unpredictable behavior or vice versa.

In the case that a delta network is within a chaotic regime, what are the bounds on its predictability? The Lyapunov time (

Predicting delta-network dynamics to five characteristic avulsion times or even just one is beyond the reach of existing state-of-the-art models such as Delft3D. Our simple estimate above for the time horizon at which a network becomes unpredictable suggests that there is substantial room for improvement in short to medium-term forecasting of delta networks before fundamental constraints on predictability take over. Predictions could be improved not only by increasing model resolution and better constraining model parameters but also by increasing our understanding of the relevant physical processes in order to improve the fundamental equations we choose for modeling delta dynamics.

Even when the dynamics are chaotic and therefore unpredictable in the long-term, they are far from random. Time series obtained from the model possess a strong degree of structure despite being nonperiodic. For example, we found that both unusually long and short waiting times between branch reoccupations are rare relative to random (*SI Appendix*, Fig. S7). Given a sufficiently long averaging window, statistics such as averages, the power spectra, and range of asymmetry values can be predicted. For instance, the sedimentation rate at any point matches the basin-wide subsidence rate when averaged over a window longer than the compensation time scale. The time scale-dependent compensation statistics we obtain from the model can be contrasted with the prediction of a compensation index of *D*). Even though chaos implies that the long-term evolution of delta networks is unpredictable, the statistical behavior is predictable (69). This is analogous to how long-term predictions of the weather are impossible, and yet one can predict the statistics of the weather (i.e., climate) on time scales well beyond the weather-predictability horizon. The compensation time scale is a measure of the time required to transition from delta “weather” to delta “climate.”

Our results show that highly complex behavior can occur even in a highly idealized system with many simplifying assumptions. In particular, we use the assumption of a fixed network geometry, implying that all avulsions are reoccupational. Although this is clearly an oversimplification of natural delta networks, a growing body of literature shows that avulsions commonly exploit previously formed channels (6, 8, 31, 70, 71). Additionally, our simulations are for a delta network consisting of just two orders of bifurcations. We have found that adding an additional order of bifurcations does not fundamentally alter the chaotic dynamics described here and produces a similar predictability horizon, although the fractal dimension is higher (*SI Appendix*, Fig. S8). Given the simplifying assumptions of our model, and the simple delta-network geometry we chose, we expect chaos to be common in delta networks in the field, which have more degrees of freedom than our model.

Our example of a simple delta network demonstrates a particular way in which chaos can arise due to coupling between processes that are not chaotic in themselves. For instance, in the Lorenz system, it is the nonlinear coupling between the state variables that gives rise to chaos. Given that most geomorphological systems are nonlinear, dissipative, and involve many coupled processes, chaos is likely common in these systems (38). Models from previous studies may be chaotic: Murray et al. (19) found that their braided river model exhibits “apparently unpredictable changes in configuration indefinitely,” despite unchanging external forcing. In Delft3D delta simulations, varying the white-noise initial condition resulted in deltas that differ in details but whose gross-scale morphology remains similar (72). These examples are suggestive of sensitive dependence on initial conditions, even though neither formally demonstrated the existence of chaos. We propose that it is useful to know whether there are fundamental limitations to the prediction of a given geomorphic system. Whenever detailed predictions are required, we suggest that varying the initial conditions within their range of uncertainty should be included as part of a Monte Carlo analysis, in addition to the more common approach of quantifying uncertainty due to parameters. While sensitive dependence on initial conditions ultimately limits our ability make accurate predictions, estimation of the Lyapunov exponent provides a clear, quantitative measure of how we expect prediction accuracy to decay.

## Conclusions

We present a simple model of a delta network consisting of an upstream bifurcation coupled to two downstream bifurcations. We find: 1) While an individual bifurcation produces periodic (nonchaotic) dynamics, two-way coupling between upstream and downstream bifurcations can produce chaos. 2) Our model produces qualitatively realistic stratigraphy and generates scale-dependent compensation statistics that are consistent with ancient and experimental delta strata. 3) The Lyapunov exponent estimated from our model implies that the loss of predictive capacity due to chaos occurs on a time scale comparable to the avulsion time. Translated to stratigraphy, this implies a loss of prediction capacity on a vertical scale of a few channel depths. 4) While chaos implies that in the long-term, delta “weather” (precise configuration) is unpredictable, delta “climate” (statistical description) can be predicted. The required time for averages to become predictable is given by the compensation time scale. 5) Based on estimates of the predictability horizon, we suggest that there is substantial room for improvement in the modeling of delta networks before fundamental limits to predictability are encountered.

## Materials and Methods

### Model Formulation.

We consider a simple distributary network with geometry as follows: an upstream branch, labeled

Each branch is described by a bed-elevation long profile,

First, we consider the evolution of these bed-elevation long profiles. The equation for mass conservation of sediment reads:

The sediment flux is calculated using a generic sediment transport formula, which reads:

The equation for the Shields stress can be written:**4**–**6** into Eq. **3**, we can see that the evolution through time of the bed elevation profile in each branch is a function of itself and the water discharge.

Next, we seek a set of nodal conditions at each bifurcation, through which the coupling between the branches of the network occur.

Firstly, the assumption of water discharge continuity can be expressed:

Next, following the approach of BRT, we allow for lateral interaction between neighboring branches occurring over a length

We assume that the water level at the entrance to the bifurcation is constant. This is written:**6**, in which the bed slope **6**, **7**, and **8** must be solved iteratively to find the discharge entering each branch and the corresponding depths. Subsequently, Eqs. **4** and **5** can be evaluated to find

Upstream cells are allowed to exchange sediment with their neighbor. This transverse sediment flux is computed:

The bed elevation of the upstream cells evolve through time according to mass conservation:

Note that in the above equations, *i* = 1 to 3. We assume that the slope **4**–**6**.

Next, we consider the downstream boundary to branches *i* = 4 to 7. Following the approach of SPV (44), we introduce the bypass fraction F. The bypass fraction is the fraction of the input sediment flux

Finally, the upstream boundary to the entire system can be specified through the water discharge, sediment flux, and grain size, or, equivalently, the width-to-depth ratio, Shields stress, and slope.

## Data Availability.

Model code data have been deposited in GitHub (https://github.com/salterg/bifurcation_network).

## Acknowledgments

This material is based upon work supported by NSF Graduate Research Fellowship Grant 00039202. We thank Rick Moeckel, Crystal Ng, and Andy Wickert for helpful discussions and two reviewers for constructive comments.

## Footnotes

↵

^{1}Present address: Division of Geological and Planetary Sciences, California Institute of Technology, Pasadena, CA 91125.- ↵
^{2}To whom correspondence may be addressed. Email: gsalter{at}caltech.edu.

Author contributions: G.S., V.R.V., and C.P. designed research; G.S. performed research; G.S., V.R.V., and C.P. analyzed data; and G.S., V.R.V., and C.P. wrote the paper.

The authors declare no competing interest.

This article is a PNAS Direct Submission.

This article contains supporting information online at https://www.pnas.org/lookup/suppl/doi:10.1073/pnas.2010416117/-/DCSupplemental.

- Copyright © 2020 the Author(s). Published by PNAS.

This open access article is distributed under Creative Commons Attribution-NonCommercial-NoDerivatives License 4.0 (CC BY-NC-ND).

## References

- ↵
- I. Overeem,
- J. P. M. Syvitski

- ↵
- Z. Tessler et al.

- ↵
- A. J. F. Hoitink et al.

- ↵
- H. Michael et al.

- ↵
- A. Tejedor et al.

- ↵
- D. A. Edmonds,
- C. Paola,
- D. C. Hoyal,
- B. A. Sheets

- ↵
- D. Mohrig,
- P. L. Heller,
- C. Paola,
- W. J. Lyons

- ↵
- ↵
- H. Jones,
- E. Hajek

- ↵
- E. Stouthamer,
- H. J. Berendsen

- ↵
- M. G. Kleinhans,
- R. I. Ferguson,
- S. N. Lane,
- R. J. Hardy

- ↵
- ↵
- V. Ganti,
- A. J. Chadwick,
- H. J. Hassenruck-Gudipati,
- M. P. Lamb

- ↵
- A. J. Chadwick,
- M. P. Lamb,
- V. Ganti

- ↵
- J. C. Roseberry,
- M. W. Schmeeckle,
- D. J. Furbish

- ↵
- B. Werner

- ↵
- L. E. Hasbargen,
- C. Paola

- ↵
- T. J. Perron,
- S. Fagherazzi

- ↵
- ↵
- D. J. Furbish

- ↵
- C. Scheidt,
- A. M. Fernandes,
- C. Paola,
- J. Caers

- ↵
- S. H. Strogatz

- ↵
- ↵
- ↵
- P. Ghilardi,
- R. Rosso

- ↵
- ↵
- J. D. Phillips

- ↵
- ↵
- A. D. Howard,
- T. R. Knutson

- ↵
- ↵
- ↵
- H. Seybold,
- J. S. Andrade,
- H. J. Herrmann

- ↵
- E. A. Hajek,
- M. A. Wolinsky

- ↵
- M. Liang,
- V. R. Voller,
- C. Paola

- ↵
- ↵
- ↵
- K. Ito

- ↵
- T. A. Cross

- R. Slingerland

- ↵
- ↵
- J. D. Phillips

- ↵
- K. Montgomery

- ↵
- ↵
- M. Bolla Pittaluga,
- R. Repetto,
- M. Tubino

- ↵
- G. Salter,
- C. Paola,
- V. R. Voller

- ↵
- G. Salter,
- V. R. Voller,
- C. Paola

- ↵
- S. Miori,
- R. Repetto,
- M. Tubino

- ↵
- M. Kleinhans,
- H. Jagers,
- E. Mosselman,
- C. Sloff

- ↵
- W. Bertoldi,
- L. Zanoni,
- S. Miori,
- R. Repetto,
- M. Tubino

- ↵
- M. Kleinhans,
- T. De Haas,
- E. Lavooi,
- B. Makaske

- ↵
- N. Gupta,
- M. G. Kleinhans,
- E. A. Addink,
- P. M. Atkinson,
- P. A. Carling

- ↵
- M. Bolla Pittaluga,
- G. Coco,
- M. G. Kleinhans

- ↵
- M. Redolfi,
- G. Zolezzi,
- M. Tubino

- ↵
- D. Edmonds,
- R. Slingerland

- ↵
- M. Redolfi,
- G. Zolezzi,
- M. Tubino

- ↵
- D. L. Turcotte

- ↵
- D. J. Jerolmack

- ↵
- M. Casartelli,
- E. Diana,
- L. Galgani,
- A. Scotti

- ↵
- E. A. Hajek,
- K. M. Straub

- ↵
- J. Allen

- ↵
- E. Hajek,
- P. Heller

- ↵
- K. M. Straub,
- C. Paola,
- D. Mohrig,
- M. A. Wolinsky,
- T. George

- ↵
- Y. Wang,
- K. M. Straub,
- E. A. Hajek

- ↵
- S. Trampush,
- E. A. Hajek,
- K. Straub,
- E. Chamberlin

- ↵
- F. Schuurman,
- M. G. Kleinhans,
- H. Middelkoop

- ↵
- M. De Vries

- ↵
- D. Hoyal,
- B. Sheets

- ↵
- M. P. Lamb,
- J. A. Nittrouer,
- D. Mohrig,
- J. Shaw

- ↵
- D. Edmonds,
- R. Slingerland,
- J. Best,
- D. Parsons,
- N. Smith

- ↵
- ↵
- A. Aslan,
- W. J. Autin,
- M. D. Blum

- ↵
- J. M. Valenza,
- D. A. Edmonds,
- T. Hwang et al.

- ↵
- D. A. Edmonds,
- R. L. Slingerland

## Citation Manager Formats

## Article Classifications

- Physical Sciences
- Earth, Atmospheric, and Planetary Sciences