Jacopo Grilli
György Barabás
Matt Michalska-Smith
Stefano Allesina

This talk was given at the annual NetSci conference in Indianapolis, IN USA on 20 June 2017. Text here will supplement slide material.

You can remove this pane (and get it back) by tapping v on your keyboard.

Some of the slides will have animations which require advancing the slide partway through these explanations. These will be indicated by "next slide" in these notes.

Slides can be advanced by tappingSpc, , PgDn, or n.

Finally, this presentation will look best in full-screen, but should be more or less intelligible on smaller screens (no guarantees about phone screens) So, whenever you're ready...

next slide

Acknowledgements

Jacopo Grilli

Gyuri Barabás

Stefano Allesina

E. L. Sander
C. A. Marcelo Serván
D. Maynard
E. Leigh
J. Aljadeff

A. Monteiro
J.A. Capitán
F. Brandl
J.M. Levine
J-F. Laslier

JEROME VALLETTE/AFP/Getty Images

Nacho Doce/Reuters

$$\sum_i x_i = 1$$

A key assumption of the model we will build is inspired by that first observation of the rainforest: that the space is entirely filled and any gap in the canopy is rapidly filled by a new individual. Because of this, the dynamics can be said to be approximately "zero-sum."

This means that the number of individuals in the community is approximately fixed and we can just model the relative species abundances $x_i$. Pictorially, this means that:

next slide

when a tree dies

next slide

a seedling is chosen at random to fill the gap and

next slide

this seedling grows into an adult

next slide

maintaining the total number of individuals in the community but potentially changing the relative abundances of the species involved.

$$\propto f_i x_i$$

$$\begin{align}\mathcal{P}(i \text{ wins against } j) &= H_{ij}\\ &= 1 - H_{ji}\end{align}$$

$$d_i x_i$$

We can model the processes at play in a more mathematically rigorous manner. For instance,
next slide
each species has a birth rate, which is proportional to its fecundity (the number of seedlings produced).
next slide
It also has a death rate. We assume that a random individual dies at each time point. Since the individual is chosen at random, the probability of any individual dying is proportional to its relative species abundance. Finally, some things happen in-between these events, which we will model as competition.
next slide
We can assume very simplified competition, where the individuals compete in pairs and the winner of any competitive bout is determined by a tournament matrix $\textbf{H}$ which encapsulates the probabilities of any species $i$ beating another species $j$. Of course, this means species $j$ has a $1-H_{ij}=H_{ji}$ probability of winning the bout.

Neutral Theory of Biodiversity (Hubbell 2001)
$$f_i = d_i = 1 \text{, } H_{ij} = H_{ji} = 0.5$$
Competitive Tournaments (Allesina & Levine 2011)
$$f_i = d_i = 1 \text{, } (H_{ij}, H_{ji}) = \begin{cases} (1,0)&\\ (0,1) \end{cases}$$
Competitive Hypertournaments
$$f_i = d_i = 1 \text{, } \text{any nonnegative } (H_{ij}, H_{ji}): H_{ij} + H_{ji} = 1$$
General Model
$$\text{any nonnegative } f_i, d_i, (H_{ij}, H_{ji}): H_{ij} + H_{ji} = 1$$

Using this simple framework, we can recapitulate some models that have already been well-studied, as well as propose some advancements. For instance, the neutral theory of biodiversity proposed by Stephen Hubbell in 2001 can be replicated in this framework by assuming
next slide
that the fecundities and death rates are equal for all species (we can consequently assume that this value is $1$ for convenience) and that the probability of any given species winning a competitive bout is equal to that of its competitor winning.
next slide
Likewise, if the probabilities of victory are exclusively $1$'s and $0$'s, that is, species $i$ always wins against species $j$ and always loses to some species $k$, we are dealing with competitive tournaments, which have been explored in an ecological context by Allesina and Levine in PNAS in 2011.
In this work, we build on these models by first allowing non-binary probabilities of victory
next slide
and second allowing full flexibility of parameterization, breaking that reliance on equal physiological parameters.
next slide
We will first consider the former and get to the latter toward the end of the presentation.

$$\propto f_i x_i$$

$$\begin{align}\mathcal{P}(i \text{ wins against } j) &= H_{ij}\\ &= 1 - H_{ji}\end{align}$$

$$d_i x_i$$

Assume $f_i = d_i = 1 \quad \forall i$

$\dfrac{\text{d} x_i }{\text{d} t} = $ $x^2_i$ $+$ $2 x_i \displaystyle\sum_{j \neq i} H_{ij} x_j$ $-$ $x_i$

$$\begin{align} \dfrac{\text{d} x_i }{\text{d} t} &= x^2_i + 2 x_i \sum_{j \neq i} H_{ij} x_j - x_i\\ &= x_i \left(2 \sum_j H_{ij} x_j - 1 \right)\\ &= \sum_j \left( H_{ij} - H_{ji}\right) x_i x_j\\ &= \sum_j P_{ij} x_i x_j \end{align}$$

$$\mathbf{H} = \left[\begin{array} 00.5 & 1 & 0 \\ 0 & 0.5 & 1 \\ 1 & 0 & 0.5 \end{array}\right] $$

$$\mathbf{P} = \left[\begin{array} 0~~~0 & ~~~1 & -1 \\ -1 & ~~~0 & ~~~1 \\ ~~~1 & -1 & ~~~0 \end{array}\right]$$

Lets look at a concrete example. Imagine three species: rock, paper, and scissors which interact in an intransitive cycle of complete dominance. The tournament matrix $\mathbf{H}$ is thus shown here, with a $1$ in the top middle, for instance, indicating that paper always wins against rock. An individual of a given species always has an equal probability of winning or losing when competing against another individual of the same species, hence the $0.5$ on the diagonal. This $\mathbf{H}$ corresponds to the following $\mathbf{P}$, which shows that, in the previously mentioned bout, the result is that rock loses an individual and paper gains one. These dynamics result is neutral cycling about an equilibrium of one third for each species.

How many species do we expect to survive?

Method

Generate many random (hyper) tournaments of $s$ species by drawing interaction strengths($H_{ij}$) at random
Run dynamics
Count the $k$ survivors
what is $\mathcal{P}(k|s)$?

$$\mathcal{P}(k|s) = \begin{cases} 0 &\text{for } k \text{ even}\\ {{s}\choose{k}} 2^{1-s} &\text{otherwise}\end{cases}$$

Result

So, we saw that in some cases, all of the species will survive and in other cases, some will go extinct. One might reasonably wonder how many species will survive on average. Put another way, if we were to open all the gates at a zoo and let the animals interact, how many species would we expect to coexist after all of the extinctions? We can approach this problem from the perspective of random matrices. In this approach, we generate a random (hyper)tournament matrix of probabilities and then run the dynamics to see how many of the original species survive. If we do this many times, we can then estimate the expected number.
next slide
When we do this, we find that this probability is $0$ for an even number of surviving species (as was also shown by Allesina & Levine 2011) and follows a binomial distribution with $p = 0.5$ otherwise. This means that, regardless of the number of initial species, the expected number that survive is half of that. Put another way, there is no upper limit for the number of surviving species, which suggests that the high biodiversity we observe in nature maybe shouldn't be so surprising after all.

One of the interesting things about the transition from strict tournaments to hyper tournaments, i.e. moving from binary probabilities of victory to real-valued probabilities, is that we can actually construct an infinite number of matrices $\mathbf{H}$ that will converge to an arbitrary equilibrium $x^\ast$. What this means is that we can fine-tune a system to allow the persistence of an even number of species.
next slide
Yet, these cases are highly degenerate and will collapse into an odd number of species with any modification of the coefficients of $\mathbf{H}$ or parameters $f_i$ or $d_i$. For this reason, I will always mark this example with an asterisk (*), denoting its special status.

Conclusions

For any matrix $\mathbf{H}$ on $s$ species

Irrespective of (positive) initial conditions, $k$ species will coexist and $s - k$ will go extinct
Those that coexist, cycle neutrally around the equilibrium point $x^\ast$
We can choose $x^\ast$ arbitrarily

For a random $\mathbf{H}$ on $s$ species

On average, $k = s/2$

Unfortunately, even this small change faces us with a dilemma. We can consider true HOI in which there is a coefficient for each unique triplet of interacting species, but then our number of parameters in the model quickly becomes intractable from an empirical standpoint. An alternative is to consider an interaction chain,
next slide
in which first two of the species compete and then the winner of this pairwise bout competes with the third species to determine the overall winner. This has the advantage of reducing the number of parameters to exactly that needed for the pairwise case, but comes with its own shortcoming: now the order of interactions matters. For instance, in this example, if we re-ordered the interactions such that the orange seedling interacts first with the purple
next slide
seedling, rather than the blue, we see that it is purple, rather than orange that advances to the second bout producing a new interaction that wasn't present in the previous chain. If blue is able to defeat purple then the final victor can also differ.

$$\dfrac{\text{d} x_i }{\text{d} t} = \sum_j P_{ij} x_i x_j$$ $$\dfrac{\text{d} x_i }{\text{d} t} = \sum_j P_{ijk} x_i x_j x_k$$

where $\mathbf{P}$ is a tensor in which $$P_{ijk} = 2 H_{ij} H_{ik} - H_{ji} H_{jk} - H_{ki} H_{kj}$$

Conclusions

The addition of higher order interactions stabilize the dynamics but do not change the equilibrium
Higher orders of higher order interactions speed up convergence to the equilibrium
This is all despite not adding any new information/parameters into the model!

Neutral Theory of Biodiversity (Hubbell 2001)
$$f_i = d_i = 1 \text{, } H_{ij} = H_{ji} = 0.5$$
Competitive Tournaments (Allesina & Levine 2011)
$$f_i = d_i = 1 \text{, } (H_{ij}, H_{ji}) = \begin{cases} (1,0)& \text{w.p. } 0.5\\ (0,1) & \text{otherwise} \end{cases}$$
Competitive Hypertournaments
$$f_i = d_i = 1 \text{, } \text{any nonnegative } (H_{ij}, H_{ji}): H_{ij} + H_{ji} = 1$$
General Model
$$\text{any nonnegative } f_i, d_i, (H_{ij}, H_{ji}): H_{ij} + H_{ji} = 1$$

What happens when these rates are allowed to vary?

In this case, we see that the rock-paper-scissors dynamics are no longer neutral cycles, but rather unstable. The cycles are amplifying and will eventually lead to extinction and monodominance. This is also seen in the second
next slide
and third
next slide
networks from the before examples.
next slide
In the case of the dynamically-fragile four-species case, we see that this modification of the physiological parameters leads to the rapid extinction of one of the species, making the number of remaining species odd, and those remaining also show the amplifying cyclical behavior. When we sample three species for each competitive bout,
next slide
we see the recovery of a stable fixed point for all cases, with the exception of the special case, where now an odd number of species coexist stably.

Conclusions

This theory relies on higher-order interactions
Leads to the stable coexistence of an arbitrary number of species

Can produce any species abundance distribution

Is robust to changes in parameters
More experimentally tractable due to parameterization of HOI as direct consequence of pairwise interactions

Such interactions chains can arise in nature due to differing scales of space or time

In conclusion, we have presented a theory which allows for the coexistence of very many species, even in the case of randomized interaction and physiological parameters. With the addition of higher-order interactions, the dynamics are stabilized, and this stability is enhanced by incorporating "higher" higher-order interactions. Remarkably, this formulation can yield stable, fixed-point dynamics leading to an arbitrarily specified equilibrium.

These changes are effected without adding additional information or parameters into the model due to our use of interaction chains of pairwise competitions rather than "true" higher order interactions. This choice both severely reduces the number of parameters in the model and makes it empirically testable.

Finally, if you see such interaction chains as a limiting assumption, this is actually a very natural interpretation when considering differing scales of space or time. For instance, two seedlings could sprout more closely to each other and only interact with a third after the first competition is already decided. Alternatively, there could be a temporal distinction in which the links of the interaction chain are separated by timescale, for instance if two seedlings sprout first and the third is delayed until after the winner of the first competition is decided.

Questions?

Slides available on Michalska-Smith.com

What happens when we add stochasticity to this framework? This plot depicts stochastic trajectories for a rock-paper-scissors system, when we vary the number of seedlings competing ($m$). Dotted lines correspond to the deterministic approximation. For $m = 1$ and equal physiological rates, one recovers a neutral model. For equal physiological rates and pairwise interactions ($m = 2$) stability is only neutral and the deterministic solution does not describe the stochastic trajectories. Densities fluctuate around the fixed point, but the amplitude of these oscillations is driven by stochasticity, rather than being fixed as in the deterministic case. For higher- order interactions ($m = 3$), densities fluctuate around the stable equilibrium of the deterministic system. When physiological rates are different, the stochastic trajectories are predicted by the deterministic equations. For $m \le 2$, extinctions occur very rapidly, and are driven by the deterministic dynamics. For $m = 3$, trajectories fluctuate around the stable equilibrium. For these simulations the total population abundance was set to $1000$.