Efficient Sampling Algorithm Research Articles

Polya tree Monte Carlo method

Markov Chain Monte Carlo (MCMC) methods have been widely used in Statistics and machine learning research. However, such methods have several limitations, including slow convergence and the inefficiency in handling multi-modal distributions. To overcome these limitations of MCMC methods, a new, efficient sampling method has been proposed and it applies to general distributions including multi-modal ones or those having complex structure. The proposed approach, called the Polya tree Monte Carlo (PTMC) method, roots in constructing a Polya tree distribution using the idea of Monte Carlo method, and then using this distribution to approximate and facilitate sampling from a target distribution that may be complex or have multiple modes. The associated convergence property of the PTMC method is established and computationally efficient sampling algorithms are developed based on the PTMC. Extensive numerical studies demonstrate the satisfactory performance of the proposed method under various settings including its superiority to the usual MCMC algorithms. The evaluation and comparison are carried out in terms of sampling efficiency, computational speed and the capacity of identifying distribution modes. Additional details about the method, proofs and simulation results are provided in the Supplementary Web Appendices online.

Automated Path Searching Reveals the Mechanism of Hydrolysis Enhancement by T4 Lysozyme Mutants

Bacteriophage T4 lysozyme (T4L) is a glycosidase that is widely applied as a natural antimicrobial agent in the food industry. Due to its wide applications and small size, T4L has been regarded as a model system for understanding protein dynamics and for large-scale protein engineering. Through structural insights from the single conformation of T4L, a series of mutations (L99A,G113A,R119P) have been introduced, which have successfully raised the fractional population of its only hydrolysis-competent excited state to 96%. However, the actual impact of these substitutions on its dynamics remains unclear, largely due to the lack of highly efficient sampling algorithms. Here, using our recently developed travelling-salesman-based automated path searching (TAPS), we located the minimum-free-energy path (MFEP) for the transition of three T4L mutants from their ground states to their excited states. All three mutants share a three-step transition: the flipping of F114, the rearrangement of α0/α1 helices, and final refinement. Remarkably, the MFEP revealed that the effects of the mutations are drastically beyond the expectations of their original design: (a) the G113A substitution not only enhances helicity but also fills the hydrophobic Cavity I and reduces the free energy barrier for flipping F114; (b) R119P barely changes the stability of the ground state but stabilizes the excited state through rarely reported polar contacts S117OG:N132ND2, E11OE1:R145NH1, and E11OE2:Q105NE2; (c) the residue W138 flips into Cavity I and further stabilizes the excited state for the triple mutant L99A,G113A,R119P. These novel insights that were unexpected in the original mutant design indicated the necessity of incorporating path searching into the workflow of rational protein engineering.

Graphical Assistant Grouped Network Autoregression Model: A Bayesian Nonparametric Recourse

Vector autoregression model is ubiquitous in classical time series data analysis. With the rapid advance of social network sites, time series data over latent graph is becoming increasingly popular. In this article, we develop a novel Bayesian grouped network autoregression model, which can simultaneously estimate group information (number of groups and group configurations) and group-wise parameters. Specifically, a graphically assisted Chinese restaurant process is incorporated under the framework of the network autoregression model to improve the statistical inference performance. An efficient Markov chain Monte Carlo sampling algorithm is used to sample from the posterior distribution. Extensive studies are conducted to evaluate the finite sample performance of our proposed methodology. Additionally, we analyze two real datasets as illustrations of the effectiveness of our approach.

Bayesian modeling and clustering for spatio-temporal areal data: An application to Italian unemployment

Spatio-temporal areal data can be seen as a collection of time series which are spatially correlated according to a specific neighboring structure. Incorporating the temporal and spatial dimension into a statistical model poses challenges regarding the underlying theoretical framework as well as the implementation of efficient computational methods. We propose to include spatio-temporal random effects using a conditional autoregressive prior, where the temporal correlation is modeled through an autoregressive mean decomposition and the spatial correlation by the precision matrix inheriting the neighboring structure. Their joint distribution constitutes a Gaussian Markov random field, whose sparse precision matrix enables the usage of efficient sampling algorithms. We cluster the areal units using a nonparametric prior, thereby learning latent partitions of the areal units. The performance of the model is assessed via an application to study regional unemployment patterns in Italy. When compared to other spatial and spatio-temporal competitors, the proposed model shows more precise estimates and the additional information obtained from the clustering allows for an extended economic interpretation of the unemployment rates of the Italian provinces.

Quantifying the Effects of Contact Tracing, Testing, and Containment Measures in the Presence of Infection Hotspots

Multiple lines of evidence strongly suggest that infection hotspots, where a single individual infects many others, play a key role in the transmission dynamics of COVID-19. However, most of the existing epidemiological models fail to capture this aspect by neither representing the sites visited by individuals explicitly nor characterizing disease transmission as a function of individual mobility patterns. In this work, we introduce a temporal point process modeling framework that specifically represents visits to the sites where individuals get in contact and infect each other. Under our model, the number of infections caused by an infectious individual naturally emerges to be overdispersed. Using an efficient sampling algorithm, we demonstrate how to estimate the transmission rate of infectious individuals at the sites they visit and in their households using Bayesian optimization (BO) and longitudinal case data. Simulations using fine-grained and publicly available demographic data and site locations from Bern, Switzerland showcase the flexibility of our framework. To facilitate research and analyses of other cities and regions, we release an open-source implementation of our framework.

A Bayesian Singular Value Decomposition Procedure for Missing Data Imputation

Missing data are common in empirical studies. Multiple imputation is a method to handle missing values by replacing them with plausible values. A common imputation method is multiple imputation with chain equations (MICE). MICE defines a series of conditional distributions to impute missing values. Although MICE is relatively easy to implement, it may not converge to a proper joint distribution. An alternative strategy is to model the variables jointly using the general location model, but this model can become complex when the number of variables increases. Both approaches require integration of prior information when there are more variables than cases. We propose a Bayesian model that is based on the singular value decomposition components of a continuous data matrix to impute missing values. The model assumes that the matrix is of low rank by applying double exponential prior distributions on the singular values. We describe an efficient sampling algorithm to estimate the model’s parameters and impute the missing data. The performance of the model is compared to current imputation methods in simulated and real datasets. Of all the methods considered and in most of the simulated and real datasets, the proposed procedure appears to be the most accurate and precise. Supplementary materials for this article are available online.

Superstatistical Wind Fields from Pointwise Atmospheric Turbulence Measurements

Accurate models of turbulent wind fields have become increasingly important in the atmospheric sciences, e.g., for the determination of spatiotemporal correlations in wind parks, the estimation of individual loads on turbine rotor and blades, or the modeling of particle-turbulence interaction in atmospheric clouds or pollutant distributions in urban settings. Because of the difficult task of resolving the fields across a broad range of scales, one oftentimes has to invoke stochastic wind field models that fulfill specific, empirically observed, properties. Whereas commonly used Gaussian random field models solely control second-order statistics (i.e., velocity correlation tensors or kinetic energy spectra), we explicitly show that our extended model emulates the effects of higher-order statistics as well. Most importantly, the empirically observed phenomenon of small-scale intermittency, which can be regarded as one of the key features of atmospheric turbulent flows, is reproduced with a very high level of accuracy and at considerably low computational cost. Our method is based on a multipoint statistical description of turbulent velocity fields that consists of a superposition of multivariate Gaussian statistics with fluctuating covariances. We propose a new and efficient sampling algorithm for this Gaussian scale mixture and demonstrate how such “superstatistical” wind fields can be constrained on a certain number of real-world measurement data points from a meteorological mast array.Received 12 April 2022Revised 26 July 2022Accepted 22 August 2022DOI:https://doi.org/10.1103/PRXEnergy.1.023006Published by the American Physical Society under the terms of the Creative Commons Attribution 4.0 International license. Further distribution of this work must maintain attribution to the author(s) and the published article's title, journal citation, and DOI.Published by the American Physical SocietyPhysics Subject Headings (PhySH)Research AreasAtmospheric scienceBoundary layersShear flowsStochastic processesStructure & turbulence of boundary layersSustainabilityTurbulenceWind energyEnergy Science & TechnologyNonlinear DynamicsFluid DynamicsStatistical Physics

A Hierarchical Random Graph Efficient Sampling Algorithm Based on Improved MCMC Algorithm

A hierarchical random graph (HRG) model combined with a maximum likelihood approach and a Markov Chain Monte Carlo algorithm can not only be used to quantitatively describe the hierarchical organization of many real networks, but also can predict missing connections in partly known networks with high accuracy. However, the computational cost is very large when hierarchical random graphs are sampled by the Markov Chain Monte Carlo algorithm (MCMC), so that the hierarchical random graphs, which can describe the characteristics of network structure, cannot be found in a reasonable time range. This seriously limits the practicability of the model. In order to overcome this defect, an improved MCMC algorithm called two-state transitions MCMC (TST-MCMC) for efficiently sampling hierarchical random graphs is proposed in this paper. On the Markov chain composed of all possible hierarchical random graphs, TST-MCMC can generate two candidate state variables during state transition and introduce a competition mechanism to filter out the worse of the two candidate state variables. In addition, the detailed balance of Markov chain can be ensured by using Metropolis–Hastings rule. By using this method, not only can the convergence speed of Markov chain be improved, but the convergence interval of Markov chain can be narrowed as well. Three example networks are employed to verify the performance of the proposed algorithm. Experimental results show that our algorithm is more feasible and more effective than the compared schemes.

Discrete Gaussian for Differential Privacy

A key tool for building differentially private systems is adding Gaussian noise to the output of a function evaluated on a sensitive dataset. Unfortunately, using a continuous distribution presents several practical challenges. First and foremost, finite computers cannot exactly represent samples from continuous distributions, and previous work has demonstrated that seemingly innocuous numerical errors can entirely destroy privacy. Moreover, when the underlying data is itself discrete (e.g., population counts), adding continuous noise makes the result less interpretable.
 With these shortcomings in mind, we introduce and analyze the discrete Gaussian in the context of differential privacy. Specifically, we theoretically and experimentally show that adding discrete Gaussian noise provides essentially the same privacy and accuracy guarantees as the addition of continuous Gaussian noise. We also present an simple and efficient algorithm for exact sampling from this distribution. This demonstrates its applicability for privately answering counting queries, or more generally, low-sensitivity integer-valued queries.

Modeling tail risks of inflation using unobserved component quantile regressions

This paper proposes methods for Bayesian inference in time-varying parameter (TVP) quantile regressions (QRs) featuring conditional heteroskedasticity. I use data augmentation schemes to render the model conditionally Gaussian and develop an efficient sampling algorithm. Regularization of the high-dimensional parameter space is achieved via dynamic shrinkage priors. The merits of the proposed approach are illustrated in a simulation study, and a simple version of TVP-QR based on an unobserved components model is applied to dynamically trace the quantiles of inflation in the United States, the United Kingdom and the euro area. In an out-of-sample forecast exercise, I find the proposed model to be competitive and perform particularly well for higher-order and tail forecasts. A detailed analysis of the resulting predictive distributions reveals that they are sometimes skewed and occasionally feature heavy tails.

Docking cyclic peptides formed by a disulfide bond through a hierarchical strategy.

Cyclization is a common strategy to enhance the therapeutic potential of peptides. Many cyclic peptide drugs have been approved for clinical use, in which the disulfide-driven cyclic peptide is one of the most prevalent categories. Molecular docking is a powerful computational method to predict the binding modes of molecules. For protein-cyclic peptide docking, a big challenge is considering the flexibility of peptides with conformers constrained by cyclization. Integrating our efficient peptide 3D conformation sampling algorithm MODPEP2.0 and knowledge-based scoring function ITScorePP, we have proposed an extended version of our hierarchical peptide docking algorithm, named HPEPDOCK2.0, to predict the binding modes of the peptide cyclized through a disulfide against a protein. Our HPEPDOCK2.0 approach was extensively evaluated on diverse test sets and compared with the state-of-the-art cyclic peptide docking program AutoDock CrankPep (ADCP). On a benchmark dataset of 18 cyclic peptide-protein complexes, HPEPDOCK2.0 obtained a native contact fraction of above 0.5 for 61% of the cases when the top prediction was considered, compared with 39% for ADCP. On a larger test set of 25 cyclic peptide-protein complexes, HPEPDOCK2.0 yielded a success rate of 44% for the top prediction, compared with 20% for ADCP. In addition, HPEPDOCK2.0 was also validated on two other test sets of 10 and 11 complexes with apo and predicted receptor structures, respectively. HPEPDOCK2.0 is computationally efficient and the average running time for docking a cyclic peptide is about 34 min on a single CPU core, compared with 496 min for ADCP. HPEPDOCK2.0 will facilitate the study of the interaction between cyclic peptides and proteins and the development of therapeutic cyclic peptide drugs. http://huanglab.phys.hust.edu.cn/hpepdock/. Supplementary data are available at Bioinformatics online.

Analysis of Professional Basketball Field Goal Attempts via a Bayesian Matrix Clustering Approach

We propose a Bayesian nonparametric matrix clustering approach to analyze the latent heterogeneity structure in the shot selection data collected from professional basketball players in the National Basketball Association (NBA). The proposed method adopts a mixture of finite mixtures framework and fully uses the spatial information via a mixture of matrix normal distribution representation. We propose an efficient Markov chain Monte Carlo algorithm for posterior sampling that allows simultaneous inference on both the number of clusters and the cluster configurations. We also establish large-sample convergence properties for the posterior distribution. The compelling empirical performance of the proposed method is demonstrated via simulation studies and an application to shot chart data from selected players in the NBAs 2017–2018 regular season. Supplementary materials for this article are available online.

Benchmarking 50-Photon Gaussian Boson Sampling on the Sunway TaihuLight

Boson sampling is expected to be an important milestone that will demonstrate quantum computational advantage (or quantum supremacy). This work establishes the benchmarking of Gaussian boson sampling (GBS) with threshold detection based on the Sunway TaihuLight supercomputer. To achieve the best performance and provide a competitive scenario for future quantum computing studies, the selected simulation algorithm is fully optimized based on a set of innovative approaches, including a parallel framework with almost perfect load balance and an instruction-level optimizing scheme based on a shortest-path-based instruction scheduling. In addition, data precision is carefully processed by an integer-instruction-based and multiple-precision fixed-point implementation, including 128- and 256-bit precison mode, which can be appropriately selected based on an adaptive precision optimizing scheme. Based on these methods, a highly efficient parallel quantum sampling algorithm is designed. The largest run enables us to obtain one Torontonian function of a <inline-formula><tex-math notation="LaTeX">$100\times 100$</tex-math></inline-formula> submatrix from 50-photon GBS within 20 hours in 128-bit precision and 2 days in 256-bit precision. To our knowledge, this was the largest quantum computing simulation based on Boson Sampling by using modern supercomputers.

Quantum mean-value approximator for hard integer-value problems

Evaluating the expectation of a quantum circuit is a classically difficult problem known as the quantum mean value problem (QMV). It is used to optimize the quantum approximate optimization algorithm and other variational quantum eigensolvers. We show that such an optimization can be improved substantially by using an approximation rather than the exact expectation. Together with efficient classical sampling algorithms, a quantum algorithm with minimal gate count can thus improve the efficiency of general integer-value problems, such as the shortest vector problem (SVP) investigated in this work.

Generic, efficient and isochronous Gaussian sampling over the integers

Gaussian sampling over the integers is one of the fundamental building blocks of lattice-based cryptography. Among the extensively used trapdoor sampling algorithms, it is ineluctable until now. Under the influence of numerous side-channel attacks, it is still challenging to construct a Gaussian sampler that is generic, efficient, and resistant to timing attacks. In this paper, our contribution is three-fold. First, we propose a secure, efficient exponential Bernoulli sampling algorithm. It can be applied to Gaussian samplers based on rejection samplings. We apply it to FALCON, a candidate of round 3 of the NIST post-quantum cryptography standardization project, and reduce its signature generation time by 13–14%. Second, we develop an isochronous Gaussian sampler based on rejection sampling. Our Algorithm can securely sample from Gaussian distributions with different standard deviations and arbitrary centers. We apply it to PALISADE (S&P 2018), an open-source lattice-based cryptography library. During the online phase of trapdoor sampling, the running time of the G-lattice sampling algorithm is reduced by 44.12% while resisting timing attacks. Third, we improve the efficiency of the COSAC sampler (PQC 2020). The new COSAC sampler is 1.46x–1.63x faster than the original and has the lowest expected number of trials among all Gaussian samplers based on rejection samplings. But it needs a more efficient algorithm sampling from the normal distribution to improve its performance.

Bayesian semi-parametric modeling of covariance matrices for multivariate longitudinal data.

The article develops marginal models for multivariate longitudinal responses. Overall, the model consists of five regression submodels, one for the mean and four for the covariance matrix, with the latter resulting by considering various matrix decompositions. The decompositions that we employ are intuitive, easy to understand, and they do not rely on any assumptions such as the presence of an ordering among the multivariate responses. The regression submodels are semi-parametric, with unknown functions represented by basis function expansions. We use spike-slap priors for the regression coefficients to achieve variable selection and function regularization, and to obtain parameter estimates that account for model uncertainty. An efficient Markov chain Monte Carlo algorithm for posterior sampling is developed. The simulation study presented investigates the gains that one may have when considering multivariate longitudinal analyses instead of univariate ones, and whether these gains can counteract the negative effects of missing data. We apply the methods on a highly unbalanced longitudinal dataset with four responses observed over a period of 20 years.

From Ties to Events in the Analysis of Interorganizational Exchange Relations.

Relational event models expand the analytical possibilities of existing statistical models for interorganizational networks by: (i) making efficient use of information contained in the sequential ordering of observed events connecting sending and receiving units; (ii) accounting for the intensity of the relation between exchange partners, and (iii) distinguishing between short- and long-term network effects. We introduce a recently developed relational event model (REM) for the analysis of continuously observed interorganizational exchange relations. The combination of efficient sampling algorithms and sender-based stratification makes the models that we present particularly useful for the analysis of very large samples of relational event data generated by interaction among heterogeneous actors. We demonstrate the empirical value of event-oriented network models in two different settings for interorganizational exchange relations-that is, high-frequency overnight transactions among European banks and patient-sharing relations within a community of Italian hospitals. We focus on patterns of direct and generalized reciprocity while accounting for more complex forms of dependence present in the data. Empirical results suggest that distinguishing between degree- and intensity-based network effects, and between short- and long-term effects is crucial to our understanding of the dynamics of interorganizational dependence and exchange relations. We discuss the general implications of these results for the analysis of social interaction data routinely collected in organizational research to examine the evolutionary dynamics of social networks within and between organizations.

Importance Sampling for a Simple Markovian Intensity Model Using Subsolutions

This article considers importance sampling for estimation of rare-event probabilities in a specific collection of Markovian jump processes used for, e.g., modeling of credit risk. Previous attempts at designing importance sampling algorithms have resulted in poor performance and the main contribution of the article is the design of efficient importance sampling algorithms using subsolutions. The dynamics of the jump processes cause the corresponding Hamilton-Jacobi equations to have an intricate state-dependence, which makes the design of efficient algorithms difficult. We provide theoretical results that quantify the performance of importance sampling algorithms in general and construct asymptotically optimal algorithms for some examples. The computational gain compared to standard Monte Carlo is illustrated by numerical examples.

An Efficient Sampling Algorithm for Difficult Tree Pairs

It is an open question whether there exists a polynomial-time algorithm for computing the rotation distances between pairs of extended ordered binary trees.The problem of computing the rotation distance between an arbitrary pair of trees, (S, T), can be efficiently reduced to the problem of computing the rotation distance between a difficult pair of trees (S', T'), where there is no known first step which is guaranteed to be the beginning of a minimal length path. Of interest, therefore, is how to sample such difficult pairs of trees of a fixed size. We show that it is possible to do so efficiently, and present such an algorithm that runs in time O(n4).

Interactive High-Resolution Simulation of Granular Material

We introduce a particle-based simulation method for granular material in interactive frame rates. We divide the simulation into two decoupled steps. In the first step, a relatively small number of particles is accurately simulated with a constraint-based method. Here, all collisions and the resulting friction between the particles are taken into account. In the second step, the small number of particles is significantly increased by an efficient sampling algorithm without creating additional artifacts. The method is particularly robust and allows relatively large time steps, which makes it well suited for real-time applications. With our method, up to 500k particles can be computed in interactive frame rates on consumer CPUs without relying on GPU support for massive parallel computing. This makes it well suited for applications where a lot of GPU power is already needed for render tasks.