Papers

We present MMbeddings, a probabilistic embedding approach that reinterprets categorical embeddings through the lens of nonlinear mixed models, effectively bridging classical statistical theory with modern deep learning. By treating embeddings as latent random effects within a variational autoencoder framework, our method substantially decreases the number of parameters – from the conventional embedding approach of cardinality $\times$ embedding dimension, which quickly becomes infeasible with large cardinalities, to a significantly smaller, cardinality-independent number determined primarily by the encoder architecture. This reduction dramatically mitigates overfitting and computational burden in high-cardinality settings. Extensive experiments on simulated and real datasets, encompassing collaborative filtering and tabular regression tasks using varied architectures, demonstrate that MMbeddings consistently outperforms traditional embeddings, underscoring its potential across diverse machine learning applications.

We consider the problem of learning the interaction strength between the nodes of a network based on dependent binary observations residing on these nodes, generated from a Markov Random Field (MRF). Since these observations can possibly be corrupted/noisy in larger networks in practice, it is important to robustly estimate the parameters of the underlying true MRF to account for such inherent contamination in observed data. However, it is well-known that classical likelihood and pseudolikelihood based approaches are highly sensitive to even a small amount of data contamination. So, in this paper, we propose a density power divergence (DPD) based robust generalization of the computationally efficient maximum pseudolikelihood (MPL) estimator of the interaction strength parameter, and derive its rate of consistency under the pure model. Along the way, we establish consistency and asymptotics for a class of general $Z$-estimators, covering our proposed DPD based estimators, under flexible assumptions that hold for a substantial class of standard models. To the best of our knowledge, these are the first central limit theorems for the class of general $Z$-estimators in such settings. Moreover, we show that the gross error sensitivities of the proposed DPD based estimators are significantly smaller than that of the MPL estimator, thereby theoretically justifying the greater (local) robustness of the former under contaminated settings. Finally, we demonstrate the superior (finite sample) performance of the DPD based variants over the traditional MPL estimator in a number of synthetically generated contaminated network datasets, and apply them to learn the network interaction strength in several real datasets from diverse domains of social science, neurobiology and genomics.

We consider the problem of estimating the trace and diagonal entries of an N-order tensor (where $N \geq 2$) under the framework where the tensor can only be accessed through tensor-vector multiplication. The aim is to estimate the tensor's diagonal entries and trace by minimizing the number of tensor-vector queries. The seminal work of Hutchinson and its extended version due to Bekas et al. give unbiased estimates of the trace and diagonal elements of a given matrix, respectively, using matrix-vector queries. However, to the best of our knowledge, no analogous results are known for estimating the trace and diagonal entries of higher-order tensors using tensor-vector queries. This paper addresses this gap and presents unbiased estimators for the trace and diagonal entries of a tensor under this model. Our proposed methods can be seen as generalizations of Hutchinson's and Bekas et al.'s estimators and reduce to their estimators when N = 2. We provide a rigorous theoretical analysis of our proposals and complement it with supporting simulations.

We study a class of degenerate diffusion generators arising in sequential testing and quickest detection problems with partial information. The observation process is driven by $k$ independent Brownian motions, while the hidden state takes $n+1$ values, with $k<n$. After transforming to posterior likelihood coordinates, we analyze Hörmander's condition both in the absence of state switching (testing) and in its presence (detection). We characterize hypoellipticity in the testing case and give explicit sufficient conditions in the detection case. We further study the stationary posterior operator, its parabolic extension, the joint observation-posterior operator, and its parabolic extension; their Hörmander conditions need not coincide. We characterize the relationships among these operators under our main structural regimes and discuss their probabilistic consequences and the regularity of the associated optimal stopping problem.

Coordinate denoising has emerged as a promising method for 3D molecular pretraining due to its theoretical connection to learning molecular force field. However, existing denoising methods rely on oversimplied molecular dynamics that assume atomic motions to be isotropic and homoscedastic. To address these limitations, we propose a novel denoising framework AniDS: Anisotropic Variational Autoencoder for 3D Molecular Denoising. AniDS introduces a structure-aware anisotropic noise generator that can produce atom-specific, full covariance matrices for Gaussian noise distributions to better reflect directional and structural variability in molecular systems. These covariances are derived from pairwise atomic interactions as anisotropic corrections to an isotropic base. Our design ensures that the resulting covariance matrices are symmetric, positive semi-definite, and SO(3)-equivariant, while providing greater capacity to model complex molecular dynamics. Extensive experiments show that AniDS outperforms prior isotropic and homoscedastic denoising models and other leading methods on the MD17 and OC22 benchmarks, achieving average relative improvements of 8.9% and 6.2% in force prediction accuracy. Our case study on a crystal and molecule structure shows that AniDS adaptively suppresses noise along the bonding direction, consistent with physicochemical principles. Our code is available at https://github.com/ZeroKnighting/AniDS.

We study low-rank matrix regression in settings where matrix-valued predictors and scalar responses are observed across multiple individuals. Rather than assuming a fully homogeneous coefficient matrices across individuals, we accommodate shared low-dimensional structure alongside individual-specific deviations. To this end, we introduce a tensor-structured homogeneity pursuit framework, wherein each coefficient matrix is represented as a product of shared low-rank subspaces and individualized low-rank loadings. We propose a scalable estimation procedure based on scaled gradient descent, and establish non-asymptotic bounds demonstrating that the proposed estimator attains improved convergence rates by leveraging shared information while preserving individual-specific signals. The framework is further extended to incorporate scaled hard thresholding for recovering sparse latent structures, with theoretical guarantees in both linear and generalized linear model settings. Our approach provides a principled middle ground between fully pooled and fully separate analyses, achieving strong theoretical performance, computational tractability, and interpretability in high-dimensional multi-individual matrix regression problems.

Alignment of large language models (LLMs) has predominantly relied on pairwise preference optimization, where annotators select the better of two responses to a prompt. While simple, this approach overlooks the opportunity to learn from richer forms of human feedback, such as multiway comparisons and top-$k$ rankings. We introduce Ranked Choice Preference Optimization (RCPO), a unified framework that bridges preference optimization with (ranked) choice modeling via maximum likelihood estimation. RCPO supports both utility-based and rank-based models, subsumes several pairwise methods (such as DPO and SimPO) as special cases, and provides principled training objectives for richer feedback formats. We instantiate this framework with two representative models (Multinomial Logit and Mallows-RMJ). Experiments on Llama-3-8B-Instruct, Gemma-2-9B-it, and Mistral-7B-Instruct across in-distribution and out-of-distribution settings show that RCPO consistently outperforms competitive baselines. RCPO shows that directly leveraging ranked preference data, combined with the right choice models, yields more effective alignment. It offers an extensible foundation for incorporating (ranked) choice modeling into LLM training.

Estimators in statistics and machine learning must typically trade off between efficiency, having low variance for a fixed target, and distributional robustness, such as multiaccuracy, or having low bias over a range of possible targets. In this paper, we consider a simple estimator, ridge boosting: starting with any initial predictor, perform a single boosting step with (kernel) ridge regression. Surprisingly, we show that ridge boosting simultaneously achieves both efficiency and distributional robustness: for target distribution shifts that lie within an RKHS unit ball, this estimator maintains low bias across all such shifts and has variance at the semiparametric efficiency bound for each target. In addition to bridging otherwise distinct research areas, this result has immediate practical value. Since ridge boosting uses only data from the source distribution, researchers can train a single model to obtain both robust and efficient estimates for multiple target estimands at the same time, eliminating the need to fit separate semiparametric efficient estimators for each target. We assess this approach through simulations and an application estimating the age profile of retirement income.

This study introduces the Spatial Health and Population Estimator (SHAPE), a spatial microsimulation framework that applies hierarchical iterative proportional fitting (IPF) to estimate two health risk behaviors and eleven health outcomes across multiple spatial scales. SHAPE was evaluated using county-level direct estimates from the Behavioral Risk Factor Surveillance System (BRFSS) and both county and census tract level data from CDC PLACES for New York (2021) and Florida (2019). Results show that SHAPE's SAEs are moderately consistent with BRFSS (average Pearson's correlation coefficient r of about 0.5), similar to CDC PLACES (average r of about 0.6), and are strongly aligned with CDC PLACES model-based estimates at both county (average r of about 0.8) and census tract (average r of about 0.7) levels. SHAPE is an open, reproducible, and transparent framework programmed in R that meets a need for accessible SAE methods in public health.

Positivity violations, which occur when some subgroups either always or never receive a treatment of interest, pose significant challenges for causal effect estimation with observational data. Recent balancing weight methods have proved to be highly effective in confounding control, however their utility is diminished in the presence of positivity violations, resulting in bias and excess variance. Approaches that deal with positivity violations, on the other hand, work by targeting a modified estimand that may be misaligned with the original research question. To address these challenges, we propose a novel balancing weights approach, which mitigates positivity violations while attempting to retain the original estimand by a targeted relaxation of the balancing constraints. Our proposed weighted estimator is consistent for the original estimand when either 1) the implied propensity score model is correct; or 2) all treatment effect modifiers are balanced to the target population. When these conditions do not hold, our estimator is consistent for a slightly modified treatment effect estimand. Furthermore, our proposed weighted estimator has reduced asymptotic variance when positivity does not hold. We evaluate our approach through applications to synthetic data, an observational study, and when transporting a treatment effect from a randomized trial.

Active-controlled trials with non-inferiority objectives are often used when effective interventions are available, but new options may offer advantages or meet public health needs. In these trials, participants are randomized to an experimental intervention or an active control. The traditional non-inferiority criterion requires that the new intervention preserve a substantial proportion of the active control effect. A key challenge is the absence of a placebo arm, which necessitates reliance on historical data to estimate the active control effect and assumptions about how well this effect applies to the target population. Another challenge arises when the active control is highly effective, as the new intervention may still be valuable even if it does not meet the traditional criterion. This has motivated alternative criteria based on sufficient efficacy relative to a hypothetical placebo. In this work, we propose a general framework for designing and evaluating non-inferiority trials that integrates all existing analytical methods and accommodates both traditional and alternative success criteria. The framework enables the systematic comparison of methods in terms of type I error, power, and robustness to misspecification of the active control effect. We illustrate its applicability in the design of a future HIV prevention trial with a highly effective active control. In this application, our framework identifies methods that provide greater efficiency and robustness than commonly used approaches and demonstrates practical advantages of the alternative non-inferiority criterion. Overall, this framework offers a comprehensive toolkit for rigorous non-inferiority trial design, supporting method selection and the evaluation of new interventions.

Generative modeling has emerged as a powerful paradigm for representation learning, but its direct applicability to challenging fields like medical imaging remains limited: mere generation, without task alignment, fails to provide a robust foundation for clinical use. We propose MAGIC-Flow, a conditional multiscale normalizing flow architecture that performs generation and classification within a single modular framework. The model is built as a hierarchy of invertible and differentiable bijections, where the Jacobian determinant factorizes across sub-transformations. We show how this ensures exact likelihood computation and stable optimization, while invertibility enables explicit visualization of sample likelihoods, providing an interpretable lens into the model's reasoning. By conditioning on class labels, MAGIC-Flow supports controllable sample synthesis and principled class-probability estimation, effectively aiding both generative and discriminative objectives. We evaluate MAGIC-Flow against top baselines using metrics for similarity, fidelity, and diversity. Across multiple datasets, it addresses generation and classification under scanner noise, and modality-specific synthesis and identification. Results show MAGIC-Flow creates realistic, diverse samples and improves classification. MAGIC-Flow is an effective strategy for generation and classification in data-limited domains, with direct benefits for privacy-preserving augmentation, robust generalization, and trustworthy medical AI.

Problem definition: Hospitals use separate pathways to identify and respond to deteriorating ward patients. At Children's Healthcare of Atlanta (CHOA), these include Watcher surveillance, Critical Care Outreach Team assessment, and transfer to the intensive care unit. When several patients need attention simultaneously, however, no common process determines who should receive which response under limited staffing and capacity. We study how to coordinate these actions to support earlier escalation without unnecessary monitoring. Methodology/results: We model ward escalation decisions using a multi-action restless multi-armed bandit model and develop Multi-action Index Learning and Optimization Policy (MILOP) as a computational solution approach. MILOP trains a neural network to generate a high-quality index policy using a decision-focused learning approach. We train and evaluate MILOP using historical CHOA data and compare it with two practical baselines. MILOP recommends an escalation before 87.4% of critical or emergent intensive-care transfers, compared with 79.5% for the stronger baseline, while using 13.5% less monitoring resources. MILOP provides a median of 19.2 hours lead time before transfer, earlier than either baseline. Moreover, 99.3% of MILOP's proactive transfer recommendations are preceded by Watcher surveillance or outreach-team assessment, demonstrating a desirable escalation of care. Managerial implications: MILOP can help leaders focus scarce clinical attention, gain more time to respond before severe deterioration, and reduce monitoring burden compared to existing practices. Broadly, hospitals may improve escalation decisions by coordinating existing pathways rather than managing each separately or using risk scores as stand-alone decision rules.

Learning mappings between functional spaces, also known as function-on-function regression, is a fundamental problem in functional data analysis with broad applications, including spatiotemporal forecasting, curve prediction, and climate modeling. Existing approaches often struggle to capture complex nonlinear relationships and/or provide reliable uncertainty quantification when data are noisy, sparse, or irregularly sampled. To address these challenges, we propose Deep Gaussian Processes for Functional Maps (DGPFM). Our method constructs a sequence of GP-based linear and nonlinear transformations directly in function space, leveraging kernel integral transforms, GP conditional means, and nonlinear activations sampled from Gaussian processes. A key insight enables a simplified and flexible implementation: under fixed evaluation locations, discrete approximations of kernel integral transforms reduce to direct functional integral transforms, allowing seamless integration of diverse transform designs. To support scalable probabilistic inference, we adopt inducing points and whitening transformations within a variational learning framework. Empirical results on both real-world and synthetic benchmark datasets demonstrate the advantages of DGPFM in terms of predictive accuracy and uncertainty calibration.

Decomposing prediction uncertainty into aleatoric (irreducible) and epistemic (reducible) components is critical for the reliable deployment of machine learning systems. While the mutual information between the response variable and model parameters is a principled measure for epistemic uncertainty, it requires access to the parameter posterior, which is computationally challenging to approximate. Consequently, practitioners often rely on probabilistic predictions from deep ensembles to quantify uncertainty, which have demonstrated strong empirical performance. However, a theoretical understanding of their success from a frequentist perspective remains limited. We address this gap by first considering a bootstrap-based estimator for epistemic uncertainty, which we prove is asymptotically correct. Next, we connect deep ensembles to the bootstrap estimator by decomposing it into data variability and training stochasticity; specifically, we show that deep ensembles capture the training stochasticity component. Through empirical studies, we show that this stochasticity component constitutes the majority of epistemic uncertainty, thereby explaining the effectiveness of deep ensembles.

Sparse variable selection improves interpretability and generalization in high-dimensional learning by selecting a small subset of informative features. Recent advances in Mixed Integer Programming (MIP) have enabled solving large-scale non-private sparse regression - known as Best Subset Selection (BSS) - with millions of variables in minutes. However, extending these algorithmic advances to the setting of Differential Privacy (DP) has remained largely unexplored. In this paper, we introduce two new pure differentially private estimators for sparse variable selection, levering modern MIP techniques. Our framework is general and applies broadly to problems like sparse regression or classification, and we provide theoretical support recovery guarantees in the case of BSS. Inspired by the exponential mechanism, we develop structured sampling procedures that efficiently explore the non-convex objective landscape, avoiding the exhaustive combinatorial search in the exponential mechanism. We complement our theoretical findings with extensive numerical experiments, using both least squares and hinge loss for our objective function, and demonstrate that our methods achieve state-of-the-art empirical support recovery, outperforming competing algorithms in settings with up to $p=10^4$.

This paper studies prediction with multiple candidate models, where the goal is to combine their outputs. This task is especially challenging in heterogeneous settings, where different models may be better suited to different inputs. We propose input adaptive Bayesian Model Averaging (IA-BMA), a Bayesian method that assigns model weights conditional on the input. IA-BMA employs an input adaptive prior, and yields a posterior distribution that adapts to each prediction, which we estimate with amortized variational inference. We derive formal guarantees for its performance, relative to any single predictor selected per input. We evaluate IABMA across regression and classification tasks, studying data from personalized cancer treatment, credit-card fraud detection, and UCI datasets. IA-BMA consistently delivers more accurate and better-calibrated predictions than both non-adaptive baselines and existing adaptive methods.

Sampling from constrained statistical distributions is a fundamental task in various fields including Bayesian statistics, computational chemistry, and statistical physics. This article considers the cases where the constrained distribution is described by an unconstrained density, as well as additional equality and/or inequality constraints, which often make the constraint set nonconvex. Existing methods for nonconvex constraint set $\Sigma \subset \mathbb{R}^d$ defined by equality or inequality constraints commonly rely on costly projection steps. Moreover, they cannot handle equality and inequality constraints simultaneously as each method only specialized in one case. In addition, rigorous and quantitative convergence guarantee is often lacking. In this paper, we introduce Overdamped Langevin with LAnding (OLLA), a new framework that can design overdamped Langevin dynamics accommodating both equality and inequality constraints. The proposed dynamics also deterministically corrects trajectories along the normal direction of the constraint surface, thus obviating the need for explicit projections. We show that, under suitable regularity conditions on the target density and $\Sigma$, OLLA converges exponentially fast in $W_2$ distance to the constrained target density $\rho_\Sigma(x) \propto \exp(-f(x))d\sigma_\Sigma$. Lastly, through experiments, we demonstrate the efficiency of OLLA compared to projection-based constrained Langevin algorithms and their slack variable variants, highlighting its favorable computational cost and reasonable empirical mixing.

The classic Mallows model is a foundational tool for modeling user preferences. However, it has limitations in capturing real-world scenarios, where users often focus only on a limited set of preferred items and are indifferent to the rest. To address this, extensions such as the top-k Mallows model have been proposed, aligning better with practical applications. In this paper, we address several challenges related to the generalized top-k Mallows model, with a focus on analyzing buyer choices. Our key contributions are: (1) a novel sampling scheme tailored to generalized top-k Mallows models, (2) an efficient algorithm for computing choice probabilities under this model, and (3) an active learning algorithm for estimating the model parameters from observed choice data. These contributions provide new tools for analysis and prediction in critical decision-making scenarios. We present a rigorous mathematical analysis for the performance of our algorithms. Furthermore, through extensive experiments on synthetic data and real-world data, we demonstrate the scalability and accuracy of our proposed methods, and we compare the predictive power of Mallows model for top-k lists compared to the simpler Multinomial Logit model.

Single-cell technologies generate high-dimensional point clouds of cells, enabling detailed characterization of complex patient states and treatment responses. Yet each patient is represented by an irregular point cloud rather than a simple vector, making it difficult to directly quantify and compare biological differences between individuals. Nonlinear methods such as kernels and neural networks achieve predictive accuracy but act as black boxes, offering little biological interpretability. To address these limitations, we adapt the Linear Optimal Transport (LOT) framework to this setting, embedding irregular point clouds into a fixed-dimensional Euclidean space while preserving distributional structure. This embedding provides a principled linear representation that preserves optimal transport geometry while enabling downstream analysis. It also forms a registration between any two patients, enabling direct comparison of their cellular distributions. Within this space, LOT enables: (i) accurate and interpretable classification of COVID-19 patient states, where classifier weights map back to specific markers and spatial regions driving predictions; and (ii) synthetic data generation for patient-derived organoids, exploiting the linearity of the LOT embedding. LOT barycenters yield averaged cellular profiles representing combined conditions or samples, supporting drug interaction testing. Together, these results establish LOT as a unified framework that bridges predictive performance, interpretability, and generative modeling. By transforming heterogeneous point clouds into structured embeddings directly traceable to the original data, LOT opens new opportunities for understanding immune variation and treatment effects in high-dimensional biological systems.

Causal discovery from observational data is a fundamental task in artificial intelligence, with far-reaching implications for decision-making, predictions, and interventions. Despite significant advances, existing methods can be broadly categorized as constraint-based or score-based approaches. Constraint-based methods offer rigorous causal discovery but are often hindered by small sample sizes, while score-based methods provide flexible optimization but typically forgo explicit conditional independence testing. This work explores a third avenue: developing differentiable $d$-separation scores, obtained through a percolation theory using soft logic. This enables the implementation of a new type of causal discovery method: gradient-based optimization of conditional independence constraints. Empirical evaluations demonstrate the robust performance of our approach in low-sample regimes, surpassing traditional constraint-based and score-based baselines on a real-world dataset. Code and data of the proposed method are publicly available at https://github$.$com/PurdueMINDS/DAGPA.

We study the problem of Offline Safe Reinforcement Learning (OSRL), where the goal is to learn a reward-maximizing policy from fixed data under a cumulative cost constraint. We propose a novel OSRL approach that frames the problem as a minimax objective and solves it by combining offline RL with online optimization algorithms. We prove the approximate optimality of this approach when integrated with an approximate offline RL oracle and no-regret online optimization. We also present a practical approximation that can be combined with any offline RL algorithm, eliminating the need for offline policy evaluation. Empirical results on the DSRL benchmark demonstrate that our method reliably enforces safety constraints under stringent cost budgets, while achieving high rewards. The code is available at https://github.com/yassineCh/O3SRL.

Safety-critical applications of machine learning require uncertainty estimates that support reliable worst-case analysis. Neural networks (NNs) provide expressive function approximation, while Gaussian processes (GPs) offer principled probabilistic uncertainty, but both face limitations in this setting: most modern neural architectures lack tractable worst-case error bounds, and Gaussian processes become computationally expensive at scale. For uncertainty to be interpretable, a central requirement is distance-awareness: that is, uncertainty increases with the distance between a test input and the nearest relevant training data. We present a general framework for worst-case distance-aware error bounds for NNs that combine dense layers with spline-based components. Our approach establishes error bounds that are both distance-aware, reflecting proximity of a test point to its nearest training data, and worst-case, providing deterministic guarantees under known Lipschitz constraints rather than probabilistic assumptions. Our algorithm, K-DAREK (Distance-Aware Error for Kurkova-Kolmogorov Networks), provides efficient and interpretable uncertainty quantification for NNs. K-DAREK is about four times faster and ten times more computationally efficient than an ensemble of KANs, 8.6 times more scalable than GP, eliminates up to 8.2% error-bound violation rate observed in DAREK, and reduces the average collision rate from 1.8% to 1.1% in the multi-agent safe control experiment. On high-dimensional real-world regression tasks (e.g., Real Estate Valuation), K-DAREK preserves distance-aware error bounds and achieves zero coverage violations, addressing the overgeneralization and inducing-point-coverage limitations exhibited by SNGP and DUE, respectively.

Selecting an appropriate evaluation metric for classifiers is crucial for model comparison, parameter optimization, and deployment decisions, yet there is no consensus on a broadly accepted evaluation paradigm explicitly aligned with Total Classification Cost (TCC) minimization. At the same time, class imbalance is often treated as a problem to be corrected per se, potentially causing misalignments with TCC minimization. To address these limitations, (i) we define Weighted Accuracy (WA), an evaluation metric for binary classifiers with a straightforward interpretation as a weighted version of accuracy and (ii) we propose a general reweighting framework for handling class imbalance in cost-sensitive scenarios, providing an alternative to resampling techniques. This framework applies to any evaluation metric or loss function that can be expressed as a linear combination of example-dependent quantities; it enables meaningful comparison of evaluation results obtained on different datasets and accounts for discrepancies between the development dataset, used for training, validation, and testing, and the target dataset, where the model will be deployed. Within this framework, we derive the conditions under which standard rebalancing techniques remain coherent with TCC minimization, and when they may instead become misleading. We prove that, under example-independent Unit Classification Costs, maximizing WA is equivalent to minimizing TCC. Finally, we analyze the robustness of WA in realistic example-dependent cost scenarios by studying its correlation with TCC across a broad range of class imbalance and cost regimes. The results show that WA maintains robust alignment with TCC across almost all examined scenarios.

Text watermarking plays a crucial role in ensuring the traceability and accountability of large language model (LLM) outputs and mitigating misuse. While promising, most existing methods assume perfect pseudorandomness. In practice, repetition in generated text induces collisions that create structured dependence, compromising Type I error control and invalidating standard analyses. We introduce a statistical framework that captures this structure through a hierarchical two-layer partition. At its core is the concept of minimal units – the smallest groups treatable as independent across units while permitting dependence within. Using minimal units, we define a non-asymptotic efficiency measure and cast watermark detection as a minimax hypothesis testing problem. Applied to Gumbel-max and inverse-transform watermarks, our framework produces closed-form optimal rules. It explains why discarding repeated statistics often improves performance and shows that within-unit dependence must be addressed unless degenerate. Both theory and experiments confirm improved detection power with rigorous Type I error control. These results provide the first principled foundation for watermark detection under imperfect pseudorandomness, offering both theoretical insight and practical guidance for reliable tracing of model outputs.

Understanding the predictions made by deep learning models remains a central challenge, especially in high-stakes applications. A promising approach is to equip models with the ability to answer counterfactual questions – hypothetical “what if?” scenarios that go beyond the observed data and provide insight into a model reasoning. In this work, we introduce the notion of causal interpretability, which formalizes when counterfactual queries can be evaluated from a specific class of models and observational data. We analyze two common model classes – blackbox and concept-based predictors – and show that neither is causally interpretable in general. To address this gap, we develop a framework for building models that are causally interpretable by design. Specifically, we derive a complete graphical criterion that determines whether a given model architecture supports a given counterfactual query. This leads to a fundamental tradeoff between causal interpretability and predictive accuracy, which we characterize by identifying the unique maximal set of features that yields an interpretable model with maximal predictive expressiveness. Experiments corroborate the theoretical findings.

Reinforcement learning (RL) promises to solve long-horizon tasks even when training data contains only short fragments of the behaviors. This experience stitching capability is often viewed as the purview of temporal difference (TD) methods. However, outside of small tabular settings, trajectories never intersect, calling into question this conventional wisdom. Moreover, the common belief is that Monte Carlo (MC) methods should not be able to recombine experience, yet it remains unclear whether function approximation could result in a form of implicit stitching. The goal of this paper is to empirically study whether the conventional wisdom about stitching actually holds in settings where function approximation is used. We empirically demonstrate that Monte Carlo (MC) methods can also achieve experience stitching. While TD methods do achieve slightly stronger capabilities than MC methods (in line with conventional wisdom), that gap is significantly smaller than the gap between small and large neural networks (even on quite simple tasks). We find that increasing critic capacity effectively reduces the generalization gap for both the MC and TD methods. These results suggest that the traditional TD inductive bias for stitching may be less necessary in the era of large models for RL and, in some cases, may offer diminishing returns. Additionally, our results suggest that stitching, a form of generalization unique to the RL setting, might be achieved not through specialized algorithms (temporal difference learning) but rather through the same recipe that has provided generalization in other machine learning settings (via scale). Project website: https://michalbortkiewicz.github.io/golden-standard/

Machine-learning methods in biochemistry commonly represent molecules as graphs of pairwise intermolecular interactions for property and structure predictions. Most methods operate on a single graph, typically the minimal free energy (MFE) structure, for low-energy ensembles (conformations) representative of structures at thermodynamic equilibrium. We introduce a thermodynamically parameterized exponential-family random graph (ERGM) embedding that models molecules as Boltzmann-weighted ensembles of interaction graphs. We evaluate this embedding on SELEX datasets, where experimental biases (e.g., PCR amplification or sequencing noise) can obscure true aptamer-ligand affinity, producing anomalous candidates whose observed abundance diverges from their actual binding strength. We show that the proposed embedding enables robust community detection and subgraph-level explanations for aptamer ligand affinity, even in the presence of biased observations. This approach may be used to identify low-abundance aptamer candidates for further experimental evaluation.

We introduce Deep Jump Gaussian Processes (DJGP), a novel method for surrogate modeling of a piecewise continuous function on a high-dimensional domain. DJGP addresses the limitations of conventional Jump Gaussian Processes (JGP) in high-dimensional input spaces by integrating region-specific, locally linear projections with JGP modeling. These projections employ region-dependent matrices to capture local low-dimensional subspace structures, making them well suited to the inherently localized modeling behavior of JGPs, a variant of local Gaussian processes. To control model complexity, we place a Gaussian Process prior on the projection matrices, allowing them to evolve smoothly across the input space. The projected inputs are then modeled with a JGP to capture piecewise continuous relationships with the response. This yields a distinctive two-layer deep learning of GP/JGP. We further develop a scalable variational inference algorithm to jointly learn the projection matrices and JGP hyperparameters. Rigorous theoretical analysis and extensive empirical studies are provided to justify the proposed approach. In particular, we derive an oracle error bound for DJGP and decompose it into four distinct sources of error, which are then linked to practical implications. Experiments on synthetic and benchmark datasets demonstrate that DJGP achieves superior predictive accuracy and more reliable uncertainty quantification compared with existing methods.

Spectral decomposition of linear operators plays a central role in many areas of machine learning and scientific computing. Recent work has explored training neural networks to approximate eigenfunctions of such operators, enabling scalable approaches to representation learning, dynamical systems, and partial differential equations (PDEs). In this paper, we revisit a classical optimization framework from the computational physics literature known as the orbital minimization method (OMM), originally proposed in the 1990s for solving eigenvalue problems in computational chemistry. We provide a simple linear-algebraic proof of the consistency of the OMM objective, and reveal connections between this method and several ideas that have appeared independently across different domains. Our primary goal is to justify its broader applicability in modern learning pipelines. We adapt this framework to train neural networks to decompose positive semidefinite operators, and demonstrate its practical advantages across a range of benchmark tasks. Our results highlight how revisiting classical numerical methods through the lens of modern theory and computation can provide not only a principled approach for deploying neural networks in numerical simulation, but also effective and scalable tools for machine learning.