key: cord-0004780-1mv1r2i0 authors: Volz, Erik title: SIR dynamics in random networks with heterogeneous connectivity date: 2007-08-01 journal: J Math Biol DOI: 10.1007/s00285-007-0116-4 sha: b02f563c662d2d521c4e4c34fdb73190b74280a0 doc_id: 4780 cord_uid: 1mv1r2i0 Random networks with specified degree distributions have been proposed as realistic models of population structure, yet the problem of dynamically modeling SIR-type epidemics in random networks remains complex. I resolve this dilemma by showing how the SIR dynamics can be modeled with a system of three nonlinear ODE’s. The method makes use of the probability generating function (PGF) formalism for representing the degree distribution of a random network and makes use of network-centric quantities such as the number of edges in a well-defined category rather than node-centric quantities such as the number of infecteds or susceptibles. The PGF provides a simple means of translating between network and node-centric variables and determining the epidemic incidence at any time. The theory also provides a simple means of tracking the evolution of the degree distribution among susceptibles or infecteds. The equations are used to demonstrate the dramatic effects that the degree distribution plays on the final size of an epidemic as well as the speed with which it spreads through the population. Power law degree distributions are observed to generate an almost immediate expansion phase yet have a smaller final size compared to homogeneous degree distributions such as the Poisson. The equations are compared to stochastic simulations, which show good agreement with the theory. Finally, the dynamic equations provide an alternative way of determining the epidemic threshold where large-scale epidemics are expected to occur, and below which epidemic behavior is limited to finite-sized outbreaks. Contact patterns constitute an important aspect of heterogeneity within a population of susceptible and infectious individuals, but it has been a difficult factor to incorporate into epidemiological models. Compartment models can capture many aspects of population heterogeneity, such as with respect to heterogeneous susceptibility and infectiousness [1, 9, 32] , however such models usually assume individuals mix homogeneously within each category. In contrast, the contact patterns responsible for the spread of many infectious diseases tend to be characterized by constant relationships marked by high levels of heterogeneity in the number of contacts per individual. An alternative approach is to model a population of susceptibles and infecteds and the contact patterns among them as a static random network [13, 18, 26, 31] . This approach has generated a new category of epidemiological models in which epidemics spread from node to node by traversing network connections [2, 8, 17, 19, 24, 25, 29, 30, 33] . Random networks with specified degree distributions have been proposed as a simple but realistic models of population structure. This case has the advantage of being well understood mathematically. The expected final size of epidemics in random networks with a given degree distribution has been solved exactly [19, 24] . The network approach has the advantage that the mathematics of stochastic branching processes [4, 15, 34] can be brought to bear on the problem. This allows for precise descriptions of the mean outbreak size well as the final size. A shortcoming of the network model is that it is difficult to describe the explicit dynamical behavior of epidemics on networks. The mean outbreak size is easy to calculate, yet the dynamic epidemic incidence, that is the number of new infecteds at a time t, has been difficult to derive. Simulation has been used in this case [12] . Heterogeneity in the number of contacts within networks makes it difficult to derive differential equations to describe the course of an epidemic. Nevertheless, several researchers [5, 7, 27, 28] have been successful modeling many of the dynamical aspects of network epidemics, particularly in the early stage where asymptotically correct equations for disease incidence are known. These solutions break down, however, when the finite size of a population becomes a significant factor. Other work has focused on pair approximation and moment closure [6, 11] . These methods produce viable approximations to SIR dynamics, but are typically high-dimensional and computationally intensive. Here we introduce an alternative strategy based on a low-dimensional system of nonlinear ordinary differential equations. This model can be used to solve for epidemic incidence at any time, from an initial infected to the final size as well as other quantities of interest. We treat the simplest possible case of the SIR dynamics with constant rate of infection and recovery. Section 2 describes the model. Several examples are given in Sect. 3, and Sect. 3.1 compares the analytical results to stochastic simulations. The networks considered here are random networks with an arbitrary degree distribution p k ( p k being the probability of a random node having degree k) [23, 26] . Nodes can be in any of three exclusive states: susceptible (S), infectious (I), or recovered (R). The dynamics are as follows. When a node is infectious, it will transmit infection to each of its neighbors independently at a constant rate r . Infectious nodes become recovered at a constant rate µ, whereupon they will no longer infect any neighbors. This will be made precise in the next section. It is desirable to determine the dynamics of the number of susceptibles and infecteds and to develop equations in terms of those quantities. This, however, turns out to be intractable due to heterogeneity in the number of contacts. The problem can be resolved by developing equations in terms of dynamic variables representing network-based quantities, for example, the number of connections to susceptible or infectious nodes at a time t. The network-and node-based quantities are defined in the next section. To bridge the divide between connection-and node-based quantities, a mathematical device known as a probability generating function (PGF) [34] is extremely useful. The PGF has many useful properties and is frequently used in probability theory and the theory of stochastic branching processes. Given a discrete probability density p k , the PGF is defined as the series: (1) The variable x in the generating function serves only as a place-holder. To illustrate the utility of this device, consider the possibility that the probability of a node being infected, say λ, is compounded geometrically according the node's degree. Then, the probability of a degree k node being susceptible is (1 − λ) k , that is, the probability of not being infected along any of k connections. If the hazard is identical for all nodes, the cumulative epidemic incidence (the fraction of nodes infectious or recovered) will be Table 1 gives a summary of the parameters used in the model. An undirected network can be defined as a graph G = {V, E} consisting of a set of vertices V corresponding to the nodes in the network, and a set of edges E with elements of unordered pairs of vertices, {a, b} where a, b ∈ V . Two vertices a, b are said to be neighbors or neighboring each other or simply connected if {a, b} ∈ E. For the purposes of this model, the terms "vertex" and "node" will often be used interchangeably. The networks considered here can be generated by the Configuration Model(CM) [22] . For our purposes, the important aspect of CM networks is that the probability of being connected to a node is proportional to the degree of the node. Denote the degree of a node v ∈ V as d v . Then given an edge {a, x} ∈ E, the probability that 296 E. Volz Table 1 Parameters and dynamic variables for the network SIR model Note that this allows multiple edges to the same node as well as loops from a node to itself, however the existence of multiple edges and loops is exceedingly rare for large sparse random networks such that results based on this case can be safely applied to networks without multiple edges. Networks of this type can be generated by a the following algorithm: 1. To each node v ∈ V assign an i.i.d. degree δ v from distribution p k 2. Generate a new set X of "half-edges" with δ v copies of node v for all nodes 3. Insure X has an even number of elements, for example, by deleting a uniform random element if odd. 4. While X is not empty, draw two elements v 1 , v 2 uniformly at random and create edge {v 1 , v 2 }. At any point in time, a vertex can be classified as susceptible, infectious, or recovered. Let S, I, and R denote the disjoint sets of vertices classified as susceptible, infectious, or recovered respectively. J = I ∪ R will denote the set of infectious or recovered nodes. S, I, and R will denote the fraction of nodes in the sets S, I, and R respectively. The cumulative epidemic incidence will be the fraction of nodes in set J . As stated in the previous section, infectious vertices a ∈ I will infect neighboring susceptible vertices b ∈ S at a constant rate r . Infectious vertices will become recovered (move to set R) at a constant rate µ. Although the network is undirected in the sense that any two neighboring vertices can transmit infection to one another, we wish to keep track of who infects who. Therefore, for each edge {a, b} ∈ E, let there be two arcs, which will be defined to be the ordered pairs (a, b) and (b, a). Let A denote the set of all arcs in the network. The first element in the ordered pair (a, b) will frequently be called the ego and the second element the alter. A XY will denote the subset of arcs such that ego ∈ X and alter ∈ Y . A X will denote the subset of arcs such that ego ∈ X . M XY = #{A XY }/#{A} will denote the fraction of arcs in the corresponding set A XY . Table 2 Network-based dynamic variables for the network SIR model • θ := The fraction of degree one nodes that remain susceptible at time t • p I := M SI /M S . The probability that an arc with a susceptible ego has an infectious alter • p S := M SS /M S . The probability that an arc with a susceptible ego has a susceptible alter Table 3 A summary of the nonlinear differential equations used to the describe the spread of a simple SIR type epidemic through a random network. The degree distribution of the network is generated by g(x)θ For example, two variables will be especially important in the derivations that follow. M SS is the fraction of arcs with a susceptible ego and a susceptible alter. M S I is the fraction of arcs with a susceptible ego and and infectious alter. M S will be the fraction of arcs with a susceptible ego and an alter of any type. Our objective is to develop a deterministic model to describe epidemic dynamics expressed with a low-dimensional system of differential equations. At first, this goal may seem incompatible with network-SIR dynamics described in the previous section. Infection spreads along links in a random network, which implies the epidemic incidence at any time as well as the final size must also be random, depending on the particular structure of a given random network. This is true, however it is possible to avoid such considerations by focusing on epidemic dynamics in the limit as population size goes to infinity. This strategy has been used in previous work to calculate the expected final size of epidemics in infinite random networks [24] expressed as a fraction of the total population size. A similar strategy is followed here by considering the fraction of nodes in sets S, I, and R, after a small fraction nodes are infected initially in a susceptible population. 1 The conclusion is the system of equations given in Table 3 in terms of the dynamic variables given in Table 2 . The dynamics predicted by these equations are compared to stochastic simulations with large but finite networks in Sect. 3.1. Consider a susceptible node ego at time t with a degree k. Then there will be a set of k arcs {(ego, alter 1 ), (ego, alter 2 ), · · · , (ego, alter k )} corresponding to ego. We will assume that for each arc (ego, alter i ) there will be a uniform probability p I = M S I /M S that alter i is infectious. Then there is an expected number kp I arcs (ego, alter) such that alter is infectious. In a time dt, an expected number rkp I dt of these will be such that the infectious alter transmits to ego. Consequently, the hazard for ego becoming infected at time t is 298 E. Volz Now let u k (t) represent the fraction of degree k nodes that remain susceptible at time t, or equivalently the probability that ego in the previous example is susceptible. Using Eq. 4, Subsequently we will use the symbol θ to denote Given θ , it is easy to determine the fraction of nodes which remain susceptible at a time t. This equation makes use of the generating function g(·) for the degree distribution which greatly simplifies this and subsequent equations. The dynamics of θ are dependent on the hazard λ 1 . Unfortunately, this does not completely specify the dynamics of θ and by extension S, which also depends on the variable p I . The derivation of the dynamics of p I follows. Our goal is to put Eq. 8 in terms of the variables θ, p S , p I and the PGF g(·). M S is easily placed in terms of these variables. M S I follows easily. Next considerṀ S . In time dt, −Ṡ nodes become infectious. Since S = g(θ ), CalculatingṀ S I requires careful consideration of the rearrangement of arcs among sets A SS and A S I as −Ṡ nodes become infected in a small time interval. Since the hazard of becoming infected is proportional to the number of arcs to an infectious alter, a susceptible node will be selected with probability proportional to the number of arcs from the node to infectious nodes. To clarify subsequent calculations, I will introduce the notation δ XY to represent the average degree of nodes in set X , selected with probability proportional to the number of arcs to nodes in set Y , not counting one arc to nodes of type Y . For example, if we select an arc (ego ∈ X, alter ∈ Y ) uniformly at random out of the set of arcs from nodes in set X to nodes in set Y (A XY ), and follow it to the node in set X , (ego), then δ XY will represent the average number of arcs (ego, alter ) not counting the arc we followed to ego. This is commonly called the "excess degree" of a node [20] . Furthermore, δ XY (Z ) will be as δ XY but counting only arcs from ego to nodes in set Z , (ego, alter ∈ Z ). To calculateṀ S I we need to first calculate δ S I , and for this it is necessary to derive the degree distribution among susceptible nodes. Considering two arcs (ego, alter 1 ) and (ego, alter 2 ) with ego ∈ S, we must suppose that the event that alter 1 ∈ X is independant of the event that alter 2 ∈ Y . 2 A consequence of this is that arcs from a susceptible ego to nodes in sets S, I, R are distributed multinomially with probabilities p S , p I , and p R = 1 − p S − p I respectively. Let d ego (X ) be the r.v. denoting the number of arcs from ego to nodes in set X . Letting c normalize the distribution, and letting the dummy variables x S , x I , and x R correspond to the number of arcs from a susceptible ego to an alter in sets S, I, R respectively, the degree distribution for susceptible nodes will be generated by 300 E. Volz Using the multinomial theorem this becomes where c = k p k θ k ( p S + p I + (1 − p S − p I )) k = g(θ ) normalizes the distribution. The degree distribution for susceptible nodes selected with probability proportional to the number of arcs to infectious nodes will be generated by the following equation. This distribution is usually called the excess degree distribution. For example, selecting a random arc in set A S I , following it to the susceptible node, and then counting the composition of arcs out of that node (not counting the one we arrived on) would give rise to this distribution. Note that this equation does not count one arc to infectious nodes (the one used to select the susceptible node). Because arcs are distributed multinomially to nodes in sets S, I, R, we have g SS (x S , x I , x R ) = g S I (x S , x I , x R ), which is easy to verify by repeating the calculation in Eq. 14. A useful property of PGF's is that the mean of the distribution they generate can be calculated by differentiating and evaluating with the dummy variables set to one [34] . Now using Eqs. 13 and 14, we have the following results. As a fraction −Ṡ nodes leave set S in time dt, the fraction of arcs between S and I, M S I is reduced by the fraction of arcs from infectious nodes to the −Ṡ newly infectious nodes (recall M S I = M I S ). Therefore M S I is reduced at rate −Ṡδ S I (I )/g (1) . Because δ S I (I ) does not count the arc along which a node was infected, M S I is also reduced at a rate r M S I to account for all arcs which have an infectious ego which transmits to the susceptible alter. And in time dt, µI nodes become recovered. The average number of arcs in A I S per infectious node is proportional to M S I /I . Then M S I is reduced at a rate µI (M S I /I ) = µM S I . The quantity M S I is also increased, as new infected nodes have links to susceptible nodes. A newly infectious node will have on average δ S I (S) arcs to susceptible nodes, so M S I is increased at a rate −Ṡδ S I (S)/g (1) . To Then applying Eqs. 16, 17, and 11 we havė Finally, it is necessary to determine the time derivative of M S . = (−r p I θg (θ ) − r p I θ 2 g (θ ))/g (1) Now applying Eqs. 3 9, 18 , and 19 to Eq. 8 we solve forṗ I in terms of the PGF and θ . This equation makes use of the variable p S which changes in time. Deriving the dynamics of this variable will complete the model. This calculation is very similar to that forṗ I .ṗ The calculation forṀ SS is very similar to that forṀ S I . Newly infected nodes have on average δ S I (S) arcs to other susceptibles, so that = −2r p I p S θ 2 g (θ )/g (1) , (22) where the factor of 2× accounts for two arcs per edge. Now applying Eqs. 9, 19, and 22 to Eq. 21, we havė The complete system of equations is summarized in Table 3 . The fraction of infectious nodes can be solved by introducing a fourth dynamic variable. The infectious class increases at a rate −Ṡ (Eq. 11) and decreases at a rate µI . Thereforeİ = r p I θg (θ ) − µI. An advantage of dynamic modeling of epidemics in networks is that the time-evolution of variables besides cumulative incidence can be calculated. Above it was shown how to calculate the degree distribution among susceptible nodes (eqn. 13). Additionally, the degree distribution among nodes which are either infectious or recovered (set J ) can be calculated by taking the complement. If a small fraction of the nodes in the network are selected uniformly at random and initially infected, we can anticipate the following initial conditions. The fraction of arcs with infectious ego will also be M I = , and since is small, there is low chance of two initial infecteds being connected. Therefore M S I ≈ M I = . θ , which can be interpreted as the fraction of degree one nodes remaining susceptible will be 1 − . We'll assume that no initial infecteds are connected to eachother, so that M Epidemic dynamics can fall into one of two qualitatively different regimes. Below a threshold in the ratio r/µ, the final size (I ∞ ) is necessarily proportional to the fraction of initial infectious nodes: I ∞ ∝ . But above this threshold, epidemics occur, and necessarily occupy a fraction of the population even as → 0. As per Eq. 4, the number of new infections in a small time interval is proportional to p I . This is in contrast to compartment models in which the number of new infections is proportional the current number of infectious. Ifṗ I (t = 0) < 0, an epidemic will necessarily die out without reaching a fraction of the population. The epidemic threshold occurs wherė Applying the initial conditions above and considering 1 giveṡ Rearranging yields the critical ratio r/µ in terms of the PGF. The epidemic threshold in Eq. 28 can also be put in terms of the the transmissibility, which is the probability that an infectious ego will transmit infection to a given alter. If ego is infectious for a duration T , the probability of transmitting to a given alter is 1 − e −r T . Integrating over an exponentially distributed duration of infectiousness T , the transmissibility τ is calculated to be Then rearranging Eq. 28 yeilds the epidemic threshold in terms of τ . This is consistent with previous results based on bond-percolation theory [24] . 4 The model has been tested on several common degree distributions: -Poisson: p k = z k e −z k! . This is generated by -Power-law. For our experiments, we utilize power-laws with exponential cutoffs κ: is the nth polylogarithm of x. This is generated by -Discrete exponential: p k = (1 − e −1/λ )e −k/λ . This is generated by Figure 1 shows the cumulative incidence for each of the degree distributions (31), (32) , and (33), with a force of infection r = .2 and recovery rate µ = .1. Initially = 10 −4 nodes are infected. The parameters of the degree distributions were chosen so that each network has an identical average degree of 3. That is, the density of connections in each network is the same. Nevertheless, there is widely different epidemic behavior due to the different degree distributions. Consistent with previous research, the degree distribution has a great impact on the final size of the epidemic [19, 24] . More importantly, the three networks exhibit widely varying dynamical behavior. The power law network experiences epidemics which accelerate very rapidly. Such epidemics enter the expansion phase (the time at which cumulative incidence increases at its maximum rate) virtually as soon as the first individual in the network is infected. In contrast, the Poisson network experiences a lag before the expansion phase of the epidemic. The discrete exponential network has behavior similar to that of the power law but not as extreme; the expansion phase is early, but not instantaneous, and the final size is not as great as the power law. These observations are consistent with the findings in [5] that the timescale of epidemics shortens with increasing contact heterogeneity. This has important implications for intervention strategies, as it is often the case that interventions are planned and implemented only after a pathogen has circulated in the population for some time. If an epidemic were to occur in the power law network, there would be little time to react before the the infection had reached a large proportion of the population. Recall from Sect. 2.4 that below the epidemic threshold τ * , only small, finite-sized outbreaks will occur. Figure 2 shows the qualitatively different dynamical behavior of outbreaks below the phase transition for networks with a Poisson distribution. Below the phase transition, the final size is always proportional to the fraction of initial infecteds . Something offered by this model and not to the author's knowledge seen previously, is an explicit calculation for how the degree distribution of susceptibles evolves over the course of the epidemic. We expect the degree distribution to become bottom-heavy, as high degree nodes are gradually weeded out of the population of susceptibles. This is indeed observed in Fig. 3 for the Poisson trial described above. Recall that the degree distribution of susceptibles is generated by the multi-variate PGF (13) . The explicit degree distribution can be retrieved from Eq. 13 by differentiation. The following gives the probability that a susceptible node has m links. For example, applying this to the Poisson PGF (Eq. (31)) gives which is simply the Poisson distribution with an adjusted parameter z × θ . Another example is illustrated in Fig. 3 , which shows the degree distribution among susceptibles for the power-law network considered above. Fig. 3 The degree distribution for susceptible nodes where the epidemic size is 50, 75, and 100% of the final size, as well as degree distribution at the beginning of the epidemic. The degree distribution for the network as a whole is a power law with exponential cutoff (Eq. 32) Simulation of SIR on networks presents two challenges: A random network must be generated with the desired degree distribution. And, the stochastic rules that govern the transmission of disease at the microscopic scale must be well-defined. In the deterministic model of the previous sections, we described dynamics using the quantities S, I, and R to represent the fraction of nodes in each category. The stochastic dynamics will be described by S n , I n , and R n which will denote the random number of susceptible, infectious, and recovered nodes respectively. Below we will give evidence that the stochastic dynamics match the deterministic dynamics in terms of time-scale and final size. Future work should provide a formal characterization of this relationship. We hypothesize that stochastic dynamics converge to the deterministic model as the population size grows to infinity, that is, lim n→∞ S n /N = S with constant fraction of initial infecteds [3] . The algorithm proposed by Molloy and Reed [22] was used to generate the random networks in these experiments. Subsequent research has shown that imperfections can arise in the networks generated by this algorithm, but such biases should be tolerably small for these purposes [21] . The simulation dynamics are as follows: -A random network of size N is generated. -A node is chosen uniformly at random from the network as an initial infected. -An infected node v will recover after an exponentially distributed random time interval t µ ∼ E x p(µ). -When a node v is infected, each arc (v, x) has a time of infection t x drawn from an exponential distribution E x p(r ). If t x < t µ , node x is infected after time t x . Otherwise x is not infected by v. This process continues until there are no more infectious nodes. Figure 4 shows the results of 450 simulations for the Poisson random network considered in the last section (z = 3) with 10 4 nodes. The black dotted line represents an independent simulation trajectory. The thick, solid line that cuts through the dense mass of simulation trajectories is the analytical trajectory based on the equations in Table 3 . The initial conditions were chosen as in the previous section using = 10 −4 . Figure 5 shows a similar series of simulations for the power law degree distribution considered in the last section. In both cases, the analytical trajectory traverses the region with the highest density of simulation trajectories. The simulation trajectories also exhibit significant variability in the time required to reach the expansion phase and final size. This is largely due to the significant impact of random events early on in the epidemic. For example, an initial infected with a low average degree, or one which takes an inordinate amount of time to infect the next infected can markedly delay the onset of the expansion phase. Figure 6 shows the median-time incidence for the exponential and Poisson networks discussed in the last section. The data points show the median time required to reach a given cumulative incidence among 450 simulation trajectories. The solid line shows the analytical trajectory based on the system of equations given in Table 3 . Intuitively, the data points are showing the path of the most central trajectory from the swarm of simulation trajectories such as in Fig. 4 . The statistical properties of SIR epidemics in random networks have been understood for some time, but the explicit dynamics have been understood mainly through simulation. This paper has addressed this shortcoming by proposing a system of nonlinear ordinary differential equations to model SIR dynamics in random networks. 308 E. Volz Table 3 It should be noted that the SI dynamics are a special case of this model (µ = 0), in which case the ultimate extent of the epidemic is simply the giant component 5 of the network. The distribution of contacts, even holding the density of contacts constant, has enormous impact on epidemic behavior. This goes beyond merely the extent of the epidemic, but as shown here, the dynamical behavior of the epidemic. In particular, the distribution of contacts plays a key role in determining the onset of the expansion phase. The distribution dynamics from Eq. 13 and shown in Fig. 3 have important implications for vaccination strategies. Previous work [14, 16] has focused on determining the critical levels of vaccination required to halt or prevent an epidemic. It is usually taken for granted that contact patterns among susceptibles are constant. Furthermore, most widespread vaccinations occur only once an epidemic is underway. Future research could be enhanced by considering optimal vaccination levels when the epidemic proceeds unhindered for variable amounts of time. It is hoped that the distribution dynamics described in this paper will find applications beyond modeling heterogeneous connectivity. The dynamic PGF approach might be used to capture other forms of heterogeneity, such as of susceptibility, mortality, and infectiousness. Infectious Diseases of Humans: Dynamics and Control Epidemic models and social networks Stochastic Epidemic Models and their Statistical Analysis Branching Processes Dynamical patterns of epidemic outbreaks in complex heterogeneous networks A versatile ODE approximation to a network model for the spread of sexually transmitted diseases Epidemic spreading in complex networks with degree correlations Halting viruses in scale-free networks Mathematical epidemiology of infectious diseases. Model building, analysis and interpretation Random Graph Dynamics Modeling dynamic and network heterogeneities in the spread of sexually transmitted diseases Modelling disease outbreaks in realistic social networks Networks of sexual contacts: implications for the pattern of spread of hiv Containging bioterrorist smallpox The Theory of Branching Processes Emergency response to a smallpox attach: the case for mass vaccination The effects of local spatial structure on epidemiological invasions The web of human sexual contacts Network theory and sars: predicting outbreak diversity Network theory and SARS: predicting outbreak diversity Uniform generation of random graphs with arbitrary degree sequences A critical point for random graphs with a given degree sequence The size of the giant component of a random graph with a given degree sequence The spread of epidemic disease on networks The Structure and Dynamics of Networks Random graph models of social networks Epidemic spreading in scale-free networks Epidemic dynamics and endemic states in complex networks Handbook of Graphs and Networks: From the Genome to the Internet. chapter Epidemics and immunization in scale-free networks Modelling development of epidemics with dynamic small-world networks Exploring complex networks On the effect of population heterogeneity on dynamics of epidemic diseases Percolation on disordered networks as a model for epidemics The author gratefully acknowledges helpful suggestions from Lauren Meyers, Shweta Bansal, Kimberly Hopkins, Steven Strogatz, Michael Marder, and the anonymous reviewers.