Markov Chain Monte Carlo Sampling Research Articles

Identifiability and convergence behavior for Markov chain Monte Carlo using multivariate probit models

Multivariate probit models have been popularly utilized to analysis multivariate ordinal data. However, the identifiable multivariate probit models entail the covariance matrix for the underlying multivariate normal variables to be a correlation matrix, which brings a rigorous task to conduct efficient statistical analysis. Parameter expansion to make the identifiable model to be non-identifiable has been inevitably explored. However, the effect of the expanded parameters on the convergence of Markov chain Monte Carlo (MCMC) is seldomly investigated; in addition, the comparison of MCMC developed based on the identifiable model and that based on the non-identifiable model is hardly ever explored, especially for data with large sample sizes. In this paper, we conduct a thorough investigation to illustrate the effect of the expanded parameters on the convergence of MCMC and compare the behavior of MCMC between the identifiable and non-identifiable models. Our investigation provides a practical guide regarding the construction of non-identifiable models and development of corresponding MCMC sampling methods. We conduct our investigation using simulation studies and present an application using data from the Russia Longitudinal Monitoring Survey-Higher School of Economics (RLMS-HSE) study.

Stochastic Nonlinear Model Updating in Structural Dynamics Using a Novel Likelihood Function within the Bayesian-MCMC Framework

The study presents a novel approach for stochastic nonlinear model updating in structural dynamics, employing a Bayesian framework integrated with Markov Chain Monte Carlo (MCMC) sampling for parameter estimation by using an approximated likelihood function. The proposed methodology is applied to both numerical and experimental cases. The paper commences by introducing Bayesian inference and its constituents: the likelihood function, prior distribution, and posterior distribution. The resonant decay method is employed to extract backbone curves, which capture the non-linear behaviour of the system. A mathematical model based on a single degree of freedom (SDOF) system is formulated, and backbone curves are obtained from time response data. Subsequently, MCMC sampling is employed to estimate the parameters using both numerical and experimental data. The obtained results demonstrate the convergence of the Markov chain, present parameter trace plots, and provide estimates of posterior distributions of updated parameters along with their uncertainties. Experimental validation is performed on a cantilever beam system equipped with permanent magnets and electromagnets. The proposed methodology demonstrates promising results in estimating parameters of stochastic non-linear dynamical systems. Through the use of the proposed likelihood functions using backbone curves, the probability distributions of both linear and non-linear parameters are simultaneously identified. Based on this view, the necessity to segregate stochastic linear and non-linear model updating is eliminated.

Real-time modeling of transient crustal deformation through the quantification of uncertainty deduced from GNSS data

We propose a new method for real-time uncertainty monitoring of earthquake and volcano source models using data from the global navigation satellite system (GNSS) and explore its application concerning observation station placement. The Geospatial Information Authority of Japan operates two main types of GNSS earth observation network system (GEONET) coordinates for crustal deformation monitoring on different time scales: post-processing analysis values and real-time GEONET analysis system for rapid deformation monitoring (REGARD). REGARD uses the Markov Chain Monte Carlo (MCMC) method termed real-time automatic uncertainty estimation of a coseismic single rectangular model using GNSS data (RUNE) for single rectangular fault model estimation to handle uncertainty. Thus far, no GNSS monitoring system can automatically detect transient crustal deformation events, such as volcanic activity and earthquake swarms, on timescales of a day or less. We extended RUNE and developed a core program for a new monitoring system for earthquake and volcanic source models and their uncertainties. Our program achieved automatic and stable MCMC utilization for rectangular fault, dike, Mogi, and spheroid models by increasing the computational speed, improving search efficiency, and adjusting hyperparameters. The program automatically determines the standard deviation of the likelihood function assuming a normal distribution with weights for each observation station. The calculation time was within 15 s for 1 × 106 samples on a standard 1U server. We assessed the reliability of the developed method using synthetic and observed GNSS data from the 2015 Sakurajima volcanic event. The results were consistent with the assumed model and previous studies and indicated an advantage in automatically quantifying uncertainty in a short computation time. Based on MCMC samples, we developed a new visualization algorithm to indicate areas on a map in which the number of observation stations should be expanded. We assessed the reliability using data from the 2023 Noto Peninsula earthquake [Mj 6.5]. The results indicate that the algorithm is helpful in studying the placement of stations. The above model extensions and their application are essential to achieve a rapid quantitative understanding of disaster events near urban areas and for utilizing this information in emergency response activities.Graphical

Flow-based nonperturbative simulation of first-order phase transitions

We present a flow-based method for simulating and calculating nucleation rates of first-order phase transitions in scalar field theory on a lattice. Motivated by recent advancements in machine learning tools, particularly normalizing flows for lattice field theory, we propose the “partitioning flow-based Markov chain Monte Carlo (PFMCMC) sampling” method to address two challenges encountered in normalizing flow applications for lattice field theory: the “mode-collapse” and “rare-event sampling” problems. Using a (2+1)-dimensional real scalar model as an example, we demonstrate the effectiveness of our PFMCMC method in modeling highly hierarchical order parameter probability distributions and simulating critical bubble configurations. These simulations are then used to facilitate the calculation of nucleation rates. We anticipate the application of this method to (3+1)-dimensional theories for studying realistic cosmological phase transitions.

Bayesian Ensemble Kalman Filter for Gaussian Mixture Models

AbstractInverse theory and data assimilation methods are commonly used in earth and environmental science studies to predict unknown variables, such as the physical properties of underground rocks, from a set of measured geophysical data, like geophysical seismic or electromagnetic data. A new Bayesian approach based on the ensemble Kalman filter using Gaussian mixture models is presented to overcome the assumption of Gaussian distribution of the unknown variables commonly used in the data assimilation literature and to generalize the algorithm to inverse problems with multimodal probability distributions. In applications of subsurface characterization, the multimodality of the unknown variables is generally due to the presence of different rock types, also known as geological facies. In the proposed method, the weights of the Gaussian mixture model represent the facies proportions, and they follow a Markov chain model. The proposed Bayesian model generates the unknown model parameters conditioned on measured data using a Markov chain Monte Carlo sampler. The validity of the method is demonstrated on a data assimilation problem where the goal is to estimate the posterior distribution of the unknown rock density from a set of repeated measurements of acoustic wave velocity measured at different times. The proposed method provides accurate estimates with efficient computational times.

Euclid preparation. LIII. LensMC, weak lensing cosmic shear measurement with forward modelling and Markov Chain Monte Carlo sampling

is a weak lensing shear measurement method developed for and Stage-IV surveys. It is based on forward modelling in order to deal with convolution by a point spread function (PSF) with comparable size to many galaxies, sampling the posterior distribution of galaxy parameters via Markov chain Monte Carlo, and marginalisation over nuisance parameters for each of the 1.5 billion galaxies observed by We quantified the scientific performance through high-fidelity images based on the Flagship simulations and emulation of the VIS images, realistic clustering with a mean surface number density of $ arcmin^ for galaxies, and $ arcmin^ for stars, and a diffraction-limited chromatic PSF with a full width at half maximum of $ $ and spatial variation across the field of view. measured objects with a density of $ arcmin^ in $4500 deg^2 $. The total shear bias was broken down into measurement (our main focus here) and selection effects (which will be addressed in future work). We found measurement multiplicative and additive biases of $m_1=(-3.6 $, $m_2=(-4.3 $, $c_1=(-1.78 $, and $c_2=(0.09 $; a large detection bias with a multiplicative component of $1.2 $ and an additive component of $-3 $; and a measurement PSF leakage of $ $ and $ $. When model bias is suppressed, the obtained measurement biases are close to requirement and largely dominated by undetected faint galaxies ($-5 $). Although significant, model bias will be straightforward to calibrate given its weak sensitivity on galaxy morphology parameters. is publicly available at

Automatic structuring of radiology reports with on-premise open-source large language models.

Structured reporting enhances comparability, readability, and content detail. Large language models (LLMs) could convert free text into structured data without disrupting radiologists' reporting workflow. This study evaluated an on-premise, privacy-preserving LLM for automatically structuring free-text radiology reports. We developed an approach to controlling the LLM output, ensuring the validity and completeness of structured reports produced by a locally hosted Llama-2-70B-chat model. A dataset with de-identified narrative chest radiograph (CXR) reports was compiled retrospectively. It included 202 English reports from a publicly available MIMIC-CXR dataset and 197 German reports from our university hospital. Senior radiologist prepared a detailed, fully structured reporting template with 48 question-answer pairs. All reports were independently structured by the LLM and two human readers. Bayesian inference (Markov chain Monte Carlo sampling) was used to estimate the distributions of Matthews correlation coefficient (MCC), with [-0.05, 0.05] as the region of practical equivalence (ROPE). The LLM generated valid structured reports in all cases, achieving an average MCC of 0.75 (94% HDI: 0.70-0.80) and F1 score of 0.70 (0.70-0.80) for English, and 0.66 (0.62-0.70) and 0.68 (0.64-0.72) for German reports, respectively. The MCC differences between LLM and humans were within ROPE for both languages: 0.01 (-0.05 to 0.07), 0.01 (-0.05 to 0.07) for English, and -0.01 (-0.07 to 0.05), 0.00 (-0.06 to 0.06) for German, indicating approximately comparable performance. Locally hosted, open-source LLMs can automatically structure free-text radiology reports with approximately human accuracy. However, the understanding of semantics varied across languages and imaging findings. Question Why has structured reporting not been widely adopted in radiology despite clear benefits and how can we improve this? Findings A locally hosted large language model successfully structured narrative reports, showing variation between languages and findings. Critical relevance Structured reporting provides many benefits, but its integration into the clinical routine is limited. Automating the extraction of structured information from radiology reports enables the capture of structured data while allowing the radiologist to maintain their reporting workflow.

Cosmology in energy-momentum squared symmetric teleparallel gravity

In this letter, we explore the f(Q,TμνTμν) gravity theory, building upon the foundations laid by the f(Q) and f(Q,T) gravity theories. Here, Q represents non-metricity and Tμν stands for the energy-momentum tensor. The proposed action encompasses an arbitrary function of both non-metricity Q and the square of the energy-momentum tensor, specifically T2=TμνTμν. We find the analytical solution for the barotropic fluid case p=ωρ for the model f(Q,TμνTμν)=Q+η(TμνTμν). We constrain parameters of the solution H(z) utilizing CC, BAO, and latest Pantheon+SH0ES samples with the help of Monte Carlo Markov Chain sampling technique along with Bayesian statistical analysis. Further, from the Om diagnostic test, we find that the assumed cosmological model favors the quintessence regime.

Bayesian Probabilistic Inference of Nonparametric Distance Distributions in DEER Spectroscopy.

Double electron-electron resonance (DEER) spectroscopy measures distance distributions between spin labels in proteins, yielding important structural and energetic information about conformational landscapes. Analysis of an experimental DEER signal in terms of a distance distribution is a nontrivial task due to the ill-posed nature of the underlying mathematical inversion problem. This work introduces a Bayesian probabilistic inference approach to analyze DEER data, assuming a nonparametric distance distribution with a Tikhonov smoothness prior. The method uses Markov Chain Monte Carlo sampling with a compositional Gibbs sampler to determine a posterior probability distribution over the entire parameter space, including the distance distribution, given an experimental data set. This posterior contains all of the information available from the data, including a full quantification of the uncertainty about the model parameters. The corresponding uncertainty about the distance distribution is visually captured via an ensemble of posterior predictive distributions. Several examples are presented to illustrate the method. Compared with bootstrapping, it performs faster and provides slightly larger uncertainty intervals.

Gradient boosting: A computationally efficient alternative to Markov chain Monte Carlo sampling for fitting large Bayesian spatio-temporal binomial regression models

Disease forecasting and surveillance often involve fitting models to a tremendous volume of historical testing data collected over space and time. Bayesian spatio-temporal regression models fit with Markov chain Monte Carlo (MCMC) methods are commonly used for such data. When the spatio-temporal support of the model is large, implementing an MCMC algorithm becomes a significant computational burden. This research proposes a computationally efficient gradient boosting algorithm for fitting a Bayesian spatio-temporal mixed effects binomial regression model. We demonstrate our method on a disease forecasting model and compare it to a computationally optimized MCMC approach. Both methods are used to produce monthly forecasts for Lyme disease, anaplasmosis, ehrlichiosis, and heartworm disease in domestic dogs for the contiguous United States. The data have a spatial support of 3108 counties and a temporal support of 108–138 months with 71–135 million test results. The proposed estimation approach is several orders of magnitude faster than the optimized MCMC algorithm, with a similar mean absolute prediction error.

Observational Constraints and Cosmographic Analysis of f(T,TG) Gravity and Cosmology

We perform observational confrontation and cosmographic analysis of f(T,TG) gravity and cosmology. This higher-order torsional gravity is based on both the torsion scalar, as well as on the teleparallel equivalent of the Gauss–Bonnet combination, and gives rise to an effective dark-energy sector which depends on the extra torsion contributions. We employ observational data from the Hubble function and supernova Type Ia Pantheon datasets, applying a Markov chain Monte Carlo sampling technique, and we provide the iso-likelihood contours, as well as the best-fit values for the parameters of the power-law model, an ansatz which is expected to be a good approximation of most realistic deviations from general relativity. Additionally, we reconstruct the effective dark-energy equation-of-state parameter, which exhibits a quintessence-like behavior, while in the future the Universe enters into the phantom regime, before it tends asymptotically to the cosmological constant value. Furthermore, we perform a detailed cosmographic analysis, examining the deceleration, jerk, snap, and lerk parameters, showing that the transition to acceleration occurs in the redshift range 0.52≤ztr≤0.89, as well as the preference of the scenario for quintessence-like behavior. Finally, we apply the Om diagnostic analysis to cross-verify the behavior of the obtained model.

Bayesian polynomial neural networks and polynomial neural ordinary differential equations.

Symbolic regression with polynomial neural networks and polynomial neural ordinary differential equations (ODEs) are two recent and powerful approaches for equation recovery of many science and engineering problems. However, these methods provide point estimates for the model parameters and are currently unable to accommodate noisy data. We address this challenge by developing and validating the following Bayesian inference methods: the Laplace approximation, Markov Chain Monte Carlo (MCMC) sampling methods, and variational inference. We have found the Laplace approximation to be the best method for this class of problems. Our work can be easily extended to the broader class of symbolic neural networks to which the polynomial neural network belongs.

Bayesian Inference for the Calibration of Progressive Damage Model of Dovetail Specimens from Laminated Composite Fan Blade

Abstract Composite fan blades are the preferred alternative for the fan stage of most advanced high bypass ratio turbofan engines. The ability of the dovetail to safely bear this load is one of the key points of the multi-level “test pyramid” approach of compliance demonstration for special airworthiness regulation. Debonding between adjacent layers is the main damage mode of laminated composite fan blades. However, there is difficulty in measuring the as manufactured interlaminar mechanical properties used in finite element models. In this study, tensile loading was applied to simulate the interacting centrifugal force and capture mixed mode damage evolution. Structural responses and material damages were calibrated with measured tensile loads through Bayesian inversion. Posterior probability distributions of maximum interface tractions and contact stresses were solved using Markov chain Monte Carlo sampler. Results indicated the two bilinear cohesive material models had a capacity of predicting empirical means of longitudinal reaction forces as that in test considering additional discrepancy term (0.035 kN and 0.96 kN respectively), while they made an significant impact on the prediction of tensile load history especially when two delamination cracks initiated and propagated. Interface elements provided a higher matching quality in predicting loading history and capturing damage mechanism in association with in-plane progressive damage analysis. This calibrated parameter set could be functioned as benchmark in numerically determining the ultimate tensile load of dovetail elements and reducing the necessary number of physical tests at elemental length level.

Parallel Sampling of Decomposable Graphs Using Markov Chains on Junction Trees

Bayesian inference for undirected graphical models is mostly restricted to the class of decomposable graphs, as they enjoy a rich set of properties making them amenable to high-dimensional problems. While parameter inference is straightforward in this setup, inferring the underlying graph is a challenge driven by the computational difficulty in exploring the space of decomposable graphs. This work makes two contributions to address this problem. First, we provide sufficient and necessary conditions for when multi-edge perturbations maintain decomposability of the graph. Using these, we characterize a simple class of partitions that efficiently classify all edge perturbations by whether they maintain decomposability. Second, we propose a novel parallel nonreversible Markov chain Monte Carlo sampler for distributions over junction tree representations of the graph. At every step, the parallel sampler executes simultaneously all edge perturbations within a partition. Through simulations, we demonstrate the efficiency of our new edge perturbation conditions and class of partitions. We find that our parallel sampler yields improved mixing properties in comparison to the single-move variate, and outperforms current state-of-the-art methods in terms of accuracy and computational efficiency. The implementation of our work is available in the Python package parallelDG. Supplementary materials for this article are available online, including a standardized description of the materials available for reproducing the work.

Multivariate moment least-squares variance estimators for reversible Markov chains

Markov chain Monte Carlo (MCMC) is a commonly used method for approximating expectations with respect to probability distributions. Uncertainty assessment for MCMC estimators is essential in practical applications. Moreover, for multivariate functions of a Markov chain, it is important to estimate not only the auto-correlation for each component but also to estimate cross-correlations, in order to better assess sample quality, improve estimates of effective sample size, and use more effective stopping rules. Berg and Song (2023) introduced the moment least squares (momentLS) estimator, a shape-constrained estimator for the autocovariance sequence from a reversible Markov chain, for univariate functions of the Markov chain. Based on this sequence estimator, they proposed an estimator of the asymptotic variance of the sample mean from MCMC samples. In this study, we propose novel autocovariance sequence and asymptotic variance estimators for Markov chain functions with multiple components, based on the univariate momentLS estimators from Berg and Song (2023). We demonstrate strong consistency of the proposed auto(cross)-covariance sequence and asymptotic variance matrix estimators. We conduct empirical comparisons of our method with other state-of-the-art approaches on simulated and real-data examples, using popular samplers including the random-walk Metropolis sampler and the No-U-Turn sampler from STAN. Supplemental materials for this article are available online.

Notes on the Practical Application of Nested Sampling: MultiNest, (Non)convergence, and Rectification

Nested sampling is a promising tool for Bayesian statistical analysis because it simultaneously performs parameter estimation and facilitates model comparison. MultiNest is one of the most popular nested sampling implementations, and has been applied to a wide variety of problems in the physical sciences. However, MultiNest results, like those of any sampling tool, can be unreliable, and accompanying convergence tests are a component of any analysis. Using analytically tractable test problems, I illustrate how MultiNest, when applied without rigorously chosen hyperparameters, (1) can produce systematically erroneous estimates of the Bayesian evidence, which are more significantly biased for problems of higher dimensionality; (2) can derive posterior estimates with errors on the order of ~100%; (3) can, particularly when sampling noisy likelihood functions, systematically underestimate posterior widths. Furthermore, I show how MultiNest, thanks to the advantageous speed at which it explores parameter space, can also be used to jump-start Markov chain Monte Carlo sampling or more rigorous nested sampling techniques, potentially accelerating more robust measurements of posterior distributions and Bayesian evidences, and overcoming the challenge of Markov chain Monte Carlo initialization.

Using Early Rejection Markov Chain Monte Carlo and Gaussian Processes to Accelerate ABC Methods

Approximate Bayesian computation (ABC) is a class of Bayesian inference algorithms that targets problems with intractable or unavailable likelihood functions. It uses synthetic data drawn from the simulation model to approximate the posterior distribution. However, ABC is computationally intensive for complex models in which simulating synthetic data is very expensive. In this article, we propose an early rejection Markov chain Monte Carlo (ejMCMC) sampler based on Gaussian processes to accelerate inference speed. We early reject samples in the first stage of the kernel using a discrepancy model, in which the discrepancy between the simulated and observed data is modeled by Gaussian process (GP). Hence, synthetic data is generated only if the parameter space is worth exploring. We demonstrate through theory, simulation experiments, and real data analysis that the new algorithm significantly improves inference efficiency compared to existing early-rejection MCMC algorithms. In addition, we employ our proposed method within an ABC sequential Monte Carlo (SMC) sampler. In our numerical experiments, we use examples of ordinary differential equations, stochastic differential equations, and delay differential equations to demonstrate the effectiveness of the proposed algorithm. We develop an R package that is available at https://github.com/caofff/ejMCMC.

Bayesian Inference and Condition Assessment Based on the Deflection of Aging Reinforced Concrete Hollow Slab Bridges

This paper presents a Bayesian inference framework for updating the structural rigidity ratio of aging hollow slab RC bridges using deflection measurements. The framework models the structural rigidity ratio as a stochastic field along the hollow RC slabs, using the Karhunen–Loeve (KL) transform to capture spatial correlation and variation. Bayesian inference is then applied using deflection data from static loading tests, supported by a finite element model (FEM) and a Kriging surrogate model to enhance computational efficiency. The posterior distribution of the structural rigidity ratio is derived using a Markov chain Monte Carlo (MCMC) sampler. The proposed method was tested on an RC bridge with hollow slabs, using deflection measurements taken before and after reinforcement. The Bayesian updates indicated increased structural rigidity ratios after reinforcement, validating the effectiveness of the reinforcement. The deflection predictions from the updated models closely matched the measurements, with the 95% confidence bounds encompassing most of the data. This demonstrates the method’s validity and robustness in capturing the structural improvements post-reinforcement.

Kinetic profile inference with outlier detection using support vector machine regression and Gaussian process regression

Abstract We propose an outlier-resilient Gaussian process regression (GPR) model supported by support vector machine regression (SVMR) for kinetic profile inference. GPR, being a non-parametric regression using Bayesian statistics, has advantages in that it imposes no constraints on profile shapes and can be readily used to integrate different kinds of diagnostics, while it is vulnerable to the presence of even a single outlier among a measured dataset. As an outlier classifier, an optimized SVMR is developed based only on the measurements. Hyper-parameters of the developed GPR model with informative prior distributions are treated in two different ways, i.e. maximum a posteriori (MAP) estimator and marginalization using a Markov Chain Monte Carlo sampler. Our SVMR-supported GPR model is applied to infer ion temperature Ti profiles using measured data from the KSTAR charge exchange spectroscopy system. The GPR-inferred Ti profiles with and without an outlier are compared and show prominent improvement when the outlier is removed by the SVMR. Ti profiles inferred with the MAP estimator and the marginalization scheme are compared. They are noticeably different when observation uncertainties are not small enough, and the marginalization scheme generally provides a smoother profile.

Predicting Colloidal Interaction Parameters from Small-Angle X-ray Scattering Curves Using Artificial Neural Networks and Markov Chain Monte Carlo Sampling.

Small-angle X-ray scattering (SAXS) is a characterization technique that allows for the study of colloidal interactions by fitting the structure factor of the SAXS profile with a selected model and closure relation. However, the applicability of this approach is constrained by the limited number of existing models that can be fitted analytically, as well as the narrow operating range for which the models are valid. In this work, we demonstrate a proof of concept for using an artificial neural network (ANN) trained on SAXS curves obtained from Monte Carlo (MC) simulations to predict values of the effective macroion valency (Z eff) and the Debye length (κ-1) for a given SAXS profile. This ANN, which was trained on 200,000 simulated SAXS curves, was able to predict values of Z eff and κ-1 for a test set containing 25,000 simulated SAXS curves, where most predicted values had errors smaller than 20%. Subsequently, an ANN was used as a surrogate model in a Markov chain Monte Carlo sampling algorithm to obtain maximum a posteriori estimates of Z eff and κ-1, as well as the associated confidence intervals and correlations between Z eff and κ-1 for an experimentally obtained SAXS profile.