Invasion fronts on graphs: the Fisher-KPP equation on homogeneous trees and Erd\H{o}s-R\'eyni graphs

We study the dynamics of the Fisher-KPP equation on the infinite homogeneous tree and Erd\H{o}s-R\'eyni random graphs. We assume initial data that is zero everywhere except at a single node. For the case of the homogeneous tree, the solution will either form a traveling front or converge pointwise to zero. This dichotomy is determined by the linear spreading speed and we compute critical values of the diffusion parameter for which the spreading speed is zero and maximal and prove that the system is linearly determined. We also study the growth of the total population in the network and identify the exponential growth rate as a function of the diffusion coefficient, $\alpha$. Finally, we make predictions for the Fisher-KPP equation on Erd\H{o}s-R\'enyi random graphs based upon the results on the homogeneous tree. When $\alpha$ is small we observe via numerical simulations that mean arrival times are linearly related to distance from the initial node and the speed of invasion is well approximated by the linear spreading speed on the tree. Furthermore, we observe that exponential growth rates of the total population on the random network can be bounded by growth rates on the homogeneous tree and provide an explanation for the sub-linear exponential growth rates that occur for small diffusion.


Introduction
We consider the Fisher-KPP equation defined on a graph and study the transition of the system from unstable to stable homogeneous states. To provide some motivation for this study we consider the dynamics of an invasive species on a transportation network. One can think of each node in the network as describing a physical location and the edges connecting these nodes as available transportation routes between these locales. When the invasive species is introduced at one location in the network we expect the species to grow and spread through the network until it resides at every node in the network that is favorable to the species. The goal would be to determine in what manner this transition occurs and to estimate the time that is required for the species to arrive at any given node in the network.
A realistic description of this process would involve complicated network topologies and a great deal of heterogeneity -both in the local dynamics of the species at each node as well as the transportation mechanisms available to move between modes. Our focus in this paper is to study a simplified version of this problem where the local dynamics are identical and prescribed by logistic growth while the movement between nodes is encapsulated by diffusion. In regards to network topologies, we concentrate mainly on the case where the network is a homogeneous tree. This restriction allows for us to make precise statements regarding the dynamics of the system and make predictions for more complicated and realistic situations. These prediction are tested on Erdős-Réyni random graphs.
To be precise, let G = (V, E) be a undirected and unweighted countable graph; while we will mainly be interested in the homogeneous tree, but leave the discussion general for the time being leave. Consider a differential equation defined on the set of vertices with local dynamics prescribed by the scalar differential equation u t = f (u). The reaction term f (u) ∈ C 1 ([0, 1]) is assumed to be of KPP type, see [13,17], and satisfies the conditions, for 0 < u < 1 f (u) < f (0)u for 0 < u < 1.
We can assume further that f (0) = 1 after rescaling the independent variable. In most of the numerical simulations contained in this paper we use the explicit example of logistic growth, f (u) = u (1 − u). Diffusion between nodes is incorporated by the graph Laplacian, ∆ G = A(G) − D(G) where A(G) is the adjacency matrix and D(G) is diagonal with ith entry equal to the degree of the ith node. Incorporating both the local dynamics and diffusion between nodes we arrive at the reaction-diffusion equation that will be the main focus of this article, The dependent variable u i ∈ R can be thought of as representing the population of a species residing at node i ∈ V .
Consider initial conditions where the population of the species is zero everywhere in the graph aside from one node where some positive concentration exists. Due to diffusion, the population will spread out into the graph. We are interested in both the pointwise and aggregate behavior of the population. By pointwise behavior, we refer to the question of whether the solution at any given node will converge to the stable steady state at one or remain near the unstable steady state at zero. We will be particularly interested in the existence/ non-existence of traveling fronts that propagate asymptotically with fixed speed and replace the unstable state at zero with the stable steady state at one. When these fronts exist, their speeds provide a measure of how long it takes for the invasive species to overtake a node as a function of the distance of that node from the node of introduction. For the aggregate dynamics, we are interested in the growth rate of the total population in the network -that is the sum of the population over all the nodes in the graph. Of particular interest is how this growth rate depends on system parameters and its relation to the existence/ non-existence of traveling fronts.
Consider first the case where G is an infinite homogeneous tree where each node has degree k + 1. Identify one node as the root, label it u 1 , and consider initial conditions wherein u 1 (0) = 1 with all other nodes initially equal to zero. For this initial data, (1.1) can be reduced to the lattice dynamical system, where u n (t) is a representative node from the set of nodes at distance n − 1 from the root.
The evolution equation for any non-root node of the tree can be expressed as In this way, the linear terms can be viewed as a competition between diffusive and advective effects. For example when α 1, the system can be viewed as a discretization of the continuum limit and the dynamics are approximated by the solution of PDE, on a semi-infinite domain with no-flux boundary conditions at the left boundary. In this scenario, the advection dominates the diffusion and localized initial data propagates to the left, i.e. up the tree, and eventually converges to zero. On the other hand, when α is small and the system is near the anti-continuum limit, the reaction terms dominate the diffusive ones and localized initial data will spread down the tree. The goal of this article is to understand the transition from spreading to extinction and, in the case of spreading, to predict the spreading speed of the solution to (1.2).
There has been a tremendous amount of effort dedicated to the study of the Fisher-KPP equation in a variety of contexts. The system was initially studied as a PDE on the real line by Fisher [13] and Kolmogorov, Petrovskii and Piscunov [17]. In the context of lattice dynamical systems, the existence of traveling waves was established in [29]. Most related to the current study is recent work in [19] where spreading speeds for the Fisher-KPP equation on the hyperbolic space H n are investigated. Homogeneous trees are often used as a model for H n and there are strong connections between the results obtained here and those in [19]. Indeed, in [19] it is shown that there exists a critical diffusion coefficient above which the solution converges uniformly to zero, while for values below that threshold compactly supported initial data form a traveling front.
The previous two decades has also witnessed an explosion of interest into the dynamics of differential equations on networks; see [21,22,24,26] and the references therein for an overview of the field. Many studies of dynamical systems on networks treat each node as an individual with discrete state. This work deals with meta-population models, where the dependent variable at each node describes the concentration of individuals; see for example [7,11,14] for studies of disease epidemics within the meta-population paradigm. The work undertaken here is related to [8]. There the Fisher-KPP equation on random graphs is studied with a focus on the growth rate of the overall population in the network. Exponential growth rates of the population growth are studied and sub-linear growth rates are observed for small values of the diffusion parameter. Finally, we also mention [18] where a bistable reaction-diffusion equation is considered and the existence of traveling and pinned, or stationary, fronts are exhibited on random networks and trees.
We now outline our main results, beginning with those that pertain to the case of the infinite homogeneous tree.
Linearly selected spreading speedsl 1 and l 2 critical diffusion coefficients A key feature of the dynamics of (1.2) is non-monotonicity of the spreading speed. We identify two two critical diffusion coefficients: α 1 (k) and α 2 (k). The spreading speed is maximal for α 1 (k) while for α 2 (k) the spreading speed is zero and marks the boundary between spreading and extinction. For any 0 < α < α 2 (k), the initial value problem (1.2) forms a traveling front spreading with the linear spreading speed. For α > α 2 (k) the solution converges pointwise exponentially fast to zero. The decay rate of the selected front is l 1 critical at α 1 and l 2 critical for α 2 .
We also consider the generalization of (1.2) to periodic trees, where the number of children at each level is periodic. In this case a similar reduction to a lattice dynamical system is possible. In analogy with the homogeneous case, we observe that the linear spreading speed for this system is zero when the selected decay rate is l 2 critical and the spreading speed has a critical point when the selected decay rate is l 1 critical.
It is interesting to contrast these dynamics with the case of the classical Fisher-KPP equation on the real line. In that case, as the diffusion parameter α is increased the solution propagates with an increased speed. This makes intuitive sense, as increasing the diffusion constant corresponds to increasing the mobility of the species being modeled. On the infinite tree, it would appear that for α > α 1 (k), increasing the mobility of the species decreases the speed at which it invades. However, we emphasize that while the solution to (1.2) slows down with increasing α > α 1 , the invasion of the species does not.
When α is small the reactive effects dominate the diffusive ones and we obtain an asymptotic expansion of the spreading speed as α(k) → 0. This expansion is independent of k to leading order. We therefore predict that dynamics of (1.2) on more general networks should also be independent of the network topology for α sufficiently small. We test this prediction on Erdős-Rényi random graphs.
Growth rate of the total population We also study the total population on the graph -both in aggregate and at each level of the tree. Our main contribution here is to identify the maximal growth rate of the total population and to locate which levels of the tree this maximal population growth is concentrated on as a function of time.
The exponential growth rate of the population in the linear system is one. In the nonlinear system, growth is saturated at u = 1 and we expect that the growth rate of the population may be slower in this case. In turns out that the critical diffusion parameter α 1 (k) demarcates the boundary between these two possibilities. For α < α 1 (k), the population growth occurs primarily at the front interface and the maximal growth rate depends on the speed and decay rate of the front. When α is larger than α 1 (k), the maximal growth rate is one and occurs in a region ahead of the front interface where the components u n (t) are near zero. The speed at which the maximal growth rate occurs can be determined in terms of the group velocity of the mode k n .
Front propagation and population growth rates on Erdős-Rényi random graphs: numerical simulations Finally, we address the potential relevance of our results to large random graphs of Erdős-Rényi type, [12]. When the number of nodes in the network is large and the average degree of each node is small it may be reasonable to approximate the dynamics on the random graph by that of a tree with k chosen appropriately. We perform numerical simulations that suggest that for small values of α, the dynamics of (1.1) are well approximated by the finite homogeneous tree. We also consider the growth rate of the total population on the graph and demonstrate that this rate can be estimated from the growth rate on the homogeneous tree. In analogy to the tree case, for small values of α the growth rate is observed to be less than the linear growth rate.
The article is organized as follows. In Section 2, we study the pointwise behavior of the solution on the homogeneous tree and prove Theorem 1. In Section 3, we study the behavior of the total population in the tree. In Section 4, we apply some of the insights gained from the homogeneous tree to Erdős-Rényi random graphs. We conclude in Section 5 with a short discussion.
2 Front propagation -the linear spreading speed and proof of Theorem 1 In this section, we study the dynamics of (1.2) and establish the spreading speed (or lack thereof) of solutions propagating down the tree. We first review the notion of linear spreading speeds. We then establish some properties of the linear spreading speed for (1.2) and identify the critical coefficients α 1 (k) and α 2 (k). Noteworthy among these properties is the independence of the spreading speed on the constant k in the asymptotic limit of α small. Finally, we establish that the linear spreading speed is the speed selected in the nonlinear system and we extend the linear spreading speed to the case of periodic trees.
The Fisher-KPP equation, whether it be posed on the real line, lattice or other more exotic domain is an example of a linearly determinate system wherein the spreading speed for the full nonlinear system is equal to the spreading speed of the system linearized about the unstable state, [28]. In this regard, a natural starting point for our analysis is system (1.2) linearized about the unstable zero state. Ignoring the equation for the root, we obtain the linearized problem, We seek separable solutions of the form e λt−γn , γ > 0, and obtain a solvability condition that is known as the dispersion relation, Simple roots of the dispersion relation relate exponentially decaying modes e −γn to their temporal growth rate e λt . Double roots, or more precisely pinched double roots, of (2.2) identify pointwise growth rates of compactly supported initial data; see for example [3,5,6,16,23,25]. Double roots differentiate between absolute and convective instabilities. In the context of (2.1), an absolute instablity refers to exponential growth of both the solution in norm (l 2 for example) as well as pointwise at each fixed node. On the other hand, a convective instablity grows in norm, but is transported away from its original location so that the solution decays pointwise at each node. A convective instability is one where the dispersion relation has simple roots with positive temporal growth rates, but all pinched double roots have negative growth rates, i.e. Re(λ) < 0.
The linear spreading speed is the speed at which compactly supported initial data spreads in the asymptotic limit as t → ∞. This speed can be thought of as the speed for which the system transitions from absolute to convective instability and is defined by the presense of a pinched double root with purely imaginary, or in our case zero, temporal growth rate. We therefore consider the dispersion relation (2.2) in a co-moving frame. Let d s (λ, γ) = d(λ − sγ, γ) denote the dispersion relation in this frame. Then the linear spreading speed for (2.2) is the speed s lin for which a γ lin simultaneously satisfies the equations, together with a pinching condition. Pinching here refers to the requirement that the roots participating in the double root can be traced to opposite sides of the imaginary axis as Re(λ) → ∞.
In order to introduce some further terminology we also mention an alternate, and equivalent, formulation of the linear spreading speed as a minimizer of the envelope velocity, The envelope velocity prescribes the speed at which an exponential solution e −γn propagates in the linear equation (2.1). For systems that obey the comparison principle, any positive solution of the linearized system places an upper bound on the speed of propagation in the nonlinear system and so minimizing over all such speeds provides an upper bound. In this way, the linear spreading speed can be alternatively defined as The equivalence of these two formulations can be observed by imposing that The quantity dλ dγ is known as the group velocity and so the linear spreading speed occurs when the group velocity equals the envelope velocity. In a frame moving at the envelope velocity, this then implies zero group velocity and hence that γ is a pinched double root with s = s env . The notion of group velocity will arise again in Section 3.
Moving forward, we will use the formulation in terms of pinched double roots and define

Properties of the linear spreading speed for (2.1)
Inspection of (2.1) reveals that explicit expressions for the linear spreading speed are not available. We compute solutions numerically and plot the linear spreading speed in Figure 1. Two features of the linear spreading speed are readily apparent: it is non-monotone as a function of α and becomes negative above some critical value. In what follows, we identify the linear spreading speed at two critical values of the diffusion parameter: the diffusion coefficient leading to the fastest spreading speed we denote α 1 (k) whereas the diffusion coefficient with zero spreading speed we denote α 2 (k). The notation is chosen since we will show that for α 1 and α 2 the decay rates of the fronts are critical in l 1 and l 2 respectively. Recall the definition of the l p norm on G, Under our assumption that the values of u at all nodes on each level of the tree are equivalent, the l p norm is expressed as The critical coefficient α 2 (k) We define α 2 as the critical diffusion parameter for which there exists a pinched double root on the imaginary axis for s = 0. To compute α 2 (k), let s = 0 and consider the second component of F set equal to zero, 2), calculated numerically as a function of α for k = 3 (red), k = 4 (black) and k = 5 (blue). Note the critical values α 2 (k) for which the spreading speed is zero and α 1 (k) where the speed is maximal. Also note that as α → 0, these spreading speeds appear to approach a common curve.
from which we compute the critical decay rate Plugging this into the dispersion relation, again with s = 0, we find Note that γ 2 is the l 2 -critical exponential decay rate, in the sense that any steeper exponential is l 2 on the tree while for any weaker exponential is not, see (2.4).
Compute the Jacobians from which the Implicit Function Theorem implies that double roots can be continued as functions s(α) and γ(α) with s (α 1 ) = 0.
Proof. Solving the second equation in F = 0 for s and plugging into the first equation we obtain an implicit equation for γ and α given by Note that G(0) = 0, lim γ→∞ G(γ) = ∞ and that G (γ) > 0. Therefore G is monotone and invertible on its range. Consequently, for every γ > 0 there exists a α for which F (α(e γ − ke −γ ), γ, α) = 0. We observe that the selected decay rate γ is a decreasing function of α and that the spreading speed α(e γ − ke −γ ) may be negative.
Applying the Implicit Function Theorem we find an expression for We observe that this zero exactly for α = α 1 . Differentiating the expression in the parenthesis we obtain This derivative changes sign exactly at α = α 2 and thus e γ − k − 1 + ke −γ is monotone decreasing for α < α 2 and the result follows.
Asymptotics for small α Consider 0 < α 1. We will derive a the leading order term in an asymptotic expansion for the linear spreading speed in this limit. Consider the implicit equation αG(γ) = 1.
In the limit as α tends to zero, it must be the case that G(γ) → ∞. Since γ > 0, G(γ) is dominated by the term γe γ and we expect Inverting, the selected decay rate is to leading order where W is the Lambert W-function. Note that s = αe γ to leading order and obtain to leading order. We remark that the linear spreading speed is independent of the degree k to leading order when α 1. This means that for α small the spreading speed of the system is close to the spreading speed of the Fisher-KPP equation on a lattice in the anti-continuum limit and the arrival time to any node in the tree is dominated by its distance from the root. This has important implications for applications to more general networks, which we will return to later.

Linear Determinacy
We have thus far focused entirely on the spreading properties of the linearized system. We now turn our attention to the nonlinear system and introduce some standard terminology. Define the invasion point for (1.2) as The selected spreading speed is then defined as We prove the following theorem which establishes that the selected spreading speed is the linear spreading speed.
Theorem 1 states that the system is linearly determinate and that the spreading speed of the nonlinear system is the same as the system linearized about the unstable state. We remark that had (1.2) been posed on Z, rather than N, then this result would be a direct consequence of [28]. Rather than adapt that proof to this context, we will construct explicit sub and super solutions that bound the spreading speed of system (1.2). Since the proof proceeds along familiar lines, we delay it until Appendix A.

The periodic tree
We now consider the generalization to an infinite rooted tree where the number of nodes per level is periodic. We again assume initial data consisting of one at the root and zero elsewhere and obtain a similar reduction to a lattice dynamical system with the evolution given by du n dt = α (u n−1 − (k n + 1)u n + k n u n+1 ) + f (u n ), n ≥ 2, where k n is periodic in n with periodic m, i.e. k n+m = k n . Our main goal is to study the linear spreading speed as a function of α and to show that the l 2 critical front is stationary while the l 1 front is a critical point of the linear spreading speed.
To this end, we linearize (2.6) about the unstable state at zero and obtain the following equation at any (non-root) node du n dt = α (u n−1 − (k n + 1)u n + k n u n+1 ) + u n .
We seek separable solutions of the form where γ > 0 and U n is periodic with period m. Since U n is periodic we focus on one segment of length m which we denote U = (U 1 , U 2 , . . . , U m ) T . Plugging in the ansatz (2.7), U must satisfy where for m ≥ 3, Note that M , after adding an appropriate multiple of the identity is positive, so it has a real principle eigenvalue λ(γ, p) and corresponding principle eigenvector with positive entries, which we denote U . In this manner, we obtain a scalar dispersion relation that relates the decay rate γ to the envelope speed s, d(γ, s) = αλ(γ) − sγ + 1 = 0.
Pinched double roots are therefore solutions of the system of equations We are interested in the linear spreading speeds corresponding to l p critical decay rates. Critical fronts in l p have decay rate Then α p and s p can be defined through the relationship F (s p , γ p , α p ) = 0 and we may express F entirely as a function of p. Set the second component of F equal to zero, from which .
Then for the first component of F to be zero it is required that Combining these two expressions we obtain an equation for the linear spreading speed as a function of p, in terms of the principle eigenvalue of M , We now express the matrix M in terms of p. We first note that m . Before stating our results, we show that M (γ p ) is similar to M (γ q ), where q is the Hölder dual of p. To do this, we compute the determinant, det(M (γ)) = (−1) m+1 k 1 . . . k m e −mγ + e mγ + ∆(k 1 , . . . , k m ), where the terms ∆(k 1 , . . . , k m ) do not depend on γ. With γ = γ p , we find Consequently, the characteristic polynomial is invariant when p is replaced with its Hölder dual q and therefore the matrices M (γ p ) and M (γ q ) are similar. Proof. Let q be the Hölder dual of p. Making use of the fact that M (p) and M (q) are similar it follows that λ(p) = λ(q) and hence that λ (p) = λ (q)q . Since q = −(q/p) 2 this rewrites as In the case p = q = 2 this reduces to λ (2) = 0. Thus s 2 = 0 as desired.
Proof. First, we compute In the case p = ∞ (1) the matrix M has row (column) sum zero, hence the vector of all ones is in its (left) null space. As this vector is positive, it is necessarily the principal (left) eigenvector, hence zero is the principal eigenvalue. With λ = 0 the expression for s p reduces to mp 1 log k l . Further s p vanishes when λ = 0.
We now present some results for periodic trees with period two, where explicit expressions can be obtained.
Example 2.4. Consider the periodic tree with m = 2, where the degrees alternate between k 1 and k 2 . We compute An explicit calculation reveals On the left, we plot α 1 as a function of k 1 with k 2 fixed to preserve the mean degree. On the right, we plot α 2 as a function of k 1 . Note that in both case the periodic heterogeneity increases the critical diffusion rates. and We plot α 1 and α 2 in Figure 2.
Note that in both cases, the heterogeneity induced by the periodic coefficients leads to increases in the critical values as compared to the case of the homogeneous tree. Take for example k 1 = 1, k 2 = 5 and compare α 2 to the critical value for the homogeneous tree with degree k = 3. In the homogeneous case, the transition occurs at approximately 1.86 whereas the periodic case has a critical value near 5.1. Heterogeneity in the tree can therefore allow for front propagation over a much larger region in parameter space.

Aggregate dynamics: growth of the total population
In Section 2, we determined the pointwise behavior of u n (t). We now turn our attention to the growth rate of the total population within the homogeneous tree. Let Note that P (t) is equivalent to the l 1 norm of the solution. We will find it convenient to express P (t) in terms of ∞ n=1 w n (t), where w n (t) is the total population at level n and is related to u n via w n (t) = k n−1 u n (t). (3.1) The evolution of w n (t) is given by the lattice dynamical system dw n dt = α (kw n−1 − (k + 1)w n + w n+1 ) + k n f (k −n w n ), n ≥ 2, This system is equivalent to (1.2) after transforming to an exponentially weighted space.
In Figure 3, we show numerical simulations of (1.2) and (3.2), with w n normalized. We observe that for small values of α, the maximal growth rate appears to occur at the front interface, while for larger values of α the maximal growth rate is achieved ahead of the front interface. This transition occurs at α = α 1 (k). We have two primary goals in this section. These are to identify the maximal growth rate of P (t) and determine the location in space-time for which the maximal growth rate is achieved. As was the case in Section 2, we begin with the system linearized about the unstable state and then consider the implications for the nonlinear system.

Linear growth rates for P (t)
The growth rate of P (t) for the linear problem can be determined from the l 1 spectrum of the linear operator α∆ G + Id. The spectrum of the graph Laplacian on homogeneous trees has been the focus of intense study and its spectrum for any l p space is known explicitly, see for example [20]. The l 1 spectrum of ∆ G is an ellipse which is contained in the left half of the complex plane aside from the point at zero. As a result, the l 1 spectrum of α∆ G + Id extends to λ = 1 and we expect a maximal growth rate of e t , modulo algebraic heat decay corrections.
Of course, we are ultimately interested in the population of the nonlinear system. It is conceivable that the nonlinearity could constrain the growth in such a way that slower growth rates for P (t) are observed. To determine whether this occurs, we require pointwise information on the growth rates of the total population.

Pointwise growth rates and spreading speeds for the total population
In a fashion analogous to the system for u n , we can linearize (3.2) about the unstable zero state and derive a dispersion relation for w n , as well as the equivalent dispersion relation in a co-moving framed s (λ, γ) =d(λ + sγ, γ). We again remark that transforming from u n to w n amounts to transforming u n to a weighted space. Any pinched double root (λ * , γ * ) of the dispersion relation d s (λ, γ) corresponds to a pinched double root of the dispersion relationd s for (λ, γ) = (λ * + s log k, γ * − log k); see the following calculation Recall that the λ * value gives the pointwise growth rate of the solution in a frame moving with speed s. The decay (growth) rate of any pinched double root for the local population u n (t) is thus increased by s log(k) in the weighted space. Thus, pointwise growth rates of w n may be positive even when growth rates for u n are negative. We have the following.
Lemma 3.1. The linear growth rate is maximized in a frame moving at speed i.e. at the u-group velocity of the mode e −n log k .
Remark 3.2. Lemma 3.1 states that the maximal growth rate of a linear problem in an exponentially weighted space occurs on a ray in space-time moving with the group velocity of the exponential weight. This is a general fact and there are several alternative avenues by which to verify the statement of Lemma 3.1. Similar results are often obtained by a asymptotic analysis of the Green's function via stationary phase. In the context of a homogeneous tree, this fact can be verified from a direct analysis of the heat kernel, see for example [9,10]. For the case of parabolic partial differential equations, this was shown in Lemma 6.3 of [16]. Finally, we mention [15], where the appearance of pointwise instabilities when transforming to weighted spaces leads to some unexpected spreading speeds in a system of coupled Fisher-KPP equations.
We now compare the speed s * to the linear spreading speed s lin and make the the following observations, • For α = α 1 , the maximal growth rate is realized at the linear spreading speed; that is, s * = s 1 • For α < α 1 we have that s * < s lin , while s * > s lin for α > α 1 .
These observations lead us to conjecture that the the population grows fastest when s * > s lin since in this regime the maximal growth rate is achieved ahead of the front interface where the nonlinearity is negligible. . The blue curves are u n (t) while the red curves depict the normalized population at each level, i.e. w n (t)/max n (w n (t)). Note that 0.2 < α 1 (3) < 0.8 < α 2 (3) < 2.2. For α = 0.2, we observe that the maximal population is concentrated at the front interface. For α = 0.8, the maximal population is concentrated ahead of the front interface. Finally, for α = 2.2 the local population at any fixed node converges to zero, but the total population grows and eventually is concentrated at the final level of the tree.

Growth rates in the nonlinear problem
We return our attention to the full nonlinear system and study the growth rate of the total population. We expect that the growth rates predicted by the linearized system to place an upper bound on the growth rates in the nonlinear system.
A key distinction is made between α < α 1 (k) and α > α 1 (k). For α below the l 1 critical value, the maximal linear growth rate at level n is achieved in a frame of reference that propagates slower than the linear spreading speed. In view of Theorem 1, this means that the invasion speed of the species exceeds s * . Furthermore, since γ lin > log k we expect the mass in the system to be added at the front interface. The dynamics of the total population can then be approximated by When γ lin > log(k) then this sum converges and we expect Therefore, we expect a maximal growth rate of e log(k)s lin t when γ lin > log(k). Note that for α 1 , it is the true that log(k)s 1 = 1.
For α > α 1 (k), the speed s * exceeds that of s lin and therefore the maximal linear growth rate is achieved ahead of the front interface. Since the solution u n (t) is pointwise close to zero in this region, we expect that the linear dynamics remain valid here and that P (t) ≈ Ce t for these parameter values. This is confirmed in numerical simulations, see Figure 4.  Figure 4: On the left, we compare predictions for the exponential growth rate of the maximum of w n (t) versus α (blue line) and the observed growth rates in direct numerical simulations (asterisks) for k = 3. On the right, we compare numerically observed spreading speeds for w n (t) (asterisks) versus linear spreading speeds determined numerically from the pinched double root criterion applied tod s (γ, λ) (blue line). Here we have taken k = 5.

Erdős-Rényi random graphs
We now turn our attention to more complicated network topologies. The analysis performed on the homogeneous trees leads us to make two qualitative predictions for more general networks. First, based upon the leading order independence of the linear spreading speed on the degree of the tree, we expect that for α asymptotically small the pointwise dynamics of the system should be well approximated by a finite tree where the only relevant feature is the distance of a node from the location of initiation. Second, based upon our analysis in Section 3 we expect sub-linear growth rates for the total population for α small. We test both these predictions in numerical simulations of Erdős-Rényi random graphs.
We first review the Erdős-Rényi random graph model [12]. Let N denote the total number of nodes in the graph and select some p ∈ (0, 1). Erdős-Rényi random graphs are constructed by assigning edges to the graph with fixed probability p. We are interested in the case where N 1 and the expected degree of each node k ER + 1 = N p is small. We select one node at random and call that node the root. We denote the solution at that node u 1 (t) and consider (1.1) with u 1 (0) = 1 and zero initial condition at every other node.
Erdős-Rényi random graphs with p small share a close relationship with trees. The root will be connected to roughly k ER +1 other nodes, which in turn will each be connected to roughly k ER other nodes. In this manner, the goal here is to employ some of the results obtained for homogeneous trees to make predictions for Erdős-Rényi random graphs.
We construct networks with size ranging from N = 60, 000 to N = 500, 000, see [2], and vary the expected degree from three (k ER = 2) to sixteen (k ER = 15). Note that for the small p values considered here, the constructed networks are disconnected with probability one. Therefore, when we discuss the network dynamics below we always restrict to the largest connected component.
Numerically, solutions are computed using explicit Euler with timestep ∆t = .0025. Changing the timestep did not significantly alter the results.

Pointwise dynamics: arrival times
In this section, we focus on the pointwise dynamics of the system and describe numerical results characterizing the spreading speed on random graphs of Erdős Réyni type.
For finite graphs the only stable steady state is the uniform state at one. We assign to each node a label i ∈ {1, 2, . . . , N }. Let d(i, j) be the shortest distance along the network between nodes i and j. Denote the node where the species is initialized as i = 1 and specify initial conditions where u 1 (0) = 1 and u j (0) = 0 for all j = 1. The primary descriptor of the dynamics in this situation is τ i , the arrival time at node i ∈ V , defined as for some threshold 0 < κ < 1.
Arrival times as a function of the distance from the initial node are plotted in Figure 5. For very small α, we observe a strong linear relationship between the arrival time and distance from the node i = 1. As α is increased this linear relationship becomes less apparent.
We compute the mean arrival time at level l, Based upon this data we find a best fit linear regression line as follows. In order to minimize the effect of initial transients or boundary effects, we consider only those M l corresponding to distances greater than two from the initial node or less than two from the furthest node. A linear regression is then computed and the reciprocal of the slope is called the numerically observed mean spreading speed. These regression lines are shown in green in the Figure 5 with green asterisks denoting M l .
Mean spreading speeds are computed and plotted in Figure 6. We compare these speeds to spreading speeds on homogeneous trees computed in Section 2. Before drawing any conclusions, some remarks regarding the relationship between random graphs and regular trees are in order. Important in any comparison of spreading speeds in the random graph and the homogeneous tree is a selection of a value of k to be used for comparison. This choice is not obvious, nor is it clear that one such k will be sufficient. For large Erdős-Rényi graphs, the average degree in the graph will converge to k ER + 1 as N → ∞. At the same time, however, if k ER is small then the Erdős-Rényi graph will be disconnected with probability one. Within the largest connected component the average degree will exceed k ER + 1. The simplest choice is then to take k to be the average degree of all nodes over which the linear regression is computed. With this choice, the observed and theoretical speeds are close, however the linear spreading speed consistently overestimates the observed mean spreading speed. This discrepancy could be due to several factors including the heterogeneity of the network. However, we also note that in the system (1.2) it takes some time for initial data to converge to a traveling front and even then the convergence is notoriously slow [4]. As a result, it is perhaps unrealistic to expect that spreading speeds measured directly in numerical simulations on the random graph to match exactly with predictions from the linear spreading speed.
We comment briefly on the dynamics for large values of α. For α sufficiently large the first arrival times cease to be monotonic with respect to distance from the original node. For even larger values of α the arrival time for most nodes on the graph is approximately constant. This is to be expected. For α sufficiently small the dominant eigenvalue of the linearization about zero is one, whose eigenvector is the constant vector 1. Projecting the initial condition onto this vector we expect leading order dynamics from which we find a homogeneous arrival time that scales with log N . This is consistent with the timescale observed for large values of α.

Aggregate dynamics: population growth rate
We now investigate the exponential growth rate of the total population for (1.1) for Erdős-Rényi random graphs. We compare these growth rates to those observed on the homogeneous tree.
The logarithm of the total population is shown in Figure 7. After an initial transient, the growth rate appears to be roughly linear corresponding to exponential growth. We measure this growth rate and compare it to similar computations in the homogeneous tree.
As was the case for the pointwise analysis, an immediate challenge lies in determining the precise value of k to use in the comparison to the homogeneous tree. Given the variation in the average degrees across levels of the network, we instead compare numerically observed growth rates in the random graphs with observed growth rates in homogeneous trees that we expect to provide bounds on the total growth rate. The results of these simulations are depicted in Figure 8.

Discussion and outlook for future work
We studied the existence and speed of invasion fronts for Fisher-KPP equations on homogeneous trees and Erdős-Réyni random graphs. For the homogeneous tree these fronts are linearly determined and propagate at the linear spreading speed. An interesting property is that these speeds are non-monotone with respect to the diffusion parameter and eventually become zero. Therefore, for values of α > α 2 the solution converges pointwise exponentially fast to zero. We find that the critical point of the linear spreading speed at α = α 1 occurs for the parameter value where the maximal linear growth rate of the total population in the network occurs at the linear spreading speed. For larger values of α the front is slowed due to some fraction of the population escaping down the graph. The increased mobility of the species allows the population to explore nodes with only small concentrations of the species where the local growth rate is maximal.
An interesting feature of the linear spreading speed is that it is independent of the degree of the nodes in the network in the asymptotic limit as α → 0. This suggested that the spreading properties observed in the homogeneous tree should be present in random graphs and explored this for the case of random networks of Erdős-Réyni type. For α asymptotically small, we observed traveling fronts and characterized their speed. It was also the case that the exponential growth rate of the population was less than the maximal growth rate of one and bounded this growth rate by the corresponding growth rate on the homogeneous tree.
Erdős-Réyni graphs are known to share a close relationship to trees. It would be interesting to extend the results here to different, more realistic, random networks such as small world [27] or scale free networks [1]. The independence of the spreading speed on the degree for small α suggests that one should continue to observe traveling fronts on these more general networks that are well approximated by the lattice Fisher-KPP equation. The situation is more complicated for larger values of α where the heterogeneity of the network can not be denied. It would be interesting to characterize whether critical diffusion rates analogous to α 1 and α 2 for these more general classes of graphs. Example 2.4 gives some indication that the effect of heterogeneity can be quite severe and so we expect that analysis of these systems would require a more general approach.
As a result, we have that the linear terms in N (ū n ) are zero and which is positive by the KPP assumption f (u) < f (0)u. Thus, the exponential is a super-solution for all n ≥ 2. Since C > 1, the super-solution is always one at the root and we have thatū is a super-solution.
We now turn our attention to the construction of sub-solutions. Let 0 < µ < 1 and consider the linear equation Exponential solutions e −γ(n−st) are obtained for γ and s satisfying the following equation For each µ, there exists s µ such that (A.1) has no real solutions for s < s µ . Let γ µ ∈ C be a complex solution. Then let φ µ (y) = e γµy + c.c., where c.c. denotes the complex conjugate and fix a < b such that φ(y) > 0 for y ∈ (a, b) and φ(a) = φ(b) = 0. We use these functions to construct sub-solutions in the following Lemma.
Lemma A.2. Fix α < α 2 (k) and let s < s lin . Let µ < 1 such that s < s µ . Then there exists T > 0, a function φ µ (y) and an * such that is a sub-solution for all 0 < < * and for all t > T .
Proof. For s < s lin , select µ < 1 so that s < s µ and consider the function φ µ in (A.2) whose support is the interval [a + st, b + st]. Let T be sufficiently large so that sT + a > 1. Now consider N (u n ) = u n − α (u n−1 − (k + 1)u n + ku n+1 ) − u n − (f (u n ) − u n ) .
For u small, there exists a C > 0 such that −Cu 2 < f (u) − u. Let be sufficiently small so that this bound holds for (A.3). Then N (u n ) < φ(n − st) (−(1 − µ) + C φ(n − st)) , and if we restrict then N (u n ) < 0 and u n is a sub-solution.
Proof. We follow the proof of Lemma A.3 and compute N (u n ). The calculation is exactly the same, aside from the root. There we calculate Since both terms are positive, we have that N (u 1 ) is positive as well and the proof is completed.
We now prove Theorem 1. The claim of pointwise convergence for α > α 2 (k) follows from Lemma A.3 and the fact that s env (γ) < 0. Consider then α < α 2 . Here, we have that the linear spreading speed is positive and obtained for some value of s lin with selected decay rate γ lin . Then since C > 1 in Lemma A.1 we have thatū n (0) ≥ u n (0). Due to the comparison principle, this relationship holds for all time and we obtain an upper bound on the location of the invasion point: κ(t) ≤ s lin (α)t. To obtain a lower bound, we consider the sub-solutions constructed in Lemma A.2. Select s < s lin (α) and consider u n from Lemma A.2. Let t 1 > T . By the maximum principle, we have that u n (t 1 ) > 0 for all n ∈ N. Since u n (t 1 ) is compactly supported, we can select > 0 sufficiently small such that u n (t 1 ) ≤ u n (t 1 ). Therefore, we have that κ(t) ≥ st. This holds for any s < s lin (α) and we find s lin (α) ≤ s sel ≤ s lin (α) and Theorem 1 is established.  Figure 6: Speed associated to the mean arrival times in numerical simulations on an Erdős Réyni graph with N = 60, 000 and expected degree k ER = 2 are shown in asterisks. The blue curve is the spreading speed predicted by the analysis in Section 2 for the homogeneous tree with k = 2.54. This value is chosen since it is one less than the mean degree of the network over those nodes with distance between 3 and 12 from the original location.   Figure 8: Numerically calculated exponential growth rate for the Erdős-Rényi graph. On the left, N = 500, 000 and observed growth rates are plotted as circles. The asterisks are the corresponding growth rates in the homogeneous tree with depth 13. The lower curve is degree k = 3 while the larger curve is degree k = 4. On the right are the same computations, but for the Erdős-Rényi graph with N = 60, 000 and for homogeneous trees with k = 2 and k = 3.