Bayesian Computing Research Articles

Complex hybridization patterns in European pond turtles (Emys orbicularis) in the Pyrenean Region

Hybrid zones are natural laboratories allowing insights in genetic processes like lineage diversification, speciation and introgression. Using large sampling, 15 microsatellite loci and a mitochondrial marker, we examined the Pyrenean contact zone of three pond turtle taxa (Emys orbicularis orbicularis, E. o. galloitalica, E. o. occidentalis). The Pyrenees are a biogeographically important region separating many lineages endemic to the Iberian Peninsula from their Western European counterparts. We discovered limited admixture, reflecting a complex biogeographic scenario. Simulations using Approximate Bayesian Computing supported that E. o. orbicularis invaded the Iberian Peninsula in the Holocene, circumventing the Pyrenees along the Mediterranean coast, and hybridized in the northern peninsula with the local coastal subspecies galloitalica, and to a lesser extent, with occidentalis. While E. o. occidentalis, and in particular E. o. orbicularis, expanded their ranges considerably during Holocene warming, E. o. galloitalica remained largely confined to its former Iberian refuge. Admixture among the three taxa is surprisingly low, and a future taxonomic investigation that includes the unstudied subspecies of E. orbicularis from North Africa, Eastern Europe and Asia has to determine whether their current status properly reflects their evolutionary divergence or whether certain taxa should be regarded as full species.

Regional population structure of the endangered Bridle Shiner (Notropis bifrenatus)

In the last 100 years, the Bridle Shiner has declined over significant parts of its range. We used mitochondrial and nuclear microsatellite markers to investigate population structure of this species. Tissue samples were collected from populations in six drainages in PA, NJ, NY and CT. One predominant haplotype was observed in the Delaware, Housatonic, Passaic, and Raritan specimens. Specimens from the Hudson and St. Lawrence drainages had a separate unique haplotype. Microsatellite loci revealed low within-population genetic diversity and pairwise population comparisons of allelic divergence showed significant genetic differences among most drainages. Population structure analysis of microsatellite loci recover the same population clusters inferred using mtDNA. Within the Delaware drainage, there was significant differentiation among areas separated by 128 km. Several alternate scenarios of population divergence and population bottlenecks were investigated using approximate Bayesian computing. These supported a scenario with a bottleneck occurring in the ancestral population of Bridle Shiner followed by recent divergence of Northern and Southern drainages.

Improving the Accuracy of Planet Occurrence Rates from Kepler Using Approximate Bayesian Computation

Abstract We present a new framework to characterize the occurrence rates of planet candidates identified by Kepler based on hierarchical Bayesian modeling, approximate Bayesian computing (ABC), and sequential importance sampling. For this study, we adopt a simple 2D grid in planet radius and orbital period as our model and apply our algorithm to estimate occurrence rates for Q1–Q16 planet candidates orbiting solar-type stars. We arrive at significantly increased planet occurrence rates for small planet candidates (R p < 1.25 R ⊕) at larger orbital periods (P > 80 day) compared to the rates estimated by the more common inverse detection efficiency method (IDEM). Our improved methodology estimates that the occurrence rate density of small planet candidates in the habitable zone of solar-type stars is per factor of 2 in planet radius and orbital period. Additionally, we observe a local minimum in the occurrence rate for strong planet candidates marginalized over orbital period between 1.5 and 2 R ⊕ that is consistent with previous studies. For future improvements, the forward modeling approach of ABC is ideally suited to incorporating multiple populations, such as planets, astrophysical false positives, and pipeline false alarms, to provide accurate planet occurrence rates and uncertainties. Furthermore, ABC provides a practical statistical framework for answering complex questions (e.g., frequency of different planetary architectures) and providing sound uncertainties, even in the face of complex selection effects, observational biases, and follow-up strategies. In summary, ABC offers a powerful tool for accurately characterizing a wide variety of astrophysical populations.

Assessing the variability in transmission of bovine tuberculosis within Spanish cattle herds

In Spain, despite years of efforts to eradicate bovine tuberculosis (bTB), the disease is still endemic, with some areas of high prevalence. In this context, the surveillance and control plans may need to be re-evaluated, and understanding the dynamics of bTB spread within Spanish herds may help to develop new strategies for reducing the time for detection of infected herds and for the elimination of bTB from the herds already infected. Here, we developed a compartmental stochastic model to simulate bTB within-herd transmission, fed it with epidemiological data from 22 herds (obtained from a previous work) and carried out parameter inference using Approximate Bayesian Computing methods We also estimated the “Within-herd transmission potential Number” (Rh), i.e. the average number of secondary cases generated by a single animal infected introduced into a totally susceptible herd, considering different scenarios depending on the frequency of controls. The median global values obtained for the transmission parameters were: for the transmission coefficient (β), 0.014 newly infected animals per infectious individual per day (i.e. 5.2 per year), for the rate at which infected individuals become infectious (α), 0.01 per day (equivalent to a latent period of 97 days), and for the rate at which infected individuals become reactive to the skin test (α1), 0.08 per day (equivalent to a period of 12 days for an infected animal to become reactive). However, the results also evidenced a great variability in the estimates of those parameters (in particular β and α) among the 22 herds. Considering a 6-month interval between tests, the mean Rh was 0.23, increasing to 0.82 with an interval of 1 year, and to 2.01 and 3.47 with testing intervals of 2 and 4 years, respectively.

Projection of Climate change impacts on the Zayandeh-Rud River Basin’s temperature by using a Bayesian Soft Computing Technique

Projection of Climate change impacts on the Zayandeh-Rud River Basin’s temperature by using a Bayesian Soft Computing Technique

A Queue Model for Reliable Forecasting of Future CPU Consumption

Statistical queuing models are popular to analyze a computer systems ability to process different types requests. A common strategy is to run stress tests by sending artificial requests to the system. The rate and sizes of the requests are varied to investigate the impact on the computer system. A challenge with such an approach is that we do not know if the artificial requests processes are realistic when the system is applied in a real setting. Motivated by this challenge, we develop a method to estimate the properties of the underlying request processes to the computer system when the system is used in a real setting. In particular we look at the problem of recovering the request patterns to a CPU processor. It turns out that this is a challenging statistical estimation problem since we do not observe the request process (rate and size of the requests) to the CPU directly, but only the average CPU usage in disjoint time intervals. In this paper we demonstrate that, quite astonishingly, we are able to recover the properties of the underlying request process (rate and sizes of the requests) by using specially constructed statistics of the observed CPU data and apply a recently developed statistical framework called Approximate Bayesian Computing. Further we apply the model to forecast future CPU consumption. Our results show that the model forecast future CPU consumption with less error than both the hidden Markov model (HMM) in (Hammer et al. 2016) and an ARIMA model. Another good property of the queue model is that we can forecast the instantaneous CPU consumption at any time point in the future, while the HMM in (Hammer et al. 2016) and time series models are limited to only forecasting the average CPU consumption in disjoint time intervals.

Lensing substructure quantification in RXJ1131-1231: a 2 keV lower bound on dark matter thermal relic mass

We study the substructure content of the strong gravitational lens RXJ1131-1231 through a forward modelling approach that relies on generating an extensive suite of realistic simulations. We use a semi-analytic merger tree prescription that allows us to stochastically generate substructure populations whose properties depend on the dark matter particle mass. These synthetic halos are then used as lenses to produce realistic mock images that have the same features, e.g. luminous arcs, quasar positions, instrumental noise and PSF, as the data. We then analyse the data and the simulations in the same way with summary statistics that are sensitive to the signal being targeted and are able to constrain models of dark matter statistically using Approximate Bayesian Computing (ABC) techniques. {In this work, we focus on the thermal relic mass estimate and fix the semi-analytic descriptions of the substructure evolution based on recent literature.} We are able, based on the HST data for RXJ1131-1231, to rule out a warm dark matter thermal relic mass below 2 keV at the 2σ confidence level.

Non parametric Bayesian analysis of the two-sample problem with censoring

ABSTRACTTesting for differences between two groups is a fundamental problem in statistics, and due to developments in Bayesian non parametrics and semiparametrics there has been renewed interest in approaches to this problem. Here we describe a new approach to developing such tests and introduce a class of such tests that take advantage of developments in Bayesian non parametric computing. This class of tests uses the connection between the Dirichlet process (DP) prior and the Wilcoxon rank sum test but extends this idea to the DP mixture prior. Here tests are developed that have appropriate frequentist sampling procedures for large samples but have the potential to outperform the usual frequentist tests. Extensions to interval and right censoring are considered and an application to a high-dimensional data set obtained from an RNA-Seq investigation demonstrates the practical utility of the method.

Bayesian Computing with INLA: A Review

The key operation in Bayesian inference is to compute high-dimensional integrals. An old approximate technique is the Laplace method or approximation, which dates back to Pierre-Simon Laplace (1774). This simple idea approximates the integrand with a second-order Taylor expansion around the mode and computes the integral analytically. By developing a nested version of this classical idea, combined with modern numerical techniques for sparse matrices, we obtain the approach of integrated nested Laplace approximations (INLA) to do approximate Bayesian inference for latent Gaussian models (LGMs). LGMs represent an important model abstraction for Bayesian inference and include a large proportion of the statistical models used today. In this review, we discuss the reasons for the success of the INLA approach, the R-INLA package, why it is so accurate, why the approximations are very quick to compute, and why LGMs make such a useful concept for Bayesian computing.

Fast and accurate Bayesian model criticism and conflict diagnostics using R‐INLA

Bayesian hierarchical models are increasingly popular for realistic modelling and analysis of complex data. This trend is accompanied by the need for flexible, general and computationally efficient methods for model criticism and conflict detection. Usually, a Bayesian hierarchical model incorporates a grouping of the individual data points, as, for example, with individuals in repeated measurement data. In such cases, the following question arises: Are any of the groups “outliers,” or in conflict with the remaining groups? Existing general approaches aiming to answer such questions tend to be extremely computationally demanding when model fitting is based on Markov chain Monte Carlo. We show how group‐level model criticism and conflict detection can be carried out quickly and accurately through integrated nested Laplace approximations (INLA). The new method is implemented as a part of the open‐source R‐INLA package for Bayesian computing (http://r-inla.org). Copyright © 2017 John Wiley & Sons, Ltd.

Beyond the Model Limit: Parameter Inference Across Scales

In multiscale studies of emergent phenomena, a common approach is to define a microscopic scale generative model and subsequently, by passing to a diffusion limit, to derive a mesoscopic or macroscopic model via a homogenization argument. Microscopic models are often inherently stochastic, while the diffusion limit model may represent distributions or mean field quantities. A canonical example of this type of limiting process is the classical diffusion equation as a limit of Brownian random motion. Being typically based on ordinary or partial differential equations formalism, large scale models are more amenable to inverse problems, where the model parameters are estimated from indirect observations. Due to the intrinsic ill-posed nature of inverse problems, the estimation of the parameters of large scale models may be itself a significant challenge. In some multiscale investigations, the quantities of primary interest are not those characterizing the large scale model, but rather the parameters of the associated microscopic model. The latter are sometimes related to our understanding of the processes that are described by the generative models. In many cases, there is no direct connection between the small scale parameters and the large scale parameters describing the aggregate effects of the underlying processes. Moreover, estimating the unknown parameters of a microscopic, stochastic model from the large scale observations poses a number of challenges, and requires the development of new methodological approaches. In the Bayesian statistical framework, stochastic forward models often lead to intractable likelihood densities, making standard statistical computations a challenge. Such models are addressed by approximate Bayesian computing (ABC). In this paper we propose an alternative approach combining Bayesian statistics and particle techniques for inference beyond the limits of model scales. The viability of the approach is demonstrated on three computed examples describing diffusion, chemical kinetics, and cell culture assay, and a connection with existing ABC methods is discussed.

Quick and energy-efficient Bayesian computing of binocular disparity using stochastic digital signals

Reconstruction of the tridimensional geometry of a visual scene using the binocular disparity information is an important issue in computer vision and mobile robotics, which can be formulated as a Bayesian inference problem. However, computation of the full disparity distribution with an advanced Bayesian model is usually an intractable problem, and proves computationally challenging even with a simple model. In this paper, we show how probabilistic hardware using distributed memory and alternate representation of data as stochastic bitstreams can solve that problem with high performance and energy efficiency. We put forward a way to express discrete probability distributions using stochastic data representations and perform Bayesian fusion using those representations, and show how that approach can be applied to diparity computation. We evaluate the system using a simulated stochastic implementation and discuss possible hardware implementations of such architectures and their potential for sensorimotor processing and robotics.

Applying ecological model evaludation: Lessons learned with the forest dynamics model Samsara2

Ecological models are increasingly used as decision-making tools and their reliability is becoming a key issue. At the same time, the sophistication of techniques for model development and analysis has given rise to a relative compartmentalization of model building and evaluation tasks. Several guidelines invite ecological modelers to follow an organized sequence of development and analysis steps and have coined the term “evaludation” for this process. The objective of this paper is to assess the feasibility and the value of a structured evaludation process, based on the working example of the Samsara2 model, a spatially explicit individual-based forest dynamics model. We implemented the six steps of model design, process level calibration, qualitative evaluation, quantitative evaluation, global sensitivity analysis, and partial recalibration using approximate Bayesian computing. We then evaluated how the evaludation process revealed model strengths and weaknesses, specified the model's conditions of use, clarified how the model works, and provided insights into forest ecosystem functioning. Finally, the efficiency/cost ratio of the process and future improvements are discussed.

A Conversation with Alan Gelfand

Alan E. Gelfand was born April 17, 1945, in the Bronx, New York. He attended public grade schools and did his undergraduate work at what was then called City College of New York (CCNY, now CUNY), excelling at mathematics. He then surprised and saddened his mother by going all the way across the country to Stanford to graduate school, where he completed his dissertation in 1969 under the direction of Professor Herbert Solomon, making him an academic grandson of Herman Rubin and Harold Hotelling. Alan then accepted a faculty position at the University of Connecticut (UConn) where he was promoted to tenured associate professor in 1975 and to full professor in 1980. A few years later he became interested in decision theory, then empirical Bayes, which eventually led to the publication of Gelfand and Smith [J. Amer. Statist. Assoc. 85 (1990) 398-409], the paper that introduced the Gibbs sampler to most statisticians and revolutionized Bayesian computing. In the mid-1990s, Alan's interests turned strongly to spatial statistics, leading to fundamental contributions in spatially-varying coefficient models, coregionalization, and spatial boundary analysis (wombling). He spent 33 years on the faculty at UConn, retiring in 2002 to become the James B. Duke Professor of Statistics and Decision Sciences at Duke University, serving as chair from 2007-2012. At Duke, he has continued his work in spatial methodology while increasing his impact in the environmental sciences. To date, he has published over 260 papers and 6 books; he has also supervised 36 Ph.D. dissertations and 10 postdocs. This interview was done just prior to a conference of his family, academic descendants, and colleagues to celebrate his 70th birthday and his contributions to statistics which took place on April 19-22, 2015 at Duke University.

On selection of statistics for approximate Bayesian computing (or the method of simulated moments)

A cross validation method for selection of statistics for Approximate Bayesian Computing, and for related estimation methods such as the Method of Simulated Moments, is presented. The method uses simulated annealing to minimize the cross validation criterion over a combinatorial search space that may contain an extremely large number of elements. A first simple example, for which optimal statistics are known from theory, shows that the method is able to select these optimal statistics out of a large set of candidate statistics. A second example of selection of statistics for a stochastic volatility model illustrates the method in a more complex case. Code to replicate the results, or to use the method for other applications, is provided at http://www.runmycode.org/companion/view/1116.

Information criteria and approximate Bayesian computing for agent-based modelling in ecology: New tools to infer on Individual-level processes

Agent or individual based models (ABMs) are computer simulation tools to represent complex systems where the unique autonomous individuals constituting the system interact, adapt or evolve. In ecology, ABMs have been applied to analyse the influence of individual variability, spatial interactions, development, behaviour, learning or genetic variability on group, population or community dynamics [1]. During the last 20 years, the approach of pattern-oriented modelling (POM) has been helping in making the use of ABMs more robust and reliable [2]. Two main elements of POM are a) ``inverse modelling'' to reduce parameter uncertainty and b) ``strong inferences'' to test how well alternative model versions changing in certain individual-level processes are at explaining multiple patterns observed at the population and/or community level. The latter has many parallels to model selection procedures in statistical modelling ( ``mutlimodel inference'' [3]). The inverse modelling borrows also concept and tools from statistics to fit unknown or uncertain parameters to real-world data. With the development of approximate Bayesian computing (ABC) in ecology [4] some propositions were made to use ABC in ABMs [5]. Bridging the two main elements of POM, an information criterion was proposed [6] that borrows concept of ABC such as Markov Chains Monte Carlo simulations without likelihood functions and multimodel inference in Bayesian context. This criterion (POMIC for pattern-oriented modelling information criterion) was applied in the original proposition paper to forest growth simulations to infer on individual tree growth function. Recently, it was applied to locust to infer on individual behaviour at hatching [7]. These kinds of cryptic processes (tree individual growth or animals' behaviour during first hours) can typically be analysed following few isolated individuals in natural or experimental conditions. However, in the context of inter-individual interactions, experimental designs to infer on how the interactions modify individual responses are hard, if not impossible to elaborate. Using population or group level data to fit ABMs and infer on the potential individual process with different ABMs version is then a way of overpassing these problems. However, the question of ABMs complexity (essential in statistical model selection procedures) needs to be addressed. The POMIC proposition includes a complexity measure and makes it therefore an interesting tool to infer on cryptic individual-level processes of complex systems. (Texte integral)

Demographic histories and genetic diversities of Fennoscandian marine and landlocked ringed seal subspecies.

Island populations are on average smaller, genetically less diverse, and at a higher risk to go extinct than mainland populations. Low genetic diversity may elevate extinction probability, but the genetic component of the risk can be affected by the mode of diversity loss, which, in turn, is connected to the demographic history of the population. Here, we examined the history of genetic erosion in three Fennoscandian ringed seal subspecies, of which one inhabits the Baltic Sea ‘mainland’ and two the ‘aquatic islands’ composed of Lake Saimaa in Finland and Lake Ladoga in Russia. Both lakes were colonized by marine seals after their formation c. 9500 years ago, but Lake Ladoga is larger and more contiguous than Lake Saimaa. All three populations suffered dramatic declines during the 20th century, but the bottleneck was particularly severe in Lake Saimaa. Data from 17 microsatellite loci and mitochondrial control-region sequences show that Saimaa ringed seals have lost most of the genetic diversity present in their Baltic ancestors, while the Ladoga population has experienced only minor reductions. Using Approximate Bayesian computing analyses, we show that the genetic uniformity of the Saimaa subspecies derives from an extended founder event and subsequent slow erosion, rather than from the recent bottleneck. This suggests that the population has persisted for nearly 10,000 years despite having low genetic variation. The relatively high diversity of the Ladoga population appears to result from a high number of initial colonizers and a high post-colonization population size, but possibly also by a shorter isolation period and/or occasional gene flow from the Baltic Sea.

Weather Derivative Risk Measures for Extreme Events

We consider pricing weather derivatives for use as protection against weather extremes by using max-stable processes to estimate risk measures. These derivatives are not currently traded on any open markets, but their use could help some institutions manage weather risks from extreme events. The central challenge is to model the dependence of payments, which increases the risk of holding multiple weather derivatives. The method described utilizes results from spatial statistics and extreme value theory to first model extremes in the weather as a max-stable process, and then simulate payments for a general collection of weather derivatives. As the joint likelihood function for max-stable processes is unavailable, we use two approaches: The first is based on the composite likelihood, and the second is based on approximate Bayesian computing (ABC). Both capture the spatial dependence of payments. To incorporate parameter uncertainty into the pricing model, we use bootstrapping with the composite likelihood approach, while the ABC method naturally incorporates parameter uncertainty. We show that the additional risk from the spatial dependence of payments can be quite substantial, and that the methods discussed can compute standard actuarial risk measures in both a frequentist and Bayesian setting.

Continued fraction enhancement of Bayesian computing

The agéd number theoretic concept of continued fractions can enhance certain Bayesian computations. The crux of this claim is due to continued fraction representations of numerically challenging special function ratios that arise in Bayesian computing. Continued fraction approximation via Lentz's Algorithm often leads to efficient and stable computation of such quantities. Copyright © 2012 John Wiley & Sons, Ltd.

Multi‐scale modeling of GMP differentiation based on single‐cell genealogies

Hematopoiesis is often pictured as a hierarchy of branching decisions, giving rise to all mature blood cell types from stepwise differentiation of a single cell, the hematopoietic stem cell. Various aspects of this process have been modeled using various experimental and theoretical techniques on different scales. Here we integrate the more common population-based approach with a single-cell resolved molecular differentiation model to study the possibility of inferring mechanistic knowledge of the differentiation process. We focus on a sub-module of hematopoiesis: differentiation of granulocyte-monocyte progenitors (GMPs) to granulocytes or monocytes. Within a branching process model, we infer the differentiation probability of GMPs from the experimentally quantified heterogeneity of colony assays under permissive conditions where both granulocytes and monocytes can emerge. We compare the predictions with the differentiation probability in genealogies determined from single-cell time-lapse microscopy. In contrast to the branching process model, we found that the differentiation probability as determined by differentiation marker onset increases with the generation of the cell within the genealogy. To study this feature from a molecular perspective, we established a stochastic toggle switch model, in which the intrinsic lineage decision is executed using two antagonistic transcription factors. We identified parameter regimes that allow for both time-dependent and time-independent differentiation probabilities. Finally, we infer parameters for which the model matches experimentally observed differentiation probabilities via approximate Bayesian computing. These parameters suggest different timescales in the dynamics of granulocyte and monocyte differentiation. Thus we provide a multi-scale picture of cell differentiation in murine GMPs, and illustrate the need for single-cell time-resolved observations of cellular decisions.