A Comparison of Bimolecular Reaction Models for Stochastic Reaction–Diffusion Systems

Stochastic reaction–diffusion models have become an important tool in studying how both noise in the chemical reaction process and the spatial movement of molecules influences the behavior of biological systems. There are two primary spatially-continuous models that have been used in recent studies: the diffusion limited reaction model of Smoluchowski, and a second approach popularized by Doi. Both models treat molecules as points undergoing Brownian motion. The former represents chemical reactions between two reactants through the use of reactive boundary conditions, with two molecules reacting instantly upon reaching a fixed separation (called the reaction-radius). The Doi model uses reaction potentials, whereby two molecules react with a fixed probability per unit time, λ, when separated by less than the reaction radius. In this work, we study the rigorous relationship between the two models. For the special case of a protein diffusing to a fixed DNA binding site, we prove that the solution to the Doi model converges to the solution of the Smoluchowski model as λ→∞, with a rigorous $O(\lambda^{-\frac{1}{2} + \epsilon})$ error bound (for any fixed ?>0). We investigate by numerical simulation, for biologically relevant parameter values, the difference between the solutions and associated reaction time statistics of the two models. As the reaction-radius is decreased, for sufficiently large but fixed values of λ, these differences are found to increase like the inverse of the binding radius.     

Interdependency and hierarchy of exact and approximate epidemic models on networks

Over the years numerous models of \(SIS\) (susceptible \(\rightarrow \) infected \(\rightarrow \) susceptible) disease dynamics unfolding on networks have been proposed. Here, we discuss the links between many of these models and how they can be viewed as more general motif-based models. We illustrate how the different models can be derived from one another and, where this is not possible, discuss extensions to established models that enables this derivation. We also derive a general result for the exact differential equations for the expected number of an arbitrary motif directly from the Kolmogorov/master equations and conclude with a comparison of the performance of the different closed systems of equations on networks of varying structure.     

When does pathogen evolution maximize the basic reproductive number in well-mixed host–pathogen systems?

Pathogen evolution towards the largest basic reproductive number, $\mathcal R _0$ , has been observed in many theoretical models, but this conclusion does not hold universally. Previous studies of host–pathogen systems have defined general conditions under which $\mathcal R _0$ maximization occurs in terms of $\mathcal R _0$ itself. However, it is unclear what constraints these conditions impose on the functional forms of pathogen related processes (e.g. transmission, recover, or mortality) and how those constraints relate to the characteristics of natural systems. Here we focus on well-mixed SIR-type host–pathogen systems and, via a synthesis of results from the literature, we present a set of sufficient mathematical conditions under which evolution maximizes $\mathcal R _0$ . Our conditions are in terms of the functional responses of the system and yield three general biological constraints on when $\mathcal R _0$ maximization will occur. First, there are no genotype-by-environment interactions. Second, the pathogen utilizes a single transmission pathway (i.e. either horizontal, vertical, or vector transmission). Third, when mortality is density dependent: (i) there is a single infectious class that individuals cannot recover from, (ii) mortality in the infectious class is entirely density dependent, and (iii) the rates of recovery, infection progression, and mortality in the exposed classes are independent of the pathogen trait. We discuss how this approach identifies the biological mechanisms that increase the dimension of the environmental feedback and prevent $\mathcal R _0$ maximization.     

On the identifiability of metabolic network models

A major problem for the identification of metabolic network models is parameter identifiability, that is, the possibility to unambiguously infer the parameter values from the data. Identifiability problems may be due to the structure of the model, in particular implicit dependencies between the parameters, or to limitations in the quantity and quality of the available data. We address the detection and resolution of identifiability problems for a class of pseudo-linear models of metabolism, so-called linlog models. Linlog models have the advantage that parameter estimation reduces to linear or orthogonal regression, which facilitates the analysis of identifiability. We develop precise definitions of structural and practical identifiability, and clarify the fundamental relations between these concepts. In addition, we use singular value decomposition to detect identifiability problems and reduce the model to an identifiable approximation by a principal component analysis approach. The criterion is adapted to real data, which are frequently scarce, incomplete, and noisy. The test of the criterion on a model with simulated data shows that it is capable of correctly identifying the principal components of the data vector. The application to a state-of-the-art dataset on central carbon metabolism in Escherichia coli yields the surprising result that only $4$ out of $31$ reactions, and $37$ out of $100$ parameters, are identifiable. This underlines the practical importance of identifiability analysis and model reduction in the modeling of large-scale metabolic networks. Although our approach has been developed in the context of linlog models, it carries over to other pseudo-linear models, such as generalized mass-action (power-law) models. Moreover, it provides useful hints for the identifiability analysis of more general classes of nonlinear models of metabolism.     

Parameter non-identifiability of the Gyllenberg–Webb ODE model

An ODE model introduced by Gyllenberg and Webb (Growth Develop Aging 53:25–33, 1989) describes tumour growth in terms of the dynamics between proliferating and quiescent cell states. The passage from one state to another and vice versa is modelled by two functions $r_o$ and $r_i$ depending on the total tumour size. As these functions do not represent any observable quantities, they have to be identified from the observations. In this paper we show that there is an infinite number of pairs ( $r_o, r_i$ ) corresponding to the same solution of the ODE system and the functions ( $r_o, r_i$ ) will be classified in terms of this equivalence. Surprisingly, the technique used for this classification permits a uniqueness proof of the solution of the ODE model in a non-Lipschitz case. The reasoning can be widened to a more general setting including an extension of the Gyllenberg–Webb model with a nonlinear birth rate. The relevance of this result is discussed in a preclinical application scenario.     

Modeling Endocrine Control of the Pituitary–Ovarian Axis: Androgenic Influence and Chaotic Dynamics

Mathematical models of the hypothalamus-pituitary-ovarian axis in women were first developed by Schlosser and Selgrade in 1999, with subsequent models of Harris-Clark et al. (Bull. Math. Biol. 65(1):157–173, 2003) and Pasteur and Selgrade (Understanding the dynamics of biological systems: lessons learned from integrative systems biology, Springer, London, pp. 38–58, 2011). These models produce periodic in-silico representation of luteinizing hormone (LH), follicle stimulating hormone (FSH), estradiol (E2), progesterone (P4), inhibin A (InhA), and inhibin B (InhB). Polycystic ovarian syndrome (PCOS), a leading cause of cycle irregularities, is seen as primarily a hyper-androgenic disorder. Therefore, including androgens into the model is necessary to produce simulations relevant to women with PCOS. Because testosterone (T) is the dominant female androgen, we focus our efforts on modeling pituitary feedback and inter-ovarian follicular growth properties as functions of circulating total T levels. Optimized parameters simultaneously simulate LH, FSH, E2, P4, InhA, and InhB levels of Welt et al. (J. Clin. Endocrinol. Metab. 84(1):105–111, 1999) and total T levels of Sinha-Hikim et al. (J. Clin. Endocrinol. Metab. 83(4):1312–1318, 1998). The resulting model is a system of 16 ordinary differential equations, with at least one stable periodic solution. Maciel et al. (J. Clin. Endocrinol. Metab. 89(11):5321–5327, 2004) hypothesized that retarded early follicle growth resulting in “stockpiling” of preantral follicles contributes to PCOS etiology. We present our investigations of this hypothesis and show that varying a follicular growth parameter produces preantral stockpiling and a period-doubling cascade resulting in apparent chaotic menstrual cycle behavior. The new model may allow investigators to study possible interventions returning acyclic patients to regular cycles and guide developments of individualized treatments for PCOS patients.     

Food and space revisited: The role of drift-feeding theory in predicting the distribution,growth, and abundance of stream salmonids

In this paper we review drift-feeding models for stream salmonids. We assess their historical development and current state, and we propose areas for future research. Drift-feeding models serve as the critical input for energetics-based habitat selection and habitat quality models, which have recently begun to see widespread use for predicting salmonid distribution, growth and abundance. We use a bibliometric approach to find drift-feeding model publications, especially those citing three landmark papers that began the quantification of drift feeding by stream fish (Fausch 1984; Hughes and Dill 1990; Hill and Grossman 1993). Subsequent drift-feeding models have largely been built upon these models. Research effort has focused on model development and applications but model testing has been neglected. To date, the only rigorous test of a drift-feeding model (Hughes et al. 2003) identified several limitations and violations of model assumptions. The most important limitation was that prey capture- and gross energy intake rates were overestimated by a factor of two, due largely to poor predictions of prey detection probabilities. Consequences of error in drift-feeding models, and consequently in the habitat selection/quality models that employ them, are greater for applications aimed at predicting growth and abundance than they are for predicting distribution. Research effort on a broad front is needed to advance both drift-feeding models and habitat selection/quality models, including: further development of drift-foraging theory, revision and testing of drift-feeding models (specifically new, functional prey detection and interception sub-models), and revision of habitat selection/quality models to incorporate spatial, temporal, and flow-dependent variation in drift concentration.     

Spiky and transition layer steady states of chemotaxis systems via global bifurcation and Helly’s compactness theorem

The most important phenomenon in chemotaxis is cell aggregation. To model this phenomenon we use spiky or transition layer (step-function-like) steady states. In the case of one spatial dimension, we carry out global bifurcation analysis on the Keller–Segel model and several variants of it, showing that positive steady states exist if the chemotactic coefficient ${\chi}$ is larger than a bifurcation value ${\bar{\chi}_1}$ which can be explicitly expressed in terms of the parameters in the models; then we use Helly’s compactness theorem to obtain the profiles of these steady states when the ratio of the chemotactic coefficient and the cell diffusion rate is large, showing that they are either spiky or have the transition layer structure. Our results provide insights on how the biological parameters affect pattern formation, and reveal the similarities and differences of some popular chemotaxis models.     

Hash: a program to accurately predict protein Hα shifts from neighboring backbone shifts

Chemical shifts provide not only peak identities for analyzing nuclear magnetic resonance (NMR) data, but also an important source of conformational information for studying protein structures. Current structural studies requiring H^α chemical shifts suffer from the following limitations. (1) For large proteins, the H^α chemical shifts can be difficult to assign using conventional NMR triple-resonance experiments, mainly due to the fast transverse relaxation rate of C^α that restricts the signal sensitivity. (2) Previous chemical shift prediction approaches either require homologous models with high sequence similarity or rely heavily on accurate backbone and side-chain structural coordinates. When neither sequence homologues nor structural coordinates are available, we must resort to other information to predict H^α chemical shifts. Predicting accurate H^α chemical shifts using other obtainable information, such as the chemical shifts of nearby backbone atoms (i.e., adjacent atoms in the sequence), can remedy the above dilemmas, and hence advance NMR-based structural studies of proteins. By specifically exploiting the dependencies on chemical shifts of nearby backbone atoms, we propose a novel machine learning algorithm, called Hash, to predict H^α chemical shifts. Hash combines a new fragment-based chemical shift search approach with a non-parametric regression model, called the generalized additive model, to effectively solve the prediction problem. We demonstrate that the chemical shifts of nearby backbone atoms provide a reliable source of information for predicting accurate H^α chemical shifts. Our testing results on different possible combinations of input data indicate that Hash has a wide rage of potential NMR applications in structural and biological studies of proteins.     

Threshold dynamics in an SEIRS model with latency and temporary immunity

A disease transmission model of SEIRS type with distributed delays in latent and temporary immune periods is discussed. With general/particular probability distributions in both of these periods, we address the threshold property of the basic reproduction number \(R_0\) and the dynamical properties of the disease-free/endemic equilibrium points present in the model. More specifically, we 1. show the dependence of \(R_0\) on the probability distribution in the latent period and the independence of \(R_0\) from the distribution of the temporary immunity, 2. prove that the disease free equilibrium is always globally asymptotically stable when \(R_0<1\) , and 3. according to the choice of probability functions in the latent and temporary immune periods, establish that the disease always persists when \(R_0>1\) and an endemic equilibrium exists with different stability properties. In particular, the endemic steady state is at least locally asymptotically stable if the probability distribution in the temporary immunity is a decreasing exponential function when the duration of the latency stage is fixed or exponentially decreasing. It may become oscillatory under certain conditions when there exists a constant delay in the temporary immunity period. Numerical simulations are given to verify the theoretical predictions.     

Approximating the dynamics of communicating cells in a diffusive medium by ODEs—homogenization with localization

Bacteria may change their behavior depending on the population density. Here we study a dynamical model in which cells of radius $R$ within a diffusive medium communicate with each other via diffusion of a signalling substance produced by the cells. The model consists of an initial boundary value problem for a parabolic PDE describing the exterior concentration $u$ of the signalling substance, coupled with $N$ ODEs for the masses $a_i$ of the substance within each cell. We show that for small $R$ the model can be approximated by a hierarchy of models, namely first a system of $N$ coupled delay ODEs, and in a second step by $N$ coupled ODEs. We give some illustrations of the dynamics of the approximate model.     

Modelling the Impact of Marine Reserves on a Population with Depensatory Dynamics

In this study, we use a spatially implicit, stage-structured model to evaluate marine reserve effectiveness for a fish population exhibiting depensatory (strong Allee) effects in its dynamics. We examine the stability and sensitivity of the equilibria of the modelled system with regards to key system parameters and find that for a reasonable set of parameters, populations can be protected from a collapse if a small percentage of the total area is set aside in reserves. Furthermore, the overall abundance of the population is predicted to achieve a maximum at a certain ratio \(A\) of reserve area to fished area, which depends heavily on the other system parameters such as the net export rate of fish from the marine reserves to the fished areas. This finding runs contrary to the contested “equivalence at best” result when comparing fishery management through traditional catch or effort control and management through marine reserves. Lastly, we analyse the problem from a bioeconomics perspective by computing the optimal harvesting policy using Pontryagin’s Maximum Principle, which suggests that the value for \(A\) which maximizes the optimal equilibrium fishery yield also maximizes population abundance when the cost per unit harvest is constant, but can increase substantially when the cost per unit harvest increases with the area being harvested.     

A dynamical systems perspective on the relationship between symbolic and non-symbolic computation

It has been claimed that connectionist (artificial neural network) models of language processing, which do not appear to employ “rules”, are doing something different in kind from classical symbol processing models, which treat “rules” as atoms (e.g., McClelland and Patterson in Trends Cogn Sci 6(11):465–472, 2002). This claim is hard to assess in the absence of careful, formal comparisons between the two approaches. This paper formally investigates the symbol-processing properties of simple dynamical systems called affine dynamical automata, which are close relatives of several recurrent connectionist models of language processing (e.g., Elman in Cogn Sci 14:179–211, 1990). In line with related work (Moore in Theor Comput Sci 201:99–136, 1998; Siegelmann in Neural networks and analog computation: beyond the Turing limit. Birkhäuser, Boston, 1999), the analysis shows that affine dynamical automata exhibit a range of symbol processing behaviors, some of which can be mirrored by various Turing machine devices, and others of which cannot be. On the assumption that the Turing machine framework is a good way to formalize the “computation” part of our understanding of classical symbol processing, this finding supports the view that there is a fundamental “incompatibility” between connectionist and classical models (see Fodor and Pylyshyn 1988; Smolensky in Behav Brain Sci 11(1):1–74, 1988; beim Graben in Mind Matter 2(2):29--51,2004b). Given the empirical successes of connectionist models, the more general, super-Turing framework is a preferable vantage point from which to consider cognitive phenomena. This vantage may give us insight into ill-formed as well as well-formed language behavior and shed light on important structural properties of learning processes.     

Threshold dynamics of an infective disease model with a fixed latent period and non-local infections

In this paper, we derive and analyze an infectious disease model containing a fixed latency and non-local infection caused by the mobility of the latent individuals in a continuous bounded domain. The model is given by a spatially non-local reaction–diffusion system carrying a discrete delay associated with the zero-flux condition on the boundary. By applying some existing abstract results in dynamical systems theory, we prove the existence of a global attractor for the model system. By appealing to the theory of monotone dynamical systems and uniform persistence, we show that the model has the global threshold dynamics which can be described either by the principal eigenvalue of a linear non-local scalar reaction diffusion equation or equivalently by the basic reproduction number ${\mathcal{R}_0}$ for the model. Such threshold dynamics predicts whether the disease will die out or persist. We identify the next generation operator, the spectral radius of which defines basic reproduction number. When all model parameters are constants, we are able to find explicitly the principal eigenvalue and ${\mathcal{R}_0}$ . In addition to computing the spectral radius of the next generation operator, we also discuss an alternative way to compute ${\mathcal{R}_0}$ .     

Stability of cluster solutions in a cooperative consumer chain model

We study a cooperative consumer chain model which consists of one producer and two consumers. It is an extension of the Schnakenberg model suggested in Gierer and Meinhardt [Kybernetik (Berlin), 12:30–39, 1972] and Schnakenberg (J Theor Biol, 81:389–400, 1979) for which there is only one producer and one consumer. In this consumer chain model there is a middle component which plays a hybrid role: it acts both as consumer and as producer. It is assumed that the producer diffuses much faster than the first consumer and the first consumer much faster than the second consumer. The system also serves as a model for a sequence of irreversible autocatalytic reactions in a container which is in contact with a well-stirred reservoir. In the small diffusion limit we construct cluster solutions in an interval which have the following properties: The spatial profile of the third component is a spike. The profile for the middle component is that of two partial spikes connected by a thin transition layer. The first component in leading order is given by a Green’s function. In this profile multiple scales are involved: The spikes for the middle component are on the small scale, the spike for the third on the very small scale, the width of the transition layer for the middle component is between the small and the very small scale. The first component acts on the large scale. To the best of our knowledge, this type of spiky pattern has never before been studied rigorously. It is shown that, if the feedrates are small enough, there exist two such patterns which differ by their amplitudes.We also study the stability properties of these cluster solutions. We use a rigorous analysis to investigate the linearized operator around cluster solutions which is based on nonlocal eigenvalue problems and rigorous asymptotic analysis. The following result is established: If the time-relaxation constants are small enough, one cluster solution is stable and the other one is unstable. The instability arises through large eigenvalues of order $O(1)$ . Further, there are small eigenvalues of order $o(1)$ which do not cause any instabilities. Our approach requires some new ideas: (i) The analysis of the large eigenvalues of order $O(1)$ leads to a novel system of nonlocal eigenvalue problems with inhomogeneous Robin boundary conditions whose stability properties have been investigated rigorously. (ii) The analysis of the small eigenvalues of order $o(1)$ needs a careful study of the interaction of two small length scales and is based on a suitable inner/outer expansion with rigorous error analysis. It is found that the order of these small eigenvalues is given by the smallest diffusion constant ${\epsilon }_2^2$ .     

Individual-Based Model for Quorum Sensing with Background Flow

Quorum sensing is a wide-spread mode of cell–cell communication among bacteria in which cells release a signalling substance at a low rate. The concentration of this substance allows the bacteria to gain information about population size or spatial confinement. We consider a model for \(N\) cells which communicate with each other via a signalling substance in a diffusive medium with a background flow. The model consists of an initial boundary value problem for a parabolic PDE describing the exterior concentration \(u\) of the signalling substance, coupled with \(N\) ODEs for the masses \(a_i\) of the substance within each cell. The cells are balls of radius \(R\) in \(\mathbb {R} ^3\) , and under some scaling assumptions we formally derive an effective system of \(N\) ODEs describing the behaviour of the cells. The reduced system is then used to study the effect of flow on communication in general, and in particular for a number of geometric configurations.     

Differential mechanical response and microstructural organization between non-human primate femoral and carotid arteries

Unique anatomic locations and physiologic functions predispose different arteries to varying mechanical responses and pathologies. However, the underlying causes of these mechanical differences are not well understood. The objective of this study was to first identify structural differences in the arterial matrix that would account for the mechanical differences between healthy femoral and carotid arteries and second to utilize these structural observations to perform a microstructurally motivated constitutive analysis. Femoral and carotid arteries were subjected to cylindrical biaxial loading and their microstructure was quantified using two-photon microscopy. The femoral arteries were found to be less compliant than the carotid arteries at physiologic loads, consistent with previous studies, despite similar extracellular compositions of collagen and elastin ( \(P> 0.05\) ). The femoral arteries exhibited significantly less circumferential dispersion of collagen fibers ( \(P< 0.05\) ), despite a similar mean fiber alignment direction as the carotid arteries. Elastin transmural distribution, in vivo axial stretch, and opening angles were also found to be distinctly different between the arteries. Lastly, we modeled the arteries’ mechanical behaviors using a microstructural-based, distributed collagen fiber constitutive model. With this approach, the material parameters of the model were solved using the experimental microstructural observations. The findings of this study support an important role for microstructural organization in arterial stiffness.