Peer Revenue Publish or profit

Papers

Fresh from arXiv, typeset properly. Highlighted words trade on the market.

Statistics Statistical Computation

stat.CO Watch Statistical Computation (sign in to use)

Sunday 20 September

Statistical Computation · Cross-listed from stat.ME

Bayesian Posterior Learning of Mixed Graphical Models

Mixed Graphical Models (MGMs) provide a flexible framework for structure learning from heterogeneous data by treating sets of both continuous and discrete. Bayesian inference for MGMs remains challenging due to the combinatorial complexity of graph space exploration and posterior computation. In this paper, we model discrete components through latent Gaussian variables and consider two likelihood specifications: a copula-based ranked likelihood, yielding the copula-MGM, and a probit formulation based on cut-off points, yielding the probit-MGM. We then propose the Mixed Graph WWA, a class of MCMC methods for posterior simulation in Bayesian MGMs. Building upon the WWA algorithm, we develop two specialized algorithms: copula-WWA for copula-MGMs and probit-WWA for probit-MGMs. Both methods exploit the latent Gaussian representations to perform posterior inference through a Gibbs sampling scheme that alternates between latent-variable augmentation and graph-structure updates. Through extensive simulation studies we demonstrate that the proposed methods achieve graph recovery accuracy comparable to or better than existing approaches, including copula-BD MCMC and probit-BD MCMC based on the Birth-Death MCMC methodology, while maintaining efficient posterior exploration and favorable effective sample size per unit computational time. We further illustrate the practical utility of our approach through an application to the PAM$50$ breast cancer gene expression dataset, where the inferred MGMs reveal meaningful dependencies between gene expression profiles and cancer subtypes. These results highlight the effectiveness of the Mixed Graph WWA method as a scalable and principled tool for Bayesian structure learning in MGMs.
Read the full abstract

Statistical Computation · Cross-listed from cs.LG

Optimal Multi-way Decision Trees for Stratified Sampling in Online Controlled Experiments

Online controlled experiments, or A/B tests, are widely used to estimate causal effects on digital platforms. A central challenge is to improve experimental sensitivity, or statistical power, without increasing the experimental sample size. Stratified sampling is a classical variance reduction technique; however, its effectiveness depends critically on how the strata are constructed. We thus propose an optimization-based stratification framework for stratified sampling using optimal multi-way decision trees. Our method, called Optimal Multi-way Stratification Trees (OMST), formulates stratification as a path-selection problem over a feature graph. The selected paths define interpretable stratification rules and are optimized using an exact variance-minimizing binary optimization formulation under continuous proportional allocation and a Neyman-type optimal allocation. We incorporate supervised optimal binning to generate outcome-relevant candidate splits for numerical features. Furthermore, we introduce reduction procedures for redundant candidate paths and assignment constraints, substantially reducing the optimization problem size. Experiments on both a real-world and a simulated dataset demonstrate that OMST achieves comparable or superior variance reduction to existing methods while maintaining shallow and interpretable stratification trees.
Read the full abstract

Saturday 19 September

Statistical Computation · Cross-listed from stat.AP

Auditing Bayesian Graph Alignment: Diagnostic Comparisons and Reference Failure

Bayesian graph alignment estimates correspondence probabilities, but convergence of an alignment-score trace need not imply accurate correspondence marginals. We audit this gap on 240 new exact graph pairs from four source families, 240 larger pairs with 20-100 vertices, and a separate 60-case exact implementation check. Under an explicit edge-flip likelihood, we compare three samplers and score, marginal, indicator, categorical, and classifier-based diagnostics. Marginal disagreement improves error discrimination over score R-hat for the exact informed sampler, but its improvement for vanilla local sampling is uncertain. Assignment-based R* and short indicator panels are competitive; no diagnostic dominates across samplers and endpoints. At larger sizes, diagnostics predict subsequent marginal changes, not posterior error, and classification performance depends on the drift threshold. Disjoint-window and held-out-chain checks attenuate but preserve positive associations. Only 22 of 240 original reference sets pass an agreement screen. On forty failure-selected cases, eightfold SMC particle escalation does not resolve disagreement, whereas additional rejuvenation helps. Longer informed runs remain unstable. An elementary feasible-alignment bound demonstrates severely unrepresentative SMC and informed-chain scores in concentrated 100-vertex cases, independently of approximate reference consensus. We also exhibit common-start chains with near-zero disagreement despite exact marginal error near .967. These results support assignment-sensitive auditing while identifying limits of finite budgets, diagnostic rankings, and reference agreement as evidence of accuracy.
Read the full abstract

Statistical Computation · Cross-listed from cs.LG

Causal Inference with Unobserved Confounding: A Mixture Learning Perspective

Unobserved confounding is a fundamental challenge in causal inference from observational data. This article develops a mixture-learning perspective, viewing latent confounders as sources of heterogeneity that induce mixture structure in observed data. Under suitable structural and identifiability assumptions, recovering the mixing distribution and component mechanisms enables estimation of interventional distributions and causal estimands. Using variants of Bernoulli mixtures as a running example, we contextualize mixture-learning techniques and their structural assumptions, and connect them to causal inference in panel-data settings, including latent factor models and synthetic interventions.We then consider high-dimensional exponential-family mixtures with dependent outcome trajectories, moving beyond counterfactual means to model counterfactual distributions. We situate this perspective relative to complementary approaches for unobserved confounding. Together, these ideas provide a bridge between mixture learning and causal inference, connecting recent advances in high-dimensional mixture learning to scalable identification and estimation of causal effects while raising new challenges for mixture learning.
Read the full abstract

Friday 18 September

Statistical Computation · Cross-listed from stat.ML

A Bayesian Vertical Federated Learning Framework for Multivariate Reduced-Rank High-Dimensional Regression

Federated learning (FL) has emerged as a leading privacy-preserving framework for collaborative machine learning across decentralized environments. While considerable progress has been made in horizontal federated learning (HFL), where data with common features is distributed across sites, vertical federated learning (VFL), where sites share observations across distinct feature sets, remains less explored. Advancing Bayesian high-dimensional multivariate reduced-rank regression methods for VFL poses unique challenges: (a) stringent privacy regulations preventing local site data sharing, and (b) fitting local regressions overlooks essential modeling aspects like inter-variable correlations. In contrast HFL allows each site to fit a comparable model independently. We present a novel Bayesian VFL framework for multivariate high-dimensional reduced-rank regression, termed BayesVFLReg, which enables precise coefficient estimation while safeguarding both feature and response privacy. Participating sites use a shared random sketching matrix to compress local variables into privacy-preserving sketches. A central server collects these sketches where Bayesian multivariate reduced-rank regression uses Gaussian scale mixture priors. For feature selection, we introduce a single-step post-processing strategy based on mixture-model clustering of the absolute posterior coefficient means to distinguish signal from noise per response variable. BayesVFLReg is computationally scalable for large, high-dimensional datasets and facilitates efficient variable selection. Theoretically, we establish sharp non-asymptotic bounds on the posterior probability that the fitted density falls within a Hellinger ball centered at the true data-generating density. Comparative simulation studies and real-world data analyses show that BayesVFLReg reliably identifies sparse feature effects, even under feature correlation.
Read the full abstract

All authors Brigham Halverson, Sharmistha Guha, Jessica Bernard, Rajarshi Guhaniyogi

Statistical Computation · Cross-listed from physics.comp-ph

Stochastic consensus dynamics for decentralized decision systems

A fundamental mechanism underlying collective phenomena in social, technological, and economic systems is decentralized decision-making, in which the behavior of individual agents follows from local interactions in the absence of central coordination. Consensus formation is a central feature of such systems, with relevance to blockchain networks, distributed artificial intelligence, autonomous multi-agent systems, distributed control, and collective decision networks. We investigate consensus formation using a stochastic consensus model on random networks and examine how initial conditions, network connectivity, and system size shape the emergence of unanimous states. We show that the stochastic dynamics strongly amplify small initial majorities, progressively suppress the competing state, and drive the system toward a predictable collective outcome. Network connectivity primarily controls the efficiency of this process: increasing connectivity accelerates the propagation of local agreement and reduces the likelihood that fluctuations reverse the initially dominant state, although these gains gradually saturate in highly connected networks. System size produces a complementary effect. Larger networks require more individual updates to reach unanimity, but they are also increasingly reliable in selecting the state favored by the initial majority. Finite-size analysis shows that the range of initial conditions associated with uncertain outcomes becomes progressively narrower as the network grows, decreasing approximately with the inverse square root of the system size. These results reveal a collective amplification mechanism by which weak initial asymmetries become increasingly decisive in large decentralized networks, and they provide a simple framework for understanding the efficiency, predictability, and reliability of consensus formation in distributed decision systems.
Read the full abstract

All authors André L. M. Vilela, Caio B. L. Silva, Kenric P. Nelson, Emilio Cobanera, Gaogao Dong

Statistical Computation · Cross-listed from math.PR

Interval-Constrained Brownian Paths: Exact Interpolation and Extrapolation

We study Brownian motion and Brownian bridge processes conditioned to remain in a fixed interval $[0,a]$, focusing on the conditional distribution at a single time. For a Brownian motion in $[0,a]$, conditioned on survival up to time $t$ we recall (and present in a self-contained form) the conditional density of its position at $t$ (exact extrapolation). For a Brownian bridge conditioned to remain in $[0,a]$ on $[0,T]$ we express the probability density at an interior time as a normalized product of killed transition densities (exact interpolation). These densities admit dual complementary series representations, which are linked via the Jacobi theta identity. Our main contribution is a unified exact sampling suite for both extrapolation and interpolation that (i) includes boundary endpoints and (ii) automatically switches between density representations to keep acceptance rates efficient across regimes. To that end, we derive simple proposal families which cover both small- and large-time regimes, with an automatic rule selecting the tighter envelope. Our procedures extend to exact simulation of discrete skeleton paths at multiple times via the Markov property.
Read the full abstract

Statistical Computation · Cross-listed from stat.ME

Optimal study interior partitioning (OSIP): An optimization approach to observational study design under fine balance constraints

Covariate matching between treatment and control groups remains a foundational methodology in observational causal inference, designed to replicate the covariate balance of randomized controlled trials. Despite substantial methodological advancements, applied researchers continue to face an inherent trade-off between achieving rigorous balance and preserving the statistical power of the matched sample. We introduce a novel framework for a matching method named Optimal Study Interior Partitioning (OSIP). OSIP frames observational matching as a joint optimization problem governed simultaneously by a maximum cardinality constraint on the number of strata, and a strict fine balance threshold on the matching carried out in each strata. OSIP is based on two-stage optimization algorithm. In both stages it considers a projection of the (high-dimensional) data into one-dimensional space by evaluating an estimated propensity score for each unit. The method is based on enforcing a partition of the input to subsets of consecutive intervals of estimated propensity scores. The novelty arises because we consider this task as a global optimization problem. In the first stage, we employ an optimal partitioning algorithm with respect to an auxiliary (one-dimensional) goal function. In the second phase, we go back to the original distance function and apply heuristics approaches locally optimizing the solution from the first stage. We compare our framework with established methodologies from the literature across well-studied empirical benchmarks. We use the standard benchmark of Lindner cardiovascular dataset to demonstrate that OSIP provides a highly efficient sweet spot solution for the corresponding multi-criteria goal. This conclusion is verified by considering also other benchmark datasets.
Read the full abstract

Statistical Computation

Neural composite likelihood estimation: simulation based inference for time series

Simulation based inference (SBI) circumvents the challenge of intractable likelihoods by using a simulator that generates data given parameter values. For instance, neural likelihood estimation (NLE) estimates the likelihood function by training a neural network to perform conditional density estimation on simulated data given corresponding parameters. However such density estimation is only feasible for relatively low dimensional data. We extend the scalability of SBI methods to a higher dimensional problem: long sequences with a complex dependency structure. We introduce Neural Composite Likelihood Estimation (NCLE). This divides the sequence into smaller, equal-sized batches. Instead of training NLE to estimate the likelihood for an entire sequence, we estimate the likelihood for each batch separately. The product of these forms an approximate composite likelihood (CL), and we perform frequentist inference using methods from the CL literature: we get a point estimate from maximising the approximate CL and obtain confidence intervals by estimating the Godambe information matrix. We demonstrate the effectiveness of NCLE with experiments on time series models.
Read the full abstract

Statistical Computation · Cross-listed from stat.ME

Model-based estimation and imputation with torus missing values

This paper addresses the problem of parameter estimation and model-based imputation for multivariate circular data lying on a p-dimensional torus in the presence of missing values. Actually, the periodic nature of the sample space invalidates conventional imputation techniques designed for Euclidean data. Then, we propose a general framework for maximum likelihood estimation under the wrapped elliptically symmetric family of distributions, with particular interest in the multivariate wrapped normal distribution, when the missing data mechanism is ignorable. The methodology leverages the conditional properties of the elliptically symmetric distributions on the unwrapped space, embedding the imputation of missing torus data into an Expectation-Maximization algorithm that treats both the wrapping coefficients and the missing entries as latent variables. Derivation of both the E and M steps is detailed and a working algorithm is discussed. Imputation methods are also taken into account. The finite-sample performance of the maximum likelihood estimator under ignorable missingness is assessed through Monte Carlo simulations under the wrapped normal specification. The methodology is also illustrated on data with artificially introduced missingness.
Read the full abstract

Statistical Computation

Scentree: a framework for generating scenario trees for multistage stochastic programming

We present scentree, an open-source Python package for constructing a scenario fan and a scenario tree for multistage stochastic programming from historical data. It combines machine learning and multivariate time series models to obtain a scenario fan that captures inter-stage dependencies in the stochastic processes. This scenario fan is subsequently transformed into a scenario tree suitable for multistage stochastic optimization, providing a flexible and extensible framework for uncertainty modeling. A key contribution is the automation of the complete workflow, including model selection, parameter estimation, scenario fan generation, and scenario tree construction. Scentree does not rely on assumptions about the underlying data distribution, reducing the statistical expertise required to produce a scenario tree. Furthermore, it is agnostic to the specific multistage stochastic problem to be solved.
Read the full abstract

Statistical Computation · Cross-listed from stat.ME

Robust Dual-Regularized Variable Selection under Outlier Contamination

Real data often contain unusual observations that can exert disproportionate effects on variable selection, especially in complex predictor settings. We propose a two-stage sparse median outer product of gradients (smOPG) method for variable selection in single index models with outlier contamination. We first estimate sparse local gradients via \(\ell_1\)-penalized local median regression and then recover the active predictor set from a rank-one sparse approximation of the resulting gradient matrix using regularized singular value decomposition. The combination of median regression and local weighting provides robustness to both response outliers and leverage points. Extensive simulations across varying dimensions and contamination mechanisms demonstrate the favorable variable selection performance of smOPG relative to existing methods. Applications to air pollution and genomic data demonstrate practical utility, while theory establishes active-set recovery without requiring selection consistency of individual local regressions.
Read the full abstract

Statistical Computation

Repulsive normalizing flow mixtures for adaptive importance sampling: reliability analysis of complex systems

Accurate rare-event estimation can be computationally expensive. Classical adaptive importance sampling (IS) schemes often rely on restrictive proposal families and can struggle under multiple failure modes. We propose FAMIS, a flow-based multiple importance sampling (MIS) framework that learns a nonuniform mixture of normalizing flow proposals for rare event estimation. The method does not require presampled failure data or prior knowledge of the number, location, or geometry of the failure modes. Instead, it adaptively learns the mixture through sequential evaluations of the limit state function. To guide training toward the failure domain, FAMIS uses a smooth rare-event surrogate and a tempered target sequence. A defensive exploration mixture improves early-stage coverage, a Rao Blackwellized update adapts the mixture weights, and a Jensen-Shannon repulsion term promotes separation and diversity among the base components. The final failure probability is computed with a deterministic-mixture MIS estimator. Numerical experiments demonstrate that FAMIS accurately approximates quasi-optimal IS densities with fewer training samples and model evaluations, providing stable variance reduction across complex reliability problems.
Read the full abstract

Thursday 17 September

Statistical Computation · Cross-listed from stat.ML

Triply-Scalable Equivariant Gaussian Process Modeling

Gaussian processes (GPs) provide principled probabilistic predictions while encoding prior knowledge, including equivariances. Yet, their use in large-scale scientific problems is limited by computational cost. Equivariant neural networks are common but typically lack the uncertainty quantification offered by GPs, which is valuable in applications such as molecular research. High-dimensional inputs and large symmetry groups further demand scalability. We establish results pertaining to the interplay of GP equivariance and conditioning and leverage them to obtain equivariant sparse GPs through suitable mean functions and covariance kernels. We instantiate this framework with a flexible class of integration-free equivariant kernels, yielding scalable and data-efficient GP inference. In particular, we introduce triply scalable equivariant Gaussian processes. We employ equivariant sparse variational Gaussian processes for $\mathrm{SO}(2)$-equivariant vector fields and molecular property prediction. Alongside the SVGP, we develop a matrix-free equivariant full-GP implementation that combines an exact Kronecker reduction with preconditioned conjugate-gradient solves, enabling fast and scalable evaluation of the full joint predictive density. We further compare different approaches for selecting inducing points in the equivariant sparse GP models. Our test cases include synthetic $\mathrm{SO}(2)$-equivariant fields as well as the prediction of electric dipole moments of N-methylformamide based on quantum chemistry simulations, achieving accurate, uncertainty-aware predictions at a fraction of the computational cost of classical GP inference.
Read the full abstract

Statistical Computation

Hankel–Christoffel–Nevai Screening for Bayesian Inverse Problems

Likelihood evaluation in Bayesian inverse problems often requires a forward-model solve. We study candidate screening using a likelihood-weighted prior moment matrix and two associated scores: a Christoffel ratio and a Nevai polynomial average. A finite pilot supplies likelihood information that is reused to rank further prior candidates. We connect finite-pilot matrix error, quantitative Legendre localization, and fixed-budget posterior-mass regret. Conditional likelihoods distinguish feature loss from polynomial and sampling errors, and a filtered Gaussian construction permits controlled feature selection and localization. In a nonlinear function-coefficient PDE, equal-budget comparisons separate methods using scalar likelihoods from surrogates using forward outputs. The Christoffel ratio outperforms the tested direct likelihood regressions, but forward-response surrogates capture more posterior mass at small retention budgets. Independent-pilot tests and reference resampling quantify two distinct sources of uncertainty. Controlled Gaussian experiments illustrate finite-filter truncation and the benefit of selecting relevant features. These results identify when moment geometry provides useful screening information, without asserting universal superiority over surrogate models or an exact posterior sampling method.
Read the full abstract

Statistical Computation · Cross-listed from stat.ME

Efficient computation of mixture confidence sequences in generalized linear models

Classical confidence intervals, when repeatedly obtained on accumulating data at different sample sizes, produce contradictory inferences with high probability. We propose a simple and efficient strategy for computing, instead, mixture confidence sequences for regression coefficients in generalized linear models under this sequential framework. Simulations demonstrate the computational convenience of our approach and the importance of drawing inferential conclusions based on these anytime-valid tools when observations become available in batches over time. The usefulness of the proposed procedure is also shown in the analysis of streaming data from the American National Automotive Sampling System.
Read the full abstract

All authors Claudia Di Caterina, Luigi Pace, Alessandra Salvan, Nicola Sartori

Statistical Computation · Cross-listed from stat.ML

Online Supervised Dimension Reduction with Random Features: Diagnostics and Computational Trade-offs

Accurate optimization of a supervised spectral objective need not produce an accurate population subspace or a better predictive representation. We investigate these distinctions for Online Kernel Supervised Principal Component Analysis (OKSPCA), which combines a centered cross-moment in finite random-feature coordinates with an Adam-style orthonormal basis update for an established objective. Fixed-map consistency, concentration and perturbation results describe the estimator and its exact subspace; same-target comparisons then assess the practical iterate separately. Across six predictive benchmarks, performance depends on the declared pipeline: replacing the tracker with the exact empirical target leaves the two regression deficits largely unchanged. Direct classification-rank models capture nearly all terminal objective energy on average, but a saved intermediate state exhibits substantial geometric deviation; a controlled sample-size study further separates empirical accuracy from population recovery. In distinct numerical-service workloads, exact on-request computation is faster in the tested classification settings, whereas Adam saves time relative to the tested full thin-SVD service for some dense wider-regression requests, alongside persistent geometric error. These diagnostics limit explanations based solely on terminal optimization accuracy and distinguish numerical cost from quality, rank coverage and freshness; they establish neither practical-tracker convergence nor predictive or deployment benefits from basis availability.
Read the full abstract

Statistical Computation · Cross-listed from stat.ME

Exponential Smoothing for Time Series of Random Objects

Time series of random objects, such as covariance matrices, probability distributions, and functional data, call for forecasting methods that do not rely on standard arithmetic operations. We introduce geodesic exponential smoothing, a generalization of exponential smoothing to time series in Hadamard spaces: the forecast level moves a fixed fraction of the way along the geodesic toward each new observation. The smoothing parameter is estimated by minimizing the average squared distance between observations and their forecasts. We further introduce an innovations mechanism under which each observation has conditional Fréchet mean equal to the current level, providing the metric-space analog of the innovations state-space model. In contrast to autoregressive models for object-valued time series, the framework involves a single scalar parameter, assumes no stationarity, and updates online in constant time per observation. Under this mechanism, we establish sample-path properties of the generative process via the quasilinearization available in Hadamard spaces, and prove almost-sure consistency of the smoothing-parameter estimator. Three real-data applications, spanning covariance-matrix, distributional, and functional time series, assess the forecasting performance of the method against structurally heavier alternatives.
Read the full abstract

All authors Takuo Matsubara, Peiwen Jiang, Wilson Ye Chen, Minh-Ngoc Tran

Statistical Computation · Cross-listed from stat.ME

A Bayesian Bi-Directional Splitting Framework for Variable Selection in Large Datasets

Modern tabular datasets are becoming increasingly large, both in the number of samples and covariates, posing significant challenges for Bayesian variable selection due to the resulting computational burden. While there is extensive literature on scaling Bayesian inference to large numbers of observations or high-dimensional covariate spaces, comparatively little work addresses scalable Bayesian variable selection when both dimensions are large simultaneously in a practical setting. This paper presents a novel Bayesian variable selection framework for efficiently analysing data with a large number of both rows and columns. The proposed framework operates via a divide-and-conquer approach, splitting data into batches along both directions and analysing each batch independently in parallel, after which the results are combined together in a two-phase consensus procedure. Experiments demonstrate the computational gains while retaining strong variable selection performance, successfully identifying relevant signals in noisy, high-dimensional settings despite the loss of information induced by data splitting. Practical guidelines are also provided, including recommendations for tuning parameter choices and effective data partitioning strategies, followed by a real application to an H3N2 influenza dataset.
Read the full abstract

Statistical Computation · Cross-listed from stat.ME

Heterogeneity-calibrated Byzantine-robust distributed composite quantile regression

We study sparse composite quantile regression (CQR) for distributed data with heterogeneous honest sites and Byzantine workers. Honest sites share a common slope but may differ in their covariate distributions, error laws, and quantile intercepts. The proposed heterogeneity-calibrated robust CQR (HC-RCQR) profiles local intercepts and calibrates scores using an approximate inverse profile Hessian. Honest workers transmit the resulting vectors, whereas Byzantine workers may send arbitrary vectors. The server updates the estimate by coordinatewise trimming and soft thresholding. A scalar example shows how unequal honest-site curvatures allow intermediate Byzantine reports to survive trimming and how ideal calibration reduces their possible effect. We also establish nonidentification of the mean of unrestricted honest-site slopes when fault identities are unknown. Under suitable conditions, we establish conditional contraction and support-recovery guarantees. The bound separates score offset, sampling fluctuation, contamination, calibration error, and the Newton remainder. Simulations and a bike-demand study examine performance under heterogeneous data and adversarial messages.
Read the full abstract

Wednesday 16 September

Statistical Computation · Cross-listed from stat.ME

Variable Selection for Infectious Disease Transmission via Penalized Pairwise Accelerated Failure Time Models

In infectious disease transmission, regression and variable selection for time-to-event outcomes are complicated by transmission-induced dependence that violates standard independence assumptions. Pairwise survival analysis addresses this dependence by modeling contact intervals in ordered pairs, defined as the time from the onset of infectiousness in the pair's first individual to their first infectious contact with the second individual. We propose a penalized pairwise accelerated failure time model for contact intervals that accommodates simultaneous covariate effects on infectiousness and susceptibility while jointly representing within-population transmission and infection from external sources. For variable selection with many candidate covariates, we develop adaptive LASSO penalized maximum likelihood estimation with a coordinate descent algorithm, treating distributional baseline parameters as unpenalized nuisance parameters. We establish selection consistency and asymptotic normality of the adaptive LASSO estimator under profile-likelihood regularity conditions in a local asymptotic regime in which information accumulates over all ordered pairs that are ever at risk of transmission during follow-up. Simulation studies show accurate variable selection and improved estimation accuracy compared with the unpenalized estimator and the LASSO. We illustrate the method using 2009 household influenza A (H1N1) surveillance data from Los Angeles County, where the adaptive LASSO produces a sparse, interpretable model and improves the precision of selected covariate effect estimates, including the effect of antiviral prophylaxis on susceptibility.
Read the full abstract

All authors Jeongjin Lee, Patrick M. Schnell, Grzegorz A. Rempala, Eben Kenah

Statistical Computation

A bridge representation of Gaussian Whittle-Matérn fields on compact metric graphs

Gaussian Whittle-Matérn fields form a flexible class of Gaussian processes on compact metric graphs, where spatial dependence is governed by the geometry and connectivity of the network through a fractional-order stochastic partial differential equation. This paper develops a new bridge representation of these fields in the case of half-integer smoothness parameters, when the fields have Markov properties. This representation decomposes the field into a finite-dimensional graph component and independent Whittle-Matérn bridge processes on the individual edges. The resulting decomposition leads to efficient likelihood evaluation, kriging prediction, and simulation methods. We show that this improves numerical stability and can greatly reduce computation time compared to previous methods. A simulation study on a Chicago street-network graph illustrates the computational efficiency of the sampling method and an application to Madrid traffic intensity data demonstrates the practical gains for likelihood-based inference and prediction. The methods are implemented in the R package MetricGraph.
Read the full abstract

Tuesday 15 September

Statistical Computation · Cross-listed from stat.ME

Online Bayesian Model Averaging with Joint Uncertainty Quantification for Models and Regression Coefficients in Binary Regression

Streaming binary response data arise in many applications, including virtual learning platforms where student responses are collected sequentially and estimated probabilities of correct responses may inform future question assignment. Logistic regression provides an interpretable framework for such analysis, but the data may support multiple plausible predictor subsets. Bayesian model averaging (BMA) accounts for such model uncertainty by averaging inference across competing models. Since repeatedly applying BMA to accumulating data can be computationally burdensome, we develop an online implementation based on renewable estimation. At each update, retained summaries and new observations are used to approximate the joint posterior of models and model-specific coefficients without revisiting historical data, enabling point and interval estimates that incorporate model uncertainty. Simulations show that online BMA closely approximates its offline counterpart in coefficient estimation, model and variable importance, prediction, and interval estimation of predictive probabilities, while substantially reducing computational cost. In a virtual learning application, we find substantial uncertainty across competing predictor subsets and large variability in credible interval widths across questions. Our method distinguishes questions with precise estimates from those with substantial uncertainty, providing information beyond point estimates for adaptive question assignment.
Read the full abstract

Statistical Computation · Cross-listed from math.NA

On the computational cost of Stochastic Gradient Langevin Dynamics

Stochastic Gradient Langevin Dynamics (SGLD) reduces the cost of Langevin-based sampling by replacing full-dataset drift evaluations with mini-batch approximations, but the resulting subsampling error may offset this computational saving. We study this trade-off for stochastic differential equations with finite-sum drifts and compare the computational cost of SGLD with that of the Euler-Maruyama (EM) method. For a prescribed mean-square accuracy $\varepsilon^2$, we derive complexity estimates that explicitly track the dependence on the dataset size $m$, mini-batch size $s$, and accuracy parameter $\varepsilon$. The resulting comparison reveals distinct parameter regimes in which either method is preferable. In particular, EM can have lower leading-order cost only in a small-data, aggressive-subsampling regime, whereas SGLD is favoured over most of the remaining parameter space. In the practically relevant regime $s \ll m$, the transition between the two methods occurs at the scale $m \asymp \varepsilon^{-1}$. We complement the theoretical analysis with numerical experiments based on a Gaussian Bayesian inference model, which examine the predicted cost regimes together with the underlying discretisation error and variance estimates.
Read the full abstract

All authors Mateusz B. Majka, Tigran Nagapetyan, Łukasz Szpruch, Yue Wu, Danqi Zhuang

Statistical Computation

Local Epochs, Averaging, and Variable Selection in Federated Lasso

Theoretical analysis and Monte Carlo simulation separate local-epoch effects from averaging and tuning effects in federated Lasso. How much local work should precede averaging when fitting a sparse regression? An orthogonal calculation shows that extra epochs can have no effect while averaging still enlarges the selected set. A correlated two-site construction gives an exact, nonmonotone limiting objective gap and its minimizing epoch count. We then compare coordinate-descent averaging, two thresholding modifications, and adapted FedDualAvg across twelve scenarios and 600 replicates. Methods share a penalty, independent validation samples, selection rules, and resource limits. FedDualAvg generally achieves smaller objective gaps but does not always recover variables better. Thresholding gains depend strongly on selection rules. Epoch effects vary with correlation, signal strength, site allocation, and the outcome measured. These results distinguish faster iteration from better optimization, prediction, and variable selection.
Read the full abstract

Statistical Computation · Cross-listed from math.ST

Exact finite-sample inference for multi-mixed fractional Brownian motion with drift

In this paper we study a linear drift perturbed by a superposition of $m$ independent fractional Brownian motions with known Hurst parameters and a common scale, observed at $N$ equidistant times. Inference for such models is usually asymptotic; we show that here it is exact. We derive the maximum likelihood estimators of the drift $\theta$ and of the scale $\alpha^{2}$ in closed form and obtain their exact finite-sample joint law: $\widehat\theta$ is Gaussian, $N\widehat\alpha^{\,2}/\alpha^{2}$ is chi-square with $N-1$ degrees of freedom, and the two are independent. As this law is free of every model parameter, we deduce Student and chi-square confidence intervals and tests of exact level for every $N\ge2$, whatever the Hurst vector. We also prove that the estimators are uniformly minimum variance unbiased with $\widehat\theta$ attaining the Cramér–Rao bound at every $N$, that both are strongly consistent and asymptotically normal, and that the drift estimators form, in law, a Brownian motion run along their own variance scale. A sharp non-asymptotic bound shows that the accuracy of the drift is governed by the length of the observation window and not by the mesh, and a Monte Carlo study confirms exact coverage, even at small sample sizes, and quantifies what is lost when the Hurst vector is misspecified.
Read the full abstract

Statistical Computation

Optimal Scaling of Langevin Proposals with Generalized Acceptance Rules

Langevin-based Markov chain Monte Carlo (MCMC) algorithms use gradient information to improve sampling, particularly in high dimensions. Classical optimal scaling theory for these algorithms has largely focused on the Metropolis-Hastings (MH) acceptance rule. However, there has been a recent surge in acceptance rules beyond MH for applications spanning differential privacy, stochastic MCMC, diffusion models, and molecular dynamics. We develop optimal scaling results for Langevin proposals employed with generalized acceptance rules belonging to a suitable class. For high-dimensional targets, we recover the usual $O(d^{-1/3})$ scaling, while different acceptance rules lead to different optimal acceptance probabilities.
Read the full abstract

Monday 14 September

Statistical Computation · Cross-listed from cs.LG

NeuralCert: certified computational discovery of extremal mathematical constructions

Neural networks are becoming popular in solving mathematical problems, but stochastic models do not provide mathematical exactness by themselves. This study introduces a discovery-to-certification framework in which high-dimensional variational trial functions are learned in a compact separable representation, spectrally diagnosed and pruned, and then certified exactly through multimodular evaluation. Exact certification makes the numerical proofs fully explicit and independently verifiable. This framework can be run on a standard personal computer. Across three extremal problems, we show that neural optimization can contribute to rigorous mathematics in three distinct ways: by discovering improved constructions, by exposing empirical invariants that lead to proofs, and by revealing optimization barriers whose geometry motivates new analytic or numerical representations. More broadly, these results suggest a path toward AI-assisted mathematics in which flexible computational discovery and exact certification become complementary components of a single rigorous workflow.
Read the full abstract

Statistical Computation

bayprior: Structured Bayesian Prior Elicitation, Conflict Diagnostics, and Regulatory Reporting

Bayesian clinical trials require a well-specified prior distribution, and regulators increasingly expect that specification to be documented, diagnosed for prior-data conflict, and shown to be robust to reasonable alternatives, yet few R tools connect the full workflow of prior construction, validation, and regulatory justification into a single package. We introduce bayprior, an R package and Shiny application that implements structured expert elicitation via quantile matching, moment matching, and the SHELF roulette method across six distribution families; linear and logarithmic expert opinion pooling with Bhattacharyya agreement diagnostics; derivation of meta-analytic-predictive (MAP) priors from historical trial data via random-effects meta-analysis; prior-data conflict assessment using the Box p-value, S-value (surprisal) surprise index, information divergence, Bhattacharyya overlap, and multivariate Mahalanobis distance; hyperparameter sensitivity analysis with tornado and influence heatmap visualisations; and robust, sceptical, and power prior alternatives. A regulatory reporting module generates self-contained HTML, PDF, or Word documents addressing expectations in the FDA's 2026 draft Bayesian methods guidance and in EMA guidance on incorporating external and historical information. The package is demonstrated on a synthetic oncology Phase II trial. bayprior is available on CRAN and includes a fully modular Shiny application for interactive use.
Read the full abstract

Statistical Computation · Cross-listed from stat.ML

Copula Adapted Directed Acyclic Graph for Cluster Representation of Biomedical Data

Diagnostic errors and mislabeling are common in biomedicine, which compromise the reliability of predictive models and data-driven outcomes. Stratifying unlabeled biomedical data based on complex relationships between features eliminates the need for data labels and overcomes the limitations of supervised learning. Traditional clustering methods assume restrictive data distributions, making them suboptimal for capturing complex dependencies in high-dimensional biomedical data. This paper introduces a novel cluster-friendly data presentation framework that integrates the non-Gaussian and non-linear feature dependence of copula models with an ensemble of causal structure discovery (CSD) methods based on Directed Acyclic Graphs (DAGs). While copulas model flexible multivariate distributions by relaxing assumptions related to multivariate normality, linear dependence, and symmetric relationships, an ensemble of DAG-based CSD methods identifies stable causal relationships between features. When clustered using K-means, the new data representation obtained by the proposed copula-adapted DAG (CopDAG) ranks first among the 12 methods in normalized clustering accuracy and adjusted Rand index across 16 biomedical datasets. Our CopDAG method predicts ground-truth class labels directly from feature relationships without data annotations and supervised learning, while also providing cluster visualizations and explainable causal structures of the biomedical data features.
Read the full abstract