  In the May-Leonard model of three cyclically competing species, we analyze
the statistics of rare events in which all three species go extinct due to
strong but rare fluctuations. These fluctuations are from the tails of the
probability distribution of species concentrations. They render a coexistence
of three populations unstable even if the coexistence is stable in the
deterministic limit. We determine the mean time to extinction (MTE) by using a
WKB-ansatz in the master equation that represents the stochastic description of
this model. This way, the calculation is reduced to a problem of classical
mechanics and amounts to solving a Hamilton-Jacobi equation with zero-energy
Hamiltonian. We solve the corresponding Hamilton's equations of motion in
six-dimensional phase space numerically by using the Iterative Action
Minimization Method. This allows to project on the optimal path to extinction,
starting from a parameter choice where the three-species coexistence-fixed
point undergoes a Hopf bifurcation and becomes stable. Specifically for our
system of three species, extinction events can be triggered along various paths
to extinction, differing in their intermediate steps. We compare our analytical
predictions with results from Gillespie simulations for two-species
extinctions, complemented by an analytical calculation of the MTE in which the
remaining third species goes extinct. From Gillespie simulations we also
analyze how the distributions of times to extinction change upon varying the
bifurcation parameter. Even within the same model and the same dynamical
regime, the MTE depends on the distance from the bifurcation point in a way
that contains the system size dependence in the exponent. It is challenging and
worthwhile to quantify how rare the rare events of extinction are.

Biological Physics