We introduce a hybrid mathematical framework for morphogenesis and regeneration motivated by bioelectric phenomena documented in highly regenerative organisms, particularly planaria. The model couples four dynamical layers on a discrete cellular network: (i) a bistable bioelectric layer, in which each cell admits stable hyperpolarized and depolarized equilibria connected through gap-junction currents; (ii) a synthetic intracellular gene regulatory network (GRN) with proliferation, differentiation, positional-identity, and regenerative-response modules; (iii) adaptive gap-junction conductances that evolve in response to electrical state, regenerative activity, and tissue identity; and (iv) a slow tissue-memory variable representing persistent cellular commitment at an epigenetic timescale. Damage is represented by a propagating wound signal on the cellular graph. The central conceptual departure from classical models is that target morphology is not prescribed externally but emerges as an attractor of the coupled multiscale dynamics, in the spirit of distributed attractor-based memory. The framework is designed to capture anatomical homeostasis, regeneration after lesion, attractor switching induced by transient electrical perturbations, regenerative thresholds, and axial polarity. The paper establishes three analytical results for reduced subsystems: single-cell bistability, absence of a Turing instability in the reduced bioelectric–regulatory subsystem, and Lyapunov descent for the pure bioelectric layer. These are complemented by a set of open mathematical questions and numerical experiments that investigate pattern nucleation, regeneration robustness, polarity reversal, and adaptive energy-landscape reshaping in the full multiscale model.
Cooperation is fundamental to biological and social systems, yet its evolutionary success depends critically on population structure. In reality, individuals often exhibit inertia, retaining current behaviors even when imitation is beneficial. The synergistic effects between heterogeneous inertia and network structure have been explored primarily through numerical simulations. In this work, we develop a mathematical framework incorporating degree-dependent inertia, wherein the tendency to maintain one's strategy is determined by node degree, and self-loops are permitted, and self-loops are permitted. Using a coalescent-theoretic approach, we derive the threshold benefit-to-cost ratio favoring cooperation in the donation game on arbitrary graphs in the weak-selection limit. Compared with the no-inertia, uniform-inertia, and inverse degree-dependent inertia cases, positive degree-dependent inertia markedly reduces this critical threshold in disassortative networks. An analytical expression for multi-star networks further reveals that, architectures with fewer hubs and more leaves promote cooperation most effectively. While prior studies emphasize that slower updating by high-degree nodes and faster updating by low-degree nodes can foster cooperation, we propose a refined perspective: cooperation is especially favored when high-degree nodes update slowly and their neighbors are predominantly low-degree individuals. This alignment of inertia and neighborhood composition provides a mechanistic explanation for the emergence of cooperation in structured populations.
Any mass action network gives rise to a parameterised family of polynomial equations whose positive solutions are the positive equilibria of the network. Here, we consider alternative systems of equations, whose solutions are in smooth, one-to-one correspondence with positive equilibria of the network, and capture degeneracy or nondegeneracy of the corresponding equilibria. The construction leads us to consider partitions of networks in a natural sense, and we explore the implications of choosing different partitions. The alternative systems are in some situations simpler than the original mass action equations, which allows us to rapidly identify various algebraic and geometric properties of the positive equilibrium set. This includes the characterisation of toricity and local toricity, bounds on the number of positive nondegenerate equilibria on stoichiometric classes, semialgebraic descriptions of the parameter regions for multistationarity, and the study of bifurcations. After discussing the construction of the alternative systems, various consequences for particular classes of networks and numerous examples are presented. We also develop additional techniques specifically for quadratic networks, the most common class of networks in applications, and use these techniques to derive strengthened results for quadratic networks.
Tuberculosis (TB) is a highly contagious chronic infectious disease that, without timely intervention, can lead to severe health consequences or even death. Improper or incomplete treatment often induces the emergence of drug-resistant TB (DR-TB), which greatly complicates disease management and intensifies its public health burden. This paper develops and analyzes a two-strain TB transmission model incorporating age structure during latency, spatial diffusion, and a treatment-induced resistance pathway. Methodologically, we establish the global existence and non-negativity of model solutions and derive explicit expressions for the basic reproduction numbers of the sensitive and resistant strains, as well as for the reproduction number associated with treatment-induced resistance. This study then extends the persistence proof method for single-strain space-age structured models to examine the dynamics of competitive exclusion and persistence between the two strains. Subsequently, we characterize the local and global stability of equilibria using spectral analysis and Lyapunov function methods. Calibrating the model with WHO data for China, we estimate that the basic reproduction number for the sensitive strain exceeds one, and owing to the presence of a treatment-induced resistance pathway, the basic reproduction number for the resistant strain displays two distinct distributions, both of which remain below one. Despite this, theoretical and numerical results demonstrate that DR-TB can persist even when its basic reproduction number is less than one or even zero. Furthermore, our projections indicate that, given the current level of TB control, China is unlikely to achieve the WHO's 2035 incidence reduction target. Despite this, significant improvement in treatment efficacy and reduction of resistance induction risk could make the goal attainable. Moreover, under comparable conditions, the elimination target appears relatively easier to achieve for DR-TB. Our findings suggest that in the absence of treatment-induced resistance, the WHO's DR-TB elimination goal could be reached approximately 2 years earlier. Notably, early increases in DR-TB cases due to improved treatment should be anticipated, underscoring that TB control efforts must not only target existing DR-TB cases but also ensure standardized treatment for drug-sensitive TB (DS-TB) infections; otherwise, treatment-induced resistance in patients will further increase the TB burden.
The relationship between two pedigree members may be summarized by their nine condensed identity coefficients (Δ _1,… ,Δ _9) , also known as Jacquard coefficients. These give the expected relative frequencies of the possible patterns of identity by descent between the alleles carried by the individuals at an autosomal locus. The collection J of all such vectors, taken over all pairs in all pedigrees, forms a subset of the probability simplex S^8 ⊂ℝ^9 . Despite decades of constant use and recurring interest in the Jacquard coefficients, remarkably little is known about the geometry of J. For instance, a long-standing open question is whether J is dense in S^8 , or, conversely, whether there exists full-dimensional unattainable regions not realizable by any pedigree. In the case of non-inbred individuals, this was famously resolved by Thompson, who identified the unattainable region bounded by a quadratic equation in the two-dimensional simplex of such relationships. In this paper we develop a general framework for constructing restrictions in the identity coefficients of any number of individuals. Applied to the pairwise case, we find several explicit quadratic inequalities in the nine Jacquard coefficients, each defining a full-dimensional unattainable region. One of these contains the known region unattainable by non-inbred relationships. To our knowledge, our results provide the first examples of full-dimensional holes in the space of pairwise relationships, arising purely from Mendelian constraints.
The complex nonlinear dynamics of tumor–immune interactions drive tumor heterogeneity and complicate treatment strategies. While bifurcation and multistability analyses of mathematical models can uncover critical system dynamics, few studies have explored bifurcation–particularly cusp bifurcation–in high-dimensional tumor–immune systems, as most existing analyses are limited to simplified two-dimensional models. In this study, we extend the tumor–immune model proposed by Anderson et al. (2024a) to incorporate combination therapy with immune checkpoint inhibitors (ICIs) and C–C chemokine receptor type 2 (CCR2) antagonists. By combining linear algebra methods, Sotomayor’s theorem, projection singularity analysis, numerical simulation, and continuation techniques, we explicitly identify transcritical and saddle-node bifurcations, conduct a systematic cusp bifurcation analysis as a classical two-parameter problem, and derive explicit quantitative conditions. Calibration against four independent murine datasets demonstrates the model’s ability to reproduce tumor growth dynamics across multiple tumor datasets. Bifurcation analysis characterizes the multiplicity of tumorous equilibria, and reveals the emergence of monostable and bistable regimes. In particular, the saddle-node bifurcation curve partitions the two-parameter space under combination therapy into monostable and bistable regions. Notably, a dose-dependent inverse correlation between ICIs and CCR2 antagonists emerges along this bifurcation boundary. Our results demonstrate that effective tumor control depends not only on treatment intensity but also on the patient’s initial tumor burden. The corresponding critical thresholds are determined by bifurcation structure and the stable manifold (or characteristic space) of the saddle point, respectively. These findings provide theoretical guidance for treatment design under different tumor states. Importantly, the identification of a stable low-tumor state–being clinically acceptable and controllable–highlights the potential for sustained tumor control as an alternative therapeutic objective to complete tumor eradication.
Cells normally combine glycolysis and oxidative phosphorylation (OXPHOS) to meet energy demands, but this balance shifts under pathological conditions. During SARS-CoV-2 infection, hypoxia, viral entry, and elevated tissue lactate alter cellular metabolism. To explore these effects, we propose a parsimonious mathematical model describing how oxygen levels, viral infiltration, and extracellular lactate jointly regulate metabolic balance through HIF-1 α protein, inside the cell, accounting for lactate’s biphasic, non-monotonic influence on glycolysis. Model simulations reveal a single steady state whose position on the glycolysis–OXPHOS phase plane depends on environmental conditions, namely, oxygen concentration, infection, and extracellular lactate. We identify four metabolic regimes, determined by sufficiency of energy production and the driving process (OXPHOS or glycolysis). Decreasing the oxygen shifts cells from OXPHOS to glycolysis dominance in both infected and non-infected states, but infected cells may become energy-deficient even with sufficient oxygen due to virus-induced mitochondrial damage. Rising extracellular lactate initially promotes glycolysis but ultimately suppresses it at high levels, pushing cells into severe energy deficit with inhibited glycolysis. Simulations of reoxygenation exhibit hysteresis: cells pass through an energy-deficient zone during hypoxia onset but return through a safer trajectory when oxygen is restored; a vulnerability is higher in infected cells. Overall, the model clarifies metabolic trajectories during viral infection, suggesting that early hypoxia is particularly dangerous and that severe acidosis can further collapse energy production. Preventing or rapidly reversing hypoxia in respiratory infection may protect cells from energy failure and limit harmful lactate accumulation.
We study mechanisms and effects of ideal free distributions (IFDs) on an ecological community consisting of n prey and m predator species, for any positive integers n and m, by considering the corresponding diffusive Lotka-Volterra system with time-periodic coefficients. We define a notion of joint IFD in a time-periodic environment, and give necessary and sufficient conditions for it to be achieved by a subcollection of prey and predator species with suitable dispersal strategies. Next, we show, via construction of a Lyapunov functional, that such dispersal strategies are evolutionarily stable, in the sense that if a subcollection of prey and predator species adopts an ideal free dispersal strategy, then the total community must converge to an IFD for large time; if a unique combination of prey-predator species adopts an ideal free strategy, then it can drive all other species to extinction. Conversely, if a combination of prey-predator species adopts a non-ideal free dispersal strategy, then it can be invaded by some suitable mutant strategies. Our results provide insight into the evolution of spatial distribution of ecological communities with predator-prey interactions.
Malaria is the most widespread and deadly parasitic disease in the world. This disease is transmitted from an infected individual to a healthy individual via the female Anopheles mosquito. Children under the age of 5 are most at risk from this infectious disease. Sub-Saharan Africa is most affected by malaria in the world. The geographical distribution and seasonality of the malaria cases observed are closely related to the climatic factors (temperature and precipitation) that influence its transmission dynamics. To control this disease, the use of insecticide treated nets and education of the population for a clean environment without mosquito breeding grounds has been advocated. To perform a mathematical analysis of malaria dynamics, a mathematical model describing malaria transmission must take into account intra-annual climate variation, compliance with measures to reduce contact with mosquitoes, and the age structure of the population. Some mathematical models developed are mainly concerned with the influence of climatic factors on transmission dynamics. Other models are being developed to assess the impact of control factors such as the use of insecticide-treated nets on the behaviour of changes in the number of malaria cases. For a more appropriate analysis, in this paper, we simultaneously assess the impact of climate, the proportion of people who comply with malaria control measures, and the age structure of the population. To take into account the temporal variation of the defined parameters as a function of temperature and precipitation, we describe the dynamics of malaria transmission using a non-autonomous system of ordinary differential equations. The mathematical analysis of the model allowed us to establish the biological properties and asymptotic behavior of the model solutions. From numerical simulations, we assessed the threshold proportion of people who comply with the rules that is necessary to control malaria dynamics, which is 0.5. Using two cities in Cameroon, Maroua and Ngaoundere with different climatic data, we highlighted the influence of climate on the evolution of the number of new malaria cases.
The LPA model is a discrete-time map that has been well-studied and well-validated with experimental data using Tribolium castaneum. We argue that the long lifespan of adults warrants the use of a continuous-time model, and thus we propose a two-dimensional system of delay differential equations to describe Tribolium dynamics. We show that the time delay allows for the periodic behavior that characterizes flour beetle populations. We illustrate the global stability of the extinction equilibrium for R 0 < 1 . Because of nonlinearities in the model, we analyze stability of the positive equilibrium in the special case when there is no adult cannibalism of pupae, exploring the biological parameter space, and finally using a quasi-steady-state argument to reduce the model. We study bifurcations numerically and, unlike previous work, find that larvae development rate and larvae mortality affect the potential for asymptotic cyclic behavior, while other parameters (adult egg production rate per day and adult mortality) influence cycles in the transient phase. This study suggests that model structure and form play an important role in finding chaos.
Oropouche virus (OROV) is a vector-borne arbovirus that has caused more than 500,000 infections across South America, Central America, and the Caribbean. Its primary vector, the midge species Culicoides paraensis, acquires the virus from wild reservoir hosts and transmits it to humans during blood feeding. As outbreaks within the Amazon basin continue to occur and confirmed cases steadily increase, a proper understanding of the interplay between climatic variables, vector behavior, and OROV transmission is needed to develop effective public health interventions. In this paper, we develop a novel mathematical model of Oropouche virus spread in Amazonas, Brazil using ordinary differential equations (ODEs). First, we derive the basic reproductive number at baseline to assess the potential for disease persistence. We then incorporate seasonal forcing and fit the model to weekly incidence data from the 2023-2024 outbreak using least squares optimization. With our model structure, we evaluated multiple candidate seasonality functions and implementations, including precipitation-based forcing functions and others intended to model general seasonality or the right-skewed incidence curve. Among these, quantitative results show that the seasonal variation in OROV transmission is best explained (in the sense of minimum l^2 error) by an exponentially decaying pulse applied to midge recruitment, which captures the incidence data with a relative l^2 error of 11.94
Prevention and control strategies play an essential role in reducing new infections of infectious diseases, while uniform and instantaneous responses to disease outbreaks are assumed in traditional epidemic models. To model the adaptive public health policies that account for hysteresis effects and heterogeneous control measures, we proposed an SIQR epidemic model with multiple two-threshold structures to describe the prevalence-dependent real-time control strategies. The formulated model incorporated the Preisach operator which effectively characterizes time-varying control strategies that depend on the historical evolution of a disease. We established the global stability of this steady-state set through an innovative analytical approach by constructing a series of Lyapunov functions corresponding to different branches of the hysteresis operator. Specifically, our results indicate that: (i) heterogeneity facilitates the convergence of system trajectories to a point within the continuum of steady-state, and (ii) greater heterogeneity leads to higher peaks of infected individuals during the evolution of an epidemic and more severe disease prevalence in the long term. Our findings reveal that the hysteresis phenomenon in control events may lead to existence of a continuum of endemic equilibria, significantly altering the effectiveness of disease suppression strategies. These results provide a theoretical foundation for designing adaptive control policies, emphasizing the importance of incorporating past intervention effects in epidemic modeling.
The COVID-19 pandemic highlighted the ability of epidemics to evolve through the emergence of successive strains of greater infectiousness, which prompted the insight that hyper-exponential growth (HEG) can arise in the development of an epidemic. The phenomenon of HEG has intrigued many researchers, because of some radical differences from exponential growth such as the finite-time singularity. However, the actual mechanism of HEG was usually hidden or needed to be added phenomenologically. Thus, there was little insight into how constraints would be triggered. In this study, we explore an SIR model with evolving parameters leading to a discrete sequence of variants. This allows us to consider the HEG phenomenon in greater depth and with greater mathematical rigour. The model yields a mechanistic description of what happens at the collapse of HEG. The model also yields closed expressions for important features, such as the critical time to the singularity, in terms of basic demographic parameters. Our analysis flags a few important issues needing further research, such as the stochastic character of the emergence of variants. Greater understanding of the HEG process will yield dividends in other fields, since modern societies exhibit HEG at several levels, such as human population growth, economic indicators and technological innovation. For this reason, more research on HEG remains imperative, especially on HEG in the presence of limited resources.
This work investigates the dynamics of solutions to a diffusive two-strain epidemic model with varying total population size. First, assuming a spatially homogeneous environment, we show that the long-term dynamics of the diffusive model mirrors that of the corresponding system of ODE epidemic model. In this setting, we establish that the competitive-exclusion principle holds. However, when the environment is spatially heterogeneous, the global dynamics of the diffusive epidemic model is more challenging. Under some biologically meaningful assumptions, we establish results concerning the existence, uniqueness, and global stability of coexistence endemic equilibrium. Our findings highlight the complex interplay between population movement and spatial heterogeneity in shaping the dynamics of multi-strain infectious diseases.
Seasonally migrating animals must navigate environments where resources shift predictably but are increasingly perturbed by climate change and human activities. Empirical work highlights the importance of cognition for these movements, yet the joint roles of perception and memory in sustaining stable seasonal migration remain poorly understood. We develop and analyze a novel PDE (partial differential equation) model that couples random dispersal with two taxis processes: perception-driven movement toward a nonlocally sensed, periodically shifting resource, and memory-driven movement guided by a spatiotemporal map of past foraging successes over seasonal time windows. We first establish global well-posedness of the system, proving existence, uniqueness, and uniform boundedness of classical solutions. Using Leray-Schauder degree theory, we then show that the model admits at least one time-periodic solution synchronized with the seasonal resource. By constructing a Lyapunov-Krasovskii functional, we further derive sufficient conditions under which this periodic migratory pattern is unique and globally asymptotically stable, revealing a key trade-off between perception- and memory-driven taxis strengths and diffusive spreading. Numerical simulations corroborate the analytical results and demonstrate how the balance of perception and memory, the precision of memory, and its match or mismatch with environmental periodicity jointly govern migration efficiency and persistence. Together, these results provide a rigorous theoretical framework linking individual-level cognitive processes to the emergence, stability, and breakdown of seasonal migration routes in changing environments.
We consider a discrete-time competition model between native and alien predator populations for a common prey. The model is described by a system of three recurrence relations governing their population dynamics, with a generic density-dependent effect for prey and a generic predation factor. In our model, the native and alien predators prey on different stages of the common prey: juvenile-specific and adult-specific predators. In this paper, we focus on the invadability of native prey-predator system to an alien predator and investigate which stage-specific alien predator could be more successful in invading the native system with a different stage-specific native predator. Our mathematical results demonstrate that the prey-predator system with an adult-specific native predator is more vulnerable to the juvenile-specific alien predator invasion compared to the native system with a juvenile-specific native predator. To illustrate the general result more clearly, we present some detailed results on a specific model with a Beverton-Holt type of density-dependent effect and a Nicholson-Bailey type of predation factor, which effectively demonstrates such the dependence of vulnerability to an alien predator invasion on the stage-specific predation.
Interspecific mutualism is a universal phenomenon that exists in biological and social populations. The generation and maintenance of such mutualism have always been a hotspot in evolutionary biology. This study enriches the research on the evolution of mutualism in structured populations from a theoretical perspective. In particular, we develop an evolutionary model of arbitrary interdependent populations using the coalescent random walk theory and employ it to derive the conditions for mutualism under weak selection. We find that evolution favors mutualism when the strength of intraspecific interactions exceeds a threshold value, which depends on the structure of multiple populations and game parameters. We thus reveal the positive role of intraspecific interactions on the evolution of mutualism, which extends prior research that considered the single effect of interspecific interactions.
Motivated by the invasion of Aedes albopictus mosquitoes and the interspecific competition between Aedes aegypti and Aedes albopictus in Florida, we formulate a two-species competition model on a compact metric graph. This model accounts for the species’ ability to inhabit and traverse along the graph’s edges. In the scenario of weak-strong competition, we prove that solutions of the competition model converge uniformly to a semi-positive equilibrium, where either Aedes albopictus survives while Aedes aegypti goes extinct, or vice versa. Whereas for weak-weak competition, the solutions converge uniformly to a positive equilibrium, enabling the coexistence of both Aedes aegypti and Aedes albopictus. We conduct numerical simulations on the two-species competition model along Route 441 in Florida, aiming to illustrate the invasion dynamics and competitive interactions of Aedes aegypti and Aedes albopictus.
We investigate periodic solutions of a genetic negative feedback loop model incorporating protein-sequestration-based repression. Stoichiometry is defined as the average ratio of the concentration of repressors to that of activators over one period of a periodic solution. We examine how stoichiometric balance plays a critical role in the emergence of oscillations. Using Hopf bifurcation analysis combined with numerical simulations, we establish the existence of periodic solutions. We further analyze the stoichiometric conditions associated with oscillation generation in the system. To this end, we approximate the stoichiometry by a quantity evaluated at the bifurcation point and justify that this quantity decreases with increasing activator concentration. Our results precisely characterize the stoichiometric range in which sustained oscillations occur, as inferred from the approximating quantity. The validity of this approximation is supported by the properties of Hopf bifurcation and is illustrated by numerical computations. Finally, we address the effects of differential degradation rates on the existence and stability of periodic solutions, as well as on stoichiometric balance.
We present a new approach to model partial differential equation (PDE) networks. The originality of our method lies in its ability to model various types of networks, including those whose nodes (domains) are non-identical, i.e. they can have various shapes and sizes. This approach makes it possible to model phenomena with possible connections both from within domains and from their boundaries, thus offering a general overview of transmission models. We propose a general model for the scalar case, demonstrating flow conservation in the network as well as the existence, uniqueness and positivity of the solution. In addition, we develop a model for the vector case, where we establish results similar to those of the scalar case, while adding an analysis on the global bounding of the solution in mass-conserving systems. We illustrate our approach with a concrete application to a predator-prey model in a network with two nodes, to which we apply all the obtained theoretical results. Finally, numerical simulations are provided to illustrate our theoretical findings.
