Bayesian Computation Methods Research Articles

Bayesian 3D X-ray Computed Tomography with a Hierarchical Prior Model for Sparsity in Haar Transform Domain.

In this paper, a hierarchical prior model based on the Haar transformation and an appropriate Bayesian computational method for X-ray CT reconstruction are presented. Given the piece-wise continuous property of the object, a multilevel Haar transformation is used to associate a sparse representation for the object. The sparse structure is enforced via a generalized Student-t distribution (), expressed as the marginal of a normal-inverse Gamma distribution. The proposed model and corresponding algorithm are designed to adapt to specific 3D data sizes and to be used in both medical and industrial Non-Destructive Testing (NDT) applications. In the proposed Bayesian method, a hierarchical structured prior model is proposed, and the parameters are iteratively estimated. The initialization of the iterative algorithm uses the parameters of the prior distributions. A novel strategy for the initialization is presented and proven experimentally. We compare the proposed method with two state-of-the-art approaches, showing that our method has better reconstruction performance when fewer projections are considered and when projections are acquired from limited angles.

Species delimitation in the presence of strong incomplete lineage sorting and hybridization: Lessons from Ophioderma (Ophiuroidea: Echinodermata)

Accurate species delimitation is essential to properly assess biodiversity, but also for management and conservation purposes. Yet, it is not always trivial to accurately define species boundaries in closely related species due to incomplete lineage sorting. Additional difficulties may be caused by hybridization, now evidenced as a frequent phenomenon. The brittle star cryptic species complex Ophioderma longicauda encompasses six mitochondrial lineages, including broadcast spawners and internal brooders, yet the actual species boundaries are unknown. Here, we combined three methods to delimit species in the Ophioderma longicauda complex and to infer its divergence history: (i) unsupervised species discovery based on multilocus genotypes; (ii) divergence time estimation using the multi-species coalescent; (iii) divergence scenario testing (including gene flow) using Approximate Bayesian Computation (ABC) methods. 30 sequence markers (transcriptome-based, mitochondrial or non-coding) for 89 O. longicauda and outgroup individuals were used. First, multivariate analyses revealed six genetic clusters, which globally corresponded to the mitochondrial lineages, yet with many exceptions, suggesting ancient hybridization events and challenging traditional mitochondrial barcoding approaches. Second, multi-species coalescent-based analyses confirmed the occurrence of six species and provided divergence time estimates, but the sole use of this method failed to accurately delimit species, highlighting the power of multilocus genotype clustering to delimit recently diverged species. Finally, Approximate Bayesian Computation showed that the most likely scenario involves hybridization between brooders and broadcasters. Our study shows that despite strong incomplete lineage sorting and past hybridization, accurate species delimitation in Ophioderma was possible using a combination of complementary methods. We propose that these methods, especially multilocus genotype clustering, may be useful to resolve other complex speciation histories.

A new POD-based approximate bayesian computation method to identify parameters for formed AHSS

Forming has significant influence on the material properties which are critical for the accuracy of simulations. However, at present, most of the existing works took the uniform raw mechanical properties of plate as the one after forming, which inevitably resulted in the inaccuracy simulation. In this study, in order to figure out the material properties of formed DP590, comprehensive experiments based on micro-indentation tests and digital image correlation (DIC) are presented. The versatile approximate Bayesian computation (ABC) inverse method is utilized to identify material parameters. Due to its difficulty of selection of summary statistics, a new flexible and feasible proper orthogonal decomposition (POD) based ABC framework is proposed. The POD can project the high dimensional observations into a much lower dimensional coefficient vector. Thus, we can make use of comparisons between simulated and observed coefficient vector of POD to determine the material parameters. With the POD, the ABC can circumvent the complex selection of summary statistics. In the computational process, the Neural Network (NN) is used to construct the nonlinear connection between material properties and its POD coefficient vector. Thus, the coefficient vector with given material properties can be determined efficiently with the well-constructed NN instead of time-consuming finite element (FE) simulations.

Bayesian Computational Methods for Sampling from the Posterior Distribution of a Bivariate Survival Model, Based on AMH Copula in the Presence of Right-Censored Data.

In this paper, we study the performance of Bayesian computational methods to estimate the parameters of a bivariate survival model based on the Ali–Mikhail–Haq copula with marginal distributions given by Weibull distributions. The estimation procedure was based on Monte Carlo Markov Chain (MCMC) algorithms. We present three version of the Metropolis–Hastings algorithm: Independent Metropolis–Hastings (IMH), Random Walk Metropolis (RWM) and Metropolis–Hastings with a natural-candidate generating density (MH). Since the creation of a good candidate generating density in IMH and RWM may be difficult, we also describe how to update a parameter of interest using the slice sampling (SS) method. A simulation study was carried out to compare the performances of the IMH, RWM and SS. A comparison was made using the sample root mean square error as an indicator of performance. Results obtained from the simulations show that the SS algorithm is an effective alternative to the IMH and RWM methods when simulating values from the posterior distribution, especially for small sample sizes. We also applied these methods to a real data set.

Variability Timescale and Spectral Index of Sgr A* in the Near Infrared: Approximate Bayesian Computation Analysis of the Variability of the Closest Supermassive Black Hole.

Sagittarius A* (Sgr A*) is the variable radio, near-infrared (NIR), and X-ray source associated with accretion onto the Galactic center black hole. We present an analysis of the most comprehensive NIR variability data set of Sgr A* to date: eight 24 hr epochs of continuous monitoring of Sgr A* at 4.5 μm with the IRAC instrument on the Spitzer Space Telescope, 93 epochs of 2.18 μm data from Naos Conica at the Very Large Telescope, and 30 epochs of 2.12 μm data from the NIRC2 camera at the Keck Observatory, in total 94,929 measurements. A new approximate Bayesian computation method for fitting the first-order structure function extracts information beyond current fast Fourier transformation (FFT) methods of power spectral density (PSD) estimation. With a combined fit of the data of all three observatories, the characteristic coherence timescale of Sgr A* is minutes (90% credible interval). The PSD has no detectable features on timescales down to 8.5 minutes (95% credible level), which is the ISCO orbital frequency for a dimensionless spin parameter a = 0.92. One light curve measured simultaneously at 2.12 and 4.5 μm during a low flux-density phase gave a spectral index α s = 1.6 ± 0.1 . This value implies that the Sgr A* NIR color becomes bluer during higher flux-density phases. The probability densities of flux densities of the combined data sets are best fit by log-normal distributions. Based on these distributions, the Sgr A* spectral energy distribution is consistent with synchrotron radiation from a non-thermal electron population from below 20 GHz through the NIR.

Estimating the joint distribution of rate parameters across multiple reactions in the absence of experimental data

Abstract A procedure for determining the joint uncertainty of Arrhenius parameters across multiple combustion reactions of interest is demonstrated. This approach is capable of constructing the joint distribution of the Arrhenius parameters arising from the uncertain measurements performed in specific target experiments without having direct access to the underlying experimental data. The method involves constructing an ensemble of hypothetical data sets with summary statistics consistent with the available information reported by the experimentalists, followed by a fitting procedure that learns the structure of the joint parameter density across reactions using this consistent hypothetical data as evidence. The procedure is formalized in a Bayesian statistical framework, employing maximum-entropy and approximate Bayesian computation methods and utilizing efficient Markov chain Monte Carlo techniques to explore data and parameter spaces in a nested algorithm. We demonstrate the application of the method in the context of experiments designed to measure the rates of selected chain reactions in the H2-O2 system and highlight the utility of this approach for revealing the critical correlations between the parameters within a single reaction and across reactions, as well as for maximizing consistency when utilizing rate parameter information in predictive combustion modeling of systems of interest.

Considerations for monitoring population trends of colonial waterbirds using the effective number of breeders and census estimates.

Detecting trends in population size fluctuations is a major focus in ecology, evolution, and conservation biology. Populations of colonial waterbirds have been monitored using demographic approaches to determine annual census size (N a). We propose the addition of genetic estimates of the effective number of breeders (N b) as indirect measures of the risk of loss of genetic diversity to improve the evaluation of demographics and increase the accuracy of trend estimates in breeding colonies. Here, we investigated which methods of the estimation of N b are more precise under conditions of moderate genetic diversity, limited sample sizes and few microsatellite loci, as often occurs with natural populations. We used the wood stork as a model species and we offered a workflow that researchers can follow for monitoring bird breeding colonies. Our approach started with simulations using five estimators of N b and the theoretical results were validated with empirical data collected from breeding colonies settled in the Brazilian Pantanal wetland. In parallel, we estimated census size using a corrected method based on counting active nests. Both in simulations and in natural populations, the approximate Bayesian computation (ABC) and sibship assignment (SA) methods yielded more precise estimates than the linkage disequilibrium, heterozygosity excess, and molecular coancestry methods. In particular, the ABC method performed best with few loci and small sample sizes, while the other estimators required larger sample sizes and at least 13 loci to not underestimate N b. Moreover, according to our N b/N a estimates (values were often ≤0.1), the wood stork colonies evaluated could be facing the loss of genetic diversity. We demonstrate that the combination of genetic and census estimates is a useful approach for monitoring natural breeding bird populations. This methodology has been recommended for populations of rare species or with a known history of population decline to support conservation efforts.

Full mitochondrial genome sequences reveal new insights about post-glacial expansion and regional phylogeographic structure in the Atlantic silverside (Menidia menidia)

One of the biggest challenges in population genetics of marine species is to uncover subtle phylogeographic patterns masked by large effective population sizes and/or high gene flow. In this paper, we use 189 full mitochondrial genome sequences from the Atlantic silverside (Menidia menidia), obtained during whole nuclear genome re-sequencing, to address this challenge. With this approach, we were able to provide a high-resolution analysis of the demographic history and current genetic structure of the species. We clearly defined a regional genetic structure that is stronger than previously thought. This structure groups the species into three regional subdivisions: (1) south of Cape Cod; (2) Gulf of Maine; and (3) Gulf of St. Lawrence. Among the regional subdivisions, our data reveal that the two northern groups show the greatest divergence despite their adjacency along the latitudinal gradient, while genetic homogeneity within the southern subdivision suggests connectivity throughout its broad geographical distribution. Furthermore, using approximate Bayesian computation (ABC) methods, we inferred that both northern populations are the result of independent colonization events from the south after the Last Glacial Maximum (LGM). Our analyses indicate that at least one of the northern populations has received two waves of colonization, one timed immediately after the LGM and the other timed after the end of the Younger Dryas glacial re-advance. Finally, we found one locus potentially under positive selection in the mitochondrial genome. The results of this study illustrate the power of full mitochondrial genome sequencing in phylogeographic research, and because the Atlantic silverside is known for its clinal phenotypic variation throughout its range, our findings have important implications for the study of local adaptation.

Massively parallel approximate Bayesian computation for estimating nanoparticle diffusion coefficients, sizes and concentrations using confocal laser scanning microscopy.

We implement a massively parallel population Monte Carlo approximate Bayesian computation (PMC-ABC) method for estimating diffusion coefficients, sizes and concentrations of diffusing nanoparticles in liquid suspension using confocal laser scanning microscopy and particle tracking. The method is based on the joint probability distribution of diffusion coefficients and the time spent by a particle inside a detection region where particles are tracked. We present freely available central processing unit (CPU) and graphics processing unit (GPU) versions of the analysis software, and we apply the method to characterize mono- and bidisperse samples of fluorescent polystyrene beads.

Can Machines Learn Respiratory Virus Epidemiology?: A Comparative Study of Likelihood-Free Methods for the Estimation of Epidemiological Dynamics.

To estimate and predict the transmission dynamics of respiratory viruses, the estimation of the basic reproduction number, R0, is essential. Recently, approximate Bayesian computation methods have been used as likelihood free methods to estimate epidemiological model parameters, particularly R0. In this paper, we explore various machine learning approaches, the multi-layer perceptron, convolutional neural network, and long-short term memory, to learn and estimate the parameters. Further, we compare the accuracy of the estimates and time requirements for machine learning and the approximate Bayesian computation methods on both simulated and real-world epidemiological data from outbreaks of influenza A(H1N1)pdm09, mumps, and measles. We find that the machine learning approaches can be verified and tested faster than the approximate Bayesian computation method, but that the approximate Bayesian computation method is more robust across different datasets.

Approximation of Bayesian Predictive p-Values with Regression ABC

In the Bayesian framework a standard approach to model criticism is to compare some function of the observed data to a reference predictive distribution. The result of the comparison can be summarized in the form of a p-value, and computation of some kinds of Bayesian predictive p-values can be challenging. The use of regression adjustment approximate Bayesian computation (ABC) methods is explored for this task. Two problems are considered. The first is approximation of distributions of prior predictive p-values for the purpose of choosing weakly informative priors in the case where the model checking statistic is expensive to compute. Here the computation is difficult because of the need to repeatedly sample from a prior predictive distribution for different values of a prior hyperparameter. The second problem considered is the calibration of posterior predictive p-values so that they are uniformly distributed under some reference distribution for the data. Computation is difficult because the calibration process requires repeated approximation of the posterior for different data sets under the reference distribution. In both these problems we argue that high accuracy in the computations is not required, which makes fast approximations such as regression adjustment ABC very useful. We illustrate our methods with several examples.

Approximate Bayesian Computation and Simulation-Based Inference for Complex Stochastic Epidemic Models

Approximate Bayesian Computation (ABC) and other simulation-based inference methods are becoming increasingly used for inference in complex systems, due to their relative ease-of-implementation. We briefly review some of the more popular variants of ABC and their application in epidemiology, before using a real-world model of HIV transmission to illustrate some of challenges when applying ABC methods to high-dimensional, computationally intensive models. We then discuss an alternative approach—history matching—that aims to address some of these issues, and conclude with a comparison between these different methodologies.

Learning Summary Statistic for Approximate Bayesian Computation via Deep Neural Network

Approximate Bayesian Computation (ABC) methods are used to approximate posterior distributions in models with unknown or computationally intractable likelihoods. Both the accuracy and computational efficiency of ABC depend on the choice of summary statistic, but outside of special cases where the optimal summary statistics are known, it is unclear which guiding principles can be used to construct effective summary statistics. In this paper we explore the possibility of automating the process of constructing summary statistics by training deep neural networks to predict the parameters from artificially generated data: the resulting summary statistics are approximately posterior means of the parameters. With minimal model-specific tuning, our method constructs summary statistics for the Ising model and the moving-average model, which match or exceed theoretically-motivated summary statistics in terms of the accuracies of the resulting posteriors.

Using Approximate Bayesian Computation by Subset Simulation for Efficient Posterior Assessment of Dynamic State-Space Model Classes

Approximate Bayesian Computation (ABC) methods have gained in popularity over the last decade because they expand the horizon of Bayesian parameter inference methods to the range of models for which an analytical formula for the likelihood function might be difficult, or even impossible, to establish. The majority of the ABC methods rely on the choice of a set of summary statistics to reduce the dimension of the data. However, as has been noted in the ABC literature, the lack of convergence guarantees induced by the absence of a vector of sufficient summary statistics that assures intermodel sufficiency over the set of competing models hinders the use of the usual ABC methods when applied to Bayesian model selection or assessment. In this paper, we present a novel ABC model selection procedure for dynamical systems based on a recently introduced multilevel Markov chain Monte Carlo method, self-regulating ABC-SubSim, and a hierarchical state-space formulation of dynamic models. We show that this formulation makes it possible to independently approximate the model evidence required for assessing the posterior probability of each of the competing models. We also show that ABC-SubSim not only provides an estimate of the model evidence as a simple by-product but also gives the posterior probability of each model as a function of the tolerance level, which allows the ABC model choices made in previous studies to be understood. We illustrate the performance of the proposed framework for ABC model updating and model class selection by applying it to two problems in Bayesian system identification: a single-degree-of-freedom bilinear hysteretic oscillator and a three-story shear building with Masing hysteresis, both of which are subject to a seismic excitation.

An approximate likelihood perspective on ABC methods

We are living in the big data era, as current technologies and networks allow for the easy and routine collection of data sets in different disciplines. Bayesian Statistics offers a flexible modeling approach which is attractive for describing the complexity of these datasets. These models often exhibit a likelihood function which is intractable due to the large sample size, high number of parameters, or functional complexity. Approximate Bayesian Computational (ABC) methods provides likelihood-free methods for performing statistical inferences with Bayesian models defined by intractable likelihood functions. The vastity of the literature on ABC methods created a need to review and relate all ABC approaches so that scientists can more readily understand and apply them for their own work. This article provides a unifying review, general representation, and classification of all ABC methods from the view of approximate likelihood theory. This clarifies how ABC methods can be characterized, related, combined, improved, and applied for future research. Possible future research in ABC is then outlined.

TreeSimGM: Simulating phylogenetic trees under general Bellman-Harris models with lineage-specific shifts of speciation and extinction in R.

Understanding macroevolutionary processes using phylogenetic trees is a challenging and complex process that draws on mathematics, computer science and biology. Given the development of complex mathematical models and the growing computational processing power, simulation tools are becoming increasingly popular.In order to simulate phylogenetic trees, most evolutionary biologists are forced to build their own algorithms or use existing tools built on different platforms and/or as standalone programmes. The absence of a simulation tool accommodating for user‐chosen model specifications limits, amongst others, model testing and pipelining with approximate Bayesian computation methods or other subsequent statistical analysis.We introduce “TreeSimGM,” an r‐package simulation tool for phylogenetic trees under a general Bellman and Harris model. This package allows the user to specify any desired probability distribution for the waiting times until speciation and extinction (e.g. age‐dependent speciation/extinction). Upon speciation, the user can specify whether one descendant species corresponds to the ancestor species inheriting its age or whether both descendant species are new species of age 0. Moreover, it is possible to scale the waiting time to speciation/extinction for newly formed species. Thus, “TreeSimGM” not only allows the user to simulate stochastic phylogenetic trees assuming several popular existing models, such as the Yule model, the constant‐rate birth–death model, and proportional to distinguishable arrangement models, but it also allows the user to formulate new models for exploration. A short explanation of the supported models and a few examples of how to use our package are presented here.As an r‐package, “TreeSimGM” allows flexible and powerful stochastic phylogenetic tree simulations. Moreover, it facilitates the pipelining of outputs or inputs with other functions in r. “TreeSimGM” contributes to the tools available to the r community in the fields of ecology and evolution, is freely available under the GPL‐2 licence and can be downloaded at https://cran.r-project.org/web/packages/TreeSimGM.

Predicting population extinction from early observations of the Lotka–Volterra system

Population extinction is one of the central themes in population biology. We propose a statistical algorithm for long-term prediction of an extinction event in the paradigmatic predator–prey model. The algorithm is based on noisy and sporadic observations of the Lotka–Volterra (LV) system at the early stages of its evolution, when the system is still very far from extinction. There are two stages in the algorithm: first, the unknown parameters (reaction rates) of the LV system are estimated using the Approximate Bayesian Computation method; then an analytical expression for the time-scale of extinction (which involves the estimated parameters) is applied to compute the probability density function of extinction time. The proposed algorithm is validated by numerical simulations for the case of a stochastic LV system specified by the birth–death rate equations. The algorithm can be seen as an initial step in the quest for long-term prediction of rare “catastrophic” events in complex stochastic dynamic systems (epidemics, host-parasite dynamics, enzyme kinetics, dynamic trading, etc.).

Bayesian approach for quantifying silica fume and water reducer impact on creep of concrete

The behaviour of concrete structures is affected by the time-dependent volume changes of concrete caused by sustained loads and moisture loss. Creep deformations have an important effect on long-term integrity of reinforced concrete structures. Incorrect or inaccurate prediction of these volume changes leads to undesirable and harmful consequences for the structures. Knowing that the Eurocode 2 does not take admixtures and additives into consideration in predicting creep deformations, the aim of this study is to apply correction coefficients to the Eurocode 2 creep formulas to take into account the effect of additives and admixtures, especially silica fume and water reducer. A large database from international laboratories and research centres is used in this study. The first step consists of comparing experimental measurements of creep deformations with creep values predicted by the Eurocode 2 formulas. The second step lies in introducing correction coefficients to the creep formulas of Eurocode 2, taking into consideration the percentage of admixtures and additives. These new correction coefficients are introduced in the Eurocode 2 formulas using an approximate Bayesian computation method based on a rejection algorithm.

Colonization history of the western corn rootworm (Diabrotica virgifera virgifera) in North America: insights from random forest ABC using microsatellite data

First described from western Kansas, USA, the western corn rootworm, Diabrotica virgifera virgifera, is one of the worst pests of maize. The species is generally thought to be of Mexican origin and to have incidentally followed the expansion of maize cultivation into North America thousands of years ago. However, this hypothesis has never been investigated formally. In this study, the genetic variability of samples collected throughout North America was analysed at 13 microsatellite marker loci to explore precisely the population genetic structure and colonization history of D. v. virgifera. In particular, we used up-to-date approximate Bayesian computation methods based on random forest algorithms to test a Mexican versus a central-USA origin of the species, and to compare various possible timings of colonization. This analysis provided strong evidence that the origin of D. v. virgifera was southern (Mexico, or even further south). Surprisingly, we also found that the expansion of the species north of its origin was recent—probably not before 1100 years ago—thus indicating it was not directly associated with the early history of maize expansion out of Mexico, a far more ancient event.

Donor-Recipient Identification in Para- and Poly-phyletic Trees Under Alternative HIV-1 Transmission Hypotheses Using Approximate Bayesian Computation.

Diversity of the founding population of Human Immunodeficiency Virus Type 1 (HIV-1) transmissions raises many important biological, clinical, and epidemiological issues. In up to 40% of sexual infections, there is clear evidence for multiple founding variants, which can influence the efficacy of putative prevention methods, and the reconstruction of epidemiologic histories. To infer who-infected-whom, and to compute the probability of alternative transmission scenarios while explicitly taking phylogenetic uncertainty into account, we created an approximate Bayesian computation (ABC) method based on a set of statistics measuring phylogenetic topology, branch lengths, and genetic diversity. We applied our method to a suspected heterosexual transmission case involving three individuals, showing a complex monophyletic-paraphyletic-polyphyletic phylogenetic topology. We detected that seven phylogenetic lineages had been transmitted between two of the individuals based on the available samples, implying that many more unsampled lineages had also been transmitted. Testing whether the lineages had been transmitted at one time or over some length of time suggested that an ongoing superinfection process over several years was most likely. While one individual was found unlinked to the other two, surprisingly, when evaluating two competing epidemiological priors, the donor of the two that did infect each other was not identified by the host root-label, and was also not the primary suspect in that transmission. This highlights that it is important to take epidemiological information into account when analyzing support for one transmission hypothesis over another, as results may be nonintuitive and sensitive to details about sampling dates relative to possible infection dates. Our study provides a formal inference framework to include information on infection and sampling times, and to investigate ancestral node-label states, transmission direction, transmitted genetic diversity, and frequency of transmission.