A Hepatitis C Virus Infection Model with Time-Varying Drug Effectiveness: Solution and Analysis

Simple models of therapy for viral diseases such as hepatitis C virus (HCV) or human immunodeficiency virus assume that, once therapy is started, the drug has a constant effectiveness. More realistic models have assumed either that the drug effectiveness depends on the drug concentration or that the effectiveness varies over time. Here a previously introduced varying-effectiveness (VE) model is studied mathematically in the context of HCV infection. We show that while the model is linear, it has no closed-form solution due to the time-varying nature of the effectiveness. We then show that the model can be transformed into a Bessel equation and derive an analytic solution in terms of modified Bessel functions, which are defined as infinite series, with time-varying arguments. Fitting the solution to data from HCV infected patients under therapy has yielded values for the parameters in the model. We show that for biologically realistic parameters, the predicted viral decay on therapy is generally biphasic and resembles that predicted by constant-effectiveness (CE) models. We introduce a general method for determining the time at which the transition between decay phases occurs based on calculating the point of maximum curvature of the viral decay curve. For the parameter regimes of interest, we also find approximate solutions for the VE model and establish the asymptotic behavior of the system. We show that the rate of second phase decay is determined by the death rate of infected cells multiplied by the maximum effectiveness of therapy, whereas the rate of first phase decline depends on multiple parameters including the rate of increase of drug effectiveness with time.     

Combining Gompertzian growth and cell population dynamics

A two-compartment model of cancer cells population dynamics proposed by Gyllenberg and Webb includes transition rates between proliferating and quiescent cells as non-specified functions of the total population, N. We define the net inter-compartmental transition rate function: Phi(N). We assume that the total cell population follows the Gompertz growth model, as it is most often empirically found and derive Phi(N). The Gyllenberg-Webb transition functions are shown to be characteristically related through Phi(N). Effectively, this leads to a hybrid model for which we find the explicit analytical solutions for proliferating and quiescent cell populations, and the relations among model parameters. Several classes of solutions are examined. Our model predicts that the number of proliferating cells may increase along with the total number of cells, but the proliferating fraction appears to be a continuously decreasing function. The net transition rate of cells is shown to retain direction from the proliferating into the quiescent compartment. The death rate parameter for quiescent cell population is shown to be a factor in determining the proliferation level for a particular Gompertz growth curve.     

A Fokker-Planck equation for growing cell populations

Closed form solutions are obtained for a Fokker-Planck model for cell growth as a function of maturation velocity and degree of maturation. For reproduction rules where daughter cells inherit their parent's maturation velocity the complete solution is derived in terms of Airy functions. For more complicated reproduction rules partial results are obtained. Emphasis is given to the relationship of these problems to time dependent linear transport theory.     

Resource competition in stage-structured populations

Two models are made to account for the dynamics of a consumer-resource system in which the consumers are divided into juveniles and adults. The resource grows logistically and a type II functional response is assumed for consumers. Resource levels determine fecundity and maturation rates in one model, and mortality rates in the other. The analysis of the models shows that the condition for establishment of consumers is that the product of per capita fecundity rate and maturation rates is higher than the product of juvenile and adult per capita decay rates at a resource level equal to its carrying capacity. This result imposes a minimal abundance of resource able to maintain the consumers. A second result shows an equilibrium stage structure, with a small instability when juveniles and adults mean saturation constants are different. The implications of these results for community dynamics are discussed.     

Distributed parameter identification for a label-structured cell population dynamics model using CFSE histogram time-series data

In this work we address the problem of the robust identification of unknown parameters of a cell population dynamics model from experimental data on the kinetics of cells labelled with a fluorescence marker defining the division age of the cell. The model is formulated by a first order hyperbolic PDE for the distribution of cells with respect to the structure variable x (or z) being the intensity level (or the log₁₀-transformed intensity level) of the marker. The parameters of the model are the rate functions of cell division, death, label decay and the label dilution factor. We develop a computational approach to the identification of the model parameters with a particular focus on the cell birth rate α(z) as a function of the marker intensity, assuming the other model parameters are scalars to be estimated. To solve the inverse problem numerically, we parameterize α(z) and apply a maximum likelihood approach. The parametrization is based on cubic Hermite splines defined on a coarse mesh with either equally spaced a priori fixed nodes or nodes to be determined in the parameter estimation procedure. Ill-posedness of the inverse problem is indicated by multiple minima. To treat the ill-posed problem, we apply Tikhonov regularization with the regularization parameter determined by the discrepancy principle. We show that the solution of the regularized parameter estimation problem is consistent with the data set with an accuracy within the noise level in the measurements.      

Optimality of Mutation and Selection in Germinal Centers

The population dynamics theory of B cells in a typical germinal center could play an important role in revealing how affinity maturation is achieved. However, the existing models encountered some conflicts with experiments. To resolve these conflicts, we present a coarse-grained model to calculate the B cell population development in affinity maturation, which allows a comprehensive analysis of its parameter space to look for optimal values of mutation rate, selection strength, and initial antibody-antigen binding level that maximize the affinity improvement. With these optimized parameters, the model is compatible with the experimental observations such as the ∼100-fold affinity improvements, the number of mutations, the hypermutation rate, and the “all or none” phenomenon. Moreover, we study the reasons behind the optimal parameters. The optimal mutation rate, in agreement with the hypermutation rate in vivo, results from a tradeoff between accumulating enough beneficial mutations and avoiding too many deleterious or lethal mutations. The optimal selection strength evolves as a balance between the need for affinity improvement and the requirement to pass the population bottleneck. These findings point to the conclusion that germinal centers have been optimized by evolution to generate strong affinity antibodies effectively and rapidly. In addition, we study the enhancement of affinity improvement due to B cell migration between germinal centers. These results could enhance our understanding of the functions of germinal centers.     

A cytokinetic model for heterogeneous mammalian cell populations. I. Cell growth and cell death

A general model is proposed for describing the growth behavior of mammalian cell populations, which features:(a) a cell cycle time distribution function with properties such that mean and variance increase with increasing population size; (b) maturation age and maturation rate functions which constrain the maturational pathways of individual cells; and (c) a death rate function, where cell death is construed as irreparable damage to a cell's reproductive apparatus. The biological implications of the model are discussed, and methods for relating the model to real cell systems by means of commonly used experimental techniques are described. The model is compared with earlier models.     

Rate of excitation propagation in a reduced Hodgkins-Huxley model. I. Rapid relaxation of the sodium current]

The question of calculating excitation propagation velocity is analyzed on the basis of the Hodgkin-Huxley model. The activation of the sodium current is assumed to be rapid as compared to the rate of potential variation. Because of slow variation of potassium activation and sodium inactivation the dynamics of these processes is assumed to be of negligible effect in the region of impulse velocity formation. By means of pieace-wise linear approximation of thus obtained voltage-current characteristics the characteristics the analytical solution of the problem was found. In two limiting cases this solution coincides with the solutions of Kolmogorov and Scott. The dependence of impulse velocity on parameters is analyzed and illustrated graphically.     

Neutrophil marrow release. A system analysis.

Neutrophil marrow egress is governed by several processes. The most important are cell maturation, functional behavior of marrow sinusoids and humoral or neuro-vascular factors. Neutrophil release cannot be observed directly but is reflected in the size, cellular composition and kinetics of the nonproliferating pool of granulocytopoiesis in bone marrow and of blood neutrophil pool. These experimentally determined parameters were used as the basis of a mathematical model study. The model describes two catenated compartments, the nonproliferating pool of granulocytopoiesis in marrow and the total blood granulocyte pool. Cell transit from one pool to the other was assumed to be age-dependent. It was expressed by a positive sloping sigmoidal function that defines the egress potential fo the cells that increases with cell maturation. During maturation granulocytopoietic cells develop intense motility which determines the morphology of the cells on smears. Relationship between cell motility and its morphology was defined by functions determining the age-dependent probabilities of cell fixation as metamyelocytes, band- and segmented forms, respectively. The parameters of this model could be so adjusted that all experimental data were matched within experimental errors. Thus, qualitative and quantitative information on neutrophil marrow egress was obtained for normal and pathological states of granulocytopoiesis.     

Transient Responses to Rapid Changes in Mean and Variance in Spiking Models

The mean input and variance of the total synaptic input to a neuron can vary independently, suggesting two distinct information channels. Here we examine the impact of rapidly varying signals, delivered via these two information conduits, on the temporal dynamics of neuronal firing rate responses. We examine the responses of model neurons to step functions in either the mean or the variance of the input current. Our results show that the temporal dynamics governing response onset depends on the choice of model. Specifically, the existence of a hard threshold introduces an instantaneous component into the response onset of a leaky-integrate-and-fire model that is not present in other models studied here. Other response features, for example a decaying oscillatory approach to a new steady-state firing rate, appear to be more universal among neuronal models. The decay time constant of this approach is a power-law function of noise magnitude over a wide range of input parameters. Understanding how specific model properties underlie these response features is important for understanding how neurons will respond to rapidly varying signals, as the temporal dynamics of the response onset and response decay to new steady-state determine what range of signal frequencies a population of neurons can respond to and faithfully encode.     

Iteroparous reproduction strategies and population dynamics

Asymptotic relationships between a class of continuous partial differential equation population models and a class of discrete matrix equations are derived for iteroparous populations. First, the governing equations are presented for the dynamics of an individual with juvenile and adult life stages. The organisms reproduce after maturation, as determined by the juvenile period, and at specific equidistant ages, which are determined by the iteroparous reproductive period. A discrete population matrix model is constructed that utilizes the reproductive information and a density-dependent mortality function. Mortality in the period between two reproductive events is assumed to be a continuous process where the death rate for the adults is a function of the number of adults and environmental conditions. The asymptotic dynamic behaviour of the discrete population model is related to the steady-state solution of the continuous-time formulation. Conclusions include that there can be a lack of convergence to the steady-state age distribution in discrete event reproduction models. The iteroparous vital ratio (the ratio between the maximal age and the reproductive period) is fundamental to determining this convergence. When the vital ratio is rational, an equivalent discrete-time model for the population can be derived whose asymptotic dynamics are periodic and when there are a finite number of founder cohorts, the number of cohorts remains finite. When the ratio is an irrational number, effectively there is convergence to the steady-state age distribution. With a finite number of founder cohorts, the number of cohorts becomes countably infinite. The matrix model is useful to clarify numerical results for population models with continuous densities as well as delta measure age distribution. The applicability in ecotoxicology of the population matrix model formulation for iteroparous populations is discussed.     

Constrained analysis of fluorescence anisotropy decay:application to experimental protein dynamics

Hydrodynamic properties as well as structural dynamics of proteins can be investigated by the well-established experimental method of fluorescence anisotropy decay. Successful use of this method depends on determination of the correct kinetic model, the extent of cross-correlation between parameters in the fitting function, and differences between the timescales of the depolarizing motions and the fluorophore's fluorescence lifetime. We have tested the utility of an independently measured steady-state anisotropy value as a constraint during data analysis to reduce parameter cross correlation and to increase the timescales over which anisotropy decay parameters can be recovered accurately for two calcium-binding proteins. Mutant rat F102W parvalbumin was used as a model system because its single tryptophan residue exhibits monoexponential fluorescence intensity and anisotropy decay kinetics. Cod parvalbumin, a protein with a single tryptophan residue that exhibits multiexponential fluorescence decay kinetics, was also examined as a more complex model. Anisotropy decays were measured for both proteins as a function of solution viscosity to vary hydrodynamic parameters. The use of the steady-state anisotropy as a constraint significantly improved the precision and accuracy of recovered parameters for both proteins, particularly for viscosities at which the protein's rotational correlation time was much longer than the fluorescence lifetime. Thus, basic hydrodynamic properties of larger biomolecules can now be determined with more precision and accuracy by fluorescence anisotropy decay.     

Catastrophic shifts in vertical distributions of phytoplankton The existence of a bifurcation set

A model of phytoplankton dynamics within a water column was analyzed with special consideration on the existence of a bifurcation set in the parameter space. We considered two resources, light and a limiting nutrient, for phytoplankton growth and assumed that the water column is separated into two layers by thermal and/or density stratification. It was shown that there exists a bifurcation set in the parameter space when the growth function meets several conditions that are general for growth functions of two essential resources. Specifically, these conditions include that a less abundant of the two resources limits the growth while the effect of the other is sufficiently small. Folded structure with two stable states separated by one unstable state appears in the catastrophe manifold when parameters move to a certain direction with a certain curvature from a point in the bifurcation set. These results suggest that occurrence of discontinuous transition between two alternative vertical patterns is possible nature of phytoplankton dynamics within a stratified water column.     

HIV model parameter estimates from interruption trial data including drug efficacy and reservoir dynamics

Mathematical models based on ordinary differential equations (ODE) have had significant impact on understanding HIV disease dynamics and optimizing patient treatment. A model that characterizes the essential disease dynamics can be used for prediction only if the model parameters are identifiable from clinical data. Most previous parameter identification studies for HIV have used sparsely sampled data from the decay phase following the introduction of therapy. In this paper, model parameters are identified from frequently sampled viral-load data taken from ten patients enrolled in the previously published AutoVac HAART interruption study, providing between 69 and 114 viral load measurements from 3-5 phases of viral decay and rebound for each patient. This dataset is considerably larger than those used in previously published parameter estimation studies. Furthermore, the measurements come from two separate experimental conditions, which allows for the direct estimation of drug efficacy and reservoir contribution rates, two parameters that cannot be identified from decay-phase data alone. A Markov-Chain Monte-Carlo method is used to estimate the model parameter values, with initial estimates obtained using nonlinear least-squares methods. The posterior distributions of the parameter estimates are reported and compared for all patients.     

Stance leg control: variation of leg parameters supports stable hopping

The spring-loaded inverted pendulum describes the planar center-of-mass dynamics of legged locomotion. This model features linear springs with constant parameters as legs. In biological systems, however, spring-like properties of limbs can change over time. Therefore, in this study, it is asked how variation of spring parameters during ground contact would affect the dynamics of the spring-mass model. Neglecting damping initially, it is found that decreasing stiffness and increasing rest length of the leg during a stance phase are required for orbitally stable hopping. With damping, stable hopping is found for a larger region of rest-length rates and stiffness rates. Here, also increasing stiffness and decreasing rest length can result in stable hopping. Within the predicted range of leg parameter variations for stable hopping, there is no need for precise parameter tuning. Since hopping gaits form a subset of the running gaits (with vanishing horizontal velocity), these results may help to improve leg design in robots and prostheses.     

A new constitutive model for multi-layered collagenous tissues

Collagenous tissues such as the aneurysmal wall or the aorta are multi-layered structures with the mean fibre alignments distinguishing one layer from another. A constitutive representation of the multiple collagen layers is not yet developed, and hence the aim of the present study. The proposed model is based on the constitutive theory of finite elasticity and is characterized by an anisotropic strain-energy function which takes the material structure into account. The passive tissue behaviour is modelled and the related mechanical response is assumed to be dominated by elastin and collagen. While elastin is modelled by the neo-Hookean material the constitutive response of collagen is assumed to be transversely isotropic for each individual layer and based on an exponential function. The proposed constitutive function is polyconvex which ensures material stability. The model has five independent material parameters, each of which has a clear physical interpretation: the initial stiffnesses of the collagen fabric in the two principal directions, the shear modulus pertaining to the non-collagenous matrix material, a parameter describing the level of nonlinearity of the collagen fabric, and the angle between the principal directions of the collagen fabric and the reference coordinate system. An extension-inflation test of the adventitia of a human femoral artery is simulated by means of the finite element method and an error function is minimized by adjusting the material parameters yielding a good agreement between the model and the experimental data.     

Mechanisms of dendritic integration underlying gain control in fly motion-sensitive interneurons

In the compensatory optomotor response of the fly the interesting phenomenon of gain control has been observed by Reichardt and colleagues (Reichardt et al., 1983): The amplitude of the response tends to saturate with increasing stimulus size, but different saturation plateaus are assumed with different velocities at which the stimulus is moving. This characteristic can already be found in the motion-sensitive large field neurons of the fly optic lobes that play a role in mediating this behavioral response (Hausen, 1982; Reichardt et al, 1983; Egelhaaf, 1985; Haag et al., 1992). To account for gain control a model was proposed involving shunting inhibition of these cells by another cell, the so-called pool cell (Reichardt et al., 1983), both cells sharing common input from an array of local motion detectors. This article describes an alternative model which only requires dendritic integration of the output signals of two types of local motion detectors with opposite polarity. The explanation of gain control relies on recent findings that these input elements are not perfectly directionally selective and that their direction selectivity is a function of pattern velocity. As a consequence, the resulting postsynaptic potential in the dendrite of the integrating cell saturates with increasing pattern size at a level between the excitatory and inhibitory reversal potentials. The exact value of saturation is then set by the activation ratio of excitatory and inhibitory input elements which in turn is a function of other stimulus parameters such as pattern velocity. Thus, the apparently complex phenomenon of gain control can be simply explained by the biophysics of dendritic integration in conjunction with the properties of the motion-sensitive input elements.     

Separable models in age-dependent population dynamics

A class of population models is considered in which the parameters such as fecundity, mortality and interaction coefficients are assumed to be age-dependent. Conditions for the existence, stability and global attractivity of steady-state and periodic solutions are derived. The dependence of these solutions on the maturation periods is analyzed. These results are applied to specific single and multiple population models. It is shown that periodic solutions cannot occur in a general class of single population age-dependent models. Conditions are derived that determine whether increasing the maturation period has a stabilizing effect. In specific cases, it is shown that any number of switches in stability can occur as the maturation period is increased. An example is given of predator-prey model where each one of these stability switches corresponds to a stable steady state losing its stability via a Hopf bifurcation to a periodic solution and regaining its stability upon further increase of the maturation period.     

Interpreting photoactivated fluorescence microscopy measurements of steady-state actin dynamics.

A continuum model describing the steady-state actin dynamics of the cytoskeleton of living cells has been developed to aid in the interpretation of photoactivated fluorescence experiments. In a simplified cell geometry, the model assumes uniform concentrations of cytosolic and cytoskeletal actin throughout the cell and no net growth of either pool. The spatiotemporal evolution of the fluorescent actin population is described by a system of two coupled linear partial-differential equations. An analytical solution is found using a Fourier-Laplace transform and important limiting cases relevant to the design of experiments are discussed. The results demonstrate that, despite being a complex function of the parameters, the fluorescence decay in photoactivated fluorescence experiments has a biphasic behavior featuring a short-term decay controlled by monomer diffusion and a long-term decay governed by the monomer exchange rate between the polymerized and unpolymerized actin pools. This biphasic behavior suggests a convenient mechanism for extracting the parameters governing the fluorescence decay from data records. These parameters include the actin monomer diffusion coefficient, filament turnover rate, and ratio of polymerized to unpolymerized actin.     

Modeling the Genetic Control of HIV-1 Dynamics After Highly Active Antiretroviral Therapy

The progression of HIV disease has been markedly slowed by the use of highly active antiretroviral therapy (HAART). However, substantial genetic variation was observed to occur among different people in the decay rate of viral loads caused by HAART. The characterization of specific genes involved in HIV dynamics is central to design personalized drugs for the prevention of this disease, but usually cannot be addressed by experimental methods alone rather than require the help of mathematical and statistical methods. A novel statistical model has been recently developed to detect genetic variants that are responsible for the shape of HAART-induced viral decay curves. This model was employed to an HIV/AIDS trial, which led to the identification of a major genetic determinant that triggers an effect on HIV dynamics. This detected major genetic determinant also affects several clinically important parameters, such as half-lives of infected cells and HIV eradication times.Key Words: Hardy-weinberg equilibrium, bi-exponential function, quantitative trait loci, HIV dynamics, functional mapping.