Elucidating relationships between P.falciparum prevalence and measures of genetic diversity with a combined genetic-epidemiological model of malaria

There is an abundance of malaria genetic data being collected from the field, yet using these data to understand the drivers of regional epidemiology remains a challenge. A key issue is the lack of models that relate parasite genetic diversity to epidemiological parameters. Classical models in population genetics characterize changes in genetic diversity in relation to demographic parameters, but fail to account for the unique features of the malaria life cycle. In contrast, epidemiological models, such as the Ross-Macdonald model, capture malaria transmission dynamics but do not consider genetics. Here, we have developed an integrated model encompassing both parasite evolution and regional epidemiology. We achieve this by combining the Ross-Macdonald model with an intra-host continuous-time Moran model, thus explicitly representing the evolution of individual parasite genomes in a traditional epidemiological framework. Implemented as a stochastic simulation, we use the model to explore relationships between measures of parasite genetic diversity and parasite prevalence, a widely-used metric of transmission intensity. First, we explore how varying parasite prevalence influences genetic diversity at equilibrium. We find that multiple genetic diversity statistics are correlated with prevalence, but the strength of the relationships depends on whether variation in prevalence is driven by host- or vector-related factors. Next, we assess the responsiveness of a variety of statistics to malaria control interventions, finding that those related to mixed infections respond quickly (∼months) whereas other statistics, such as nucleotide diversity, may take decades to respond. These findings provide insights into the opportunities and challenges associated with using genetic data to monitor malaria epidemiology.     

An agent-based modelling approach to estimate error in gyrodactylid population growth

Comparative studies of gyrodactylid monogeneans on different host species or strains rely upon the observation of growth on individual fish maintained within a common environment, summarised using maximum likelihood statistical approaches. Here we describe an agent-based model of gyrodactylid population growth, which we use to evaluate errors due to stochastic reproductive variation in such experimental studies. Parameters for the model use available fecundity and mortality data derived from previously published life tables of Gyrodactylus salaris, and use a new data set of fecundity and mortality statistics for this species on the Neva stock of Atlantic salmon, Salmo salar. Mortality data were analysed using a mark-recapture analysis software package, allowing maximum-likelihood estimation of daily survivorship and mortality. We consistently found that a constant age-specific mortality schedule was most appropriate for G. salaris in experimental datasets, with a daily survivorship of 0.84 at 13°C. This, however, gave unrealistically low population growth rates when used as parameters in the model, and a schedule of constantly increasing mortality was chosen as the best compromise for the model. The model also predicted a realistic age structure for the simulated populations, with 0.32 of the population not yet having given birth for the first time (pre-first birth). The model demonstrated that the population growth rate can be a useful parameter for comparing gyrodactylid populations when these are larger than 20-30 individuals, but that stochastic error rendered the parameter unusable in smaller populations. It also showed that the declining parasite population growth rate typically observed during the course of G. salaris infections cannot be explained through stochastic error and must therefore have a biological basis. Finally, the study showed that most gyrodactylid-host studies of this type are too small to detect subtle differences in local adaptation of gyrodactylid monogeneans between fish stocks.     

An Alternative to Moment Closure

Moment closure methods are widely used to analyze mathematical models. They are specifically geared toward derivation of approximations of moments of stochastic models, and of similar quantities in other models. The methods possess several weaknesses: Conditions for validity of the approximations are not known, magnitudes of approximation errors are not easily evaluated, spurious solutions can be generated that require large efforts to eliminate, and expressions for the approximations are in many cases too complex to be useful. We describe an alternative method that provides improvements in these regards. The new method leads to asymptotic approximations of the first few cumulants that are explicit in the model’s parameters. We analyze the univariate stochastic logistic Verhulst model and a bivariate stochastic epidemic SIR model with the new method. Errors that were made in early applications of moment closure to the Verhulst model are explained and corrected.     

A model for estimating total parasite load in falciparum malaria patients

We describe an age-structured mathematical model of the malaria parasite life cycle that uses clinical observations of peripheral parasitaemia to estimate population dynamics of sequestered parasites, which are hidden from the clinical investigator. First, the model was tested on parasite populations cultured in vitro, and was found to account for approximately 72% of the variation in that sub-population of parasites that would have been sequestered in vivo. Next, the model was applied to patients undergoing antimalarial therapy. Using individual data sets we found that although the model fitted the peripheral parasite curves very well, unique solutions for the fit could not be obtained; therefore, robust estimates of sequestered parasite dynamics remained unavailable. We conclude that even given detailed data on individual parasitaemia, estimates of sequestered numbers may be difficult to obtain. However, if data on individuals undergoing similar therapy are collected at equal time intervals, some of these problems may be overcome by estimating specific parameters over groups of patients. In this manner we estimated sequestered parasite density in a group of patients sampled at identical time points following antimalarial treatment. Using this approach we found significant relationships between changes in parasite density, age structure and temperature that were not apparent from the analysis of peripheral parasitaemia only.     

Host–parasite determinants of parasite population structure: lessons from bats and mites on the importance of time

By definition, parasitic organisms are strongly dependant on their hosts, and for a great majority, this dependence includes host-to-host transmission. Constraints imposed by the host's spatial distribution and demography, in combination with those of the parasite, can lead to a metapopulation structure, where parasite populations are highly stochastic (i.e. prone to frequent extinctions and re-colonizations) and where drift becomes a major force shaping standing genetic variation. This, in turn, will directly affect the observed population structure, along with the ability of the parasite to adapt (or co-adapt) to its host. However, only a specific consideration of temporal dynamics can reveal the extent to which drift shapes parasite population structure; this is rarely taken into account in population genetic studies of parasitic organisms. The study by Bruyndonckx et al. in this issue of Molecular Ecology does just this and, in doing so, illustrates how a comparison of host–parasite co-structures in light of temporal dynamics can be particularly informative for understanding the ecological and evolutionary constraints imposed by the host. More specifically, the authors examine spatial and temporal population genetic data of a parasitic mite Spinturnix bechsteini that exclusively exploits the Bechstein's bat Myotis bechsteinii and consider these data in relation to host–parasite life histories and the population structure of the host.     

Estimating the variation in S phase duration from flow cytometric histograms

A stochastic model for interpreting BrdUrd DNA FCM-derived data is proposed. The model is based on branching processes and describes the progression of the DNA distribution of BrdUrd-labelled cells through the cell cycle. With the main focus on estimating the S phase duration and its variation, the DNA replication rate is modelled by a piecewise linear function, while assuming a gamma distribution for the S phase duration. Estimation of model parameters was carried out using maximum likelihood for data from two different cell lines. The results provided quite a good fit to the data, suggesting that stochastic models may be a valuable tool for analysing this kind of data.     

Technical and biological variance structure in mRNA-Seq data: life in the real world

ABSTRACT: BACKGROUND: mRNA expression data from next generation sequencing platforms is obtained in the form of counts per gene or exon. Counts have classically been assumed to follow a Poisson distribution in which the variance is equal to the mean. The Negative Binomial distribution which allows for over-dispersion, i.e., for the variance to be greater than the mean, is commonly used to model count data as well. RESULTS: In mRNA-Seq data from 25 subjects, we found technical variation to generally follow a Poisson distribution as has been reported previously and biological variability was over-dispersed relative to the Poisson model. The mean-variance relationship across all genes was quadratic, in keeping with a Negative Binomial (NB) distribution. Over-dispersed Poisson and NB distributional assumptions demonstrated marked improvements in goodness-of-fit (GOF) over the standard Poisson model assumptions, but with evidence of over-fitting in some genes. Modeling of experimental effects improved GOF for high variance genes but increased the over-fitting problem. CONCLUSIONS: These conclusions will guide development of analytical strategies for accurate modeling of variance structure in these data and sample size determination which in turn will aid in the identification of true biological signals that inform our understanding of biological systems.     

Distribution and abundance of parasite nematodes: ecological specialisation,phylogenetic constraint or simply epidemiology?

We investigate the patterns of abundance‐spatial occupancy relationships of adult parasite nematodes in mammal host populations (828 populations of nematodes from 66 different species of terrestrial mammals). A positive relationship between mean parasite abundance and host occupancy, i.e. prevalence, is found which suggests that local abundance is linked to spatial distribution across species. Moreover, the frequency distribution of the parasite prevalence is bimodal, which is consistent with a core‐satellite species distribution. In addition, a strong positive relationship between the abundance (log‐transformed) and its variance (log‐transformed) is observed, the distribution of worm abundance being lognormally distributed when abundance values have been corrected for host body size.
Hanski et al. proposed three distinct hypotheses, which might account for the positive relationship between abundance and prevalence in free and associated organisms: 1) ecological specialisation, 2) sampling artefact, and 3) metapopulation dynamics. In addition, Gaston and co‐workers listed five additional hypotheses. Four solutions were not applicable to our parasitological data due to the lack of relevant information in most host‐parasite studies. The fifth hypothesis, i.e. the confounded effects exerted by common history on observed patterns of parasite distributions, was considered using a phylogeny‐based comparison method. Testing the four possible hypotheses, we obtained the following results: 1) the variation of parasite distribution across host species is not due to phylogenetic confounding effects; 2) the positive relationship between mean abundance and prevalence of nematodes may not result from an ecological specialisation, i.e. host specificity, of these parasites; 3) both a positive abundance‐prevalence relationship and a negative coefficient of variation of abundance‐prevalence relationship are likely to occur which corroborates the sampling model developed by Hanski et al. We argue that demographic explanations may be of particular importance to explain the patterns of bimodality of prevalence when testing Monte‐Carlo simulations using epidemiological modelling frameworks, and when considering empirical findings. We conclude that both the bimodal distribution of parasite prevalence and the mean‐variance power function simply result from demographic and stochastic patterns (highlighted by the sampling model), which present compelling evidence that nematode parasite species might adjust their spatial distribution and burden in mammal hosts for simple epidemiological reasons.     

The beta-binomial model for host specificity among organisms in trophic interactions

In this paper, we present a new stochastic model where the host specificity among organisms in trophic interactions in a community, say parasite-host interactions, is estimated by a beta-binomial model. The expected proportion of the host species in a community that a given parasite species is utilizing is modeled as a realization from an inhomogeneous Poisson process, where the rate of this process is assumed to be proportional to a beta probability distribution. The observed number of host species utilized by the parasites is then binomially distributed with the number of trials equaling the number of different host species in the sample. When the degree of polyphagy is estimated by the parameters of the beta-binomial model, quantities like community host specificity and the expected total number of parasite species that will utilize the host species in the community can be predicted as functions of the number of host species available. The predictions can then be applied in analysis of, e.g., symbiotic interactions among organisms, local species richness, and community structure.     

Gain-of-function mutations reveal expanded intermediate states and a sequential action of two gates in MscL

The tension-driven gating transition in the large mechanosensitive channel MscL proceeds through detectable states of intermediate conductance. Gain-of-function (GOF) mutants with polar or charged substitutions in the main hydrophobic gate display altered patterns of subconducting states, providing valuable information about gating intermediates. Here we present thermodynamic analysis of several GOF mutants to clarify the nature and position of low-conducting conformations in the transition pathway. Unlike wild-type (WT) MscL, which predominantly occupies the closed and fully open states with very brief substates, the mild V23T GOF mutant frequently visits a multitude of short-lived subconducting states. Severe mutants V23D and G22N open in sequence: closed (C) --> low-conducting substate (S) --> open (O), with the first subtransition occurring at lower tensions. Analyses of equilibrium state occupancies as functions of membrane tension show that the C-->S subtransition in WT MscL is associated with only a minor conductance increment, but the largest in-plane expansion and free energy change. The GOF substitutions strongly affect the first subtransition by reducing area ((Delta)A) and energy ((Delta)E) changes between C and S states commensurably with the severity of mutation. GOF mutants also exhibited a considerably larger (Delta)E associated with the second (S-->O) subtransition, but a (Delta)A similar to WT. The area changes indicate that closed conformations of GOF mutants are physically preexpanded. The tension dependencies of rate constants for channel closure (k(off)) predict different positions of rate-limiting barriers on the energy-area profiles for WT and GOF MscL. The data support the two-gate mechanism in which the first subtransition (C-->S) can be viewed as opening of the central (M1) gate, resulting in an expanded water-filled "leaky" conformation. Strong facilitation of this step by polar GOF substitutions suggests that separation of M1 helices associated with hydration of the pore in WT MscL is the major energetic barrier for opening. Mutants with a stabilized S1 gate demonstrate impeded transitions from low-conducting substates to the fully open state, whereas extensions of S1-M1 linkers result in a much higher probability of reverse O-->S transitions. These data strongly suggest that the bulk of conductance gain in the second subtransition (S-->O) occurs through the opening of the NH2-terminal (S1) gate and the linkers are coupling elements between the M1 and S1 gates.     

Modelling variability in lymphatic filariasis: macrofilarial dynamics in the Brugia pahangi--cat model.

A striking feature of lymphatic filariasis is the considerable heterogeneity in infection burden observed between hosts, which greatly complicates the analysis of the population dynamics of the disease. Here, we describe the first application of the moment closure equation approach to model the sources and the impact of this heterogeneity for macrofilarial population dynamics. The analysis is based on the closest laboratory equivalent of the life cycle and immunology of infection in humans--cats chronically infected with the filarial nematode Brugia pahangi. Two sets of long-term experiments are analysed: hosts given either single primary infections or given repeat infections. We begin by quantifying changes in the mean and aggregation of adult parasites (inversely measured by the negative binomial parameter, kappa in cohorts of hosts using generalized linear models. We then apply simple stochastic models to interpret observed patterns. The models and empirical data indicate that parasite aggregation tracks the decline in the mean burden with host age in primary infections. Conversely, in repeat infections, aggregation increases as the worm burden declines with experience of infection. The results show that the primary infection variability is consistent with heterogeneities in parasite survival between hosts. By contrast, the models indicate that the reduction in parasite variability with time in repeat infections is most likely due to the ''filtering'' effect of a strong, acquired immune response, which gradually acts to remove the initial variability generated by heterogeneities in larval mortality. We discuss this result in terms of the homogenizing effect of host immunity-driven density-dependence on macrofilarial burden in older hosts.     

Influence of Vectors’ Risk-Spreading Strategies and Environmental Stochasticity on the Epidemiology and Evolution of Vector-Borne Diseases: The Example of Chagas’ Disease

Insects are known to display strategies that spread the risk of encountering unfavorable conditions, thereby decreasing the extinction probability of genetic lineages in unpredictable environments. To what extent these strategies influence the epidemiology and evolution of vector-borne diseases in stochastic environments is largely unknown. In triatomines, the vectors of the parasite Trypanosoma cruzi, the etiological agent of Chagas’ disease, juvenile development time varies between individuals and such variation most likely decreases the extinction risk of vector populations in stochastic environments. We developed a simplified multi-stage vector-borne SI epidemiological model to investigate how vector risk-spreading strategies and environmental stochasticity influence the prevalence and evolution of a parasite. This model is based on available knowledge on triatomine biodemography, but its conceptual outcomes apply, to a certain extent, to other vector-borne diseases. Model comparisons between deterministic and stochastic settings led to the conclusion that environmental stochasticity, vector risk-spreading strategies (in particular an increase in the length and variability of development time) and their interaction have drastic consequences on vector population dynamics, disease prevalence, and the relative short-term evolution of parasite virulence. Our work shows that stochastic environments and associated risk-spreading strategies can increase the prevalence of vector-borne diseases and favor the invasion of more virulent parasite strains on relatively short evolutionary timescales. This study raises new questions and challenges in a context of increasingly unpredictable environmental variations as a result of global climate change and human interventions such as habitat destruction or vector control.     

Parameter inference for stochastic biochemical models from perturbation experiments parallelised at the single cell level

Understanding and characterising biochemical processes inside single cells requires experimental platforms that allow one to perturb and observe the dynamics of such processes as well as computational methods to build and parameterise models from the collected data. Recent progress with experimental platforms and optogenetics has made it possible to expose each cell in an experiment to an individualised input and automatically record cellular responses over days with fine time resolution. However, methods to infer parameters of stochastic kinetic models from single-cell longitudinal data have generally been developed under the assumption that experimental data is sparse and that responses of cells to at most a few different input perturbations can be observed. Here, we investigate and compare different approaches for calculating parameter likelihoods of single-cell longitudinal data based on approximations of the chemical master equation (CME) with a particular focus on coupling the linear noise approximation (LNA) or moment closure methods to a Kalman filter. We show that, as long as cells are measured sufficiently frequently, coupling the LNA to a Kalman filter allows one to accurately approximate likelihoods and to infer model parameters from data even in cases where the LNA provides poor approximations of the CME. Furthermore, the computational cost of filtering-based iterative likelihood evaluation scales advantageously in the number of measurement times and different input perturbations and is thus ideally suited for data obtained from modern experimental platforms. To demonstrate the practical usefulness of these results, we perform an experiment in which single cells, equipped with an optogenetic gene expression system, are exposed to various different light-input sequences and measured at several hundred time points and use parameter inference based on iterative likelihood evaluation to parameterise a stochastic model of the system.     

A mixed model for two-state Markov processes under panel observation

Many chronic medical conditions can be meaningfully characterized in terms of a two-state stochastic process. Here we consider the problem in which subjects make transitions among two such states in continuous time but are only observed at discrete, irregularly spaced time points that are possibly unique to each subject. Data arising from such an observation scheme are called panel data, and methods for related analyses are typically based on Markov assumptions. The purpose of this article is to present a conditionally Markov model that accommodates subject-to-subject variation in the model parameters by the introduction of random effects. We focus on a particular random effects formulation that generates a closed-form expression for the marginal likelihood. The methodology is illustrated by application to a data set from a parasitic field infection survey.     

Beating irregularity of single pacemaker cells isolated from the rabbit sinoatrial node.

Single pacemaker heart cells discharge irregularly. Data on fluctuations in interbeat interval of single pacemaker cells isolated from the rabbit sinoatrial node are presented. The coefficient of variation of the interbeat interval is quite small, approximately 2%, even though the coefficient of variation of diastolic depolarization rate is approximately 15%. It has been hypothesized that random fluctuations in interbeat interval arise from the stochastic behavior of the membrane ionic channels. To test this hypothesis, we constructed a single channel model of a single pacemaker cell isolated from the rabbit sinoatrial node, i.e., a model into which the stochastic open-close kinetics of the individual membrane ionic channels are incorporated. Single channel conductances as well as single channel open and closed lifetimes are based on experimental data from whole cell and single channel experiments that have been published in the past decade. Fluctuations in action potential parameters of the model cell are compared with those observed experimentally. It is concluded that fluctuations in interbeat interval of single sinoatrial node pacemaker cells indeed are due to the stochastic open-close kinetics of the membrane ionic channels.     

Different genomic relationship matrices for single-step analysis using phenotypic,pedigree and genomic information

Background

The incorporation of genomic coefficients into the numerator relationship matrix allows estimation of breeding values using all phenotypic, pedigree and genomic information simultaneously. In such a single-step procedure, genomic and pedigree-based relationships have to be compatible. As there are many options to create genomic relationships, there is a question of which is optimal and what the effects of deviations from optimality are.

Methods

Data of litter size (total number born per litter) for 338,346 sows were analyzed. Illumina PorcineSNP60 BeadChip genotypes were available for 1,989. Analyses were carried out with the complete data set and with a subset of genotyped animals and three generations pedigree (5,090 animals). A single-trait animal model was used to estimate variance components and breeding values. Genomic relationship matrices were constructed using allele frequencies equal to 0.5 (G05), equal to the average minor allele frequency (GMF), or equal to observed frequencies (GOF). A genomic matrix considering random ascertainment of allele frequencies was also used (GOF*). A normalized matrix (GN) was obtained to have average diagonal coefficients equal to 1. The genomic matrices were combined with the numerator relationship matrix creating H matrices.

Results

In G05 and GMF, both diagonal and off-diagonal elements were on average greater than the pedigree-based coefficients. In GOF and GOF*, the average diagonal elements were smaller than pedigree-based coefficients. The mean of off-diagonal coefficients was zero in GOF and GOF*. Choices of G with average diagonal coefficients different from 1 led to greater estimates of additive variance in the smaller data set. The correlation between EBV and genomic EBV (n = 1,989) were: 0.79 using G05, 0.79 using GMF, 0.78 using GOF, 0.79 using GOF*, and 0.78 using GN. Accuracies calculated by inversion increased with all genomic matrices. The accuracies of genomic-assisted EBV were inflated in all cases except when GN was used.

Conclusions

Parameter estimates may be biased if the genomic relationship coefficients are in a different scale than pedigree-based coefficients. A reasonable scaling may be obtained by using observed allele frequencies and re-scaling the genomic relationship matrix to obtain average diagonal elements of 1.     

Estimation of Parameters for Macroparasite Population Evolution Using Approximate Bayesian Computation

Summary We estimate the parameters of a stochastic process model for a macroparasite population within a host using approximate Bayesian computation (ABC). The immunity of the host is an unobserved model variable and only mature macroparasites at sacrifice of the host are counted. With very limited data, process rates are inferred reasonably precisely. Modeling involves a three variable Markov process for which the observed data likelihood is computationally intractable. ABC methods are particularly useful when the likelihood is analytically or computationally intractable. The ABC algorithm we present is based on sequential Monte Carlo, is adaptive in nature, and overcomes some drawbacks of previous approaches to ABC. The algorithm is validated on a test example involving simulated data from an autologistic model before being used to infer parameters of the Markov process model for experimental data. The fitted model explains the observed extra‐binomial variation in terms of a zero‐one immunity variable, which has a short‐lived presence in the host.     

A model for nematodiasis in New Zealand lambs.

A strategic model is described for the epidemiology of mixed nematode infections in New Zealand lambs. The model successfully reproduces known patterns of parasite epidemiology and production loss in lambs under currently implemented control strategies. The variation in model output during sensitivity analysis was within acceptable limits defined by field data. Model output was most sensitive to variation in parameters affecting survival and migration of the free-living stages and host resistance to infection, suggesting that these factors are most influential in regulating parasite populations. It is intended to use the model to focus research on key aspects of nematode epidemiology and control and, following the incorporation of appropriate genetic mechanisms, anthelmintic resistance.     

Evolution in the real world: stochastic variation and the determinants of fitness in Carlina vulgaris

Empirical studies of life histories often ignore stochastic variation, despite theoretical demonstrations of its potential impact on life-history evolution. Here we use a novel approach to explore the effects of stochastic variation on life-history evolution and estimate the selection pressures operating on the monocarpic perennial Carlina vulgaris, in which flowering may be delayed by up to eight years. The approach is novel in that we use modern theoretical techniques to estimate selection pressures and the fitness landscape from a fully parameterised individual-based model. These approaches take into account temporal variation in demographic rates and density dependence. Analysis of 16 years' data revealed significant temporal variation in growth, mortality, and recruitment in our study population. Flowering was strongly size dependent and, unusually for such a species, also age dependent. Individual-based models of the flowering strategy, parameterized using field data, consistently underestimated the size at flowering, when temporal variation in demographic rates was ignored. In contrast, models that incorporated temporal variation in growth, mortality, and recruitment predicted sizes at flowering not significantly different from those observed in the field. Temporal variation in mortality, which had the largest effect on the flowering strategy, selected for increased size at flowering. An analytical approximation is presented to explain this result, extending the "1-year look-ahead criterion" presented in Rees et al. (2000). A fitness landscape generated by following the fate of rare mutant invaders with a broad range of alternative flowering strategies demonstrated that the observed parameters were adaptive. However, the fitness landscape reveals that approximately equal fitness is achieved by a broad range of strategies, providing a mechanism for the maintenance of genetic variation. To understand how the different parameters that defined our models determine the fitness of rare mutants, we numerically estimated the elasticities and sensitivities of mutant fitness. This demonstrated strong selection on a number of the parameters. Elasticities and sensitivities estimated in constant and random environments were significantly positively correlated, and both were negatively related to the standard error of the parameter. This last result is surprising and, we argue, reflects the genetic and phenotypic responses to selection.