Papers

What can a finished text reveal about the process that produced it? Drawing on Turing's imitation game and statistical decision theory, this article examines the limits of AI-text detection as an inference from a completed object to an unobserved history. For two known source distributions with equal prior probabilities, a standard identity gives the minimum average classification error as half their probability overlap. This overlap is one minus their total variation distance. Perfect detection therefore requires nonoverlapping distributions; useful discrimination does not. In practice, the problem is harder: human and machine writing form changing families of distributions, posterior probabilities depend on base rates, and “AI-written” becomes ambiguous when people and software contribute to the same text. Turing's interrogator can ask another question; a detector restricted to finished prose cannot. Distinguishability, source attribution, authentication, and compliance are different inferential tasks. The paper's own documented human-AI provenance shows what a production history can reveal that a binary label cannot.

On evolving graphs, node classifiers must satisfy two key requirements: inductive generalization to newly arriving nodes under distribution shift and calibrated uncertainty for safety-sensitive applications. Standard graph neural networks (GNNs) are typically trained once and address neither requirement. We adapt the Bayesian last-layer (BLL) model by placing random last-layer parameters on top of a deterministic GNN encoder for uncertainty quantification. The categorical softmax likelihood required for classification breaks Gaussian conjugacy, so neither the training posterior nor the test-time streaming update has a closed-form solution. To address both challenges, we introduce a variational Bayesian last-layer (VBLL) objective that jointly trains the encoder and an approximate last-layer posterior by maximizing an evidence lower bound with a Monte Carlo expected log-likelihood. At test time, we freeze the encoder and apply an online Laplace update to the last-layer posterior. This update corresponds to a power-prior Bayesian model with exponential forgetting and a Kullback-Leibler anchor to the training posterior. Across five node-classification benchmarks under distribution shift, online GVBLL is the only method to achieve the best accuracy and negative log-likelihood on every dataset. It improves accuracy by up to 17 percentage points on Cora and 14 percentage points on ogbn-arxiv over the strongest non-GVBLL baseline, while remaining competitive in calibration with MC Dropout, Deep Ensembles, Temperature Scaling, and Gaussian-process classifiers.

This work addresses decentralized online Riemannian optimization on Hadamard manifolds. Prior work under geodesic convexity (g-convexity) may require curvature information in the optimization analysis, typically through a finite lower bound on the sectional curvature. Curvature may also enter the step size or contraction factor of tangent-space Riemannian consensus schemes. In this work, we relax the curvature dependence for a narrower class of horospherical convex (h-convex) functions. We study Distributed Riemannian Online Gradient Descent (D-ROGD), which combines local Riemannian h-subgradient updates with an implicit Fréchet-mean consensus. For h-convex and strongly h-convex local objectives, we establish $O(\sqrt{T})$ and $O(\log T)$ static regret, respectively, matching the corresponding Euclidean rates with respect to $T$, with network dependence governed solely by the spectral gap. To our knowledge, these are the first curvature-independent regret guarantees for decentralized online optimization on Hadamard manifolds. Experiments on hyperbolic embeddings corroborate the predicted rates, with no observable degradation due to curvature.

Network analysis for multivariate time series is popular in many fields, from neuroscience to seismology. The inverse spectral density is a common choice for time series network analysis due to its representation of the frequency domain correlation between two variables after removing the best linear predictor of all other variables. In many applications, the goal is to study how these networks change across different conditions. For example, in neuroscience, one might be interested in how the brain connectivity network changes before and after stimulation. Towards this goal, we develop an inference framework based on a direct estimate of the difference in two high-dimensional inverse spectral densities. We develop a new Gaussian approximation error bound for any de-biased D-trace estimation procedure which is then leveraged to both inform optimal window sizes of Welch's estimators of the spectral density and establish asymptotic normality of our de-biased D-trace estimator. Moreover, we develop an efficient algorithm based on a generalized D-trace estimation procedure to overcome the computational complexity of high-dimensional inference. The method is illustrated on synthetic data experiments and on experiments with electroencephalography data.

Long records of rainfall contain critical information about the temporal organization of precipitation within a year, including wet and dry spells. Detecting when this temporal pattern changes is central to climate-impact assessment and environmental monitoring.A common data structure features observations from m locations over T years. For location i and year t, the record is a vector of daily measurements $X_{i,t} = (X_{i,t,1}, ..., X_{i,t,D})^\top$, with D = 365. The goal is to identify an unknown year $t_0$ when the temporal rainfall pattern changes across the region. Because changes in the intra-annual pattern affect the frequency-domain behavior of the daily sequence, we represent each yearly curve through a spectral feature capturing the strength of temporal oscillation at a fixed frequency.Classical changepoint methods often rely on parametric likelihoods or Gaussian approximations, which are hard to justify for skewed and heavy-tailed rainfall data. Furthermore, many procedures only return a point estimate. A confidence set for the change year is crucial in environmental applications to quantify temporal uncertainty and distinguish sharp transitions from weak evidence.This paper develops a distribution-free framework for detecting a synchronized structural change in yearly rainfall patterns across independently monitored locations. For each location-year pair, we compute a fixed-frequency nonparametric spectral density estimate. These features form a $T \times m$ data matrix. By exploiting the exchangeability of yearly spectral features before and after the true changepoint, alongside spatial independence, we construct finite-sample valid confidence sets for the unknown change year.

AI-based item generation and NLP-based prediction of item parameters are producing item banks that are substantially larger, sparser, and more frequently updated than conventional banks. Hierarchical Bayesian item response theory (IRT) is a natural calibration framework for such banks, but the common practice of refitting the entire accumulated response history at each update is costly and can exceed available memory. We describe consensus calibration, a divide-and-conquer procedure that calibrates each time period independently and reconstructs the pooled posterior in two layers. First, the posterior draws of each period are mapped to a common metric by a robust characteristic-curve linking (Haebara) that is solved separately for each draw, which propagates the uncertainty of the linking transformation into the linked posteriors. Second, the linked item posteriors are combined as a product of Gaussian densities from which the population prior contributed by each period is removed and a single prior—obtained by consensus across the per-period population posteriors—is reinstated. The correction targets the posterior dispersion, not only its location. As evidence for consensus calibration, we compare it to a pooled single-run analysis on a large operational assessment in terms of item-parameter recovery, an uncertainty-by-exposure diagnostic, and the ability distributions.

We develop a causal framework for matrix completion under missing not at random (MNAR) data. Drawing on synthetic controls from the econometric panel data literature, our approach relaxes two assumptions common in MNAR matrix completion: positivity and independence of observation indicators. Unlike traditional panel data models, which often require prescribed block-sparse geometries, our framework accommodates flexible, heterogeneous observation patterns through target-specific local information structures. We propose synthetic nearest neighbors (SNN), a local synthetic-controls-inspired estimator, and establish finite-sample entrywise error bounds and consistency for mean recovery under suitable conditions. We further derive asymptotic normality under heteroskedastic noise and develop feasible entrywise inference. To estimate entry-specific noise variances, we apply the same local principle to squared outcomes, obtaining consistency under bounded noise and asymptotic unbiasedness under general subgaussian noise. Simulation studies corroborate the theoretical findings across a range of missingness designs and observation patterns.

Objectives in scientific machine learning are often prescribed as a sum of several terms, such as the residual, boundary, initial, and data losses of a physics-informed neural network. In the pullback-corrected scalar auxiliary variable (PB–SAV) method, one scalar tracks the shifted objective while the component gradients build a positive semidefinite curvature correction of rank at most the number of components. We carry that correction into an optimizer with momentum and an adaptive mobility, applying it to the gradient and the stored momentum in a single implicit solve. A mobility that is nonincreasing in the Loewner order yields an exact modified energy law, covering Euclidean and AMSGrad-type choices; the corresponding identity for momentum appended after the solve carries a cross term of indefinite sign. For a fixed mobility we give a necessary and sufficient condition for local stability at a stationary point, depending on the Hessian minus twice the correction, and show that it also gives local geometric convergence for every scalar relaxation sequence. The implicit solve reduces to a dense system whose order is the number of components. In the forward Burgers comparison, four components reduce the mean tail objective by 64.7% and the final solution error by 50.2% relative to one component at the same learning rate and momentum settings.

Statistical validation of synthetic multivariate data requires assessing whether a generator preserves the joint dependence structure of the target population without merely reproducing observed records. We develop a model-agnostic framework based on full conditional distributions. For each coordinate, we normalize the conditional probability assigned to the observed value by the largest conditional probability available in the same record context; averaging this quantity yields a one-sided MAP-alignment statistic that can be estimated using a conditional model fitted on held-out real data. The mathematical contribution is twofold: under strict positivity and compatibility, the complete normalized conditional profile identifies the joint distribution, and its integrated L1 difference defines a metric on finite-state generative processes; we also establish consistency and finite-sample concentration for the corresponding empirical estimators. Because high conditional alignment alone can arise from copying or concentration on conditional modes, we pair it with nearest-real similarity as a separate record-level novelty diagnostic. We evaluate the framework on NSHAP health and aging data, influenza B genomic surveillance, and 34 General Social Survey waves. In GSS, the Large Science Model matched the original-data control in mean conditional alignment while retaining substantial novelty, indicating preservation of conditional structure without row reuse. In influenza B, a Chow-Liu generator matched the control alignment but had almost no novelty, revealing near-reproduction of observed records. The framework therefore distinguishes three statistically different failure modes: loss of dependence, record reuse, and mode concentration, and provides a principled basis for validating synthetic health, surveillance, and population data.

We study switching regret in adversarial multi-armed bandits, where the learner competes with an arm sequence that changes at most $S$ times. When $S$ is known, an optimal expected regret of $\widetilde{\mathcal{O}}(\sqrt{(S+1)KT})$ is obtainable [Auer et al., 2002]. However, when $S$ is unknown, Marinov and Zimmert [2021] show that this guarantee is impossible under an adaptive adversary. In this paper, we show that a single algorithm achieves $\widetilde{\mathcal{O}}(\sqrt{(S+1)KT})$ expected regret for every $S$ against an oblivious adversary, resolving an open problem of Auer et al. [2019b]. Our algorithm combines a fixed-share learner initialized with a small learning rate and dyadic-interval subroutines that search for local improvements using randomized learning rates and implicit exploration. Importantly, a non-uniform prior favors following the main learner, keeping the cost of maintaining many subroutines small. When the subroutines accumulate sufficient improvement over the main learner, its learning rate doubles, allowing adaptation to the unknown number of comparator switches $S$.

The basic question in perturbation analysis of Markov chains is how small changes in their transition kernels affect their stationary distributions. Classical perturbation bounds typically require the kernel error to be much smaller than $1/\tau$, where $\tau$ is a mixing or relaxation time. Although this scaling is sharp for Markov chains in general, we investigate a general "square-rooting" phenomenon in which one-step errors of order roughly $1/\sqrt{\tau}$ can be sufficient for local updates. We proved a form of this phenomenon in Lin, Liu and Smith (2025) under strong assumptions. Here we substantially weaken these assumptions, and prove this phenomenon occurs using three distinct approaches. First, block factorization applies under structural assumptions on the stationary measures of both chains. Second, approximate block-update arguments extend the result to statistically relevant Markov chain Monte Carlo (MCMC) settings, where structural guarantees are available for the exact posterior and its associated sampler, but not for the perturbed posterior. Third, we use direct calculations for a class of models with hard constraints where neither general result is directly available. We illustrate these results in three MCMC settings and show how they directly inform the tuning of approximate MCMC algorithms.

Detection and attribution (D&A) analyses provide a statistical framework for quantifying the contribution of external forcings to observed climate change. Optimal fingerprinting, the primary approach for D&A, is formulated as an errors-in-variables regression with a high-dimensional covariance structure. Reliable inference is challenging because covariance matrices must be estimated from a limited number of control simulations, and climate models may exhibit variability patterns that differ from those of the observed climate system. We develop a spiked optimal fingerprinting framework that exploits the spiked structure of dominant climate variability modes while providing stable covariance estimation in high-dimensional settings. The proposed method develops bias-corrected spiked covariance estimation for constructing adaptive weight matrices and incorporates variability inflation adjustments to account for model–observation differences in internal variability. We further develop a valid uncertainty quantification procedure for the scaling factor estimators and a residual consistency test for assessing model adequacy. Numerical studies demonstrate improved estimation accuracy and uncertainty quantification compared with existing practical approaches. Applied to annual mean near-surface temperature data, the proposed framework produces shorter and better-calibrated confidence intervals and leads to different detection and attribution conclusions in several regions, providing new insights into the interpretation of climate attribution results.

We study generative bootstrap processes obtained by resampling from fitted generative distributions. We establish necessary and sufficient conditions for their conditional weak convergence to the $P$-Brownian bridge, formulated in terms of finite-dimensional conditional weak convergence and conditional asymptotic equicontinuity. We further provide general sufficient conditions and a Giné-Zinn-type characterization relating this convergence to the $P$-Donsker property. As applications, we verify these conditions for prominent classes of generative models, including triangular normalizing flows, flow matching, score-based diffusion models, and Wasserstein generative adversarial networks.

Reconstructing stationary space-time Gaussian processes at unobserved locations is costly when many sites share regular temporal records. We develop a spectral likelihood combining a temporal discrete Fourier transform (DFT) with a Hilbert-space Gaussian process (HSGP) representation of frequency-specific spatial covariance. We derive the exact covariance of the finite-record DFT coefficients and use a Whittle likelihood that approximates distinct frequencies as independent spatial problems. At each temporal frequency, HSGP approximates spatial covariance by evaluating the sampled spectral multiplier at retained Laplacian eigenfrequencies. For models specified by a joint spectral density whose half-spectrum lacks a convenient closed form, this construction avoids repeated Fourier inversion. The fixed spatial basis also permits cached feature projections to be reused in likelihood fitting and held-site reconstruction, with approximation accuracy depending on domain extension and basis size. In a simulation study of such a model, HSGP achieved reconstruction accuracy comparable to a high-accuracy quadrature reference while reducing mean fitting time by 69%. Additional applications to wind-field reconstruction tasks examine the effects of basis rank, temporal record length, and separability, and demonstrate accurate inference on held-out test sites when sufficiently rich bases are used. Taken together, the results indicate that computational savings are attainable when the spectral multiplier can be evaluated directly, numerical spatial inversion is costly, and an adequate basis has rank lower than the number of fitting sites.

As compute scales, models evolve, and training algorithms advance, our ability to explain the increasingly powerful AI systems they enable is eroding. To help safeguard interpretability, we introduce a specialized training algorithm (MACCHIATO) that jointly constructs (i) an explicitly structured $\operatorname{ReLU}$-MLP from partial truth-table observations and (ii) an explicit Boolean circuit over signed literals with $\{\operatorname{AND},\operatorname{OR},\operatorname{XOR}\}$ gates certifying what its subnetworks compute and how they compose. Intuitively, we iteratively project the residuals of a Boolean function onto low-dimensional $\{\operatorname{AND},\operatorname{OR},\operatorname{XOR}\}$-circuit classes and exactly compile the resulting circuit into a $\operatorname{ReLU}$-MLP; we combine $\operatorname{ReLU}$-MLP circuit compilation, ESPRESSO logic minimization, and influence-based variable selection. Roughly speaking, our interpretability certificate is complemented by a statistical guarantee: under the theorem's influence-recovery conditions, if each of the $m$ stage-wise residuals depends on at most $\log_2(B)$ bits, a sample-splitting variant of our algorithm trained on $T$ observations returns a six-layer $\operatorname{ReLU}$-MLP (counting the input layer) of width $\mathcal{O}(mB)$ with truth-table error $\mathcal{O}\bigl(\sqrt{m(B+\log(m/\delta))/T}\bigr)$. On synthetic random-junta tasks, our networks outperform depth- and hidden-width-matched Adam-trained MLPs in several data-sparse or projection-aligned regimes, while the trained ReLU-MLPs are stronger in others. Moreover, in our explicit PyEDA truth-table implementation, the iterative procedure completes in regimes where flat ambient-dimensional ESPRESSO exceeds the three-hour computational budget.

We study the mixing time of a systematic-scan analogue of Kac's walk that was proposed as a fast surrogate for Haar-distributed orthogonal matrices in randomized high-dimensional algorithms and was conjectured to approach Haar measure after only logarithmically many sweeps. We show that this conjectured speed-up does not occur for convergence of the full matrix law to Haar measure in Frobenius Wasserstein distance. At fixed normalized accuracy, the mixing time lies between order $n/\log n$ and order $n$ sweeps; at fixed absolute Frobenius accuracy, the corresponding bounds are between order $n$ and order $n\log n$. More strongly, below the scale $n/\log n$, the normalized Wasserstein distance remains asymptotically at its extremal value. We also show that the output law is singular with respect to Haar measure for fewer than $n/2$ sweeps. Thus the sampler may provide effective application-specific randomization without exhibiting the much faster full-Haar mixing.

In this paper, we study a stochastic Susceptible-Infected-Removed (SIR) model where the infection and the recovery rates depend on individual covariates for susceptibility and infectiousness of the infector and the infectee. Such models allow explicit nonlinearity in the incidence term. They are also important from a practical perspective, as they allow for the incorporation of individual heterogeneity into the epidemic process. Statistical estimates for crucial epidemiological parameters, such as the basic reproduction number, herd immunity threshold, could be vastly different, and even biased, when the population heterogeneity is ignored in the mathematical model. We describe our epidemic model as an Interacting Particle System (IPS) of Stochastic Differential Equations (SDEs) driven by Poisson Random Measures. Our main mathematical contributions are a Functional Law of Large Numbers (FLLN), which approximates the empirical random measure of the IPS by means of a deterministic measure-valued function, and the propagation of chaos phenomenon, which establishes asymptotic independence of the particles as the population size goes to infinity with an explicit construction of McKean–Vlasov type Kac's “nonlinear process”. We also briefly mention how the propagation of chaos phenomenon leads to a product-form likelihood function, which forms the basis of the so-called Dynamic Survival Analysis (DSA) method for parameter inference based on sparse data.

Multi-objective Bayesian optimisation (MOBO) is a sample-efficient approach for optimising expensive black-box functions with multiple objectives. In MOBO, the goal is to adequately approximate the Pareto front; that is, to obtain a high-quality solution set with 1) good convergence (closeness to the Pareto front) and 2) good diversity (spread across the Pareto front). Existing MOBO methods typically aim to accomplish these two tasks simultaneously, i.e., driving the search towards the Pareto front while maintaining a diverse set of nondominated solutions, such that the solutions, ideally, can gradually approach the entire front. When sufficient search budgets are available, this approach is effective. However, considering both convergence and diversity throughout the search is not easy and requires careful design. Under very tight budgets, there may not be enough solutions generated to be able to simultaneously approach the entire Pareto front. To address this issue, this paper proposes a converge-then-diversify (CTD) approach that decouples convergence and diversity into two stages. In the first stage, CTD focuses on convergence, aiming to quickly drive the search toward a single point on the Pareto front. In the second stage, CTD focuses on diversity, aiming to spread solutions across the front. We present two simple instantiations of CTD by using widely adopted acquisition functions in the area. Experimental results show that, across all 446 pairwise comparisons, CTD statistically outperforms state-of-the-art methods in 72.9% of the cases, performs equivalently in 21.1%, and is statistically worse in only 6.1%, with the advantage being particularly evident in settings with very tight evaluation budgets or in high-dimensional problems.

Probabilistic forecasting is central to decision-making under uncertainty, yet its methodological landscape has become increasingly fragmented across temporal and spatiotemporal forecasting, statistical modeling, machine learning, and deep generative modeling. This survey develops a unified perspective by organizing probabilistic forecasting methods according to where and how uncertainty is introduced into the forecasting pipeline. Our taxonomy connects model-agnostic approaches including ensembles and distribution-free calibration, with model-intrinsic approaches spanning Bayesian modeling, parametric predictive distributions, distributional regression, and modern generative models, and further examines the emerging role of time series foundation models. Beyond methodological synthesis, we identify the assumptions, computational demands, and forms of uncertainty represented by different paradigms, and translate these distinctions into data-driven and domain-specific guidance for method selection. We complement the survey with a cross-paradigm empirical study on univariate, multivariate, and spatiotemporal forecasting tasks. The results reveal that no single uncertainty-quantification paradigm dominates across settings. Calibration, sharpness, predictive accuracy, and computational efficiency can lead to substantially different model preferences, while expressive generative models and zero-shot foundation models exhibit markedly different accuracy-efficiency trade-offs. Lastly, we identify unresolved challenges surrounding uncertainty in evolving dependency structures, physics-informed predictive distributions, forecasting extreme events, handling count-valued, directional, and continuous-time series, and the development of unified software resources. Our survey provides both a conceptual framework and a practical roadmap for probabilistic forecasting research.

Stochastic gradient descent (SGD) admits diffusion approximations that replace the complicated randomness of stochastic gradients by Gaussian noise, providing a powerful tool for understanding its dynamics and long-time behavior. We investigate whether an analogous approximation principle holds for optimization over probability measures, where the objective is a functional defined on the Wasserstein space P2. The nonlinear geometry and infinite-dimensional nature of P2 prevent a direct extension of the classical Euclidean theory. Using Lions differentiability, we lift the problem to a linear Hilbert space, where higher-order differential calculus becomes available. We then construct a Gaussian random-field approximation whose velocity field matches the mean and covariance of the original stochastic gradient. By exploiting this moment matching through higher-order Taylor expansions, we show that the Gaussian approximation captures the SGD dynamics with second-order weak accuracy. Our result provides a rigorous foundation for replacing sample-driven randomness by analytically tractable Gaussian fluctuations in stochastic optimization over probability measures.

Statistical models of ecological and evolutionary dynamics often include latent variables that are either continuous (e.g., average body size) or discrete (e.g., numerical abundance). Mixed-type hierarchical models containing both are typically fitted using Markov chain Monte Carlo (MCMC), which can be prohibitively slow for large models. Here, we introduce an alternative in the R package RTMB that automates the sequential reduction of small groups of related discrete variables, allowing them to be efficiently marginalized. Sequential reduction is combined with automatic differentiation and the Laplace approximation to estimate parameters and predict both continuous and discrete variables. We demonstrate speed and flexibility using demographic examples (occupancy, dynamic occupancy, N-mixture, and open dynamic N-mixture models), benchmarking RTMB against JAGS and unmarked. We then develop two novel applications. The first is a multi-site open N-mixture model with a spatial latent variable governing site-specific initial abundance and recruitment, which shows that continuous Gaussian Markov random fields can be estimated jointly with discrete abundance dynamics in under a minute. The second is phylogenetic trait imputation for a published data set of female Liolaemus lizards, where we jointly impute a binary trait (viviparity), estimate its state-switching rates, and estimate its effect on a continuous trait (body size) during ancestral state reconstruction. This indicates that phylogenetic comparative methods can estimate linkages among discrete and continuous traits. We envision that intuitive and efficient specification of mixed-type models will allow more expressive representation of ecological and evolutionary dynamics.

Computing matrix inverses (or, e.g., left inverses) using known algorithms is in general computationally prohibitive in high dimension, such as in an operational weather forecasting context. In this article, a swift, scalable algorithm to calculate the left inverse (through singular values decompositions techniques) of a square root of a horizontally localized ensemble covariance matrix is introduced. In this algorithm, calculations involving model variables located in different grid vertical columns can be performed separately. This is the aspect that makes this algorithm scalable by allowing for parallel computing. The algorithm is developed in the case when the number of ensemble members equals the number of model variables in each grid vertical column (this number is of order $10^3$ for current high-resolution models). The algorithm is tested in a synthetic experiment on a three-dimensional grid, and it is found to have very good accuracy. Additionally, a version of the algorithm outputting the left-multiplication of a generic vector by that above-mentioned left inverse is also introduced. For operational applications, if executed on a single processor, that version of the algorithm is already $9$ orders of magnitude faster than a traditional algorithm, and, with multiple, parallel processors, it scales until the number of processors reaches the order of $10^6$. In addition, an application of the algorithm is pointed out: performing the transform from model variables to an alternative set of variables, defined herein as standardized variables, being climatologically uncorrelated and each having unit climatological error variance. Furthermore, it is pointed out that the standardized variables allow constructing a new Hybrid covariance model, which is positive definite and climatologically unbiased, and it features a geographically dependent climatological error covariance matrix.

When providing forecasted probabilities with a predictive model, the ideal model offers perfect calibration: the true probability of the outcome (i.e., the probability that $Y=1$) exactly matches the forecasted probability $f(X)$. In practice, models inevitably exhibit calibration error, and it is therefore important to be able to measure this miscalibration to assess a model's reliability. The Expected Calibration Error (ECE) is the most widely used measure of miscalibration, but is known to be impossible to estimate the ECE with guaranteed accuracy in an assumption-free setting. In this work, we propose an alternative measure, the rankECE, that is based on comparing points with neighboring values of the predicted probability $f(X)$. Our theoretical guarantees and empirical results establish that rankECE provides a better proxy for ECE as compared to binned approximations to ECE, which are the most commonly-used approximations in practice.

Deep neural networks are expressive enough to contain worst-case targets that can be evaluated in polynomial time but cannot be learned in polynomial time by gradient descent. For practical tasks they nonetheless learn well, raising the question of what non-generic structure of real-world targets enables this. Existing surrogate models cannot pose this question because they either lack hard-to-learn targets entirely (deep linear networks) or cannot evaluate such targets efficiently (kernel methods, infinite-width limits). We study tree tensor networks (TTNs), a model class that generalizes deep linear networks and Tucker decompositions. We show they embed arbitrary read-once Boolean formulas, and thus contain polynomial-size targets that cannot be learned by gradient descent in polynomial time under the same mechanism as neural networks. Despite this, we prove that their loss landscapes are conditionally benign for every realizable target: every local minimum that is minimum-norm is global. Thus, surprisingly, bad local minima are not what distinguishes between typical and worst-case problems in TTNs. Instead, learning difficulty in TTNs can arise from high-order degenerate saddle points, which we show are caused by rank-deficiency. This is explored through a case study of the parity function, illustrating the potential for TTNs to relate landscape geometry to computational hardness.

Transfer learning explores how to leverage knowledge from various tasks or domains (sources) to enhance predictive performance in related tasks or domains (targets). Typically, transfer learning research is segmented into several isolated sub-areas (such as domain generalisation, domain adaptation, or multi-domain learning), each making distinct assumptions about target data availability, namely how much data and how many labels are available at training time. However, in many real-world applications, data availability is not fixed but evolves over time, as instances and labels are progressively collected from a new domain. Each of the classical settings then describes only a snapshot of a trajectory that a deployed system must traverse in full. We formalise this trajectory as a transfer learning problem in its own right, Transfer Learning for Evolving Domains (TrED), specified by a data availability process fixed by the environment, a learning protocol that the method is free to choose, and an evaluation criterion that scores the whole trajectory of models rather than a single one. Within this formalism, the classical settings are recovered as regimes that a learner may pass through, rather than as separate problems that TrED concatenates. We then examine the transfer learning literature to identify mechanisms that are promising building blocks for a solution, and find that most methods are tailored to a single regime and that even the strongest existing candidates do not yet optimise the whole trajectory. We argue that TrED is a well-posed and unsolved problem, and an important direction for future research.

We develop a unified framework for $c$-sample homogeneity testing on the circle. The proposed class of $c$-sample Sobolev tests generalizes the two-sample Sobolev tests based on uniform scores and encompasses the existing multisample tests as particular cases. We further embed this class into a broader aggregation framework, where the samplewise Sobolev components are combined through general merging functions, including average-type, maximum-type, and interpolating families. Within this class, we introduce the first Anderson-Darling-type homogeneity test for circular data, and we further propose two new tests designed to detect multimodal departures from homogeneity, constructed from softmax and Poisson kernels. We derive the asymptotic null distribution of the class, prove its consistency against a broad family of fixed alternatives, and obtain the asymptotic distribution under shift-type local alternatives. The tests are distribution-free and therefore do not require resampling. A comprehensive simulation study demonstrates the strong power of the Anderson-Darling-type test compared with Cramér-von Mises-type competitors across several scenarios, and the effectiveness of the softmax and Poisson tests under several alternatives. The testing toolbox is applied to analyze the nursing patterns of polar bears.

Increasing wildfire activities in the United States have generated far-reaching effects on air pollution and human health. While prescribed fires are used as a preventive strategy to reduce wildfire risks, the health impacts of their smoke remain uncertain. Small-scale case studies indicate that prescribed fires generally produce lower pollutant emissions than uncontrolled wildfires, potentially benefiting air quality. However, the extent to which this translates into reduced population-level smoke exposure remains uncertain due to a lack of robust, large-scale quantitative evidence. Moreover, the factors driving differences in air quality impacts between prescribed fires and wildfires were underexplored. Using fire data in California, Florida, and Georgia (2006-2020) and a causal forest approach, we quantified and compared the impacts of prescribed fires and wildfires on air quality, assessed their differences in population-level smoke exposure, and identified key factors influencing their smoke exposure levels. On a per square-kilometer basis, prescribed fires resulted in over six times less population-level smoke exposure compared to wildfires. Additionally, prescribed fires were associated with lower smoke exposure in regions with higher relative humidity and precipitation. This work could enable the broad risk-benefit assessment to guide the judicious use of prescribed fires and optimize public health outcomes while mitigating wildfire threats.

Characterizing complex Multi-Input Multi-Output (MIMO) systems presents two issues: feedback controllers mask fault variance, rendering single-variable monitoring ineffective, and repeated testing across multiple conditions inflates false positive rates due to the multiple comparisons problem. This study proposes a multivariate statistical framework to resolve these limitations. We extend a statistical library that identifies system dynamics via Frequency Response Functions (FRFs) by transforming them into time-domain Pseudo-Impulse Responses (PIRs). This functional representation captures the complete dynamic signature of the system, offering a richer diagnostic profile than traditional static scalar metrics. The framework, originally developed for SISO systems, is extended to the MIMO case by introducing supervectors that aggregate multiple PIRs, evaluated using Permutational Multivariate Analysis of Variance (PERMANOVA). This non-parametric approach handles the high dimensionality and complex correlation structures inherent in MIMO functional data. We demonstrate that the PIR Supervector is more powerful than single-variable analysis, and the MIMO PERMANOVA approach outperforms traditional SISO anomaly detection. By distinguishing between normal operation, external thermal disturbance effectively compensated within the investigated operating range, and severe parametric faults that standard univariate methods miss, the framework provides a single, rigorous metric for specific fault diagnosis ($p < 0.001$).

We consider regression problems in which the marginal distribution of the covariates is informative for estimation and variable selection, rather than merely auxiliary. Motivated by random-design, high-dimensional, and latent-effect settings, we propose a general design-assisted regression framework in which the estimating criterion depends on both the conditional model for $Y \mid \bfX$ and structured features of the covariate distribution. The framework identifies two roles of design information: stabilizing weak design directions through quadratic regularization and correcting latent-effect distortion through nuisance augmentation. We establish oracle properties for the resulting estimator, separate the effects of stochastic error, shrinkage, and approximation, and compare it with a benchmark sparse procedure that ignores design information. These results show that the proposed framework improves estimation while preserving first-order prediction performance. Numerical studies and two real-data applications illustrate the practical impact of incorporating design information.

Classical convergence guarantees for stochastic gradient methods typically assume Lipschitz-smooth objectives and finite-variance gradient noise, both frequently violated in practice. In contrast, we study nonconvex stochastic optimization under the joint relaxation of these assumptions: objectives with $(L,s)$-Hölder continuous gradients, $s\in(0,1]$, and gradient noise satisfying only a bounded $\alpha$-th moment condition for $\alpha\in(1,2]$. We establish three convergence results. Firstly, that standard SGD converges at rate $O(T^{-s/(1+s)})$ whenever $\alpha\ge1+s$, extending the classical nonconvex SGD rate to heavy-tailed noise and Hölder smoothness simultaneously. Secondly, we analyze $\delta$-regularized gradient clipping ($\delta$-GClip), a provable trainer of wide and deep nets, and establish a stationarity rate of $O(T^{-2s(\alpha-1)/[(1+s)(2\alpha-1)]})$ under the same condition. Thirdly, we analyze standard gradient clipping (G-Clip) and show that it recovers the above rate for $\alpha\ge1+s$ while in the very heavy-tailed regime $\alpha<1+s$, it has a convergence rate $O(T^{-2s(\alpha-1)/[(\alpha-1)+s(2\alpha-1)]})$ — the first convergence guarantee in this regime for any stochastic gradient based method.