Mean and quasideterministic equivalence for linear stochastic dynamics

In linear, stochastic dynamics it is shown that the quasideterministic population size is equivalent to the mean population size. The quasideterministic dynamics are defined by the conditional infinitesimal mean of the process. The stochastic component of the dynamics includes both Gaussian and Poisson white noise, with amplitude coefficients proportional to the population size. Generalizations are given for nonautonomous coefficients and for distributed Poisson jump amplitudes. A counter example--an exactly integrable nonlinear jump model--shows that the equivalence result does not hold for nonlinear stochastic dynamics.     

Optimal harvesting of a logistic population in an environment with stochastic jumps

Dynamic programming is employed to examine the effects of large, sudden changes in population size on the optimal harvest strategy of an exploited resource population. These changes are either adverse or favorable and are assumed to occur at times of events of a Poisson process. The amplitude of these jumps is assumed to be density independent. In between the jumps the population is assumed to grow logistically. The Bellman equation for the optimal discounted present value is solved numerically and the optimal feedback control computed for the random jump model. The results are compared to the corresponding results for the quasi-deterministic approximation. In addition, the sensitivity of the results to the discount rate, the total jump rate and the quadratic cost factor is investigated. The optimal results are most strongly sensitive to the rate of stochastic jumps and to the quadratic cost factor to a lesser extent when the deterministic bioeconomic parameters are taken from aggregate antarctic pelagic whaling data.Research supported in part by the National Science Foundation under grants MCS 81-01698 and MCS 83-00562.     

Varieties of stochastic model: a comparative study of the Gompertz effect

A comparative study is made of various models for the Gompertz phenomenon, which is a form of growth rate limitation in population dynamics. Deterministic, Markov birth-death, diffusion and stochastic differential equation models are studied, with a view to assessing their advantages and limitations.     

Stochastic models for the Trojan Y-Chromosome eradication strategy of an invasive species

The Trojan Y-Chromosome (TYC) strategy, an autocidal genetic biocontrol method, has been proposed to eliminate invasive alien species. In this work, we develop a Markov jump process model for this strategy, and we verify that there is a positive probability for wild-type females going extinct within a finite time. Moreover, when sex-reversed Trojan females are introduced at a constant population size, we formulate a stochastic differential equation (SDE) model as an approximation to the proposed Markov jump process model. Using the SDE model, we investigate the probability distribution and expectation of the extinction time of wild-type females by solving Kolmogorov equations associated with these statistics. The results indicate how the probability distribution and expectation of the extinction time are shaped by the initial conditions and the model parameters.     

A comparison of persistence-time estimation for discrete and continuous stochastic population models that include demographic and environmental variability

A discrete-time Markov chain model, a continuous-time Markov chain model, and a stochastic differential equation model are compared for a population experiencing demographic and environmental variability. It is assumed that the environment produces random changes in the per capita birth and death rates, which are independent from the inherent random (demographic) variations in the number of births and deaths for any time interval. An existence and uniqueness result is proved for the stochastic differential equation system. Similarities between the models are demonstrated analytically and computational results are provided to show that estimated persistence times for the three stochastic models are generally in good agreement when the models satisfy certain consistency conditions.     

The impact of short term synaptic depression and stochastic vesicle dynamics on neuronal variability

Neuronal variability plays a central role in neural coding and impacts the dynamics of neuronal networks. Unreliability of synaptic transmission is a major source of neural variability: synaptic neurotransmitter vesicles are released probabilistically in response to presynaptic action potentials and are recovered stochastically in time. The dynamics of this process of vesicle release and recovery interacts with variability in the arrival times of presynaptic spikes to shape the variability of the postsynaptic response. We use continuous time Markov chain methods to analyze a model of short term synaptic depression with stochastic vesicle dynamics coupled with three different models of presynaptic spiking: one model in which the timing of presynaptic action potentials are modeled as a Poisson process, one in which action potentials occur more regularly than a Poisson process (sub-Poisson) and one in which action potentials occur more irregularly (super-Poisson). We use this analysis to investigate how variability in a presynaptic spike train is transformed by short term depression and stochastic vesicle dynamics to determine the variability of the postsynaptic response. We find that sub-Poisson presynaptic spiking increases the average rate at which vesicles are released, that the number of vesicles released over a time window is more variable for smaller time windows than larger time windows and that fast presynaptic spiking gives rise to Poisson-like variability of the postsynaptic response even when presynaptic spike times are non-Poisson. Our results complement and extend previously reported theoretical results and provide possible explanations for some trends observed in recorded data.     

Gene expression dynamics in randomly varying environments

A simple model of gene regulation in response to stochastically changing environmental conditions is developed and analyzed. The model consists of a differential equation driven by a continuous time 2-state Markov process. The density function of the resulting process converges to a beta distribution. We show that the moments converge to their stationary values exponentially in time. Simulations of a two-stage process where protein production depends on mRNA concentrations are also presented demonstrating that protein concentration tracks the environment whenever the rate of protein turnover is larger than the rate of environmental change. Single-celled organisms are therefore expected to have relatively high mRNA and protein turnover rates for genes that respond to environmental fluctuations.     

Neural rate equations for bursting dynamics derived from conductance-based equations

A method of obtaining rate equations from conductance-based equations is developed and applied to fast-spiking and bursting neocortical neurons. It involves splitting systems of conductance-based equations into fast and slow subsystems, and averaging the effects of fast terms that drive the slowly varying quantities by showing that their average is closely proportional to the firing rate. The dependence of the firing rate on the injected current is then approximated in the analysis. The resulting behavior of the slow variables is then substituted back into the fast equations, with the further approximation of replacing the fast voltages in these terms by effective values. For bursting neurons the method yields two coupled limit-cycle oscillators: a self-exciting oscillator for the slow variables that commences limit-cycle oscillations at a critical current and modulates a fast spike-generating oscillator, thereby leading to slowly modulated bursts with a group of spikes in each burst. The dynamics of these coupled oscillators are then verified against those of the conductance-based equations. Finally, it is shown how to place the results in a form suitable for use in mean-field equations for neural population dynamics.     

Synchronization of an Excitatory Integrate-and-Fire Neural Network

In this paper, we study the influence of the coupling strength on the synchronization behavior of a population of leaky integrate-and-fire neurons that is self-excitatory with a population density approach. Each neuron of the population is assumed to be stochastically driven by an independent Poisson spike train and the synaptic interaction between neurons is modeled by a potential jump at the reception of an action potential. Neglecting the synaptic delay, we will establish that for a strong enough connectivity between neurons, the solution of the partial differential equation which describes the population density function must blow up in finite time. Furthermore, we will give a mathematical estimate on the average connection per neuron to ensure the occurrence of a burst. Interpreting the blow up of the solution as the presence of a Dirac mass in the firing rate of the population, we will relate the blow up of the solution to the occurrence of the synchronization of neurons. Fully stochastic simulations of a finite size network of leaky integrate-and-fire neurons are performed to illustrate our theoretical results.     

Bayesian estimators for conditional hazard functions

This article introduces a new Bayesian approach to the analysis of right-censored survival data. The hazard rate of interest is modeled as a product of conditionally independent stochastic processes corresponding to (1) a baseline hazard function and (2) a regression function representing the temporal influence of the covariates. These processes jump at times that form a time-homogeneous Poisson process and have a pairwise dependency structure for adjacent values. The two processes are assumed to be conditionally independent given their jump times. Features of the posterior distribution, such as the mean covariate effects and survival probabilities (conditional on the covariate), are evaluated using the Metropolis-Hastings-Green algorithm. We illustrate our methodology by an application to nasopharynx cancer survival data.     

Anti-phase,asymmetric and aperiodic oscillations in excitable cells—I. Coupled bursters

I seek to explain phenomena observed in simulations of populations of gap junction-coupled bursting cells by studying the dynamics of identical pairs. I use a simplified model for pancreatic β-cells and decompose the system into fast (spike-generating) and slow subsystems to show how bifurcations of the fast subsystem affect bursting behavior. When coupling is weak, the spikes are not in phase but rather are anti-phase, asymmetric or quasi-periodic. These solutions all support bursting with smaller amplitude spikes than the in-phase case, leading to increased burst period. A key geometrical feature underlying this is that the in-phase periodic solution branch terminates in a homoclinic orbit. The same mechanism also provides a model for bursting as an emergent property of populations; cells which are not intrinsic bursters can burst when coupled. This phenomenon is enhanced when symmetry is broken by making the cells differ in a parameter.     

The role of regulated mRNA stability in establishing bicoid morphogen gradient in Drosophila embryonic development

The Bicoid morphogen is amongst the earliest triggers of differential spatial pattern of gene expression and subsequent cell fate determination in the embryonic development of Drosophila. This maternally deposited morphogen is thought to diffuse in the embryo, establishing a concentration gradient which is sensed by downstream genes. In most model based analyses of this process, the translation of the bicoid mRNA is thought to take place at a fixed rate from the anterior pole of the embryo and a supply of the resulting protein at a constant rate is assumed. Is this process of morphogen generation a passive one as assumed in the modelling literature so far, or would available data support an alternate hypothesis that the stability of the mRNA is regulated by active processes? We introduce a model in which the stability of the maternal mRNA is regulated by being held constant for a length of time, followed by rapid degradation. With this more realistic model of the source, we have analysed three computational models of spatial morphogen propagation along the anterior-posterior axis: (a) passive diffusion modelled as a deterministic differential equation, (b) diffusion enhanced by a cytoplasmic flow term; and (c) diffusion modelled by stochastic simulation of the corresponding chemical reactions. Parameter estimation on these models by matching to publicly available data on spatio-temporal Bicoid profiles suggests strong support for regulated stability over either a constant supply rate or one where the maternal mRNA is permitted to degrade in a passive manner.     

Stochastic models for virus and immune system dynamics

New stochastic models are developed for the dynamics of a viral infection and an immune response during the early stages of infection. The stochastic models are derived based on the dynamics of deterministic models. The simplest deterministic model is a well-known system of ordinary differential equations which consists of three populations: uninfected cells, actively infected cells, and virus particles. This basic model is extended to include some factors of the immune response related to Human Immunodeficiency Virus-1 (HIV-1) infection. For the deterministic models, the basic reproduction number, R₀, is calculated and it is shown that if R₀<1, the disease-free equilibrium is locally asymptotically stable and is globally asymptotically stable in some special cases. The new stochastic models are systems of stochastic differential equations (SDEs) and continuous-time Markov chain (CTMC) models that account for the variability in cellular reproduction and death, the infection process, the immune system activation, and viral reproduction. Two viral release strategies are considered: budding and bursting. The CTMC model is used to estimate the probability of virus extinction during the early stages of infection. Numerical simulations are carried out using parameter values applicable to HIV-1 dynamics. The stochastic models provide new insights, distinct from the basic deterministic models. For the case R₀>1, the deterministic models predict the viral infection persists in the host. But for the stochastic models, there is a positive probability of viral extinction. It is shown that the probability of a successful invasion depends on the initial viral dose, whether the immune system is activated, and whether the release strategy is bursting or budding.     

Transitions in Smoking Behaviour and the Design of Cessation Schemes

The intake of nicotine by smoking cigarettes is modelled by a dynamical system of differential equations. The variables are the internal level of nicotine and the level of craving. The model is based on the dynamics of neural receptors and the way they enhance craving. Lighting of a cigarette is parametrised by a time-dependent Poisson process. The nicotine intake rate is assumed to be proportional with the parameter of this stochastic process. The effect of craving is damped by a control mechanism in which awareness of the risks of smoking and societal measures play a role. Fluctuations in this damping may cause transitions from smoking to non-smoking and vice versa. With the use of Monte Carlo simulation the effect of abrupt and gradual cessation therapies are evaluated. Combination of the two in a mixed scheme yields a therapy with a duration that can be set at wish.     

Solving the chemical master equation for monomolecular reaction systems analytically

The stochastic dynamics of a well-stirred mixture of molecular species interacting through different biochemical reactions can be accurately modelled by the chemical master equation (CME). Research in the biology and scientific computing community has concentrated mostly on the development of numerical techniques to approximate the solution of the CME via many realizations of the associated Markov jump process. The domain of exact and/or efficient methods for directly solving the CME is still widely open, which is due to its large dimension that grows exponentially with the number of molecular species involved. In this article, we present an exact solution formula of the CME for arbitrary initial conditions in the case where the underlying system is governed by monomolecular reactions. The solution can be expressed in terms of the convolution of multinomial and product Poisson distributions with time-dependent parameters evolving according to the traditional reaction-rate equations. This very structured representation allows to deduce easily many properties of the solution. The model class includes many interesting examples. For more complex reaction systems, our results can be seen as a first step towards the construction of new numerical integrators, because solutions to the monomolecular case provide promising ansatz functions for Galerkin-type methods.