Routing functions for parameter space decomposition to describe stability landscapes of ecological models

Elizabeth Gross; Jonathan D. Hauenstein; Joseph Cummings; Kyle J.-M. Dahlin

REVIEW 3 major objections 5 minor 1 cited by

Routing functions for parameter space decomposition to describe stability landscapes of ecological models

T0 review · 3 major / 5 minor · reviewed 2026-08-16 · deepseek-v4-flash

Challenge this review Re-run · record.json Download PDF Read on arXiv ↗ Review JSON

Pith's one-line read Routing functions decompose ecological parameter space into regions with a fixed set of stable steady states.

desk verdict A genuinely useful application of routing functions to ecological stability landscapes, but the limit-cycle claim is too strong and the unshown 31-factor boundary polynomial makes the coral maps conditional. read the letter →

arxiv 2505.00128 v1 pith:SZTAC4G7 submitted 2025-04-30 q-bio.PE math.AG

Joseph Cummings , Kyle J.-M. Dahlin , Elizabeth Gross , Jonathan D. Hauenstein This is my paper · ORCID

classification q-bio.PEmath.AG MSC 92D4037N2514P1034D20

keywords stabilitylandscaperoutingfunctionrealalgebraicgeometrysemi-algebraicsetsRouth-Hurwitzboundarymultistabilitycoral-bacteriasymbiosislimitcycles

verification ladder T0 review T1 audit T2 compute T3 formal

The pith

A machine-rendered reading of the paper's core claim, the machinery that carries it, and where it could break.

The reading

The paper aims to turn the stability question for ecological ODE models into a symbolic geometry problem. It shows that outside a boundary built from singular, Routh-Hurwitz, and coordinate conditions, each connected region of parameter space has a constant number and type of stable steady states, and routing functions compute those regions. Applied to a coral-bacteria symbiosis model with several rates fixed, the method maps where coexistence, pathogen exclusion, and bistability are stable, and it identifies a region with no stable steady state, implying sustained oscillations. If the framework generalizes as claimed, ecologists could read collapse and recovery thresholds directly off such maps.

What carries the argument

The key object is a routing function: for boundary polynomials $g_1,\dots,g_m$, the rational map $r_c(a)=g_1(a)\cdots g_m(a)/(1+\sum_i(a_i-c_i)^2)^D$ with $D$ large enough, whose nondegenerate critical points (routing points) and gradient-ascent paths form a graph whose connected components match the connected components of $\mathbb{R}^k_{>0}\setminus B$. The boundary $B$ is assembled from singular, Routh-Hurwitz, and coordinate boundaries, each computed as an elimination ideal; this mechanism turns the stability question into a finite graph computation.

What would settle it

Recompute the elimination ideal for the coral model with an independent symbolic computation and compare the 31 irreducible factors; then at representative parameter values inside each computed region, numerically integrate the ODEs and check that the predicted stable steady states are present and that in the white region no steady state is stable and orbits do not converge to equilibrium.

Watch

Extended reading notes

Core claim

The central claim is that routing functions—rational functions built from the product of the boundary polynomials—provide a computational route to the exact connected components of parameter space on which steady-state feasibility and stability are constant. The paper relies on a theorem stating that the graph of routing points and gradient connections is in bijection with the connected components of the complement of the total boundary, and it instantiates the method for the classical two-species competition-colonization model and for a three-species coral-bacteria model. In the coral model, after fixing death rates and one colonization rate, the total boundary is a principal ideal with 31 irreducible factors, and varying the remaining two parameters yields stability landscapes labeled by the stable equilibria, including bistable regions and a white region where no equilibrium is stable, so the dynamics there must be non-stationary.

Load-bearing premise

The load-bearing premise is that the computed 31-factor polynomial really is the complete total boundary for the coral model; if any boundary piece is missing or extraneous, every connected region and its stability label, including the claimed no-steady-state region, could change.

Editorial extensions

If this is right

The same boundary construction can be applied to any polynomial or rational ODE model, giving a complete stability-landscape map rather than a sampled approximation.
For the coral-bacteria model with the chosen fixed parameters, crossing any computed boundary changes the set of stable states, so the maps identify thresholds for pathogen exclusion versus coexistence.
The white region with no stable steady state implies sustained oscillations or more complex attractors, so the method can detect limit cycles without integrating trajectories.
Bistable regions mean that state perturbations, not just parameter changes, can switch the system between coral-only and coexistence, which is directly relevant to whether a coral population recovers after a shock.

Reading between the lines

Editorial extensions of the paper, not claims the author makes directly.

I would expect the same routing-function decomposition to work when initial conditions are treated as parameters, which would turn questions of basin size and separatrix geometry into the same graph computation; the paper notes this possibility but does not develop it.
Because the no-steady-state region is identified by algebraic elimination rather than by simulation, the method offers a certificate-style indication of oscillations for parameter sets, something bifurcation continuation typically only suggests.
A natural stress test is to perturb the fixed rates ($d=\gamma_y=\gamma_z=1$, $\beta_y=5$); the qualitative landscape, including the position of the limit-cycle window, may be sensitive to those choices, so the biological conclusions are conditional on the chosen slice.

Share X Bluesky LinkedIn Reddit HN

Signed reviews

No signed human review yet.

Editorial analysis

A structured set of objections, weighed in public.

Desk editor's note, referee report, and a circularity audit.

Desk Editor's Note

read the letter

Rough take: this is a faithful application of routing functions to ecology, and the coral-bacteria stability maps are new. The exposition of the three boundary types is clear, and the Levins-Culver reproduction is a good sanity check. But the abstract overclaims limit cycles, and the coral results rest on a boundary polynomial g that is never shown or certified, so the maps should be read as conditional on that computation being complete. What is actually new: the systematic package of singular, Routh-Hurwitz, and coordinate boundaries applied to a three-species ecological model, plus the first stability landscape decomposition for the tripartite symbiosis model. The 16 slices, including the bistability regions and the white no-stable-state region, are not in the cited literature. That is real output. Where it gets soft: Section 4.1.1 says a region with no stable steady states must have a limit cycle. In three dimensions that does not follow without ruling out chaotic or other non-periodic attractors. The conclusion uses the safer phrase "limit cycles or more complex attractors," so the abstract's "regions supporting limit cycles" is not what the computation established. Fix the wording. The bigger issue is the unshown total boundary. The routing function graph is only as good as the hypersurface arrangement. If any of the 31 irreducible factors is missing or spurious, connected components merge or split and every region label can change. The paper gives degree counts but does not display g, list the elimination ideals for the five primary components, or provide a certificate. The GitHub repository is mentioned but no URL or commit hash is given. So a referee cannot verify the load-bearing step. This is fixable by making the code and data available, but until then the coral landscape is conditional, not final. The citation pattern is fine: [7] is the source of routing functions, published with independently available software, and the authors are explicit about the lineage. The self-citation is appropriate. No invented entities or fitted constants. Who this is for: mathematical biologists who want to map stability regions in polynomial or rational ODE models, and applied algebraic geometers looking for new application domains. It is not a theoretical advance beyond [7], but it is a well-motivated application with useful results. Recommendation: this deserves a serious referee. Conditional accept, with the expectation that the authors supply the missing boundary data and correct the limit-cycle language.

Referee Report

3 major / 5 minor

Summary. This paper applies the routing-function method of [7] to ecological ODE models, decomposing the positive parameter space minus a boundary hypersurface arrangement into connected components on which the number and type of stable steady states is constant. The boundaries considered are the singular, Routh-Hurwitz, and coordinate boundaries, each obtained via elimination ideals. The method is first demonstrated on the Levins-Culver competition-colonization model, reproducing the known coexistence-stability condition. It is then applied to a coral-bacteria model with parameters fixed to d = gamma_y = gamma_z = 1 and beta_y = 5; after computing a claimed total boundary with 31 irreducible factors, the authors slice at b = 2 and b = 0.5 and at eight values of beta_z, labeling each two-dimensional component by the stable steady states found at a routing point. The abstract claims that the method reveals regions supporting limit cycles, while Section 4.1.1 states that a white region with no stable steady states must have a limit cycle; the conclusion more cautiously says 'limit cycles or more complex attractors.'

Significance. The proposed pipeline is a worthwhile use of recent computational real algebraic geometry in mathematical biology: it gives a concrete, algorithmic way to obtain stability landscapes, and the Levins-Culver section is a clean sanity check in which the known condition beta_z > beta_y (beta_y + gamma_z - gamma_y)/gamma_y and beta_y > gamma_y is recovered without case-by-case analysis. The paper also provides a useful pedagogical summary of the three boundary types and of routing functions, and it names concrete software (HypersurfaceRegions.jl, HomotopyContinuation.jl, Macaulay2) that makes the computations reproducible in principle. However, the coral-landscape results are conditional on a boundary computation that is not verifiable from the manuscript, and one headline inference (no stable steady state implies a limit cycle) is not justified for a three-dimensional system. With added certificates and corrected language, the contribution would be useful to the q-bio audience.

major comments (3)

[Section 4, total boundary computation] The coral analysis rests entirely on the claim that eliminating x, y, z from the equilibrium ideal plus the Routh-Hurwitz conditions yields a principal ideal <g> with 31 irreducible factors, whose degrees are listed but whose factors are not shown. This is load-bearing for Figures 3 and 4 and for Theorem 2.1: if g omits any hypersurface across which the number or type of stable steady states changes, connected components of the complement can merge regions with different labels; a spurious factor would create artificial boundaries. The manuscript does not display g, does not list the five primary-component elimination ideals, and gives no certificate (e.g., a Macaulay2 script and transcript, or a repository URL with a commit hash). Please supply the elimination data in a supplement or make the repository identifiable, or explicitly mark the coral landscape results as conditional on this computation.
[Section 4.1.1 and Abstract] The white region in Figure 3(e) is described as a region where 'none of the steady states are stable and the model must have a limit cycle.' In a three-dimensional autonomous system, absence of stable steady states does not imply existence of a limit cycle; chaotic attractors, heteroclinic cycles, or unbounded dynamics are not excluded by the routing-function computation, which, even assuming the boundary is exact, only certifies that no stable steady state exists in that component and does not identify the attractor. The conclusion's phrase 'limit cycles or more complex attractors' is the defensible version. Please revise the abstract and Section 4.1.1 accordingly, or supply an additional argument (e.g., boundedness plus a Poincare-Bendixson-type analysis) that actually proves a periodic orbit.
[Section 4, slice decomposition] The manuscript does not state how the 31-factor boundary arrangement is restricted to each fixed (b, beta_z) slice before the two-dimensional routing function is constructed, and it does not list the active factors for the 16 slices. This information is needed to reproduce Figures 3 and 4 and to check that slice-dependent factors (for example, factors that become constant on a slice) are handled correctly. Please describe the restriction procedure and provide the restricted boundary data for each slice, either in the text or in a supplement.

minor comments (5)

[Section 2.1] The symbol beta_y is overloaded: it denotes a Routh polynomial in the stability criterion and also a model parameter later in the paper. Please rename the Routh polynomials to avoid confusion.
[Section 3] The primary decomposition notation 'If = T4i=1 J(i)f' is ambiguous; it should be an intersection, e.g., If = intersection_{i=1}^4 J_f^(i).
[Various] There are several typos and infelicities: 'paramter landscape' in Section 2, 'for for' before 'the Levins-Culver model' in Section 4, and 'the large the number of parameters' in Section 4. These should be corrected.
[Section 4] The statement that the computations 'can be found in a GitHub repository' should include a URL, a persistent identifier, and a commit hash or version tag so that the reported results are actually reproducible.
[Section 4.2] The phrase 'as beta_z crosses four' and the description of an 'abrupt' landscape change are based on two sampled slices on either side of beta_z = 4; please phrase these observations as being observed at the sampled slices rather than as a proven characterization for all intermediate values.

Circularity Check

0 steps flagged · score 0.0 of 10

No circularity: the routing-function theorem is imported from published prior work, and every landscape label is computed from elimination ideals derived from the model equations, not fitted to the conclusion.

full rationale

The paper's derivation chain is self-contained relative to the model equations. The stability boundaries are constructed by elimination from the equilibrium ideal together with the Jacobian determinant, Routh-Hurwitz polynomials, and coordinate conditions, as in Section 2.1, with no parameter fitted to data. The routing-function graph method is imported from the published reference [7], whose theorem statements are parameter-free and do not presuppose the present models or results; although two of the present authors are co-authors of [7], the citation is real, independently available evidence and not a circular justification. The Levins-Culver section reproduces a known stability condition, serving as an external benchmark. In the coral case study, each connected component is labeled by computing the steady states at one routing point per component, which is a valid sampling procedure if the boundary computation is complete; the main concern is that the 31-factor boundary polynomial g is not displayed and no certificate or repository URL is given, a reproducibility and verification issue rather than a circular one. The abstract's claim of 'regions supporting limit cycles' overstates what the computation shows, since the conclusion correctly says 'limit cycles or more complex attractors', but that is an inference-strength concern, not a case of the result being equivalent to its inputs by construction. No step in the paper reduces to its own inputs by definition, and no fitted parameter is renamed as a prediction.

Assumptions & free parameters 4 free parameters · 6 assumptions · 0 invented entities

The paper introduces no new physical entities. It relies on several imported mathematical tools and on two unverified computational assumptions: that the 31-factor boundary polynomial is correct and complete, and that numerical routing computations are accurate. The no-stable-state-to-limit-cycle inference is an additional unsupported premise in a 3D system.

free parameters (4)

d = gamma_y = gamma_z = 1 = 1
Set equal to one in Section 4 to reduce the parameter count before boundary computation; all coral landscape maps depend on this choice.
beta_y = 5 = 5
Parasite colonization rate fixed in Section 4 to make Groebner basis computations finish; the stability maps are conditional on this value.
b = 2 or 0.5 = 2, 0.5
Intrinsic coral birth rate chosen to represent viable and non-viable coral regimes (b above and below mortality d=1).
beta_z slice values = 2.5, 3, 3.5, 3.9, 4.1, 5.9, 6.1, 10
Mutualist colonization rate slices; values 3.9, 4.1, 5.9, 6.1 chosen because beta_z-4 and beta_z-6 are boundary factors when b=2. Results are slice-specific, not a full 4D landscape.

assumptions (6)

domain assumption Generic finiteness: for a generic parameter value a*, f(x;a*)=0 has finitely many nonsingular complex solutions, and for all a the system has exactly d solutions counted with multiplicity.
Stated at the start of Section 2.1; without it the boundary description by singular, Routh-Hurwitz, and coordinate hypersurfaces is incomplete.
standard math Routh-Hurwitz conditions characterize local asymptotic stability.
Used throughout Section 2.1 and Section 4 to define the stability boundary; standard result.
standard math Routing function graph theorem (Theorem 2.1 from [7]) correctly computes connected components of the complement of a hypersurface arrangement.
The paper relies on [7, Thm. 4.4] and Algorithm 1; this is imported background, not proved here.
ad hoc to paper Numerical computations (HomotopyContinuation.jl, HypersurfaceRegions.jl) return correct routing points and connections.
No certificates or exact verification are provided for the 16 slices; errors would change the landscape.
ad hoc to paper No-stable-steady-state implies limit cycle.
Used in Section 4.1.1 to assert the white region 'must have a limit cycle.' The system is three-dimensional, where a flow with no stable steady states can have chaotic or other non-equilibrium attractors; the conclusion later hedges to 'limit cycles or more complex attractors.'
ad hoc to paper Boundary polynomial g is complete.
Section 4 states the elimination ideal is principal with 31 irreducible factors, but the polynomial is not shown; completeness of g is assumed for the component decomposition.

how reviews work

0 comments

Cite this review

Pith. "Pith review of Routing functions for parameter space decomposition to describe stability landscapes of ecological models." pith.science (2026). https://pith.science/paper/SZTAC4G7

@misc{pith2026250500128,
  author       = {Pith},
  title        = {Pith review of: Routing functions for parameter space decomposition to describe stability landscapes of ecological models},
  year         = {2026},
  howpublished = {\url{https://pith.science/paper/SZTAC4G7}},
  note         = {Machine review of arXiv:2505.00128}
}

read the original abstract

Changes in environmental or system parameters often drive major biological transitions, including ecosystem collapse, disease outbreaks, and tumor development. Analyzing the stability of steady states in dynamical systems provides critical insight into these transitions. This paper introduces an algebraic framework for analyzing the stability landscapes of ecological models defined by systems of first-order autonomous ordinary differential equations with polynomial or rational rate functions. Using tools from real algebraic geometry, we characterize parameter regions associated with steady-state feasibility and stability via three key boundaries: singular, stability (Routh-Hurwitz), and coordinate boundaries. With these boundaries in mind, we employ routing functions to compute the connected components of parameter space in which the number and type of stable steady states remain constant, revealing the stability landscape of these ecological models. As case studies, we revisit the classical Levins-Culver competition-colonization model and a recent model of coral-bacteria symbioses. In the latter, our method uncovers complex stability regimes, including regions supporting limit cycles, that are inaccessible via traditional techniques. These results demonstrate the potential of our approach to inform ecological theory and intervention strategies in systems with nonlinear interactions and multiple stable states.

Figures

Figures reproduced from arXiv: 2505.00128 by the authors.

**Figure 1.** Computing the 6 regions in R 2 >0 \B where the boundaries are the two ellipses (thick blue) and coordinate axes. Routing points (index 0 in blue, index 1 in red) and gradient paths (thin blue) are also plotted demonstrating 6 connected components. Given a point a ∈ R k >0 \ B which is not a routing point, the output of Algorithm 1 can be used to determine which connected component of R k >0 \ B contains a. For this,… view at source ↗

**Figure 2.** The coral-bacteria symbioses compartmental model represented by system ( [PITH_FULL_IMAGE:figures/full_fig_p013_2.png] view at source ↗

**Figure 3.** Stability landscapes for increasing values of [PITH_FULL_IMAGE:figures/full_fig_p016_3.png] view at source ↗

Figures from the paper (1 more)

**Figure 4.** Figure 4: Stability landscapes for increasing values of [PITH_FULL_IMAGE:figures/full_fig_p016_4.png]

Discussion (0). Continue with ORCID to comment.

Forward citations

Cited by 1 Pith paper

Reviewed papers in the Pith corpus that reference this work. Sorted by Pith novelty score. Full citation record

Strata of Ecological Coexistence via Grassmannians
math.AG 2025-08 conditional novelty 6.0 of 10

A Grassmannian and oriented matroid sign-completion algorithm proves three n=4 Lotka-Volterra interaction networks admit no feasible-stable coexistence.

Reference graph

Works this paper leans on

19 extracted references · 18 canonical work pages · cited by 1 Pith paper

[7]

Cummings, J

J. Cummings, J. D. Hauenstein, H. Hong, and C. D. Smyth. Smooth connectivity in real algebraic varieties. Numerical Algorithms, Oct 2024

work page 2024
[1]

Amarasekare

P. Amarasekare. Competitive coexistence in spatially structured environments: a syn- thesis. Ecology Letters, 6(12):1109–1122, Dec. 2003

work page 2003
[2]

E. A. Bernal, J. D. Hauenstein, D. Mehta, M. H. Regan, and T. Tang. Machine learning the real discriminant locus. Journal of Symbolic Computation , 115:409–426, 2023. 18

work page 2023
[3]

Breiding, B

P. Breiding, B. Sturmfels, and K. Wang. Computing arrangements of hypersurfaces, 2024

work page 2024
[4]

Breiding and S

P. Breiding and S. Timme. HomotopyContinuation.jl: A package for homotopy con- tinuation in Julia. In J. H. Davenport, M. Kauers, G. Labahn, and J. Urban, editors, Mathematical Software – ICMS 2018 , pages 458–465, Cham, 2018. Springer Interna- tional Publishing

work page 2018
[5]

Conradi, E

C. Conradi, E. Feliu, M. Mincheva, and C. Wiuf. Identifying parameter regions for multistationarity. PLoS Computational Biology, 13(10):e1005751, 2017

work page 2017
[6]

D. A. Cox, J. Little, and D. O’Shea. Ideals, varieties, and algorithms . Undergradu- ate Texts in Mathematics. Springer, Cham, fourth edition, 2015. An introduction to computational algebraic geometry and commutative algebra

work page 2015
[8]

Dorf and R

R. Dorf and R. Bishop. Modern Control Systems. Pearson, 2011

work page 2011

Show all 19 references

[9]

Feliu and C

E. Feliu and C. Wiuf. Simplifying biochemical models with intermediate species.Journal of the Royal Society Interface , 10(87):20130484, 2013

2013
[10]

I. M. Gelfand, M. M. Kapranov, and A. V. Zelevinsky. Discriminants, resultants and multidimensional determinants. Modern Birkh¨ auser Classics. Birkh¨ auser Boston, Inc., Boston, MA, 2008. Reprint of the 1994 edition

2008
[11]

T. L. Gibbs, K. Dahlin, J. Brennan, C. B. Silveira, and L. C. McManus. Coexistence of bacteria with a competition-colonization tradeoff on a dynamic coral host. bioRxiv, page 2024.09.15.612558, Sept. 2024

2024
[12]

D. R. Grayson and M. E. Stillman. Macaulay2, Version 1.20, 2022. http://www.math.uiuc.edu/Macaulay2/

2022
[13]

H. Hong. Connectivity in semi-algebraic sets. In 2010 12th International Symposium on Symbolic and Numeric Algorithms for Scientific Computing , pages 4–7, 2010

2010
[14]

H. Hong, J. Rohal, M. Safey El Din, and E. Schost. Connectivity in semi-algebraic sets I. arXiv:2011.02162, 2020

2011 arXiv
[15]

A. Hurwitz. Ueber die Bedingungen, unter welchen eine Gleichung nur Wurzeln mit negativen reellen Theilen besitzt. Math. Ann., 46(2):273–284, 1895
[16]

Joshi and A

B. Joshi and A. Shiu. A survey of methods for deciding whether a reaction network is multistationary. Mathematical Modelling of Natural Phenomena , 10(5):47–67, 2015

2015
[17]

Lazard and F

D. Lazard and F. Rouillier. Solving parametric polynomial systems. Journal of Symbolic Computation, 42(6):636–667, 2007. 19

2007
[18]

Levins and D

R. Levins and D. Culver. Regional coexistence of species and competition between rare species. Proceedings of the National Academy of Sciences of the United States of America, 68(6):1246–1248, June 1971

1971
[19]

D. W. Yu and H. B. Wilson. The competition-colonization trade-off is dead; long live the competition-colonization trade-off. The American Naturalist , 158(1):49–63, July 2001. 20

2001

Pith tools

Reviewed August 16, 2026 · model on record in the stance chip above.

[1] [7]

Cummings, J

J. Cummings, J. D. Hauenstein, H. Hong, and C. D. Smyth. Smooth connectivity in real algebraic varieties. Numerical Algorithms, Oct 2024

work page 2024

[2] [1]

Amarasekare

P. Amarasekare. Competitive coexistence in spatially structured environments: a syn- thesis. Ecology Letters, 6(12):1109–1122, Dec. 2003

work page 2003

[3] [2]

E. A. Bernal, J. D. Hauenstein, D. Mehta, M. H. Regan, and T. Tang. Machine learning the real discriminant locus. Journal of Symbolic Computation , 115:409–426, 2023. 18

work page 2023

[4] [3]

Breiding, B

P. Breiding, B. Sturmfels, and K. Wang. Computing arrangements of hypersurfaces, 2024

work page 2024

[5] [4]

Breiding and S

P. Breiding and S. Timme. HomotopyContinuation.jl: A package for homotopy con- tinuation in Julia. In J. H. Davenport, M. Kauers, G. Labahn, and J. Urban, editors, Mathematical Software – ICMS 2018 , pages 458–465, Cham, 2018. Springer Interna- tional Publishing

work page 2018

[6] [5]

Conradi, E

C. Conradi, E. Feliu, M. Mincheva, and C. Wiuf. Identifying parameter regions for multistationarity. PLoS Computational Biology, 13(10):e1005751, 2017

work page 2017

[7] [6]

D. A. Cox, J. Little, and D. O’Shea. Ideals, varieties, and algorithms . Undergradu- ate Texts in Mathematics. Springer, Cham, fourth edition, 2015. An introduction to computational algebraic geometry and commutative algebra

work page 2015

[8] [8]

Dorf and R

R. Dorf and R. Bishop. Modern Control Systems. Pearson, 2011

work page 2011