Exact fundamental solutions of linear parabolic equations with spatially varying coefficients

The exact general solution is obtained to a linear second order ordinary differential equation which has quite general coefficients depending on an arbitrary function of the independent variable. From this, the exact fundamental solution is derived for the corresponding linear parabolic partial differential equation with coefficients depending on the single space coordinate. In a special case this latter equation reduces to one of the Fokker-Planck type. These coefficients are then generalised and the appropriate fundamental solution is obtained. Extensions are given to linear parabolic equations in two andn space dimensions. The paper provides a collection of basic examples which illustrate and develop the theory for the generation of the exact fundamental solutions. Reduction to, and the corresponding fundamental solutions of the Fokker-Planck equations is presented, where appropriate.     

Toxicokinetics of Metals in Terrestrial Invertebrates: Making Things Straight with the One-Compartment Principle

In this analysis, we first performed a critical review of one-compartment models used to describe metal toxicokinetics in invertebrates and found mathematical or conceptual errors in almost all published studies. In some publications, the models used do not represent the exact solution of the underlying one-compartment differential equations; others use unrealistic assumptions about constant background metal concentration and/or zero metal concentration in uncontaminated medium. Herein we present exact solutions of two differential-equation models, one describing simple two-stage toxicokinetics (metal toxicokinetic follows the experimental phases: the uptake phase and the decontamination phase) and another that can be applied for more complex three-stage patterns (toxicokinetic pattern does not follow two phases determined by an experimenter). Using two case studies for carabids exposed via food, based on previously published data, we discuss and compare our models to those originally used to analyze the data. Our conclusion is that when metal toxicokinetic follows a one-compartment model, the exact solution of a set of differential equations should be used. The proposed models allow assimilation and elimination rates to change between toxicokinetic stages, and the three-stage model is flexible enough to fit patterns that are more complex than the classic two-stage model can handle.     

Zooplankton mortality and the dynamical behaviour of plankton population models

We investigate the dynamical behaviour of a simple plankton population model, which explicitly simulates the concentrations of nutrient, phytoplankton and zooplankton in the oceanic mixed layer. The model consists of three coupled ordinary differential equations. We use analytical and numerical techniques, focusing on the existence and nature of steady states and unforced oscillations (limit cycles) of the system. The oscillations arise from Hopf bifurcations, which are traced as each parameter in the model is varied across a realistic range. The resulting bifurcation diagrams are compared with those from our previouswork, where zooplankton mortality was simulated by a quadratic function—here we use a linear function, to represent alternative ecological assumptions. Oscillations occur across broader ranges of parameters for the linear mortality function than for the quadratic one, although the two sets of bifurcation diagrams show similar qualitative characteristics. The choice of zooplankton mortality function, or closure term, is an area of current interest in the modelling community, and we relate our results to simulations of other models.     

Analytical and Simulation Results for Stochastic Fitzhugh-Nagumo Neurons and Neural Networks

An analytical approach is presented for determining the response of a neuron or of the activity in a network of connected neurons, represented by systems of nonlinear ordinary stochastic differential equations—the Fitzhugh-Nagumo system with Gaussian white noise current. For a single neuron, five equations hold for the first- and second-order central moments of the voltage and recovery variables. From this system we obtain, under certain assumptions, five differential equations for the means, variances, and covariance of the two components. One may use these quantities to estimate the probability that a neuron is emitting an action potential at any given time. The differential equations are solved by numerical methods. We also perform simulations on the stochastic Fitzugh-Nagumo system and compare the results with those obtained from the differential equations for both sustained and intermittent deterministic current inputs withsuperimposed noise. For intermittent currents, which mimic synaptic input, the agreement between the analytical and simulation results for the moments is excellent. For sustained input, the analytical approximations perform well for small noise as there is excellent agreement for the moments. In addition, the probability that a neuron is spiking as obtained from the empirical distribution of the potential in the simulations gives a result almost identical to that obtained using the analytical approach. However, when there is sustained large-amplitude noise, the analytical method is only accurate for short time intervals. Using the simulation method, we study the distribution of the interspike interval directly from simulated sample paths. We confirm that noise extends the range of input currents over which (nonperiodic) spike trains may exist and investigate the dependence of such firing on the magnitude of the mean input current and the noise amplitude. For networks we find the differential equations for the means, variances, and covariances of the voltage and recovery variables and show how solving them leads to an expression for the probability that a given neuron, or given set of neurons, is firing at time t. Using such expressions one may implement dynamical rules for changing synaptic strengths directly without sampling. The present analytical method applies equally well to temporally nonhomogeneous input currents and is expected to be useful for computational studies of information processing in various nervous system centers.     

Coagulation in cell suspensions: extensions of the von Smoluchowski model

Recently, a long-range interactive force between erythrocytes has been proposed (Rowlands et al., 1981, 1982a,b) based on an apparent increase in the rate of aggregation of erythrocytes in an aqueous suspension over that predicted by one model of Brownian aggregation. Here, we examine the assumptions underlying this model and propose modifications compatible with the biological constraints on the model. The refined model is represented as a series of coupled differential equations representing the change in particle number density as a function of time. Numerical solution of these equations is consistent with the absence of an intercellular force.     

Analysis of three- and four-compartment models forin vivo radioligand-neuroreceptor interaction

Currently applied three-copartment models for analyzing kinetic data derived fromin vivo positron emission tomographic (PET) studies of radioligand-neuroreceptor interactions require assumptions which may not be strictly valid. Such assumptions include very rapid kinetics for nonspecific binding and the absence of multiple specific receptors or subtypes. Computer simulations, based on an exact analytical solution of the relevant differential equations, indicate the numerical errors that can arise when the assumptions are invalid. We propose a fourcompartment model which requires fewer assumptions. A simple relationship is derived for expressing the microscopic rate constants of either the three- or four-compartment model as explicit functions of the experimentally-observed macroscopic rate constants. This could eliminate the need for time-consuming, iterative, non-linear, curve-fitting approaches and numerical integration. The usefulness of the four-compartment model is limited, however, by the sensitivity and temporal resolution of current PET imaging devices.     

Mathematical analysis of a model for the renal concentrating mechanism

A system of ordinary differential equations, designed to model the counterflow system in the renal medulla, is studied. An existence theorem for solutions of the model equations is obtained. An exact solution of the system is obtained in the limiting case of infinite water permeability. If there is diffusion in the core, evaluation of the exact solution leads to multiple stable solutions of the model equations. One solution has a large concentration ratio, which tends to a finite asymptotic limit as the pump strength tends to infinity.     

Steady state kinetics of consecutive enzyme catalysed reactions involving single substrates: procedures for the interpretation of coupled assays

Sets of differential rate equations are written describing a linear sequence of reactions occurring in solution each catalysed by a control enzyme or one of the Michaelis-Menten type. It is shown that the solutions of these equations may be formulated as a set of Maclaurin polynomials, expressing the concentration of each reactant and of final product as a function of time. From arrays of such polynomials, general expressions are induced for the first non-zero term of the series. These are used to formulate a procedure (illustrated with an example simulated by numerical integration) by which results of coupled enzymic assays may be analysed in terms of maximal velocities and apparent Michaelis constants: correlation is made with other established methods for conducting coupled assays. The present procedure assumes a steady state of enzyme-substrate complexes but not of intermediate reactants.     

Pair formation models with maturation period

The standard model for pair formation is generalized to include a maturation period. This model in the form of three coupled delay equations is a special case of the general age-structured model for a two-sex population. The exact conditions for the existence of an exponential (persistent) two-sex solution are derived. It is shown that this solution is unique and locally stable. In order to achieve these results the theory of homogeneous differential equations is extended to a class of homogeneous delay equations.     

Complexity in a prey-predator model with prey refuge and diffusion

Prey-predator interaction is one of the most commonly observed relationships in ecosystem. In the study of prey-predator models, it is frequently assumed that the changes in population densities are only time-dependent and the dynamics is generally represented by coupled nonlinear ordinary differential equations. In natural system, however, either prey or predator or both move from one place to another for various reasons. In such a case, their dynamic interaction depends both on time and space and requires coupled nonlinear partial differential equations for its dynamic representation. It is also well documented that prey refuges affect the interaction between prey and predator significantly. In this paper, we studied the dynamics of a diffusive prey-predator interaction with prey refuge and type III response function. We have considered both one and two dimensional diffusivity in the model system and presented different stability results under the assumptions that one or both species may be mobile or sedentary. Our results showed that the system may exhibit different spatiotemporal (non-Turing) patterns, like spiral waves, patchy structures, spot pattern, or even spatiotemporal chaos depending on the refuge availability and diffusion rate of species. Another interesting finding was that the dynamic complexity in a prey-predator model increases in case of mobile predator and sedentary prey compare to mobile prey and sedentary predator while refuge availability is varied.     

New Moment Closures Based on A Priori Distributions with Applications to Epidemic Dynamics

Recently, research that focuses on the rigorous understanding of the relation between simulation and/or exact models on graphs and approximate counterparts has gained lots of momentum. This includes revisiting the performance of classic pairwise models with closures at the level of pairs and/or triples as well as effective-degree-type models and those based on the probability generating function formalism. In this paper, for a fully connected graph and the simple SIS (susceptible-infected-susceptible) epidemic model, a novel closure is introduced. This is done via using the equations for the moments of the distribution describing the number of infecteds at all times combined with the empirical observations that this is well described/approximated by a binomial distribution with time dependent parameters. This assumption allows us to express higher order moments in terms of lower order ones and this leads to a new closure. The significant feature of the new closure is that the difference of the exact system, given by the Kolmogorov equations, from the solution of the newly defined approximate system is of order 1/N(2). This is in contrast with the O(1/N) difference corresponding to the approximate system obtained via the classic triple closure. The fully connected nature of the graph also allows us to interpret pairwise equations in terms of the moments and thus treat closures and the two approximate models within the same framework. Finally, the applicability and limitations of the new methodology is discussed in detail.     

Models of central pattern generators for quadruped locomotion I. Primary gaits

In this paper we continue the analysis of a network of symmetrically coupled cells modeling central pattern generators for quadruped locomotion proposed by Golubitsky, Stewart, Buono, and Collins. By a cell we mean a system of ordinary differential equations and by a coupled cell system we mean a network of identical cells with coupling terms. We have three main results in this paper. First, we show that the proposed network is the simplest one modeling the common quadruped gaits of walk, trot, and pace. In doing so we prove a general theorem classifying spatio-temporal symmetries of periodic solutions to equivariant systems of differential equations. We also specialize this theorem to coupled cell systems. Second, this paper focuses on primary gaits; that is, gaits that are modeled by output signals from the central pattern generator where each cell emits the same waveform along with exact phase shifts between cells. Our previous work showed that the network is capable of producing six primary gaits. Here, we show that under mild assumptions on the cells and the coupling of the network, primary gaits can be produced from Hopf bifurcation by varying only coupling strengths of the network. Third, we discuss the stability of primary gaits and exhibit these solutions by performing numerical simulations using the dimensionless Morris-Lecar equations for the cell dynamics.     

Exact results for fixation probability of bithermal evolutionary graphs

One of the most fundamental concepts of evolutionary dynamics is the “fixation” probability, i.e. the probability that a mutant spreads through the whole population. Most natural communities are geographically structured into habitats exchanging individuals among each other and can be modeled by an evolutionary graph (EG), where directed links weight the probability for the offspring of one individual to replace another individual in the community. EGs have recently spurred huge interest, as it has been shown that some topology can amplify or suppress the effect of beneficial mutations. Very few exact analytical results however are known for EGs. In this article we show that the use of a new technique, the fixed point of probability generating function, allows us to compute the exact fixation probability for a large subset of bithermal graphs. We also show by numerical simulations that the computed solution holds for all bithermal graphs. Moreover, the analytical solution allows us to clarify the opposing consequences of birth–death versus death–birth processes as amplifier or suppressor of beneficial mutations for the same bithermal topology.     

Theoretical analysis of compartmented coupling in linear enzyme systems

Exact equations which describe the kinetic patterns of enzyme/enzyme complexes, when compartmented coupling occurs between them, are presented. Compartmented coupling refers to the creation of a local environment in which the concentration of an intermediate, shared by two enzymes, is higher than its solution concentration. This results in a higher coupling enzyme activity, a condition reflected in a shorter transition time for the system. In this paper, equations are presented which allow experimenters to quantitate the effect of compartmented coupling in terms of changes in the apparent Km and Vmax values. The equations presented in this paper are more exact than those previously derived since they do not incorporate first order assumptions before derivation.     

Model reduction using a posteriori analysis

Mathematical models in biology and physiology are often represented by large systems of non-linear ordinary differential equations. In many cases, an observed behaviour may be written as a linear functional of the solution of this system of equations. A technique is presented in this study for automatically identifying key terms in the system of equations that are responsible for a given linear functional of the solution. This technique is underpinned by ideas drawn from a posteriori error analysis. This concept has been used in finite element analysis to identify regions of the computational domain and components of the solution where a fine computational mesh should be used to ensure accuracy of the numerical solution. We use this concept to identify regions of the computational domain and components of the solution where accurate representation of the mathematical model is required for accuracy of the functional of interest. The technique presented is demonstrated by application to a model problem, and then to automatically deduce known results from a cell-level cardiac electrophysiology model.     

Stability analysis for a peri-implant osseointegration model

We investigate stability of the solution of a set of partial differential equations, which is used to model a peri-implant osseointegration process. For certain parameter values, the solution has a ‘wave-like’ profile, which appears in the distribution of osteogenic cells, osteoblasts, growth factor and bone matrix. This ‘wave-like’ profile contradicts experimental observations. In our study we investigate the conditions, under which such profile appears in the solution. Those conditions are determined in terms of model parameters, by means of linear stability analysis, carried out at one of the constant solutions of the simplified system. The stability analysis was carried out for the reduced system of PDE’s, of which we prove, that it is equivalent to the original system of equations, with respect to the stability properties of constant solutions. The conclusions, derived from the linear stability analysis, are extended for the case of large perturbations. If the constant solution is unstable, then the solution of the system never converges to this constant solution. The analytical results are validated with finite element simulations. The simulations show, that stability of the constant solution could determine the behavior of the solution of the whole system, if certain initial conditions are considered.     

On emmetropization

The growth of the ocular elements, such as axial length, corneal power, and lens power, is represented by a set of first order, linear, coupled differential equations. Since the rate of growth and the size of the ocular elements vary from one individual to the next, the parameters which arise in these differential equations are assumed to be Gaussian random variables. These assumptions provide a basis for understanding the controlled growth of the eye. They predict that the distribution of refractive error is normal at birth and leptokurtic in adulthood. They also predict that the variance of refractive error will decrease with age provided the second order moments of the random variables obey some rather general conditions.     

Mathematical properties of pump-leak models of cell volume control and electrolyte balance

Homeostatic control of cell volume and intracellular electrolyte content is a fundamental problem in physiology and is central to the functioning of epithelial systems. These physiological processes are modeled using pump-leak models, a system of differential algebraic equations that describes the balance of ions and water flowing across the cell membrane. Despite their widespread use, very little is known about their mathematical properties. Here, we establish analytical results on the existence and stability of steady states for a general class of pump-leak models. We treat two cases. When the ion channel currents have a linear current-voltage relationship, we show that there is at most one steady state, and that the steady state is globally asymptotically stable. If there are no steady states, the cell volume tends to infinity with time. When minimal assumptions are placed on the properties of ion channel currents, we show that there is an asymptotically stable steady state so long as the pump current is not too large. The key analytical tool is a free energy relation satisfied by a general class of pump-leak models, which can be used as a Lyapunov function to study stability.     

Exact epidemic models on graphs using graph-automorphism driven lumping

The dynamics of disease transmission strongly depends on the properties of the population contact network. Pair-approximation models and individual-based network simulation have been used extensively to model contact networks with non-trivial properties. In this paper, using a continuous time Markov chain, we start from the exact formulation of a simple epidemic model on an arbitrary contact network and rigorously derive and prove some known results that were previously mainly justified based on some biological hypotheses. The main result of the paper is the illustration of the link between graph automorphisms and the process of lumping whereby the number of equations in a system of linear differential equations can be significantly reduced. The main advantage of lumping is that the simplified lumped system is not an approximation of the original system but rather an exact version of this. For a special class of graphs, we show how the lumped system can be obtained by using graph automorphisms. Finally, we discuss the advantages and possible applications of exact epidemic models and lumping.