Papers

Large language models (LLMs) have been rapidly improving in long-context tasks, powered by Chain-of-Thought (CoT) reasoning. However, the internal mechanisms underlying this improvement remain unclear. We investigate these mechanisms through a needle-in-a-haystack (NIAH) counting task, where an LLM is asked to count the number of records dispersed in a long text. Across twelve model comparison groups, Thinking (or reasoning) improves counting accuracy over Non-thinking, with pronounced gains at larger counts. This motivates our mechanistic analysis, which identifies two contrasting mechanisms: (i) broad retrieval, where Non-thinking models broadly attend to multiple needles; (ii) targeted retrieval, where Thinking models use enumeration in CoT traces to successively retrieve needles. Targeted retrieval concentrates attention on individual needles and is accompanied by more compact internal representations. Moreover, causal intervention analysis suggests that Thinking models use the CoT trace to maintain and update an internal counter as needles are successively retrieved, even without explicit numbering. In small controlled experiments, both retrieval mechanisms and counter states emerge under standard autoregressive training. Together, our results connect long-context retrieval with representation geometry of counting, supporting a state-tracking account of CoT reasoning.

Regional seismic assessments often assign structural properties by building archetype, limiting their representation of building-specific responses. This study develops a modular physics-regularized framework for identifying equivalent single-degree-of-freedom Bouc–Wen models from seismic input–response histories. A differentiable forward solver enables joint estimation of hysteretic parameters and viscous damping through an objective combining response matching with independently weighted dynamic-equilibrium and hysteretic state-evolution residuals. Before regional application, synthetic systems and four-story steel-frame hybrid simulation data are used to evaluate identification and prediction under unseen inputs, showing improved parameter estimation relative to a genetic algorithm and improved predictive performance with physics regularization, respectively. The framework is then applied to a numerical benchmark of 14,280 buildings in Milpitas, California, using models identified from simulated responses to a smaller earthquake to predict responses under larger scenarios. With the inventory, ground-motion inputs, and damage-and-loss models held fixed, archetypal and response-informed portfolios yield distinct damage-fraction distributions and substantially different structural repair costs. For the earthquake scenario with moment magnitude $M_w=6.8$, the response-informed portfolio yields a mean structural repair cost nearly twice that of the archetypal portfolio, despite nearly identical mean damage fractions. These results highlight the influence of structural model assignment on regional damage and loss estimates and the importance of incorporating building-specific response information into regional assessment.

In semantic segmentation, a recent line of RankSEG methods directly optimizes Dice/IoU scores at inference time, improving alignment with evaluation metrics without modifying model training. Despite its theoretical and empirical success, RankSEG relies on the restrictive Conditional Independence Assumption (CIA), which ignores crucial label correlations and therefore degrades performance in ambiguous or low-contrast scenarios. However, accounting for full label dependence is computationally prohibitive, requiring $\mathcal{O}(d^3)$ time. To address this, we replace the CIA with a Spatially Localized Dependence (SLD) structure that captures local label correlations while keeping the dependence model tractable. We further overcome the remaining computational bottleneck via a Reciprocal Moment Approximation coupled with a novel fixed-point optimization strategy that eliminates exhaustive search. The proposed algorithm achieves a highly practical $\mathcal{O}(d \log d)$ complexity and consistently outperforms conventional argmax and CIA-based RankSEG across diverse segmentation benchmarks. Improvements are significant in low-contrast or small-object scenarios, where label dependence offers valuable signals complementary to image information for accurate segmentation. The code of experiments is available at https://github.com/ZixunWang/RankSEG-DEP.

Flow Matching (FM) learns a velocity field whose ODE transports a simple source distribution to a target law. Existing finite-sample theory largely treats ambient-space regularity or data supported exactly on low-dimensional sets. We study linear FM under a noisy latent-generator model, where a low-dimensional Hölder map is perturbed by nondegenerate ambient Gaussian noise, so the target law is full-dimensional despite its latent structure. We construct a spatially regular ReLU velocity class and establish non-asymptotic high-probability approximation and estimation bounds whose leading sample-size exponent is governed by the latent dimension rather than the ambient dimension, with ambient and noise dependence kept explicit. Fixed positive target noise keeps the interpolation nondegenerate over the full time interval. The same spatial regularity propagates the learned velocity error through the transport ODE, yielding a corresponding Wasserstein convergence guarantee. These results show that exact low-dimensional support is not necessary for Flow Matching to retain latent-dimensional statistical behavior.

Any-dimensional machine learning models, such as graph neural networks (GNNs), can be naturally trained and evaluated on inputs of different sizes and dimensions. Inspired by the GNN transferability literature, we show mathematical conditions under which a partial differential equation (PDE) learning-based solver can be trained in small dimensions and directly applied to solve a higher dimensional PDE in a zero-shot fashion. These conditions are based on symmetries in both the partial differential equation and the initial data. When the equations satisfy the symmetries but the data does not, which is the case for many PDEs arising from physics, we show that our theory gives a principled way of warm-starting low-dimensional PDE solvers for higher dimensional PDEs. We apply this method on the heat equation, Burgers' equation, and the compressible Navier–Stokes equations, improving the performance in both zero-shot and typical training regimes on high dimensional data. For example, we train a surrogate model on 2D Navier–Stokes data and achieve better results on 3D test data than a baseline surrogate model trained on 3D data, while only using 12$\%$ of the flops and 20$\%$ of the total data size.

For $\beta \in \mathbb{R}$, we define the partition function $$ Z_\beta(X) = \sum_{i=1}^N \exp(\beta X_i) $$ of a centered Gaussian random vector $X=(X_1, \ldots, X_N)$ with $\mathbb{E}[X_i^2]= 1$. For $q \in \mathbb{R}$, we prove a complete phase transition for the generalized $L_q$-means of $Z_\beta(X)$ at $q=1$. The centered Gaussian vector with covariance matrix $\Delta_N$ obtained from a regular simplex configuration of unit vectors maximizes these means for $q<1$, and minimizes them for $q>1$. The two regimes are governed by different principles. For $q > 1$ {and $\beta\ne0$}, the moment functional is globally strongly convex on the entire set of correlation matrices, with an explicit modulus of convexity and a quantitative centroid-shape stability estimate. For $q<1$, the strategy is different. We prove a universal comparison for log-concave profiles of reverse Brascamp–Lieb type. Specializing this result to the Gumbel profile yields a Laplace-transform comparison between the partition functions.

Augmenting small concurrent studies with external or historical cohorts is attractive in drug development, where enrollment is slow, follow-up is expensive, and closely related trial or real-world data are often already available. Bayesian dynamic borrowing (BDB) provides a principled framework for adaptively controlling the influence of external data, but classical implementations often depend on hand-specified priors and MCMC-based inference, which can be computationally expensive and not generalizable. In this work, we study amortized neural posterior estimation (NPE) as a flexible alternative. A single network is pretrained on simulated current/external dataset pairs spanning covariate shift, outcome drift, and joint non-exchangeability, and then returns an approximate posterior for a scalar current-study target in a single forward pass. Through simulation studies, we find that NPE is most useful under outcome drift and joint mismatch: in the harder outcome-drift regimes, it gives up to about five-fold lower absolute bias than the best classical baseline and keeps Type I error close to nominal. After pretraining, posterior summaries are obtained in about 8 ms per dataset, roughly $10^3\times$ faster than MCMC-based borrowing baselines in our timing experiment. We further analyze Alzheimer's Disease Neuroimaging Initiative (ADNI) data and show that, when mild cognitive impairment outcomes differ across cohorts, the NPE formulation recovers the later-cohort risk level in this example without claiming greater precision. Code is available at https://github.com/ChinHungScott/NPE-for-Bayesian-Dynamic-Borrowing-MLHC-.

We study the last iterate of standard importance-weighted temporal-difference (TD) learning from a single trajectory of a finite Markov decision process. The trajectory follows a behavior policy $b$ with an irreducible, possibly periodic state chain; a target policy $\pi$ covered by $b$ is evaluated with action importance ratios bounded by $\varrho$. Write $H=(1-\gamma)^{-1}$, with $\mu_{\min}$ the minimum stationary probability and $\mathfrak h_\star$ the largest expected hitting time of the reverse behavior state chain. We prove that the last iterate achieves sup-norm error at most $\varepsilon$ with high probability after $\widetilde O(\varrho H^3/(\mu_{\min}\varepsilon^2)+\mathfrak h_\star)$ transitions for $0<\varepsilon\le1$. The statistical term retains the cubic effective-horizon dependence of synchronous TD, while the initialization cost is additive in $\mathfrak h_\star$. Setting $b=\pi$ and $\varrho=1$ recovers the on-policy rate. For aperiodic behavior chains, $\mathfrak h_\star=\widetilde O(t_{\rm mix}/\mu_{\min})$ recovers a mixing-time bound. The proof combines an anchored local Poisson equation in reverse time to control fluctuations without a mixing-time factor, a compensation argument for successive hitting costs, and an expected cumulative propagation mass bound for decreasing step sizes. Given the policy ratios and their bound, an empirical state-count rule attains the same rate without prior knowledge of model parameters. The framework covers nonnegative affine importance-weighted prediction updates, including successor features, transition-dependent discounting, and coupled finite-mode predictions. A perturbation estimate yields a time-uniform guarantee for online return-variance estimation. A three-state, two-action construction gives matching minimax lower bounds for both terms up to logarithmic factors on the considered parameter classes.

This paper develops a framework for estimation and inference on the volumes of sets that are projections of critical function sets, focusing particularly on the convex body beneath the optimal receiver operating characteristic (ROC) surface. Specifically, we propose a volume calculation method that first uses an Aumann expectation representation and then applies Minkowski mixed volumes. Using this framework, we show that the population volume under the ROC surface (VUS) is proportional to the expectation of a symmetric U-statistic kernel. We then propose a double/debiased machine learning estimator of the VUS, derive its asymptotic properties, and develop an inference procedure. Further applications of this framework include an analysis of the feasible error set across pre-defined groups and a natural generalization of the Gini coefficient for measuring inequality.

Learn-then-differentiate (LTD) estimates gradients by fitting a model to simulation outputs and differentiating it. We develop a unified framework explaining what LTD differentiates and how accurately it estimates gradients. For models with a weighted representation, LTD differentiates a learned representation of the underlying probability measure. We then show how accuracy guarantees for fitted models translate into guarantees for gradients and higher-order derivatives, with rates approaching the standard Monte Carlo rate under suitable smoothness conditions. The framework recovers established results for kernel regression, local polynomial regression, and kernel ridge regression, and yields further guarantees for multiple kernel learning and smooth neural networks. These results provide a common foundation for understanding and analyzing LTD across learning methods.

Contrastive learning is a successful paradigm for learning $d$-dimensional geometric representations from a collection of “anchor–positive–negative” triplets $(i,j^{+},k^{-})$, indicating that “item $i$ is closer to $j$ than to $k$.” Despite its success, understanding why contrastive learning leads to representations of high generalization quality—beyond the often pessimistic predictions from PAC-learning—remains a central question. Recently, \citet*{alon2024optimal} proved that, for PAC-learning $d$-dimensional Euclidean representations of $n$-point datasets, $\Theta(\min(nd, n^2))$ triplets are necessary and sufficient, while they posed as an open question whether their VC dimension bounds for the more realistic setting of contrastive learning with a margin can be improved. For a margin parameter $\alpha >0$, a triplet $(i,j^{+},k^{-})_{\alpha}$ is satisfied by the embedding $\phi:[n]\rightarrow \mathbb{R}^{d}$, if $\|\phi(i)-\phi(k)\|_2>(1+\alpha)\cdot\|\phi(i)-\phi(j)\|_2$. In this work, we resolve their question by proving that the VC dimension of contrastive learning under any margin $\alpha\in(0,1)$ is in fact $O(n/\alpha^2)$, improving on the previous bound of $O(n\log(n)/\alpha^2)$. We also establish that the bounds are optimal up to constant factors, by providing a matching lower bound of $\Omega(\frac{n}{\alpha^2})$ (the previously known lower bound was $\Omega(\frac{n}{\alpha})$), for $\alpha\geq \max(n^{-1/2},d^{-1/2})$.

Mental health support hotlines log help-seeking call attempts in real time, yet limited service capacity allows only a small proportion of calls to be answered and assessed. The time of every attempt is recorded, whereas issue types and marks such as caller demographics and crisis severity are observed only for assessed calls. Caller anonymity also prevents linking repeated attempts to the same individual. Nevertheless, service planning and crisis monitoring require estimating unassessed demand, recovering temporal patterns in call content, and quantifying how much the partially observed marks improve these estimates. We develop a Joint Marked Dynamic Factor Model (JM-DFM) in which a shared low-dimensional latent state drives attempt intensity, issue composition, and mixed-type mark distributions. We establish identifiability, prove consistency of the proposed sieve estimator, and show that incorporating marks improves the asymptotic estimation precision. The proposed framework is then applied to records of 68,371 call attempts to China's national 12356 mental health support hotline. We estimate that about 35 unassessed attempts per day involve suicidal ideation, representing a higher proportion than among assessed calls. Incorporating the marks narrows the latent-state uncertainty bands by 18% to 49%, most during periods with few issue-labeled calls.

Policy evaluation and learning aim to assess and optimize treatment assignment rules based on individual characteristics. A fundamental challenge arises when outcomes are truncated by death, rendering conventional policy utilities undefined. To address this challenge, we study policy evaluation and learning under truncation by death within the principal stratification framework. We propose the survivor average utility and subgroup survival rates for evaluating treatment policies. Under treatment ignorability and monotonicity, the subgroup survival rates are point identified, whereas the survivor average utility is only partially identified through sharp bounds. We then develop semiparametrically efficient estimators for the subgroup survival rates and hybrid estimators for the survivor average utility bounds. Building on this evaluation framework, we formulate a constrained minimax optimization problem for policy learning that minimizes worst-case regret relative to a benchmark policy while requiring the learned policy to achieve a survival rate no lower than that of the benchmark. We show that the optimal policy admits a threshold representation based on two priority scores. We develop an estimation procedure for the optimal policy and establish regret and feasibility guarantees. Simulation studies and an application to the MIMIC-III clinical dataset demonstrate the practical performance of the proposed methods.

The Bayesian optimal phase II (BOP2) design accommodates simple and complex endpoints but does not address biomarker-guided enrichment or multiplicity control across pooled and subgroup claims. We propose BOP2-ENR, which monitors biomarker-negative (NG) and biomarker-positive (PG) cohorts separately, permits interim enrichment to PG, and calibrates pooled and PG-restricted claims. Strong family-wise error rate (FWER) control is defined over a prespecified monotone-activity parameter space. Under the mixed configuration in which NG is inactive and PG is active, the probability of a false joint claim is bounded by the probability that NG survives futility monitoring. This exact single-cohort bound does not involve the final efficacy cutoff and also yields a ceiling on joint-claim power. For efficacy-toxicity monitoring, subgroup inactivity is defined by a union null and safety is established separately within each subgroup. Across the scenarios considered and binary, co-primary and efficacy-toxicity endpoints, all selected designs met the prespecified exact-bound and Monte Carlo verification criteria; in sensitivity analyses over odds ratios between component outcomes assumed independent in calibration, point estimates of the error-control quantities remained below the nominal level, although power was association-sensitive. The bound was nearly attained as PG efficacy approached the boundary of the nuisance space, whereas disabling enrichment produced mixed-configuration error rates near one. BOP2-ENR therefore extends the BOP2 framework to adaptive enrichment while clarifying that the NG futility gate governs the mixed-configuration component of strong error control and the attainable power of a joint claim.

We study the problem of sampling from a continuous density $\pi\propto \exp(-V)$ on $\mathbb R^d$, where $V\in C^2(\mathbb R^d)$ has a $\beta$-Lipschitz gradient and $\pi$ satisfies a logarithmic Sobolev inequality with constant $\alpha^{-1}$, and write $\kappa = \beta/\alpha$. We introduce the Gaussian cloud sampler, which achieves total variation accuracy $\varepsilon$ using $\widetilde O(\kappa d^{1/5}\,\text{polylog}(1/\varepsilon))$ gradient queries in expectation. The algorithm uses first-order rejection sampling (FORS) to correct the law of smoothed Picard HMC trajectories. To do so, we represent the iterates of the ideal Picard iteration via Gaussian clouds, whose centers are never evaluated, and we develop a suite of likelihood correction gadgets which only use samples from this indirect cloud representation.

Generative AI has dramatically accelerated the rate at which new treatments—from novel pharmaceuticals to online marketing campaigns—can be conceived and deployed. As a result, modern experimentation platforms often run continuously, with treatments added as they are ready and removed when they underperform. We formalize this "Always-On" experimental setting, in which treatments can be dynamically generated, added to, and removed from a running experiment, and study the statistical problem of deciding whether to accept or reject each treatment while controlling for the false discovery rate. We develop sequential tests that achieve time-uniform Type-I error control under arbitrary stopping times and "predictable" treatment schedules. Our approach builds on the testing-by-betting framework: we construct test supermartingales for testing the average treatment effect of each treatment, and show that the construction of these test supermartingales is growth-rate optimal in an almost-sure sense.

Batting average factors exactly as BA = c times f, where c = (AB - SO)/AB is the rate of avoiding a strikeout and f = H/(AB - SO) is the rate at which non-strikeout at-bats become hits. Using the 256 major-league hitters with at least 300 at-bats in 2025, we show that the two factors behave very differently. Strikeout avoidance is highly repeatable, with a median year-to-year correlation of 0.86 over 19 consecutive-season pairs from 2004 to 2025. The finishing rate is not: its median correlation is 0.44, and batting average itself (0.44) is no more repeatable than its noisier factor. A bivariate logistic-normal random-effects model fit to the 2025 season estimates the correlation between the two talents at -0.68 (95% profile interval -0.85 to -0.50), far stronger than the raw correlation of -0.44, and implies single-season reliabilities of 0.91 for c and 0.37 for f. The model also yields a closed-form bivariate shrinkage estimator in which a hitter's strikeout rate informs the estimate of his finishing rate. Applied to decades of American and National League data, the decomposition revisits Gould's explanation for the disappearance of the .400 hitter. Since the dead-ball era, the talent variance of strikeout avoidance has grown roughly 4.6-fold while that of finishing has halved, and the correlation between them has moved from near zero to about -0.53. That emerging trade-off, rather than a general narrowing of talent, accounts for the reduced spread of modern batting averages.

Imputation-based causal estimation is typically viewed as relying exclusively on an outcome model, in contrast to augmented inverse-probability weighting, whose consistency is protected by fitting two nuisance models. This paper argues that this view can be misleading by highlighting a hidden dual weighting structure in least-squares sieve regression imputation. Although only outcome regressions are explicitly fitted, the resulting imputation estimator admits an exact weighting representation whose induced weights balance every function in the sieve space and the corresponding population weighting functions are the $L^2$ projections of the inverse propensity scores onto the same sieve space. This projection structure yields an implicit form of double robustness and, under standard sieve approximation and growth conditions, asymptotic linearity with the efficient influence function. Thus, weighting, covariate balance, double robustness, and semiparametric efficiency can all emerge from imputation alone through the geometry of least-squares projection.

We develop a three-stage adaptive procedure for principal component analysis (PCA) when high-dimensional observations are collected sequentially and additional sampling incurs a cost. The procedure balances PCA compression loss against sampling cost while selecting the retained dimension through a prescribed explained-variance criterion. Starting from a pilot sample, an intermediate stage updates the PCA quantities before determining the final sample size, thereby avoiding reliance on unknown population eigenvalues. Under suitable regularity conditions, we establish both first- and second-order efficiency relative to the population oracle. Comparison with the corresponding two-stage rule shows that the additional recalibration yields sharper second-order control and reduces the influence of the pilot stage on the final sampling decision. The theory allows the ambient dimension to exceed the sample size under appropriate covariance and spectral conditions. Simulation studies demonstrate the strong finite-sample performance of the procedure across increasing dimensions and several dense covariance structures. As a real-data application, we conduct a retrospective study of gene-expression data from 32 cancer-type cohorts in The Cancer Genome Atlas, illustrating both cost-effective early stopping and settings in which additional observations are recommended.

Modern experiments increasingly produce rich post-treatment objects—such as text, images, and video—that may provide early surrogate signals for outcomes that are delayed, expensive, or observed only for a subset of participants. A pretrained encoder can turn such an object into a strong generic embedding, but the embedding is not trained for the estimation tasks required by a causal estimator and may emphasize information irrelevant to the causal target. We propose SIRNET, a Surrogate-Informed Representation adaptation framework for neural-NETwork-based causal estimation in the presence of unstructured post-treatment surrogates. Starting from pretrained surrogate embeddings, SIRNET learns a shared adapter jointly for outcome prediction, treatment propensity, and arm-specific surrogate marginalization, beyond the generic pretraining objective. The resulting covariate and surrogate representations provide a common interface for these nuisance tasks. The final treatment-effect estimate and confidence interval are obtained from cross-fitted orthogonal scores that attain the semiparametric efficiency bound under the stated nuisance-rate conditions. We also validate SIRNET's performance in synthetic and semisynthetic experiments, as well as the Anti-Amyloid in Asymptomatic Alzheimer's Disease (A4) trial, using frozen pretrained brain-MRI embeddings as surrogates for a long-term outcome.

On-policy distillation (OPD) learns from teacher feedback on student-generated responses and has shown promise in reducing forgetting relative to supervised fine-tuning (SFT). However, its benefits and fragility remain incompletely understood. We study sequential distillation from multiple teachers, where the student minimizes its average divergence from the teachers. Forward Kullback–Leibler (KL) divergence yields a weighted arithmetic mixture, while reverse KL yields a normalized weighted geometric aggregate. We develop algorithms that learn these targets under off-policy and on-policy feedback, respectively, establishing logarithmic regret bounds in the tabular setting and extending the analysis to function approximation. By analyzing these aggregation targets, we identify mechanisms that help explain both the benefits and fragility of OPD. Relative to forward KL, reverse KL can better retain a confident expert's preferences under uninformative feedback, but is more sensitive to teachers that assign very low probabilities to correct responses. Its token-level conditionals also reveal a dependence on continuation distributions that can favor incorrect prefixes over long horizons.

We study multi-armed bandits (MAB) with multiple optimal arms, motivated by the fact that many practical decision making problems admit multiple correct answers. For $K$-armed bandits with $A$ optimal arms, we first provide a sharper analysis of previous sub-sampling algorithms (De Heide et al., 2021; Zhu and Nowak, 2020), establishing a $\tilde{O}\Big(\frac{K-A}{\sqrt{KA}}\sqrt{T} \Big)$ minimax regret, where $T$ is the total number of interactions and $\tilde O(\cdot)$ drops all constant and logarithmic factors, improving the previous $\tilde{O}(\sqrt{KT/A})$ regret. We then provide a matching lower bound up to logarithmic factors, indicating that our established rate is nearly minimax-optimal. We further show that the knowledge of $A$ up to $\tilde{O}(1)$ factors is necessary to achieve near-optimal regret, as near-optimal algorithms for one number of optimal arms must incur substantially larger regret than optimal regret for a smaller number. Overall, our results provide a comprehensive minimax characterization of $K$-armed bandits with $A$ over the entire range of $1 \leq A \leq K-1$.

Correlated predictors can support competing sparse explanations with similar predictions, making joint uncertainty about variable inclusion difficult to capture with mean-field approximations. We develop an adaptive fitting procedure for mixtures of product distributions in Gaussian regression with a point-mass spike-and-slab prior. It minimizes reverse Kullback-Leibler divergence directly on inclusion indicators and active coefficients, jointly refining component parameters and weights as the mixture grows. This avoids an additional divergence penalty on unused latent coefficients under independent augmentation. Our analysis relates approximation accuracy to mixture size, support coverage and dependence within supports, and establishes contraction, selection consistency and a Bernstein-von Mises approximation under explicit conditions on the prior, posterior concentration and variational error. On all 250 simulated datasets with exact posterior references, mixtures reduce errors in inclusion probabilities, grouped support probabilities and coefficient covariance relative to multistart mean field. Comparisons at fixed mixture size and common initialization favor direct joint refinement over augmented or restricted refinement in posterior divergence. Complete stagewise fitting can nevertheless be more accurate near collinearity. The results support direct joint refinement for posterior approximation while showing that local gains do not ensure superiority of the full adaptive search.

High-frequency athlete-tracking data can be summarized as quantile cubes (Thomas and Hannig, Journal of Quantitative Analysis in Sports, 2026): counts of time spent in bins of velocity, acceleration and movement angle. In our application each athlete-match session has 100 such bins, with counts from under a hundred to several thousand. Whether a bin's negative binomial distribution is close to normal depends on its shape, the product of mean and dispersion, not on the count alone. We develop a hybrid generalized fiducial framework that uses exact negative binomial structural equations where the estimated shape is small and moment-matched normal approximations where it is large, providing uncertainty quantification for hundreds of regression and dispersion parameters without a prior. Only the normal components contribute to the fiducial Jacobian, which is block diagonal in closed form, and we prove a Bernstein-von Mises theorem for the correctly specified hybrid model. Simulations show near-nominal coverage, agreement with a flat-prior Bayesian analysis, and that the partition should be based on the estimated shape rather than on the mean or the observed responses. In 216 sessions from 17 professional women's soccer athletes, the estimated position and age effects concentrate at the velocity extremes and highest accelerations, but an athlete-level permutation check shows these patterns cannot be distinguished from label noise.

Global Dust Storms (GDSs) dominate the Martian climate, yet the mechanisms of dust lifting, transport, and deposition remain poorly constrained because in situ observations are limited: landers are stationary, and their Radioisotope Thermoelectric Generators (RTGs) contaminate local measurements. Our solution is a swarm of lightweight, wind-propelled spherical rovers – Dust Rovers (DRs) – that follow natural wind patterns to collect high-spatial-resolution in situ dust data during storms. Here, we quantify the dispersion and communication connectivity of this robotic architecture. Using 625 wind states calculated using the Mars Climate Database (MCD) for the MY34 storm at the Curiosity site, we build a bivariate-normal model of the near-surface wind velocity and characterise its diurnal structure. Rover trajectories and swarm dispersion are propagated with three integrators: the Euler and fourth-order Runge-Kutta schemes for the deterministic dynamics, and a stochastic Ornstein-Uhlenbeck formulation that simulates turbulent gust forcing. Dispersion is controlled by the unmeasured gust correlation time: with a static rolling threshold enforced, the white-gust case (tau_g = Delta t = 1 s) provides a two-rover separation below 1 m over an hour (0.44 km over a sol), whereas correlation times of 10-60 s yield 6-25 m over an hour (1.4-3.3 km over a sol). We demonstrate that a 20-rover swarm needs only a 1.6-4.1 km per-rover range to remain connected over a sol, and stays connected across the plausible envelope – for gust amplitudes up to twice the nominal value with tau_g <= 60 s, and at the nominal amplitude up to tau_g ~ 150 s.

When can we learn from biased samples? We study regression when outcomes are observed only after passing through selection filters that depend on both covariates and outcomes themselves, a ubiquitous challenge spanning clinical trials with patient dropout, labor markets with self-selection, and auctions with strategic entry. Ignoring such selection yields systematically biased conclusions with real-world consequences. This challenge has a long history in econometrics and statistics, starting with Heckman's seminal two-stage model and followed by numerous generalizations. While these works provide various sufficient conditions for identification, a complete characterization of when such regression is possible has remained elusive. In this work, we provide a characterization for when regression is possible in the presence of sample selection bias. Our results establish the minimal assumptions required on the functional forms of selection processes under which regression remains possible, which are particularly relevant in modern settings where selection mechanisms are increasingly complex and opaque. As a corollary of our characterization, we show that there are settings where the regression function can be identified even when the selection filter itself cannot. This observation already goes beyond the “estimate selection filter, then debias regression” paradigm that is followed by virtually all existing approaches. Under natural strengthenings of our identification conditions, we also establish finite-sample estimation guarantees with explicit convergence rates and provide oracle-efficient algorithms. This yields the first general-purpose estimation method for this broad class of selection problems. Finally, we explore the implications of our results for several well-studied econometric settings with complex selection mechanisms such as auctions with entry costs and labor markets.

Tax fraud remains a central challenge for public revenue authorities worldwide, imposing fiscal losses estimated to reach up to 1 trillion euros annually in the EU alone. Fraudulent firms intentionally manipulate reported figures to conceal their activity. Business networks can provide complementary information and reveal fraud that covariates alone may miss. We propose a Bayesian positive-unlabeled (PU) classification framework that combines firm-level covariates with multilayer network information. As a policy relevant case study, we apply the framework to firms in a regulated segment of the Greek energy market, a sector exposed to excise tax evasion. Confirmed fraud labels exist only for audited cases, while the remaining observations are unlabeled rather than verified compliant. We develop a Bayesian Gaussian process (GP) classifier integrating covariates with multilayer network information through a Product of Experts (PoE) construction and incorporating a nondetection probability for missed positives. We derive identification bounds for the nondetection rate under an anchor condition and show that, for a fixed latent risk function, a common nondetection rate affects calibration but not ranking. Simulations show strong ranking performance relative to PU and network-based alternatives, while illustrating the difficulty of estimating the nondetection rate. In real audit data, the method identifies high risk firms with posterior uncertainty, estimates the number of undetected fraudulent cases, and identifies whether risk is driven by covariates, networks, or both. Of the six highest ranked firms, the two that had already been reinspected independently were both confirmed as noncompliant, providing an independent validation of these cases. The remaining four were selected by the tax authority for follow-up inspection.

Tree-based scan statistics (TBSSs) are machine learning methods for disproportionality analyses. They simultaneously scan for thousands of hierarchically related health outcomes to detect potential signals of harm from drugs and vaccines, controlling for multiplicity. TBSSs have been extensively used to mine insurance claims databases to detect potential drug adverse events. Current TBSS implementations do not allow comparative safety evaluations with time-to-event outcomes. Explicitly accounting for event timing in analyses can improve power to detect signals compared to methods that use only event counts. We propose three TBSS methods that explicitly leverage event timing to detect potentially harmful effects of drugs. The first assumes proportional hazard rates for each node of the outcome hierarchy and uses a permutation scheme for inference. The second builds on exponential survival models for the terminal nodes of the hierarchy, assuming constant hazard rates at each node, and uses a parametric bootstrap for inference. The third uses robust asymptotic approximations of the hazard rates in connection with an approximate parametric bootstrap. We compare the proposed methods with standard event count based TBSSs in simulation scenarios. Finally, we present results from a database study comparing two glucose-lowering medications among adults with type 2 diabetes.

Defining a weighted mean over probability measures under probability metrics is a central tool in probabilistic machine learning. Under the Wasserstein metric, these are called Wasserstein barycenters. While most approaches compute barycenters for a fixed weight vector, approximating the whole family of barycenters over the simplex, which we call the Wasserstein simplex, remains underexplored. We refer to this problem as Universal Barycenter Approximation, and propose BaryFM, a flow matching model transporting the marginal measures into any barycenter in the Wasserstein simplex. Once trained, the network can draw samples from measures in the Wasserstein simplex through an ordinary differential equation. We validate our method on 4 downstream tasks: domain adaptation, generalization, Bayesian posterior aggregation and algorithmic fairness. BaryFM achieves the best average rank among 15 competing methods across 10 domain adaptation benchmarks, matching or surpassing non-universal solvers.

We consider estimating the one-step-ahead conditional distribution of a multivariate stochastic process. Many existing approaches rely on assumptions such as finite-range memory, sparsity, or additivity, which can be poorly suited to processes with long-range nonlinear interactions. However, without such structural assumptions, nonparametric estimation is challenging due to the curse of dimensionality. To address this challenge, we introduce a new estimation approach based on the predictive states of a process, possibly with infinite-range memory. We show that our estimator achieves fast convergence rates when the past history can be compressed into a low-dimensional statistic that is sufficient for predicting the future. Specifically, we show that the statistical complexity of the estimation problem is determined by the intrinsic dimension of the predictive state space. We establish guarantees for an instantiation of our method based on deep neural network estimators, and we support these theoretical results with experiments.