Reversible Jump Markov Chain Monte Research Articles

A Bayesian hybrid method for the analysis of generalized linear models with missing-not-at-random covariates

Missing data handling is one of the main problems in modelling, particularly if the missingness is of type missing-not-at-random (MNAR) where missingness occurs due to the actual value of the observation. The focus of the current article is generalized linear modelling of fully observed binary response variables depending on at least one MNAR covariate. For the traditional analysis of such models, an individual model for the probability of missingness is assumed and incorporated in the model framework. However, this probability model is untestable, as the missingness of MNAR data depend on their actual values that would have been observed otherwise. In this article, we consider creating a model space that consist of all possible and plausible models for probability of missingness and develop a hybrid method in which a reversible jump Markov chain Monte Carlo (RJMCMC) algorithm is combined with Bayesian Model Averaging (BMA). RJMCMC is adopted to obtain posterior estimates of model parameters as well as probability of each model in the model space. BMA is used to synthesize coefficient estimates from all models in the model space while accounting for model uncertainties. Through a validation study with a simulated data set and a real data application, the performance of the proposed methodology is found to be satisfactory in accuracy and efficiency of estimates.

Sequential model identification with reversible jump ensemble data assimilation method

In data assimilation (DA) schemes, the form representing the processes in the evolution models are pre-determined except some parameters to be estimated. In some applications, such as the contaminant solute transport model and the gas reservoir model, the modes in the equations within the evolution model cannot be predetermined from the outset and may change with the time. We propose a framework of sequential DA method named Reversible Jump Ensemble Filter (RJEnF) to identify the governing modes of the evolution model over time. The main idea is to introduce the Reversible Jump Markov Chain Monte Carlo (RJMCMC) method to the DA schemes to fit the situation where the modes of the evolution model are unknown and the dimension of the parameters is changing. Our framework allows us to identify the modes in the evolution model and their changes, as well as estimate the parameters and states of the dynamic system. Numerical experiments are conducted and the results show that our framework can effectively identify the underlying evolution models and increase the predictive accuracy of DA methods.

Split-and-merge model selection of mixtures of Gaussian processes with RJMCMC

The mixture of Gaussian processes is a powerful statistical learning model that can be effectively applied to curve clustering and prediction. However, the corresponding model selection problem, that is, selecting an appropriate number of components in the mixture, is rather difficult to solve. In our previous work, we established the split-and-merge automatic model selection algorithm for mixtures of Gaussian processes along the output space under the framework of Reversible Jump Markov Chain Monte Carlo (RJMCMC), which can not only determine the number of actual Gaussian processes but also dynamically adjust the Gaussian process components to avoid dependence on parameter initialization and initial partitioning of the dataset during the parameter learning on a given dataset. In this study, we propose two algorithms: Penalized Likelihood RJMCMC and Penalized Prior RJMCMC. The former integrates a penalized term into the likelihood, while the latter incorporates a penalized term into the prior and operates within the full Bayesian inference framework, both aiming to focus more sharply on determining the number of components in the convergence process. Furthermore, we prove the geometric ergodicity of the RJMCMC algorithm for the mixture of Gaussian processes model, ensuring convergence of the posterior distribution with sufficient iterations. The experimental results further demonstrate the robustness of our PP-RJMCMC algorithm in model selection, showing superior performance compared to traditional approaches in curve classification and clustering. Additionally, the prediction performance is comparable to the EM algorithm. Although not directly explored in this study, the RJMCMC results can be used to initialize the EM algorithm, which could potentially improve prediction accuracy and accelerate computation.

Bayesian Estimation of the Semiparametric Spatial Lag Model

This paper proposes a semiparametric spatial lag model and develops a Bayesian estimation method for this model. In the estimation of the model, the paper combines Reversible Jump Markov Chain Monte Carlo (RJMCMC) algorithm, random walk Metropolis sampler, and Gibbs sampling techniques to sample all the parameters. The paper conducts numerical simulations to validate the proposed Bayesian estimation theory using a numerical example. The simulation results demonstrate satisfactory estimation performance of the parameter part and the fitting performance of the nonparametric function under different spatial weight matrix settings. Furthermore, the paper applies the constructed model and its estimation method to an empirical study on the relationship between economic growth and carbon emissions in China, illustrating the practical application value of the theoretical results.

A deep learning-based surrogate model for trans-dimensional inversion of discrete fracture networks

Fractures and their geometrical patterns are usually required to analyze the mechanical and flow properties of porous media in the subsurface. Fracture characterization is therefore regarded of crucial importance for optimizing production management or achieving maximum storage capacity. In this research, we propose to invert the fracture networks under the Bayesian framework for the uncertainty quantification. In particular, the number of fractures in the modelling system is treated as unknown, leading to a trans-dimensional inverse problem, and the reversible jump Markov chain Monte Carlo algorithm is applied to sample the model space with possible model moves proposed in the sampling process. A deep learning network is further applied as a surrogate model in the sampling process for increasing the computational efficiency, instead of using the physical forward simulator. We apply the proposed methodology to estimate the spatial distribution of fracture networks based on the head measurements from the steady-state flow simulation. The prior distributions of fracture parameters such as position, orientation and length are described using the discrete fracture networks approach that is deeply rooted in stochastic modelling. Due to the high non-uniqueness, the correct spatial distribution of fracture patterns cannot be successfully recovered in this case study, even a good match between observed and simulated head data is reached. More analysis could be performed in the future with the production historical data or more informative priors.

Bayesian Structure Learning for Climate Model Evaluation

AbstractA Bayesian structure learning approach is employed to compare and contrast interactions between the major climate teleconnections over the recent past as revealed in reanalyses and climate model simulations from leading Meteorological Centers. In a previous study, the authors demonstrated a general framework using homogeneous Dynamic Bayesian Network models constructed from reanalyzed time series of empirical climate indices to compare probabilistic graphical models. Reversible jump Markov Chain Monte Carlo is used to provide uncertainty quantification for selecting the respective network structures. The incorporation of confidence measures in structural features provided by the Bayesian approach is key to yielding informative measures of the differences between products if network‐based approaches are to be used for model evaluation, particularly as point estimates alone may understate the relevant uncertainties. Here we compare models fitted from the NCEP/NCAR and JRA‐55 reanalyses and Coupled Model Intercomparison Project version 5 (CMIP5) historical simulations in terms of associations for which there is high posterior confidence. Examination of differences in the posterior probabilities assigned to edges of the directed acyclic graph provides a quantitative summary of departures in the CMIP5 models from reanalyses. In general terms the climate model simulations are in better agreement with reanalyses where tropical processes dominate, and autocorrelation time scales are long. Seasonal effects are shown to be important when examining tropical‐extratropical interactions with the greatest discrepancies and largest uncertainties present for the Southern Hemisphere teleconnections.

Stochastic inversion of fracture networks using the reversible jump Markov chain Monte Carlo algorithm

Characterization of fracture networks is essential in the production optimization or storage calculation for enhanced geothermal systems or geologic carbon storage. A novel inversion approach is proposed for estimating the fracture networks in this research. The discrete fracture network technique is adopted to probabilistically describe various fracture parameters such as trace length, midpoint position or azimuthal angle. The reversible jump Markov chain Monte Carlo algorithm is applied to explore the target posterior distribution of model parameters of differing dimensionality, in which the number of fractures is assumed to be unknown. More specifically, the birth-death strategy is utilized to perturb the fracture number iteratively in the sampling process. The proposed methodology is applied with two different types of observational datasets, namely the head records from steady-state flow simulation and the acoustic impedance obtained from seismic inversion. The sampling results can successfully recover the fracture geometry in the observed domain, and the number of fractures in the system can be retrieved as well. Benchmarked on multiple Markov chain trials, the technique of parallel tempering can greatly improve the convergence efficiency and increase the diversity of sampled posterior models, through the random swapping of model states across the whole temperature ladder.

Generalized Bayesian Additive Regression Trees Models: Beyond Conditional Conjugacy

Bayesian additive regression trees have seen increased interest in recent years due to their ability to combine machine learning techniques with principled uncertainty quantification. The Bayesian backfitting algorithm used to fit BART models, however, limits their application to a small class of models for which conditional conjugacy exists. In this article, we greatly expand the domain of applicability of BART to arbitrary generalized BART models by introducing a very simple, tuning-parameter-free, reversible jump Markov chain Monte Carlo algorithm. Our algorithm requires only that the user be able to compute the likelihood and (optionally) its gradient and Fisher information. The potential applications are very broad; we consider examples in survival analysis, structured heteroskedastic regression, and gamma shape regression.

Adaptive tempered reversible jump algorithm for Bayesian curve fitting

Bayesian curve fitting plays an important role in inverse problems, and is often addressed using the reversible jump Markov chain Monte Carlo (RJMCMC) algorithm. However, this algorithm can be computationally inefficient without appropriately tuned proposals. As a remedy, we present an adaptive RJMCMC algorithm for the curve fitting problems by extending the adaptive Metropolis sampler from a fixed-dimensional to a trans-dimensional case. In this presented algorithm, both the size and orientation of the proposal function can be automatically adjusted in the sampling process. Specifically, the curve fitting setting allows for the approximation of the posterior covariance of the a priori unknown function on a representative grid of points. This approximation facilitates the definition of efficient proposals. In addition, we introduce an auxiliary-tempered version of this algorithm via non-reversible parallel tempering. To evaluate the algorithms, we conduct numerical tests involving a series of controlled experiments. The results demonstrate that the adaptive algorithms exhibit significantly higher efficiency compared to the conventional ones. Even in cases where the posterior distribution is highly complex, leading to ineffective convergence in the auxiliary-tempered conventional RJMCMC, the proposed auxiliary-tempered adaptive RJMCMC performs satisfactorily. Furthermore, we present a realistic inverse example to test the algorithms. The successful application of the adaptive algorithm distinguishes it again from the conventional one that fails to converge effectively even after millions of iterations.

Anisotropic Reversible‐Jump McMC Shear‐Velocity Tomography of the Eastern Alpine Crust

AbstractThe eastern Alpine crust has been shaped by the continental collision of the European and Adriatic plates beginning at 35 Ma and was affected by a major reorganization after 20 Ma. To better understand how the eastern Alpine surface structures link with deep seated processes, we analyze the depth‐dependent seismic anisotropy based on Rayleigh wave propagation. Ambient noise recordings are evaluated to extract Rayleigh wave phase dispersion measurements. These are inverted in a two step approach for the azimuthally anisotropic shear velocity structure. Both steps are performed with a reversible jump Markov chain Monte Carlo (rj‐McMC) approach that estimates data errors and propagates the modeled uncertainties from the phase velocity maps into the depth inversion. A two layer structure of azimuthal anisotropy is imaged in the Alpine crust, with an orogen‐parallel upper crust and approximately orogen‐perpendicular layer in the lower crust and the uppermost mantle. In the upper layer, the anisotropy tends to follow major fault lines and may thus be an apparent, structurally driven anisotropy. The main foliation and fold axis orientations might contribute to the anisotropy. In the lower crust, the N‐S orientation of the fast axis is mostly confined to regions north of the Periadriatic Fault and may be related to European subduction. Outside the orogen, no clearly layered structure is identified. The anisotropy pattern in the northern Alpine foreland is found to be similar compared to SKS studies which is an indication of very homogeneous fast axis directions throughout the crust and the upper mantle.

Exploring the Relationship Between Technical and Comfort Factors in Designing Efficient Buildings: A Statistical Analysis of a Real-World Dataset

This study investigates the three factors that contribute to designing efficient buildings, namely technical solutions, facade systems, and occupant requirements, through the use of a real-world dataset consisting of 49 efficient buildings from various locations across the globe. Each factor comprises distinct elements that are essential in achieving building efficiency. Statistical methods, including correlation and Kruskal-Wallis methods, as well as advanced statistical methods such as the reversible jump Markov chain Monte Carlo method, were employed to estimate parameters that represent the conditional dependence between the elements of each factor. The undirected graphs were generated for each factor based on the conditional depence between the elements of the factor which is shown by a link. Through the analysis of these graphs, designers can enhance their comprehension of the correlation between the various elements of each factor, which can ultimately result in improved building efficiency. This, in turn, may lead to a decrease in air pollution and energy consumption while enhancing human comfort.

A Bayesian survival treed hazards model using latent Gaussian processes.

Survival models are used to analyze time-to-event data in a variety of disciplines. Proportional hazard models provide interpretable parameter estimates, but proportional hazard assumptions are not always appropriate. Non-parametric models are more flexible but often lack a clear inferential framework. We propose a Bayesian treed hazards partition model that is both flexible and inferential. Inference is obtained through the posterior tree structure and flexibility is preserved by modeling the log-hazard function in each partition using a latent Gaussian process. An efficient reversible jump Markov chain Monte Carlo algorithm is accomplished by marginalizing the parameters in each partition element via a Laplace approximation. Consistency properties for the estimator are established. The method can be used to help determine subgroups as well as prognostic and/or predictive biomarkers in time-to-event data. The method is compared with some existing methods on simulated data and a liver cirrhosis dataset.

Artificial Bee Colony Algorithm with Adaptive Parameter Space Dimension: A Promising Tool for Geophysical Electromagnetic Induction Inversion

Frequency-domain electromagnetic induction (FDEMI) methods are frequently used in non-invasive, area-wise mapping of the subsurface electromagnetic soil properties. A crucial part of data analysis is the geophysical inversion of the data, resulting in either conductivity and/or magnetic susceptibility subsurface distributions. We present a novel 1D stochastic optimization approach that combines dimension-adapting reversible jump Markov chain Monte Carlo (MCMC) with artificial bee colony (ABC) optimization for geophysical inversion, with specific application to frequency-domain electromagnetic induction (FDEMI) data. Several solution models of simplified model geometry and a variable number of model knots, which are found by the inversion method, are used to create re-sampled resulting average models. We present synthetic test inversions using conductivity models based on 14 direct-push (DP) EC logs from Greece, Italy, and Germany, as well as field data applications using multi-coil FDEMI devices from three sites in Azerbaijan and Germany. These examples show that the method can effectively lead to solutions that resemble the known DP input models or image reasonable stratigraphic and archaeological features in the field data. Neighboring 1D solutions on field data examples show high coherence along profiles even though each 1D inversion is independently handled. The computational effort for one 1D inversion is less than 120,000 forward calculations, which is much less than usually needed in MCMC inversions, whereas the resulting models show more plausible solutions due to the dimension-adapting properties of the inversion method.

TIME SERIES WITH MULTIPLE CHANGE POINTS AND CENSORED OBSERVATIONS

This article examines a Bayesian model for a nonstationary time series with an unknown number of change points and censored observations. Each segment is assumed to be an autoregressive process with order one. To estimate the number and locations of change points, we use the reversible jump Markov chain Monte Carlo (RJMCMC) algorithm. The censored problem is solved by imputing the censored values from a multivariate normal distribution based on the observed part. A numerical example shows that the estimates of the number of change points and their localizations have little bias. Additionally, the estimates are robust to the censoring percentage.

Bayesian posterior density estimation reveals degeneracy in three-dimensional multiple emitter localization

Single-molecule localization microscopy requires sparse activation of emitters to circumvent the diffraction limit. In densely labeled or thick samples, overlap of emitter images is inevitable. Single-molecule localization of these samples results in a biased parameter estimate with a wrong model of the number of emitters. On the other hand, multiple emitter fitting suffers from point spread function degeneracy, which increases model and parameter uncertainty. To better estimate the model, parameters and uncertainties, a three-dimensional Bayesian multiple emitter fitting algorithm was constructed using Reversible Jump Markov Chain Monte Carlo. It reconstructs the posterior density of both the model and the parameters, namely the three-dimensional position and photon intensity, of overlapping emitters. The ability of the algorithm to separate two emitters at varying distance was evaluated using an astigmatic point spread function. We found that for astigmatic imaging, the posterior distribution of the emitter positions is multimodal when emitters are within two times the in-focus standard deviation of the point spread function. This multimodality describes the ambiguity in position that astigmatism introduces in localization microscopy. Biplane imaging was also tested, proving capable of separating emitters up to 0.75 times the in-focus standard deviation of the point spread function while staying free of multimodality. The posteriors seen in astigmatic and biplane imaging demonstrate how the algorithm can identify point spread function degeneracy and evaluate imaging techniques for three-dimensional multiple-emitter fitting performance.

Trans-Dimensional Bayesian Inversion of Rayleigh Wave Dispersion Curves for Shallow Surface Two-Dimensional Imaging

In this study, we achieve 2-D spatial structure directly from dispersion measurements using the reversible jump Markov chain Monte Carlo (rj-McMC) algorithm with a fully 2-D model parameterization. New parsimonious parameterization with low dimensional representation is introduced, which enable us to accurately reproduce the wave velocity field and the associated uncertainties. Field data test demonstrates that the proposed method estimates directly the shear-wave velocity profile with better lateral resolution. The associated uncertainties provide directly interpretable uncertainty on structures of interest.

Adaptive MCMC for Bayesian Variable Selection in Generalised Linear Models and Survival Models.

Developing an efficient computational scheme for high-dimensional Bayesian variable selection in generalised linear models and survival models has always been a challenging problem due to the absence of closed-form solutions to the marginal likelihood. The Reversible Jump Markov Chain Monte Carlo (RJMCMC) approach can be employed to jointly sample models and coefficients, but the effective design of the trans-dimensional jumps of RJMCMC can be challenging, making it hard to implement. Alternatively, the marginal likelihood can be derived conditional on latent variables using a data-augmentation scheme (e.g., Pólya-gamma data augmentation for logistic regression) or using other estimation methods. However, suitable data-augmentation schemes are not available for every generalised linear model and survival model, and estimating the marginal likelihood using a Laplace approximation or a correlated pseudo-marginal method can be computationally expensive. In this paper, three main contributions are presented. Firstly, we present an extended Point-wise implementation of Adaptive Random Neighbourhood Informed proposal (PARNI) to efficiently sample models directly from the marginal posterior distributions of generalised linear models and survival models. Secondly, in light of the recently proposed approximate Laplace approximation, we describe an efficient and accurate estimation method for marginal likelihood that involves adaptive parameters. Additionally, we describe a new method to adapt the algorithmic tuning parameters of the PARNI proposal by replacing Rao-Blackwellised estimates with the combination of a warm-start estimate and the ergodic average. We present numerous numerical results from simulated data and eight high-dimensional genetic mapping data-sets to showcase the efficiency of the novel PARNI proposal compared with the baseline add-delete-swap proposal.

Fracture Network Characterization in Reservoirs by Joint Inversion of Microseismicity and Thermal Breakthrough Data: Method Development and Verification

AbstractThe spatial configuration of fractures often regulates flow and transport in the subsurface. However, the characterization of fracture networks is a challenging task, especially in deep reservoirs, because here only a limited number of boreholes is available to perform downhole logs and cross‐well hydraulic and tracer tests. In this study, we develop a joint inversion procedure to infer fracture number, geometry, and aperture based on Bayesian principles, utilizing microseismic (MS) events and thermal breakthrough data; both are typically available in enhanced geothermal systems and unconventional gas/oil reservoirs. A basic discrete fracture geometry, including orientations, lengths and positions of fractures, is generated from an MS events cloud by minimizing the distance between fractures and MS events. This geometric configuration, together with the aperture, is then further adjusted, based on the reversible jump Markov chain Monte Carlo algorithm, to minimize the misfit between observed and calculated thermal breakthrough curves by iterative forward flow and heat transport modeling. The inversion model is applied to two synthetic test cases. Following the sensitivity analysis of temperatures to fracture parameters, and the robustness analysis on the model performance under supportive data featuring good and poor quality, it is confirmed that the probabilistic joint inversion procedure approximates the fracture geometry and aperture well, and errors in predicting temperatures based on realizations of fracture networks are well below 5%. The methodology provides a new way to characterize fracture networks in the subsurface, without restrictions on predetermined fracture sizes or the number of fractures.

IsoFrog: a reversible jump Markov Chain Monte Carlofeature selection-based method for predicting isoform functions.

A single gene may yield several isoforms with different functions through alternative splicing. Continuous efforts are devoted to developing machine-learning methods to predict isoform functions. However, existing methods do not consider the relevance of each feature to specific functions and ignore the noise caused by the irrelevant features. In this case, we hypothesize that constructing a feature selection framework to extract the function-relevant features might help improve the model accuracy in isoform function prediction. In this article, we present a feature selection-based approach named IsoFrog to predict isoform functions. First, IsoFrog adopts a reversible jump Markov Chain Monte Carlo (RJMCMC)-based feature selection framework to assess the feature importance to gene functions. Second, a sequential feature selection procedure is applied to select a subset of function-relevant features. This strategy screens the relevant features for the specific function while eliminating irrelevant ones, improving the effectiveness of the input features. Then, the selected features are input into our proposed method modified domain-invariant partial least squares, which prioritizes the most likely positive isoform for each positive MIG and utilizes diPLS for isoform function prediction. Tested on three datasets, our method achieves superior performance over six state-of-the-art methods, and the RJMCMC-based feature selection framework outperforms three classic feature selection methods. We expect this proposed methodology will promote the identification of isoform functions and further inspire the development of new methods. IsoFrog is freely available at https://github.com/genemine/IsoFrog.

SBFEM and Bayesian inference for efficient multiple flaw detection in structures

Bayesian inference is a powerful technique for damage/flaw detection in critical structures. This paper explores the application of Bayesian inference to identify the flaws (voids/inclusions) and its parameters, assuming no information about the flaws is known a priori. Multiple flaws are represented by a parameter vector that contains the locations and the geometric information. In the inverse problem framework, the Scaled Boundary Finite Element Method (SBFEM) with quadtree decomposition is used to solve the forward problem. The flaw parameters are statistically quantified using the Bayesian inference. The likelihood function and prior information update the joint distribution of the number of flaws and their corresponding flaw parameters. The sampling to estimate the posterior distribution of the flaw parameters is based on the trans-dimensional Reversible Jump Markov Chain Monte Carlo (RJMCMC). Also, the impact of additive Gaussian noise on the observed data is investigated. Numerical examples include identifying multiple voids of different shapes such as circle and ellipse and identifying voids and inclusions, using the proposed approach with input sensor data at different noise levels.