A non-Markovian model of rill erosion

A non-Mark o vian mo del of rill erosion Mic hael Damron ∗ C.L. Win ter † Decem b er 2008 Abstract W e in tro duce a new mo del for rill erosion. W e start with a net work similar to that in the Dynamical Discrete W eb** and instantiate a dy- namics whic h makes the pro cess highly non-Mark ovian. The b ehavior of no des in the streams is similar to the b eha vior of Poly a urns with time-dep endent input. In this pap er we use a combination of rigorous argumen ts and sim ulation results to show that the mo del exhibits man y prop erties of rill erosion; in particular, no des whic h are deep er in the net work tend to switch less quickly . 1 In tro duction 1.1 Reinforcemen t and Rill Erosion Sto c hastic pro cesses with reinforcemen t are inheren tly non-Marko vian and there- fore ma y mo del some real phenomena more accurately than can their Mark o vian coun terparts. Reinforcement is a mechanism that provides a bias to a system, making it more likely to o ccup y states the more often those states are visited. Some well-studied examples include v ariations on the urn of P´ oly a (the original in tro duced in [4] and this and subsequen t mo dels studied, for example, in [1] and [9]) and reinforced random walks ([3, 15]). The infinite memory e xhibited in these examples can force a system to sp end most (or almost all) of its time in a small subset of its state space. Man y natural phenomena exhibit similar b eha vior; for instance, the ov erall pattern of erosion on a hillslop e is relativ ely stable once it is established, although small details of the pattern may c hange frequen tly and catastrophes that p ermanen tly alter it may o ccasionally o ccur. ∗ Courant Institute of Mathematical Sciences, 251 Mercer St., New Y ork, NY 10012, USA. Email: damron@cims.nyu.edu; Research supp orted in part by NSF grant number OISE- 0730136. † National Center for Atmospheric Research, 1850 T able Mesa Dr., Boulder, CO 80305, USA. Email: lwinter@ucar.edu 0 MSC2000: Primary: 60K35, 82D99 Secondary: 60G09 0 ** See arXiv:0704.2706 and arXiv:math/0702542 0 Keywords: erosion; rill erosion; P´ olya urn; exchangeabilit y; dynamical discrete web 1 W e inv estigate a discrete time, infinite-memory random pro cess defined on the no des and edges of an oriented diagonal lattice (Figure 1) that we prop ose as a simple model of hillslope erosion. The lattice starts out smooth in the sense that it has no edges initially , but it sprouts edges ev erywhere the instan t the pro cess starts, muc h as rain can start soil erosion ev erywhere on a hillslop e at once. Edges ma y connect an interior no de to t wo, one, or neither of the tw o no des directly abov e it. Exactly one edge descends from each in terior no de, and it p oin ts either left or right. At every no de and at every time step a simple tw o parameter reinforcing la w, based on the entire history of the net work ab o v e a giv en in terior no de, randomly determines the direction of the nodes descending edge and then is up dated. Ob vious mo difications of these statements apply to no des at the top or b ottom (if one exists) of the lattice. The curren t pattern of connections among nodes represents the present state of the process, and the patterns stabilit y – measured by the tendency of the same state, or one similar to it, to o ccur on subsequent iterations of the pro cess – represen ts the patterns strength as a memory . The degree of reinforcement is set b y tuning tw o parameters, r and α . At any given moment the current pattern is a collection of dendritic netw orks that app ears similar to drainage netw orks found in nature; indeed, lattice mo dels hav e often b een used to inv estigate the morphology of natural drainage netw orks (e.g. [17]). W e fo cus on the surficial dynamics of rill net works [10], rather than their morphology . Put in terms of erosion, we are more interested in the pro cess of erosion than we are in the result. The analogy b et ween our mo del and erosion, specifically rill erosion, is straigh tforward: r can b e interpreted as a rainfall rate (or equiv alently , as the rate of sediment generation) and α − 1 as the resistance of soil to erosion, while the reinforcement dynamics corresp ond to the ov erland flow of water and sedi- men t down a hill. Rills are small, ephemeral c hannels that transp ort sediment do wn hillslop es when it rains [18]. They form when rainfall and runoff dislo dge particles from the soil surface and transp ort them along flow paths gov erned b y v ariations in the surface roughness of soils and the soil’s ability to resist erosion. Flow depths in rills are typically on the order of a few centimeters or less, while the longest channels in rill netw orks can b e several meters long. Pro- cesses affecting rill erosion take place ov er timescales ranging from milliseconds to hours. The top ology of rill netw orks is relativ ely unstable when compared to larger scale natural drainage systems (of which rills may b e a part) like gulley systems and riv er basins. Rill net wor ks are most unstable at their tops where b ound- aries b et ween rills and inter-rill areas are not well defined and shift often, but connectivit y can change downhill as well, usually at a slow er rate than uphill. Some rills grow throughout a rainfall ev ent, others are filled by sediment and disapp ear, still others alternate. A detailed description of rill erosion 1) m ust accoun t for complicated in teractions among rainfall, soil prop erties, and top og- raph y [19], and 2) often dep ends on obtaining a set of physical parameters that are difficult to measure. Despite the high degree of complexity of rill erosion at small scales, at macro- 2 scopic scales it is principally determined by particle detachmen t, o verland flow, and sediment transp ort ([16]). In turn, each of flow, detachmen t, and transp ort dep ends critically on the rate of rainfall and the soils resistance to erosion. It is not completely surprising that our simple tw o parameter model exhibits some imp ortan t elements of the macroscopic b eha vior of rills formation. In fact, sim- ilar to rill erosion, each no de in the mo del netw ork switc hes direction infinitely man y times but the switching rate depends on p osition up or do wn hill. F ur- thermore, flo ods that carry unusually large amounts of w ater and catastrophes that significantly alter the flow pattern o ccur o ccasionally in the mo del, as they do in nature. 1.2 Definition of the Mo del Consider the vertices in the ev en sub-lattice of Z 2 whic h hav e second co ordinate non-p ositiv e. That is, the set Z 2 ev en = { ( x, y ) ∈ Z 2 : x + y is even and y ≤ 0 } and edges E 2 ev en = { < ( x, y ) , ( x + 1 , y − 1) > : ( x, y ) ∈ Z 2 ev en } ∪{ < ( x, y ) , ( x − 1 , y − 1) > : ( x, y ) ∈ Z 2 ev en } Let v = ( x, y ) b e a no de with left parent w 1 = ( x − 1 , y + 1), right parent w 2 = ( x + 1 , y + 1) (for those no des with s econd co ordinate 0, paren ts will not exist), and with left child ( x − 1 , y − 1), righ t c hild ( x + 1 , y − 1). W e will use the term depth k to refer to those no des with y co ordinate equal to 1 − k . Con versely , for any no de v the term depth( v ) will denote the numerical v alue of the depth of v . First w e describe the algorithm for the b eha vior of v heuristically . A t the end of the 0 th second, v receives I v (0) = r units of rain (but does nothing else). During the n th second, for n ≥ 1, the following sequence o ccurs: 1. v flips a coin, heads-biased with probabilit y P L v ( n ), whic h reflects T L v ( n ), the total ”sedimen t” load v has sent to its left child b y time n . 2. If this coin shows heads (tails), v sends its current input of sediment I v ( n − 1) to its left (righ t) c hild. v adds this num b er to the total sedimen t, T L v ( n ) ( T R v ( n )), it has sen t to the left (right) for all time. 3. v receives sedimen t load from 0,1, or 2 paren ts and receiv es r units of rain. Call the sum of these t wo I v ( n ). Increment time and return to step 1. The evolution of the no de’s b eha vior dep ends on t wo parameters: the rainfall rate r > 0 and a term α − 1 > 0 that resists c hange. T o make this precise, we mak e several definitions. W e start by initializing v ariables. F or eac h v ∈ Z 2 ev en , let 3 T L v (0) = T R v (0) = T v (0) = 0 and I v (0) = r , P L v (0) = P R v (0) = 1 / 2 F or eac h n ≥ 1, v ∈ Z 2 ev en , we define a Bernoulli v ariable, D L v ( n ), the biased coin, that is conditionally indep endent from v ertex to vertex giv en the v ariables { D L v ( i ) : v ∈ Z 2 ev en , i < n } , with parameter P L v ( n − 1). Next, let T L v ( n ) = I v ( n − 1) D L v ( n ) + T L v ( n − 1) T v ( n ) = I v ( n − 1) + T v ( n − 1) , T R v ( n ) = T v ( n ) − T L v ( n ) W e create the bias for the next coin: P L v ( n ) = T L v ( n ) + α T v ( n ) + 2 α = T L v ( n ) r + η T v ( n ) r + 2 η (1) where η = α/r compares the effect of the rain to the system’s inheren t resistance to c hange. In this paper, we shall alwa ys tak e r = 1 so that η = α . Last we define the input I v ( n ) = ( r depth( v ) = 1 I w 1 ( n − 1)(1 − D L w 1 ( n )) + I w 2 ( n − 1) D L w 2 ( n ) + r otherwise and the filtration F n = σ ( { D L v ( k ) : v ∈ Z 2 ev en , k = 1 , ..., n } ) See Figure 1 for an illustration of the pro cess at no de v . Denote by d v = ( d v (1) , d v (2) , ... ) the sequence of directions that no de v c ho oses (for example (L,R,L,...)). A t the end of time t , after all no des ha ve sen t loads to their children, we may up date certain edge v ariables. Define a sequence of edge configurations { ω n } n ≥ 0 , where for eac h n , ω n is a map from E 2 ev en → { 0 , 1 } , using the following rule. If the no de v = ( x, y ) has d v ( n ) = L then let ω n ( < x − 1 , y − 1 > ) = 1 , ω n ( < x + 1 , y − 1 > ) = 0 If, on the other hand, d v ( n ) = R then let ω n ( < x − 1 , y − 1 > ) = 0 , ω n ( < x + 1 , y − 1 > ) = 1 See Figure 1 for a realization of ω n . W e say that no des v and w are connected at time n if there exists a path of distinct adjacent edges e 1 , .., e m with ω n ( e i ) = 1 for all i so that e 1 connects v 4 Figure 1: Left: Input and output b eha vior at no de v . The darkened line seg- men ts indicate paths of sediment flo w. The light line represents a p oten tial flow path. Righ t: A small piece of a visualization of ω n for some n . Only edges e with ω n ( e ) = 1 are solid. to one of its children and e m connects w to one of its parents. Denote by C v ,n the set of vertices which are connected to v at time n and define the backw ard (uphill) comp onen t of v = ( x, y ) at time n by C + v ,n := C v ,n ∩ { ( x 0 , y 0 ) : y 0 ≥ y } Finally , let ω − 1 n = { e : ω n ( e ) = 1 } . 1.3 Regimes for η The parameter η pla ys an imp ortan t role in the b ehavior of the mo del. F or a fixed no de v (at depth k ) we hav e that for all n ≥ 1 lim η → 0 P ( d v ( n ) = L | d v (0) = R ) = 0 This indicates that when η is small the no de v chooses a direction at time 0 and has a high probabilit y of sticking to this direction for most v alues of n ≥ 1. Since this is true for each no de v , the evolution of { ω n } is somewhat simple. In the limit as η → 0, eac h no de picks a direction and stays with that direction for all time. That is, for each n ≥ 1, and for each finite subset E ⊂ E 2 ev en , lim η → 0 P ( ω − 1 0 (1) ∩ E = ω − 1 n (1) ∩ E ) = 1 (2) and the dynamics has no effect on the configuration in an y finite subset of Z 2 ev en . The configurations at any moment are the same as those in the discrete web ([2], [7]). In the other direction, as η → ∞ eac h node ”forgets” its history . That is, for each node v , the conditional probability giv en F n that it chooses left at time n + 1 is giv en in (1), and the limit of this quantit y is 1 / 2. By symmetry , 5 P ( d v ( n + 1) = L ) = 1 / 2 = P ( d v ( n + 1) = R ) and so lim η →∞ [ P ( d v ( n + 1) = L ) − P ( d v ( n + 1) = L |F n )] = 0 . Therefore the v ariables in any finite subset of { d v ( n ) : n ≥ 0 } conv erge in distribution to i.i.d. Bernoulli(1/2) v ariables. F urthermore, the v ariables in an y finite subset of { d v ( n ) : n ≥ 0 , v ∈ Z 2 ev en } conv erge in distribution to i.i.d. Bernoulli(1/2) v ariables. Intuitiv ely this holds b ecause distinct no des only interact with each other through their input and output loads, and b oth of these are ev en tually dominated b y large η . These statemen ts indicate that when η is large, the dynamics of our erosion mo del are similar to those in a netw ork in whic h each no de flips a fair coin at eac h time n , independently from site to site and from time to time, to determine in which direction to send sediment. Th us the configurations { ω n } resem ble those tak en from the dynamical discrete w eb ([11], [7]). Giv en the relation b oth these extreme cases hav e to the v ariables { ω n } , it is natural to view the present mo del (with 0 < η < ∞ ) as an in terp olation betw een the discrete w eb and the dynamical discrete web. Indeed, for each fixed n , the distribution of ω n is the same as that in the case of η = 0 at time n = 0 or that in the case of η → ∞ at any time n . (In b oth cases, all directions are c hosen by indep enden t fair coin flips). As we shall see in section 3.2.2, the mo del with 0 < η < ∞ can b e lik ened to the case η → ∞ in the follo wing wa y . Each lev el k is asso ciated to a measure θ k (defined in eq. (8)). Eac h no de at level k samples (non-indep enden tly) from this measure a v alue p v . F or an y sequence n 1 , n 2 , ..., n m of times and for any sequence x 1 , ..., x m of elemen ts from the set { L, R } , lim T →∞ P ( d v ( n 1 + T ) = x 1 , ..., d v ( n m + T ) = x m ) = p N L v (1 − p v ) N R where N L ( N R ) is the num b er of i for which d i = L . Because of this fact, we ma y view the mo del (for large time) as one in which each node fixes a Bernoulli parameter p v and flips a p v -biased coin independently each second n (but not indep enden tly from site to site) to determine the direction in which to mov e sedimen t. 1.4 Outline of the Paper In section 2 w e discuss the (relatively simple) b eha vior of no des at depth 1. Since these no des receive constant input load ov er time, we can use the w ell kno wn mo del of P´ oly a’s Urn to analyze their output. In section 3 we discuss the more complicated b eha vior of no des at depth at least 2. Here we make use of results of P emantle [14] for the time-dep enden t P´ oly a Urn. W e lo ok more closely at prop erties of the input load, of the output load, and of the dynamics of these lo wer-depth no des. 6 2 T op Lev el Since our top lev el no des are equiv alent to the mo del of P´ oly a’s urn, we recall basic facts of P´ oly a’s mo del Start with an urn containing R 0 red balls and B 0 blac k balls and draw one ball from the urn. Return this ball to the urn, along with another ball of the same color. After this round there are R 1 red balls in the urn and B 1 blac k balls in the urn, with either R 0 = R 1 or B 0 = B 1 . Rep eat this pro cess infinitely many times, creating sequences { R n } n ≥ 0 and { B n } n ≥ 0 so that for eac h n , P ( R n +1 − R n = 1 | R n , B n ) = R n R n + B n It is well kno wn that the fraction F R n = R n R n + B n has an almost sure limit and that this limit is distributed as β ( R 0 , B 0 ) (see e.g. [8]). Let v be a no de at depth 1. At the b eginning of each second, v receiv es an amoun t of sediment equal to 1 and this input load amoun t do es not c hange with time. The no de sends this load either to the right or left, dep ending on the bias rule in (1). W e are in terested in the fraction of total load the no de sends left (righ t) up to time n . T o this end, define the load fraction LF L v ( n ) = T L v ( n ) T v ( n ) , LF R v ( n ) = T R v ( n ) T v ( n ) , n ≥ 1 Theorem 2.1. The quantities LF L v ( n ) and LF R v ( n ) have limits as n → ∞ . These limits ar e r andom: they ar e distribute d as β ( η , η ) . Pr o of. W e will indicate the proof only for the case LF L v ( n ). An easy calculation sho ws that P L v ( n ) is a martingale w.r.t. {F n } and, since it is bounded for all n , it has an almost sure limit. Solving for the limiting distribution is similar to solving for the related quantit y in the standard P´ oly a urn model. See, for instance, [5]. This gives lim n →∞ P L v ( n ) = lim n →∞ T L v ( n ) + η T v ( n ) + 2 η = lim n →∞ T L v ( n ) T v ( n ) + η T v ( n ) 1 + 2 η T v ( n ) = lim n →∞ T L v ( n ) T v ( n ) = lim n →∞ LF L v ( n ) (3) b ecause η is constant w.r.t. n and T v ( n ) → ∞ . Note that the limiting distribution in Theorem 2.1 is supp orted on [0,1] and has no atoms. This implies that with probabilit y 1, the no de v switc hes states (L,R) infinitely often and that neither of these states is transient. This is quite unlik e the ”stic king” associated to the dynamics in the η → 0 limit (refer to (2)). The distribution from the ab o ve theorem for differen t v alues of η is pictured in Figure 2. F or 0 < η < 1 the limiting load fraction has a bimodal distribution, and for η > 1 the distribution is unimo dal, symmetric ab out 1 2 . This means 7 η = 1 2 η = 1 η = 2 Figure 2: Asymptotic distributions for LF L v ( n ): η = . 5 , 1 , 2 resp ectively . that when η is small, each no de is likely to hav e a relatively strong preference for one direction and that when η is large, each no de is lik ely to fa vor L and R somewhat equally . The case η = 1 giv es a uniform distribution. Here v is equally lik ely to hav e a strong directional preference as it is not to. 3 Lo w er Lev els The simplicity of b eha vior at the top level comes from the fact that each no de has an input load which is constant w.r.t. time. This is not true at low er levels. Eac h no de has an input load whose magnitude is non-trivially time dep enden t. T o make this more apparent, isolate an arbitrary no de v at depth k . If at time t = n , v is not connected to either of its parents in ω n , then its input load is 1 unit (coming only from rain). If, on the other hand, v is connected to at least one of its parents, then its input load will b e strictly greater than 1 unit. Therefore, the geometry of the connected comp onen ts of ω n determines the b eha vior of each no de. This relationship is complex for at least tw o reasons. First, not only do es the geometry of the net work influence no de b eha vior, the no de b ehavior in turn determines the future geometry of the netw ork. In this sense, our system generates its own randomness. Second, the method by which this randomness arises inv olves propagation. The geometry of no des at depth k − l at time m affects the b eha vior of no des at depth k at time n if and only if m = n − l . In other w ords, it takes l seconds for the output load from depth k − l to reach no des at depth k . In spite of these complications, we set out to analyze these lo wer level nodes. The no de v has an input load sequence I v = ( I v (1) , I v (2) , ... ), left out- put load sequence T L v = ( T L v (1) , T L v (2) , ... ), and output direction sequence d v = ( d v (1) , d v (2) , ... ). W e are interested in analyzing the nature of this in- put sequence, the nature of the output sequence, and the relationship b et ween the t wo. 3.1 Input Load Figure 3 shows a histogram of input load v alues for all no des at (a) depth 5, (b) depth 7, and (c) depth 8 at t = 300s with η = 1 (the precise v alue 8 Ro w 5 Ro w 7 Ro w 8 Figure 3: Load distribution for k = 5 , 7 , 8 resp ectively . of η do es not matter, as a consequence of Theorem 3.1). The sim ulation w as conducted with perio dic boundary conditions, with 10 6 no des per ro w, and with 10 rows. Therefore, the histogram for depth k at time n = 300s should closely appro ximate the probability mass function of the distribution of the input load for depth k at time n = 300s. One notices a few things. First, the supp ort of the distribution at depth k is integers in the interv al [1 , 1 2 k ( k + 1)]. Next, the mass function appears to decrease from load v alue 1 to a lo cal minim um at k − 1, to increase for a bit to a lo cal maximum, and then to decrease to the edge of its supp ort. Ab out 1 / 4 of no des are at the heads of rills, while the fraction of rills starting short of the top increases with depth. The ”bump” in the load distribution to the right of the v alue k − 1 appears to trav el to the right as depth increases. Lo oking at Figure 3, it is tempting to guess that the load distribution at a given lev el is a mixture of a distribution for loads that start at the top and one for loads that do not. Last, the different mass functions hav e sev eral common v alues. F or example, the probabilities for load v alues 1 to 4 are the same in each figure, and the probabilities for load v alues 1 to 6 are the same in the cen ter and right figures. W e presen t three structural theorems regarding the load distribution. The first gives basic information needed to make calculations, and the second giv es us the v alue of the first moment of the distribution. The third discusses a limiting measure for the family of loads that do not originate at the top. Because of the simplicit y of the first theorem, we state it without pro of. Theorem 3.1. Basic pr op erties of the lo ad distribution. L et n 0 b e a fixe d time and let v b e a no de at depth k . a. Al l r andom variables d v ( n 0 ) ar e i.i.d. with pr ob ability 1 2 of b eing L or R. b. The distribution of I v ( n 0 ) is later al ly tr anslation invariant (i.e. along the x-axis) and is invariant in time for n 0 ≥ k . c. The distribution of I v ( n 0 ) is e qual to the distribution of | C + w,n 0 | for any no de w with depth e qual to min ( n 0 , k ) . Ther efor e I v ( n 0 ) takes values in [1 , n 0 ( n 0 +1) 2 ] . 9 Theorem 3.2. L et v b e a no de at depth k . The me an of the lo ad distribution is E ( I v ( n )) = ( n n ≤ k k n > k Pr o of. W e prov e by induction on k . F or k = 1, the statement is trivial, so consider k > 1. Since the distribution of I v ( n ) is constant for n ≥ k , we assume n ≤ k . Let N v ,k − 1 b e the num b er of nodes at lev el k − 1 which send sediment to v at the end of time n − 1. This v ariable takes v alues in { 0 , 1 , 2 } with probabilities { 1 / 4 , 1 / 2 , 1 / 4 } , resp ectiv ely . Call w 1 ( w 2 ) the left (righ t) parent of v . E ( I v ( n )) = 2 X i =1 E ( I v ( n ) | N v ,k − 1 = i ) P ( N v ,k − 1 = i ) = 1 / 4 + 1 / 2 [1 + E ( I w 1 ( n − 1))] + 1 / 4 [1 + E ( I w 1 ( n − 1) + I w 2 ( n − 1))] = 1 + E ( I w 2 ( n − 1)) = 1 + ( n − 1) = n where to go from the second line to the third line, w e use the fact that the v ariables I w 1 ( n − 1) and I w 2 ( n − 1) hav e the same distribution (see b. under Theorem 3.1). Theorem 3.1 lets us use geometric prop erties of clusters of a static netw ork ( ω n 0 ) to study something which is dynamic: the load at time n at no de v . That load may hav e come from a pathw ay that no longer even exists at time n . W e further exploit this relationship, but to do this we must consider the concept of the dual web, defined in, for example, [7], and of whose definition w e remind the reader. Consider the o dd sublattice Z 2 odd = { ( x, y ) ∈ Z 2 : x + y o dd and y ≤ 1 } F or an y no de v ∗ = ( x ∗ , y ∗ ) ∈ Z 2 odd w e call the no de ( x ∗ + 1 , y ∗ + 1) the righ t c hild of v ∗ and we call the no de ( x ∗ − 1 , y ∗ + 1) the left child of v ∗ . Similarly , w e call the no de ( x ∗ + 1 , y ∗ − 1) the righ t parent of v ∗ and we c all the node ( x ∗ − 1 , y ∗ − 1) the left paren t of v ∗ . W e define the set E 2 odd in the obvious w ay . The set of configurations { ω n : n ≥ 0 } induces a set of configurations { ω ∗ n : n ≥ 0 } ⊂ { 0 , 1 } E 2 odd b y the following rule. If, in the configuration ω n , a no de v = ( x, y ) is connected to its left child, we form a connection b etw een the no de v ∗ = ( x, y − 1) and its right child in the configuration ω ∗ n b y setting the image under ω ∗ n of the edge in E 2 odd b et w een v ∗ and its righ t child to 1, and the image of the edge b etw een v ∗ and its left c hild to 0. If, on the other hand, v is connected to its right child in ω n then we set the image of the edge from v ∗ to its left c hild under ω ∗ n to 1 and the image of the edge from v ∗ to its righ t c hild to 0. See Figure 4 and notice that we construct clusters in ω ∗ n so that no o ccupied edges in ω ∗ n cross an y o ccupied edges in ω n . 10 Figure 4: Portion of an realization of the erosion netw ork, along with its dual w eb. The solid lines indicate paths of sediment flow and the dotted lines show paths of the dual (courtesy of [6]). The up ward paths in ω ∗ n no w resemble the down ward paths in ω n . That is, the up ward path starting at the no de v ∗ is a simple symmetric random w alk whic h is killed at depth 1. Random walks starting at different no des are indep enden t until they meet, at whic h p oint they coalesce into one random w alk. (This is similar to the coalescing random walks picture of the discrete w eb, describ ed in [2], [11], [7].) There is an obvious physical interpretation for the paths in the dual web. F or an y tw o adjacent paths in the configuration ω n , there is a path in ω ∗ n separating them. If the paths in ω n represen t rills or drains, the paths in ω ∗ n represen t the divides or ridges b et ween them. Just as divides b et ween rills do not cross rills, paths in ω ∗ n do not cross paths in ω n . W e now characterize the load distributions for our model. F or any no de v = ( x, y ) (with depth k ), let v ∗ L = ( x − 1 , y ) and let v ∗ R = ( x + 1 , y ). Consider the set of edges in the dual lattice con tained in the paths emanating from the v ertices v ∗ R and v ∗ L in ω ∗ n un til either (a) they meet at some vertex w ∗ or (b) they reach a depth of 1. The set of no des in Z 2 ev en in the in terior of this set of edges is exactly the bac kward cluster of v in the configuration ω n . W e now mak e some definitions so that we can work with this load distri- bution. Let { X L i : i ≥ 2 } and { X R i : i ≥ 2 } b e indep enden t sets of ran- dom v ariables (also indep endent of eac h other) which take the v alues 1 and -1 eac h with probability 1 2 . F or i ≥ 2 let Y i = 1 2 ( X R i − X L i ) and for i ≥ 1, let W i = 1 + Y 2 + ... + Y i . Consider the stopping time τ = min { n : W n = 0 } Up un til the stopping time τ , the random v ariable W i represen ts the width 11 of the bac kward cluster of the node v in the real lattice (only v alleys and not separating ridges), where we only consider nodes in this cluster whose depths are b etw een k − i + 1 and k . Therefore the total num b er of no des in this partial cluster should b e L i := W 1 + ... + W i No w we can mak e an equiv alent definition of the distribution of the load I v ( n ) by saying that for each fixed n , it is the same as the distribution of the random v ariable L k ( n ) := L min( τ ,n,k ) (4) This v ariable is essen tially a discrete integral of the symmetric random walk { W i : i ≥ 1 } . Theorem 3.3. L et v b e a no de at depth k v and let w b e a no de at depth k w ≥ k v . F or any n ≥ k v and for any l < k v we have P ( I v ( n ) = l ) = P ( I w ( n ) = l ) Ther efor e the limit lim k v →∞ I v ( k v ) (5) exists in distribution. This limit is a.s. finite but has infinite me an. Pr o of. On the ev ent τ ≥ k v , i.e. the load originated from the top, I v ( n ) = L min( τ ,n,k v ) = L k v ≥ k v > l and I w ( n ) = L min( τ ,n,k w ) ≥ L min( τ ,n,k v ) > l Hence, w e need only consider τ < k v . P ( I v ( n ) = l ) = P ( I v ( n ) = l, τ < k v ) = P ( L min( τ ,n ) = l, τ < k v ) = P ( L min( τ ,n,k w ) = l, τ < k v ) = P ( I w ( n ) = l ) (6) The random v ariable L k v ( n ) is constan t for n ≥ k v , so P ( I v ( k v ) = l ) = P ( I v ( k w ) = l ) = P ( I w ( k w ) = l ) where in the last equalit y we use (6). Consequently , for an y fixed l , the limit lim k v →∞ P ( I v ( k v ) = l ) 12 exists. By the definition (4), a random v ariable with this limiting distribution is L ∞ := lim k →∞ L min( τ ,k,k ) = L τ Since τ ≤ L τ ≤ τ ( τ + 1) 2 the third statemen t of the theorem will follo w if we sho w that τ is a.s. finite and has infinite mean. But since the increments { Y i } of the random walk { W i } hav e mean zero, the w alk is recurrent. In addition, it is a standard result that the en trance time of the set { 0 } has infinite mean. This completes the pro of. 3.2 Dynamics No w we in vestigate some asp ects of the effect of η on the stability of config- urations ov er time. As noted, the dynamics creates an interpolation b et ween the discrete web and the discrete dynamical web. The evolution of the system mirrors some asp ects of rill erosion, one b eing that no des through which a large amoun t of water passes at time n 0 ha ve a non-trivial probabilit y to c hannel a large amount of w ater at any time n 1 > n 0 . The degree to whic h this is true dep ends on the parameter η , as we will see. 3.2.1 Load Correlation W e start our analysis b y insp ecting simulation results. F or any tw o p ositiv e in tegers M , N , let V M ,N b e an en umeration of the M N nodes in the b o x [0 , M − 1] × [ − N , − 1] and define the load correlation co efficien t at time n (for n ≥ N ) b y K M ,N ( n ) = P v ∈ V M,N I 0 v ( N ) I 0 v ( n ) q ( P v ∈ V M,N I 0 v ( N ) 2 )( P v ∈ V M,N I 0 v ( n ) 2 ) where I 0 v ( n ) = I v ( n ) − E ( I v ( n )). This quan tity is only defined for n ≥ N b ecause t wo load v ectors for a box of depth N are in some sense incomparable if they are tak en at times n 0 , n 1 with n 0 < N ≤ n 1 . F or example, a no de at depth n only has a maximum p ossible load of n ( n +1) 2 at time n < N , whereas its maxim um p ossible load is N ( N +1) 2 for n ≥ N . In Figure 5 we ha ve graphed simulation results for a netw ork 49 no des deep and 49 no des wide. The x-axis represents v alues of the parameter η , as it v aries from 0 to 5. The y-axis represents v alues of a time av eraged correlation co efficien t, namely the quan tity 1 N 0 − N N 0 X n = N K M ,N ( n ) 13 Figure 5: Averaged load correlation (v ertical axis) versus η (horizontal axis). Data p oin ts are taken at η -in terv als of 0.1 un til η = 2, and then at interv als of 1. for M = N = 49 and N 0 = 200. F urthermore, we a v eraged this v alue o ver 6 indep enden t trials. This quantit y is mean t to approximate v alues of K M ,N ( n ) for M , N , n large. The time av eraging seemed necessary b ecause of fluctuations, most lik ely due to finite size conditions, in the quantit y K M ,N ( n ). W e can see in Figure 5 that the co efficien t approac hes 1 as η → 0. This mak es sense b ecause, as remarked in section 1.3, the η → 0 limit of the dynam- ics (in any fixed b o x) is the same as the dynamics (or rather non-dynamics) of the discrete web. Therefore the load vector for this box should be similar (if not the same) at any tw o times. As η increases, the correlation co efficien t decreases and app ears to approac h 0. Indeed, additional simulations giv e the following data: for η = 10 , 100 , 1000 , 10000, the coefficients were . 2391 , . 0896 , . 0189 , and . 0101 , . F rom the discussion of the η → ∞ limit giv en in section 1.3, the cor- relation co efficien t should approach that computed from tw o load vectors from indep enden t realizations of the discrete w eb. 3.2.2 de Finetti Measures Whereas we can compare no des at the top lev el to standard P´ oly a urns, w e can compare low er lev el no des to time-dep endent input [14] or r andom input [13] P´ oly a urns. W e start with an urn with R 0 red balls and B 0 blac k balls, as b efore, but we also hav e a time-dep enden t (or random) input sequence I = ( I 0 , I 1 , ... ). A t time t = n w e dra w a ball from the urn and w e return it to the urn along with I n balls of the same color. Notice that this pro cess with I = (1 , 1 , ... ) is just the standard P´ oly a urn. T o analyze these low er level no des, we will also mak e use of a fundamental result in the theory of exc hangeable v ariables. 14 Definition 3.4. { 0 , 1 } -v alued v ariables X 1 , X 2 , ... are exchangeable if for any x 1 , ..., x m ∈ { 0 , 1 } and for any p erm utation σ of m elements we ha ve P ( X 1 = x 1 , ..., X m = x m ) = P ( X σ (1) = x 1 , ..., X σ ( m ) = x m ) Theorem 3.5 (de Finetti) . L et (Ω , F , P ) b e a pr ob ability sp ac e and supp ose that { X n } n ≥ 0 ar e exchange able { 0 , 1 } -value d r andom variables define d on Ω . Then ther e exists a r andom variable F on Ω so that c onditione d on F , the r andom variables X n ar e indep endent Bernoul li with p ar ameter F . It is easy to v erify that if depth( v ) = 1, then the { 0 , 1 } -v alued v ariables { D L v ( n ) } n ≥ 0 are exchangeable. In our case, the v ariable F from Theorem 3.5 is actually p v := lim n →∞ P L v ( n ) . (7) Th us if we know the asymptotic fraction of left choices for a no de, then our no de is just flipping indep enden t coins each second with the same bias. A t lo wer levels, the v ariables { D L v ( n ) } are not exchangeable. How ever, they are asymptotically exc hangeable. W e use the definition of Kingman [12]. Definition 3.6. { 0 , 1 } -v alued random v ariables X 1 , X 2 , ... are called asymp- totically exc hangeable if there exists a sequence Y 1 , Y 2 , ... of exc hangeable random v ariables so that for each x 1 , ..., x m ∈ { 0 , 1 } , lim N →∞ P ( X 1+ N = x 1 , ..., X m + N = x m ) = P ( Y 1 = x 1 , ..., Y m = x m ) In the language of Theorem 3.5, let F b e the random v ariable asso ciated with the exchangeable v ariables { X n } . W e call F the de Finetti measure for the sequence { Y n } . Let v b e a no de with depth k ≥ 1. Theorem 3.7. The variables { P L v ( n ) } n ≥ 0 form a b ounde d martingale se quenc e w.r.t. F n . Ther efor e they have an almost sur e limit p v . Pr o of. Similar to the pro of of [13, Theorem 2.1] R emark 3.8 . Using the same equations which produce (3), the limit in Theorem 3.7 is the same as the limit of the v ariables { LF L v ( n ) } n ≥ 0 . F or any num b er 0 ≤ p ≤ 1, define the measure Q p on the set { 0 , 1 } b y Q p ( { 0 } ) = 1 − p , Q p ( { 1 } ) = p Let { v 1 , ..., v r } b e a finite set of vertices. F or any vector of real num b ers ( p 1 , ..., p r ), each b etw een 0 and 1, define the pro duct measure Q ~ p on v ectors in { 0 , 1 } r to b e the pro duct measure Q r i =1 Q p i . 15 Theorem 3.9. F or fixe d v , the variables { D L v ( n ) : n ≥ 1 } ar e asymptot- ic al ly exchange able with de Finetti me asur e e qual to the distribution of p v . F urthermor e, let ~ v = ( v 1 , ..., v r ) b e a ve ctor of vertic es and for e ach n , let D L ~ v ( n ) = ( D L v 1 ( n ) , ..., D v r ( n )) . If ~ d 1 , ..., ~ d s ar e ve ctors in { 0 , 1 } r , with pr ob a- bility one, lim T →∞ P ( D L ~ v (1 + T ) = ~ d 1 , ..., D L ~ v ( s + T ) = ~ d s ) = E ( s Y i =1 Q ~ p ( ~ d i )) wher e ~ p = ( p v 1 , ..., p v r ) . Pr o of. Similar to the pro of of [13, Theorem 2.2]. Because of lateral translation in v ariance, the de Finetti measure for v de- p ends only on the depth k . In light of this, we define θ k = de Finetti measure for row k (8) With this framework we will be able to study the switching rate of eac h no de v once w e hav e the following lemma. Lemma 3.10. F or any no de v , almost sur ely, lim n →∞ 1 n n X i =1 D L v ( i ) = p v (9) Pr o of. Similar to the pro of of [13, Theorem 2.3]. W e no w define the switching function s v for n ≥ 2 b y s v ( n ) = D L v ( n )(1 − D L v ( n − 1)) + D L v ( n − 1)(1 − D L v ( n )). Define the switching rate S v ( n ) to b e the time a verage of s v , that is S v ( n ) = 1 n − 1 n X i =2 s v ( i ) Theorem 3.11. The n → ∞ limit of S v ( n ) exists a.s. lim n →∞ S v ( n ) = 2 p v (1 − p v ) (10) Pr o of. The pro of is similar to the pro of of [13, Theorem 2.3]. Let d n = D L v ( n ) and p n = P L v ( n ). A straightforw ard calculation gives E ( d n +1 |F n ) = p n , E ( p n +1 |F n ) = p n (11) No w, lim n →∞ 1 n n X i =2 s v ( i ) = lim n →∞ 1 n n X i =2 ( d i (1 − d i − 1 ) + d i − 1 (1 − d i )) 16 = lim n →∞ 1 n n X i =2 ( d i + d i − 1 ) − 2 lim n →∞ 1 n n X i =2 d i d i − 1 = 2 p v − 2 lim n →∞ 1 n n X i =2 d i d i − 1 b y Lemma 3.10. W e must show that the last limit ab o v e is a.s. equal to p 2 v . Define M n = P n i =2 ( d i − p n ) d i − 1 . M n is a martingale with resp ect to F n : E ( M n +1 |F n ) = E ( n X i =2 ( d i d i − 1 ) − n X i =2 ( p n +1 d i − 1 ) + d n +1 d n − p n +1 d n |F n ) = n X i =2 ( d i d i − 1 ) − n X i =2 ( p n d i − 1 ) + p n d n − p n d n = M n where w e use b oth equations in (11). Note also lim n →∞ 1 n n X i =2 p n d i − 1 = p v lim n →∞ 1 n n X i =2 d i − 1 whic h is p 2 v , by Lemma 3.10. Therefore it suffices to show that with probabilit y one, M n n → 0 . (12) By summing the series 1 n ( M 1 + ... + M n ) = n X i =2 M i i  i n  b y parts, it can b e sho wn that (12) will follow once w e show that ∞ X i =2 M i +1 − M i i con verges. T o this end, define M 0 n = P n − 1 i =2 M i +1 − M i i . W e lea ve the reader to v erify that M 0 n is a martingale. Using L 2 -orthogonalit y of martingale differences, E ( M 0 n ) 2 = E ( n − 1 X i =2  M i +1 − M i i  2 ) = E ( n − 1 X i =2  ( d i − p n )( d i − 1 ) i  2 ) ≤ n X i =2 1 i 2 < ∞ Therefore M 0 n is an L 2 b ounded martingale and conv erges a.s. This completes the pro of. 17 If p v ∈ (0 , 1) then neither of the choices L or R are transient for v . This prompts the question of whether or not the de Finetti measures θ k ha ve atoms at 0 or 1. F or any fixed k , the answ er is no. Theorem 3.12. F or e ach k ≥ 1 , the me asur e θ k has no atoms. Pr o of. In [14, Theorem 4] it is sho wn that a time-dep endent input P´ olya urn’s de Finetti measure cannot ha v e atoms if there is a C so that I v ( n ) ≤ C for all n . F or each realization of the dynamics and for each v , we hav e I v ( n ) ≤ k v ( k v +1) 2 for all n . The result follows. Corollary 3.13. Each no de v has a nonzer o asymptotic switching r ate. Ther e- for e, for e ach v , the states L and R ar e r e curr ent. Pr o of. This is a direct consequence of Lemma 3.10 and Theorem 3.12. Corollary 3.14. With pr ob ability one, for e ach no de v , the variable I v ( n ) takes e ach value in [1 , k v ( k v +1) 2 ] for infinitely many values of n . Pr o of. Let v 1 , ..., v m b e the m = k v ( k v +1) 2 − 1 no des ab o ve v whic h can send sedimen t to v and let ~ d 1 , ..., ~ d v k ∈ { 0 , 1 } m . Let N ≥ 1 and write D L ~ v ( n ) for the v ector ( D L v 1 ( n ) , ..., D L v m ( n )). lim N →∞ lim T →∞ P ( D L ~ v ( T + j k v +1) = ~ d 1 , ..., D L ~ v ( T +( j +1) k v ) = ~ d k v for some 1 ≤ j ≤ N ) equals zero almost surely , by Theorem 3.9 and Theorem 3.12. But this probabil- it y dominates the probability of the even t { D L ~ v ( n + 1) = ~ d 1 , ..., D L ~ v ( n + v k + 1) = ~ d v k for infinitely man y n } c . The result follows. Here an in teresting picture of our net work emerges. On the one hand we may view the system as an infinite lattice (the lo wer half plane), where eac h no de is a random input P´ oly a urn. The output of the urns at depth k at time n b ecomes the input of the urns at depth k + 1 at time n + 1. On the other hand, as remark ed in section 1.3, we may first sample (non-indep enden tly) v alues { p v : v ∈ Z 2 ev en } from the de Finetti measures { θ k : k ≥ 1 } to create an infinite arra y . As time n approac hes infinity , the b eha vior of the system approaches the b ehavior of the same net work in whic h each no de v chooses to send its current load left with probabilit y p v and right with probability 1 − p v , independently at each second. Therefore this picture is of a netw ork of tw o v ariables, a realization of v alues p v from the de Finetti measures, and realization of dynamics whic h coincides with the dynamics of a muc h simpler net work. This second net w ork is an obvious generalization of the Dynamical Discrete W eb. Figure 6 shows histograms for the de Finetti measures θ k for k = 2 , 5 , 9 and for v alues of η = . 5 , 1 , 2. One sees that the measures b ecome more biased as k increases (for fixed η ). In other w ords, the mass of θ k is concentrated on domains closer to 0 and 1 than is the mass of θ k − 1 . This would seem to imply that the 18 η = 1 2 , ro w 2 η = 1 2 , ro w 5 η = 1 2 , ro w 9 η = 1, row 2 η = 1, row 5 η = 1, row 9 η = 2, row 2 η = 2, row 5 η = 2, row 9 Figure 6: Definetti measures for k = 2 , 5 , 9 (from left to right) and η = . 5 , 1 , 2 (from top to b ottom). exp ected asymptotic switching rate of a no de at level k (whic h is 2 p v (1 − p v )) m ust decrease with k . Similarly , if k is fixed and η decreases to 0, it seems that the exp ected switching rate should decrease. Figure 7 represents data given by simulations conducted with an erosion net work with width 10 5 , depth 50, and η = . 1 , 1 , or 10. The sim ulation ran for n = 1000 steps and at the end, switch rates for eac h no de in the netw ork were computed. In each row, each no de’s rate was av eraged. Since tw o no des v 1 and v 2 with the same depth k hav e indep enden t b eha vior as long as they are at least a distance of 2 k apart, the ergo dic theorem gives that, as the netw ork size ap- proac hes infinity , the resulting a verage should resem ble the expected switch rate for a row. The ab o ve results results were av eraged by row ov er 3 indep enden t trials. Finally , the data w ere plotted by row. Not only do the av erage switch rates appear to decrease as k increases, there appears to be a non-trivial (i.e. non-zero and η dep endent) limit for the switch rate. This indicates that for at least some v alues of η , the limit of the measures θ k (if it exists) is most likely not equal to 1 2 ( δ 0 + δ 1 ). 19 Figure 7: Average switch rate (vertical axis) versus ro w (horizon tal axis). The v alues of η are 10 (top curve), 1 (middle curv e), .1 (b ottom curv e). 3.2.3 Catastrophes Next, we study the follo wing situation. Supp ose a no de v at a large depth k (for this section w e assume the depth is at least 2) starts with a small input load and k eeps a relatively small input load until a m uch later time. Then v ’s load c hanges dramatically . If this new load is sufficiently large, it could bring v ’s de Finetti measure muc h closer to 1 2 ( δ 0 + δ 1 ). This analysis is from the point of view of the no de v , whereas the analysis of the last half of the section will b e from the p oin t of view of the parent. Let A v ( n ) = T v ( n − 1) n , n ≥ 1 Definition 3.15. F or an y n ≥ 1, define the flo o d ratio F v ( n ) = I v ( n ) A v ( n ) . F or c ≥ 1, we say that a flo o d of order c o ccurs at time n if F v ( n ) ≥ c . R emark 3.16 . Since I v ( n ) , A v ( n ) ∈ [1 , 1 2 ( k ( k + 1))], we hav e 2 k ( k + 1) ≤ F v ( n ) ≤ k ( k + 1) 2 (13) Prop osition 3.17. F or any v , A v ( n ) has a limit a.s. Ther efor e, lim inf n →∞ F v ( n ) < lim sup n →∞ F v ( n ) Pr o of. W e show the first statement b y induction on depth( k ). Clearly this is true if depth( v ) = 1. Otherwise, let w 1 b e the left parent of v and assume that 20 for all no des w 0 with depth equal to that of w 1 , A w 0 ( n ) has a limit. Let N R w 1 ( n ) b e the num b er of i ≤ n such that D R w 1 ( i ) = 0. By Lemmas 3.7 and 3.10, lim n →∞ T R w 1 ( n ) n + 1 = lim n →∞ T R w 1 ( n ) T w 1 ( n ) ( n + 1) T w 1 ( n ) = (1 − p w 1 ) lim n →∞ A w 1 ( n + 1) exists. The same argument shows that, if w 2 is the righ t parent of v , lim n →∞ T L w 2 ( n ) n +1 exists. Therefore, lim n →∞ T v ( n − 1) n = lim n →∞ n − 1 + T R w 1 ( n − 2) + T L w 2 ( n − 2) n exists. F or the second statemen t of the prop osition, we use Lemma 3.14 to see that lim inf n →∞ F v ( n ) = 1 lim n →∞ A v ( n ) < k ( k + 1) 2 lim n →∞ A v ( n ) = lim sup n →∞ F v ( n ) R emark 3.18 . The ab o ve pro of shows also that lim n →∞ T R v ( n ) n +1 exists and equals (1 − p v ) lim n →∞ A v ( n ). R emark 3.19 . A simple extension of the ab o ve pro of using Theorem 3.12 sho ws that for an y lev el k , the distribution of the time-a verage of the load for lev el k has no atoms. F rom the previous prop osition we know that if v has an asymptotic a verage input load A v , then v will ha ve infinitely man y floo ds of order k ( k +1) 2 A v +  for an y  > 0. This num b er can be quite large if A v is small. W e study numerically the rates of these large flo ods as they relate to b oth depth( v ) and η . F or a fixed no de v , define the random measure µ v ,n = n X i =2 δ F v ( i ) Let { v i } i ≥ 1 b e an enumeration of the no des with the same depth as that of v . By the ergo dic theorem, the av erage 1 M M X i =1 µ v i ,n → E ( µ v ,n ) as M → ∞ . See Figure 8. The graphs come from a sim ulation run for 10 5 seconds on a netw ork with a width of 10 4 no des and a depth of 50 no des. W e graphed the densit y function for the measure 1 10 8 P 10000 i =1  µ v i , 10 5 − µ v i , 9000  in an attempt to approximate 1 n − 9000 E ( µ v ,n ) for n large and for v with depth 5, 20, and 50. The reason we subtracted µ v i , 9000 is to decrease the effect of small times, during whic h the ratio F v ( n ) is lik ely to b e an in teger. 21 η = 1 10 , ro w 5 η = 1 10 , ro w 20 η = 1 10 , ro w 50 η = 1, row 5 η = 1, row 20 η = 1, row 50 η = 10, row 5 η = 10, row 20 η = 10, row 50 Figure 8: Histograms approximating the measure E ( µ v ,n ) for n large. As η → 0 with a fixed row or as the ro w increases with fixed η , each measure seems to concentrate its mass at 0 and 1. In other words, the measure of any in terv al whic h do es not include either of these tw o points app ears to approac h 0. In addition, T able 1 sho ws that as η → 0 with a fixed row or as the ro w increases with fixed η , the exp ected fraction of time during which a large flo o d o ccurs (ratio ab ov e 5) increases. Since the time av erage of flo od ratios approaches 1 as n → ∞ , the ab o ve facts indicate a trend that as η → 0 or as the row increases, the time v ariance of measures increases, giving more p ossible v ariabilit y of the flo od ratios. Although this v ariance seems to increase, v alues near 1 (on the x-axis) sho w that the fraction of time flo od ratios sp end near 1 increases. In spirit, this is in accordance with previous results, as we explain. Figure 7 shows that as the ro w increases, the expected switc hing rate of a no de decreases. It is reasonable to believe that the same conclusion holds if the depth is fixed but η decreases. Therefore the netw ork prefers to b e more static in these circum- stances and w e w ould expect a node to receiv e a relativ ely constan t load, forcing flo od ratios to b e near 1. W e now c hange fo cus to the paren t. Mak e the following definition. Definition 3.20. F or any n ≥ 2, define the righ t catastrophe ratio 22 η Ro w 5 Ro w 20 Ro w 50 .1 .00035 .00205 .00308 1 .00136 .01235 .01789 10 .00024 .02356 .03711 T able 1: Exp ected fraction of time during which a large flo od (ratio at least 5) o ccurs. C R v ( n ) = I v ( n ) A R v ( n ) whenev er D L v ( n ) = 0 Here, A R v ( n ) = T R v ( n − 1) N R v ( n − 1) , where N R v ( n ) is the num b er of i ≤ n such that D R v ( i ) = 0. F or c ≥ 1, we say that a righ t catastrophe of order c occurs at time n if C R v ( n ) ≥ c . Make similar definitions for left catastrophe ratio and left catastrophe. R emark 3.21 . F rom a similar argument to that used in Prop osition 3.17, we ha ve lim inf n →∞ C R v ( n ) < lim sup n →∞ C R v ( n ) W e inv estigate the relationship betw een floo ds and catastrophes. Figure 9 sho ws the exp ected fraction of right catastrophes from the left paren t whic h result in a flo od of at least the same magnitude. The sim ulation was p erformed with a netw ork of width 1000, depth 50, η equal to either .1, 1, or 10, and for a duration of 10 5 seconds. The calculation of fractions w as only made b etw een 9000s and 10000s and we only consider catastrophes with ratio at least 1. It is clear from the figure that as the row increases, this exp ectation decreases. As η decreases for a fixed row, the exp ectation also decreases. As n → ∞ , A R v ( n ) = ( n − 1) T R v ( n − 1) ( n − 1) N R v ( n − 1) → lim n →∞ A v ( n ) b y Remark 3.18 and Lemma 3.10. Therefore, we may use Remark 3.19 to show that almost surely for large n a right catastrophe of order c o ccurs for the no de v at time n if v has a flo od of order c at time n − 1. In other words, whenever a node receives a flo od of order c at a large time, it has either a left or right catastrophe of the same order at the next second. No w w e ma y in terpret the probability that a no de has a flo od given that its left paren t has a righ t catastrophe as the probability that a paren t’s catastrophe incites a catastrophe in the child. This w ould b e a step of a p ossible catastro- phe c asc ade . The simulation results indicate that cascades b ecome less presen t at low er levels (on av erage) but that they should never cease to exist. Two questions naturally arise. Given that a node has a righ t catastrophe of order c , ho w far does its catastrophe cascade trav el? A t eac h step of the cascade, the relev ant (righ t or left) catastrophe ratio will generally increase. Indee d, a child 23 Figure 9: Exp ected fraction of right catastrophes from left paren t w hic h result in a flo od of at least the same magnitude. The v alues of eta are 10 (top), 1 (middle), and .1 (b ottom). The data are plotted b y row. no de may ev en receiv e a catastrophe from b oth paren ts. Ho w large do es this ratio b ecome in a typical cascade? 4 Conclusion W e hav e shown that the erosion model exhibits many prop erties of rill erosion. Eac h no de c ho oses a random initial direction (right or left) in which to send sedimen t and further such choices b ecome biased at a rate largely determined b y the parameter eta. This is similar to the metho d by which rills are cut into a hillslop e. As more water and sedimen t flo ws through a rill, a c hannel is cut deep er, giving reinforcemen t to the path, making it more lik ely to carry sedimen t in the future. Though the dynamics manifests itself through reinforcement, no fixed no de can b ecome fully biased (i.e. ha ve a de Finetti measure equal to a sum of tw o delta masses). That is, since each no de has a non-trivial asymptotic switc hing rate, sediment flow emerging from it will take b oth a left and right path a p ositiv e fraction of time. This rate of switching app ears to decrease as w e mov e further down the hill. There are a num ber of questions which deserv e careful analysis. Do the measures θ k ha ve a limit? If so, one w ould expect the limit to be the de Finetti measure asso ciated with the ”infinity pro cess”. T o define this process, we start with a lattice of no des whic h extends infinitely far in b oth p ositiv e and negativ e y directions. Since the b eha vior of a no de v at time n in the present mo del dep ends only on the no des in the n − 1 levels ab o ve it we may consider the input to the no de v at time n in the infinity mo del to b e a function of the 24 output of this finite num ber of ancestors. In the same wa y w e ha ve analyzed in this pap er, it is p ossible to show that a de Finetti measure θ ∞ for this pro cess exists and that θ ∞ = lim n →∞ θ n ( n ) , (14) where the term inside the limit is the measure given by θ n ( n )([ a, b ]) = P ( P L v ( n ) ∈ [ a, b ]) for depth( v ) = n Do es the measure θ ∞ ha ve atoms for some v alues of η ? If so, is there a critical η ∗ so that for 0 < η < η ∗ , θ ∞ has atoms? If the limit of θ k (assuming it exists) is not the same as θ ∞ , do es this limit hav e atoms and is there a critical η asso ciated with it? Ac knowledgmen ts . W e thank Ev erett Springer for careful reading and commen ts. W e also thank Charles Newman, Daniel Stein, Emmanuel Sc hertzer, and Itai Benjamini for very helpful discussions. Researc h at NCAR was funded b y the National Science F oundation. References [1] Athrey a, K.: On a c haracteristic prop ert y of P´ oly a’s urn. Stud. Sci. Math. Hung. 4 , 31–35 (1969) [2] Arratia, R.: Coalescing Brownian motions and the voter mo del on Z , Un- published partial manuscript (1981), av ailable from rarratia@math.usc.edu. [3] Diaconis, P .: Recent progress on de Finetti’s notion of exc hangeability . In J. Bernardo, M. de Gro ot, D. Lindley , and A. Smith, editors, Bayesian Statistics , pp. 111–125. Oxford Univ ersity Press, Oxford (1988) [4] Eggenberger, F., P´ oly a, G.: ¨ Ub er die Statistik V erketter V org¨ ange. Z. A ngew. Math. Me ch. 3 , 279–289 (1923) [5] F eller, W. Intr o duction to Pr ob ability The ory and Its Applic ations, V ol. 2, Wiley , New Y ork. [6] F on tes, L ., Isopi, M., Newman, C.M., Ra vishank ar, K.: Coarsening, nu- cleation and the marked Bro wnian w eb. Ann. Inst. H. Poinc ar e. 42 37–60 (2006) [7] F on tes, L., Newman, C., Ra vishank ar, K., Schertzer, E.: The Dynamical Discrete W eb. Arxiv: math.PR/07042706 (2007) [8] F reedman, D.: Bernard F riedman’s urn. Ann. Math. Statist. , 36 956–970 (1965) 25 [9] F riedman, B.: A simple urn mo del. Commun. Pur e Appl. Math. 2 , 59–70 (1949) [10] Glo c k, W.S.: The developmen t of drainage system – a synoptic view. Ge- o gr aph. R ev. , 21 , 475 (1931) [11] Howitt, C., W arren, J.: Dynamics for the Brownian web and the erosion flo w. Arxiv: math.PR/0702542 . [12] Kingman, J.: Uses of exc hangeability . Annals of Pr ob ability 2 183–197 (1978) [13] May , C., Paganoni, A., Secchi, P .: On a tw o-color generalized P´ olya urn. International Journal of Statistics 63 , 115–134 (2005) [14] Peman tle, R.: A time-dependent version of P´ oly a’s urn. Journal of The o- r etic al Pr ob ability 3 , 627–637 (1990) [15] Peman tle, R.: A survey of random pro cesses with reinforcement. Pr ob ability Surveys 4 1–79 (2007) [16] Quinton, John N.: Erosion and sedimen t transp ort, in Envir onmental mo d- el ling : finding simplicity in c omplexity , ed. By John W ainwrigh t and Mark Mulligan, Wiley , Hob oken, NJ, 2004. [17] Ro dr ´ ıguez-Iturb e, I., Rinaldo, A.: F r actal R iver Basins: Chanc e and Self- Or ganization , Cam bridge Universit y Press, Cam bridge, UK, 1997. [18] T o y , T.J., F oster, G.R., Renard, K.G.: Soil Er osion: Pr o c esses, Pr e diction, Me asur ement, and Contr ol , John Wiley and Sons, N.Y., 2002. [19] Willgo ose, G. R., Bras, R. L., and Ro driguez-Iturbe, I.: 1991a. A Ph ysically Based Coupled Netw ork Growth and Hillslop e Evolution Mo del: 1 Theory . WRR, 27(7), 1671–1684. 26