Papers

We propose a new approach for analyzing Time-Dependent Perturbed Angular Distribution data. In this, a likelihood is constructed with the help of conditional probabilities of single events and binning of time data is avoided. This likelihood is used in a Bayesian framework to calculate the posterior probability density function for the $g$ factor of interest. This approach is compared to more traditional approaches that use binned data by analyzing simulated datasets. In many cases the resulting posterior probability densities are observed to be multimodal and the results can thus often not be summarized with a single Gaussian approximation. We find that results from the new approach are more reliable for low-statistics datasets.

Ranks depend on the observations used for comparison. Reusing those observations can add association between forecast and outcome rank contrasts even when the evaluated forecast and outcome stay fixed. We characterize assignments that preserve association between specified population-rank contrasts, including forecast rank minus baseline rank compared with outcome rank minus baseline rank. Conditional on independent training, whole trajectories are sampled independently from a common law, with unrestricted dependence within each trajectory. The expected score separates into its target and an explicit interaction between map pairs. When reassigning references, we keep the learned maps, reference law and coefficient row sums fixed. Zero weighted reference overlap for every map pair is necessary and sufficient for preservation uniformly over permitted maps and laws. An unbiased three-trajectory kernel estimates the interaction; independent evaluation and validation provide finite-sample lower bounds. Complete U-statistics estimate the same target directly when all draws can be recombined. Sharp shared-baseline ranges, including ties, tighten both constructions. In a Beijing air-quality archive, interaction accounts for 91.4% to 94.5% of seven learned forecasts' expected shared scores under the empirical archive law. Separate results address category-fitting error and temporal feedback. In a matched category-control null experiment, distinct references reduce rejections from 402 to 41 out of 1,000 panels.

Quantifying the heterogeneity is an important issue in meta-analysis. The conventional heterogeneity statistic, $I^2$, heavily depends on study sample sizes and tends to approach one as the sample sizes increase. To overcome this limitation, a sample size independent measure, ${\rm ICC}_{\rm MA}$, has recently been developed to quantify the heterogeneity, yet this method is restricted to meta-analysis with continuous outcomes. In practice, however, binary outcomes are more widely encountered, with the log odds ratio (lnOR) frequently adopted as an effect size to compare the effects of two treatments in terms of odds. To fill the gap, we propose a new heterogeneity measure for meta-analysis with binary outcomes by introducing a latent continuous variable model. Under this framework, the proposed measure admits a closed-form expression within the random-effects meta-analysis model with lnOR as the effect size. We define this sample size independent measure as ${\rm ICC}_{\rm MA}^{\rm OR}$, and moreover propose a new estimator, $I^2_A$, to evaluate this measure. Simulation studies demonstrate that our new estimator is nearly unbiased and aligns closely with its counterpart applied for standardized mean differences when the continuous variable is observable. A real data application further demonstrates that the proposed heterogeneity statistic for binary outcomes and its counterpart for continuous outcomes provide interpretable and coherent estimates across different outcome types.

Graph neural networks are widely used for transductive node classification, with accuracy typically measured on randomly drawn train/validation/test splits. Reported accuracy has been shown to shift substantially across different random splits of the same dataset, making published comparisons between architectures unreliable. The classical remedy in non-graph settings is stratified $k$-fold cross-validation, which ensures each test fold reflects the full class distribution of the dataset. We argue that class stratification alone is insufficient for graphs: nodes are not isolated but connected, and folds that differ in their distribution of local neighbourhood homophily expose the model to systematically different relational conditions that directly affect message-passing behaviour. The resulting cross-fold variation reflects the homophily composition of each split, inflating reported variance beyond what model behaviour alone would produce. To address this, we propose \hp{}, a topology-aware stratification procedure that treats node homophily as the primary stratification axis, aligning folds with respect to local relational consistency alongside the class marginal that standard stratification already controls. Stratifying on homophily alone does not guarantee class balance, so \hp{} incorporates class label as a secondary axis, preserving class representativeness as a natural consequence of the procedure. We evaluate \hp{} on a broad benchmark suite comprising 15 node-classification datasets spanning the full homophily spectrum and 7 GNN architectures. \hp{} achieves a mean stability rank of 1.49 compared to 2.31 for random $k$-fold, achieving the lowest mean stability rank on 13 of 15 datasets while preserving class balance close to class-stratified splits and substantially better than random. We argue that homophily-aware split construction merits broader adoption for GNN evaluation.

We study federated soft clustering over federated learning (FL) networks of devices that each hold a private local dataset and fit a personalized Gaussian mixture model (GMM). Generalized total variation minimization (GTVMin) couples the local maximum likelihood problems through a graph regularizer that penalizes a discrepancy between the models of connected nodes. The choice of discrepancy measure is a key design decision: we compare a squared Euclidean distance between model parameters, which requires component matching, with two measures that compare the local model distributions directly and hence need no matching: a Monte-Carlo approximated Kullback-Leibler (KL) divergence and a closed-form maximum mean discrepancy (MMD). All three resulting GTVMin instances are optimized by synchronous projected gradient updates; for the smooth MMD instance we provide a convergence guarantee to stationary points. We characterize their computational cost and evaluate their robustness to data heterogeneity.

Visual signals require compact yet sufficient representations for robust downstream prediction. Convolutional sparse coding (CSC) provides an explicit mechanism for suppressing redundant components while preserving signal content, but its sparsity coefficient is often manually selected during training. We propose a training-adaptive convolutional sparse coding framework for robust visual signal representation. Specifically, we unfold the CSC optimization with the Fast Iterative Shrinkage-Thresholding Algorithm (FISTA) and treat the sparsity coefficient as a differentiable variable jointly learned with the network parameters. From the information bottleneck perspective, this coefficient controls the trade-off between information retention and compression: the sparsity term promotes compact representations, while the reconstruction term together with task loss preserves task-relevant signal content. We further introduce a label-free post-training strategy that adjusts the compression strength for corrupted inputs with the main network parameters fixed. Experiments on CIFAR and ImageNet demonstrate competitive clean-data recognition and greatly improved robustness under different input perturbations.

Randomized trials with limited concurrent control information—including but not limited to unequal-randomization settings—can offer ethical and practical advantages, especially in pediatric and rare diseases, but they often lose power because fewer control patients are available for direct comparison. Borrowing information from external controls may improve efficiency, but can also inflate the Type I error rate when the external and trial populations are not sufficiently comparable. We propose RECaST-Surv, a Bayesian transfer-learning framework for time-to-event outcomes that extends the RECaST method to survival settings. The method learns a structural survival model from external control data and calibrates it to the concurrent control arm of the current trial through a Cauchy random effect. To improve frequentist operating characteristics, we further develop a bootstrap-based procedure to calibrate the testing rule for Type I error control. RECaST-Surv can accommodate multiple external datasets and requires only summary-level information from external sources. Simulation studies show that the method maintains near-nominal Type I error across challenging settings while improving power by roughly 10%–12% over standard analyses. In an Amyotrophic Lateral Sclerosis trial emulation, RECaST-Surv increased power from 82.8% to 95.7% relative to the no-borrowing RCT analysis, while maintaining acceptable error control.

Kernel methods, and Gaussian Processes (GPs) in particular, require a Hilbertian distance measure—one whose square is conditionally negative definite (CND)—to guarantee positive semi-definiteness (PSD) of the kernel matrix; a condition that fails for many natural input spaces, including smooth manifolds and spaces of probability distributions. We propose the Sparse Landmark Embedding (SLE) kernel, which eliminates this requirement entirely. Each input is embedded into a sparse feature vector via compactly supported bump functions centered at all |D| training points; applying any standard PSD kernel in this embedding space yields a kernel that is provably PSD for arbitrary distance measures. The compact support automatically controls embedding sparsity, keeping kernel matrices well-conditioned and computationally tractable despite the high ambient dimension. We provide theoretical guarantees on PSD, sparsity, stability, and universal approximation, and demonstrate, using geodesic and Wasserstein distances, that the SLE kernel matches or substantially exceeds domain-specific baselines in both predictive accuracy and uncertainty quantification.

In the canonical Difference-in-Differences design, the control group's post-treatment change serves as an imputation of the treated group's counterfactual change in the same period, an imputation justified by parallel trends. However, differences in group composition can produce between-group differences in how outcomes would evolve over time, rendering this imputation vulnerable to confounding. An alternative imputation – such as one based on the treated group's pre-treatment change – avoids such between-group confounding but introduces the risk of confounding from within-group temporal shifts. Ideally, both imputations, each vulnerable to different sources of confounding, would have concordant values, thereby yielding the same causal conclusions. When the imputations are discordant, conclusions under parallel trends hinge more critically on that assumption since alternative imputations would point to different results. Yet in these scenarios, existing pretrends-based sensitivity analyses can show low sensitivity because they ignore post-treatment deviations from pretrends in the control group. This paper therefore proposes a discordance-based sensitivity model in which parallel pretrends are necessary but not sufficient for low sensitivity. I formally justify this model in terms of the expected distance between the ATT under parallel trends and under alternative assumptions, weighted by the joint plausibility of those assumptions. I then provide a decision-theoretic rationale for benchmarking violations of parallel trends using the worst-case discordance between the parallel trends imputation and alternative imputations. Finally, I apply both pretrends- and discordance-based sensitivity models to assess how a labor supply shock influenced electoral support for apartheid-era policies in South Africa, showing how the two approaches yield different results.

The test error of a model plotted against its number of parameters $d$ falls, peaks when the model can just fit the training data, and falls again, exhibiting the double descent phenomenon. We explain the phenomenon with statistical mechanics. The training trajectory of a stochastic gradient-based method is a particle wandering over the energy landscape of the training loss at an induced temperature $T$, and a run that has equilibrated visits every parameter vector of a given training loss equally often, the fundamental postulate of statistical mechanics, with probability given by the Boltzmann distribution. Because training starts at an initial point and has only finite time to diffuse, it carries an effective weight decay, which makes every parameter a quadratic degree of freedom. The equipartition theorem then distributes the energy among the $d$ degrees of freedom in shares of $T/2$, so at a fixed training loss adding parameters lowers the temperature and drives the Boltzmann distribution toward the stationary path. Finally, adding parameters can only lower the $L^2$ norm of the stationary path, so a solution sampled at fixed loss is less likely to be large with increasing $d$, effectively increasing weight regularization.

In physics-informed machine learning, a target function $u^*$ is learned from noisy value observations $y_i=u^*(x_i)+ \varepsilon_i$, together with differential information, given either by noisy observations $d_j=(Du^*)(z_j)+\xi_j$ or by a known physical constraint $Du^*=v$. We consider the setting where $D$ is a linear differential operator and analyze a physics-informed kernel estimator $\hat u$ combining $n$ value observations and $m$ differential observations. In this context, we ask how much can differential information improve predictions, and how does this improvement depend quantitatively on $n$, $m$, and $D$. We prove finite-sample bounds, supported by numerical simulations, revealing a two-regime structure for the prediction error. When $m$ is limited, the rate depends jointly on $n$ and $m$; when $m$ exceeds a problem-dependent threshold, the rate saturates and matches the oracle rate obtained when the perfect constraint $D \hat u = Du^*$ is imposed. Examples are discussed for Sobolev spaces which are reproducing kernel Hilbert spaces and include partial Laplacian constraints on the torus and gradient observations on bounded domains. These examples illustrate the range of possible learning rate improvements — from the standard nonparametric $n^{-1/4}$ to the parametric rate $n^{-1/2}$. Finally, we derive physically consistent rates in a stronger norm that jointly controls the errors in $\hat u$ and $D\hat u$.

Many datasets arise as multivariate time series observed at irregular time intervals, particularly in Earth observation and astronomical measurements. Classical time-series models assume that time is discrete and observations occur at equally spaced intervals, limiting their ability to capture dynamics when observation gaps vary. Existing approaches to irregular sampling rely on continuous-time formulations, which assume that observation intervals are sufficiently small. We address this limitation by introducing a discrete-time framework that accommodates irregular observation gaps while preserving multivariate dependence structures. We introduce a framework for the analysis of irregularly observed multivariate time series based on hypercomplex autoregressive processes. The framework embeds multivariate observations into a hypercomplex algebraic structure, allowing multiple variables and their interactions to be represented within a single mathematical entity while incorporating irregular observation gaps. The resulting process admits a structured matrix representation, which enables the model to be expressed as a state-space system. This formulation provides an estimation methodology in which model parameters can be inferred using Kalman filtering techniques. The state-space representation allows recursive estimation and prediction in the presence of irregular sampling intervals. The performance of the methodology is illustrated through data and applications to remote sensing and astronomical datasets, where observations occur at nonuniform time intervals. The results demonstrate that the approach captures dynamics and cross-variable interactions difficult to model using traditional time-series techniques. The proposed framework provides a methodology for analyzing irregularly sampled multivariate time series and opens possibilities for modeling observational data across disciplines.

Randomized controlled trials provide internally valid estimates of treatment effects, but their results may not directly apply to broader target populations due to differences in baseline covariate distributions, adherence to treatment or variations in outcome mechanisms. Under standard transport and time-to-event identifiability assumptions, we develop nonparametric, debiased machine learning estimators for transporting and generalizing causal survival treatment effect differences from a source population to a target population in discrete time. We derive the efficient influence functions for the transport and generalization survival difference estimands and propose cross-fitted one-step estimators that are doubly robust and achieve semiparametric efficiency bounds under weak regularity conditions. We further introduce estimators that exploit known effect modifier subsets through an additive parameterization of the survival function, reducing the dimensionality of the reweighting and yielding smaller or equal asymptotic variance. We establish asymptotic normality, double robustness, and rates of convergence for all proposed estimators. Finite-sample properties are illustrated through Monte Carlo simulations under flexible and misspecified nuisance estimation scenarios. We apply the methods to data from the Women's Health Initiative to estimate the effect of hormone therapy on coronary heart disease across trial and observational populations.

Generating signals on graphs requires permutation-equivariant models that exhibit stability with respect to relative structural perturbations. While recent graph-aware Schrödinger bridge models incorporate topology information directly into their reference dynamics, it is unclear how perturbations of the graph propagate through these dynamics and affect the resulting generated distributions. In this paper, we analyze the structural stability of graph-aware continuous-time generative models whose drift combines a graph filter with a learned graph neural network. We derive explicit Wasserstein stability bounds that quantify the effect of relative graph perturbations on the generated distributions. Motivated by these bounds, we introduce a principled framework for designing stable graph filters that preserve the smoothing behavior of graph heat diffusion, while boosting structural stability. Experiments on synthetic and fMRI signals show our stable filters enhance structural robustness while matching or exceeding the generative quality of the heat equation baseline.

Antibody lead optimization calls for a small, bounded set of edits to an existing candidate: substitutions, but also insertions and deletions. Edit-based generative models are the only ones that allocate such an edit budget without fixing the edit positions, the edit count, or the output length in advance. However, the existing approaches Edit Flows and EvoFlows did not release code or complete training specifications. Here, we show that both methods follow the same underlying process – edits firing one at a time, at learned rates, in continuous time – the pure-jump case of generator matching over finite sequences. With EditJumps we introduce the first open implementation of this framework, with a single generalist antibody editor trained on 1.66M Observed Antibody Space homolog pairs to propose homolog-like variants of a seed sequence, editing unseen leads zero-shot, without the per-family retraining original approaches require. Replicating this system from scratch exposes why open code is essential for generative biology: reconciling published edit distributions required reverse-engineering an undocumented rate-scaling hyperparameter that dictates realized mutation counts. Moreover, we show that published evaluation metrics are highly sensitive to reference sample size, frequently flipping method rankings. We release our full codebase, automated test suite, and configurations at: https://github.com/VisiumCH/editjumps

Parametric copulas capable of capturing full-range tail dependence in both lower and upper tails are useful in many applications, including as pair-copula building blocks in vine copula models. However, it is challenging to create such copulas that have closed-form density functions, which limits their practical applicability in complex modeling scenarios. To the best of our knowledge, all existing full-range tail dependence copulas do not have closed-form density functions. In this paper, we propose a new Archimedean copula family, called the full-range Archimedean copula of type I (FRA1), that captures full-range tail dependence in both tails and admits a closed-form density function. The copula has two parameters, one controlling the strength of dependence in the lower tail and the other controlling that in the upper tail. We provide results on dependence properties of the copula, and discuss how to effectively estimate the parameters based on maximum likelihood. We build a vine copula model only using the FRA1 copula, and demonstrate that the model can achieve comparable performance of an ordinary vine copula model while avoiding the need of pair-copula model selection. The new copula is also applied to data from the Medical Expenditure Panel Survey, demonstrating that an FRA1 vine model can effectively capture diverse temporal and cross-sectional dependence patterns in medical expenditures. Finally, the new copula is implemented in CopulaOne, an R package for full-range tail dependence copulas.

We develop a class of linear state space models for matrix-valued time series data where the state is a latent matrix normal process. We derive matrix versions of the Kalman filter, log-likelihood, and smoother enabling estimation of the latent state matrix as well as the model's parameters. To conduct Bayesian inference, we provide algorithms that draw from the joint posterior distribution of the latent state matrices conditional on the observed data and parameters. We apply these methods to a large panel of U.S. macroeconomic time series across the 50 U.S. states. The proposed framework accommodates mixed-frequency data, heteroskedasticity, and outliers within a unified matrix-valued structure. Empirically, we find that a small number of latent factors captures the joint dynamics across states and variables, providing a parsimonious and scalable approach to modeling high-dimensional macroeconomic systems.

Time-to-event outcomes are central in oncology and rare diseases, where treatment effects are often summarized by differences in survival curves or Restricted Mean Survival Time (RMST). In real-world data, estimating these causal effects relies on the absence of unobserved confounding, an assumption that is rarely satisfied. We develop a sensitivity analysis framework for causal treatment effects with survival outcomes under the Marginal Sensitivity Model (MSM). We introduce doubly valid and doubly sharp (DVDS) bounds for differences in survival functions and RMST, extending recent DVDS results to the time-to-event setting while accounting for informative censoring. In practice, our method yields tighter bounds and improved computational efficiency compared to a previous approach from the literature, on simulated and real data. For tractability, we assume independence between censoring and unobserved confounding, a limit that should be addressed in future works.

How many compensatory seats are needed to ensure (full) proportionality in two-tier electoral systems with multi-member constituencies? This question is answered by estimating the discrepancy between seats and votes rigorously, probabilistically and empirically. In 1919, Pólya showed that the Jefferson/D'Hondt method favours large parties. This seat surplus accumulates across constituencies. In some countries (e.g., Portugal, Poland, Spain and Turkey), this accumulation directly leads to quantifiable disproportionality in the final parliament. Elsewhere compensatory seats counteract this accumulated seat surplus. Using a logic-first approach, closed formulas for the required number of compensatory seats are derived given regional and national apportionment methods, the number of relevant parties $n$ and the size of the largest party $p$. If there are $c$ constituencies, it is derived logically that a two-tier system built on local D'Hondt apportionment in multi-member constituencies needs around $\frac{c}{2}(n-\frac{1}{p})$ compensatory seats to ensure proportionality. In words, that is half a compensatory seat per constituency times a number that increases with the number of parties and the size of the largest party. The validity of the derived formulas is analysed empirically using the extensive recent CLEA dataset. Many quantities do not depend much on the total parliament sizes, which means that the variance can be reduced by artificially varying the house sizes. Portugal, Poland, Spain and Turkey all employ a one-tier D'Hondt system and the analysis shows that if they were to add compensatory seats to the system, reasonable numbers would be 40, 75, 75 and 120, assuming that no other properties of the systems change.

This paper provides a self-contained proof of the rank of the pathlifting Jacobian of a DAG ReLU network by performing an induction on the network's number of hidden nodes. In fact, the induction is elementary, and the key recipe is to consider the skeleton matrix of the network, a sparse matrix encoding the network paths, and transform the representation of one of its hidden neurons into an output node. The proof relies on intermediate propositions which link the pathlifting, its Jacobian, the network parameters, and its skeleton matrix, which, on top of permitting to conclude on the rank of the pathlifting Jacobian, also provide a way to compute it without backpropagation and whose computation cost is super efficient in practice compare to usual backpropagation. The paper is provided with a Python module that implements the different propositions of the paper for feed forward networks and is used to experimentally quantifies the computational gain of computing the pathlifting Jacobian with the proposed theory.

We consider the problem of noisy gradient-free minimization of the k-th order partial derivative of a $\beta$-Hölder function supported on a d-dimensional cube. We show that T ^{($\beta$+d+k)/(2$\beta$+d)} log(T )^{(\beta-k)/(2\beta+d)} is a non-asymptotic minimax rate of the T step cumulative regret for all $\beta$ \ge 0. In the special case k = 0, our results cover the problem of noisy gradient-free minimization of $\beta$-Hölder functions, closing the existing gap between the known upper and lower bounds. We show that a minimizer of a suitably chosen local polynomial estimator is rate-optimal. The minimax optimal upper bound is achieved under the passive design, that is, when the query points are i.i.d. Thus, there is no advantage in considering sequential designs when it is only known that f is a $\beta$-Hölder function with no additional property. We propose an algorithm feasible in polynomial time that constructs a proxy of the minimizer of the local polynomial estimator. The procedure requires computing the estimator on auxiliary random points. The resulting polynomial time algorithm matches the lower bound.

Sign-based methods reduce communication costs in distributed environments, but aggregating local signs can introduce bias when data are heterogeneous. As a result, existing sign-based variance reduction methods fail to obtain the optimal convergence rates. In this paper, we solve this problem and obtain optimal rates for both nonconvex stochastic and finite-sum optimization. We first give a counterexample showing that majority voting can fail to approach stationary points even with exact local gradients. Motivated by this limitation, we propose tracking the global gradient at the server through unbiased compression of recursive gradient increments. As a result, we can obtain the convergence rates of $O(\sqrt{d/K}+\sqrt d (a/(nK))^{1/3})$ for the $\ell_1$-norm and $O(\sqrt{a/K}+\sqrt a/(nK)^{1/3})$ for the $\ell_2$-norm. Here, $K$ is the iteration number, $n$ is the number of workers, $d$ is the dimension, and $a=1+\omega$, with $\omega$ denoting the compressor's relative variance. For finite-sum problems with $M$ components, we combine periodic exact gradient refreshes with compressed component-gradient differences. The resulting total sample complexities are $O(M+d\sqrt{aM}\epsilon^{-2})$ and $O(M+a\sqrt M\ epsilon^{-2})$ for $\ell_1$ and $\ell_2$ gradient norms at most $\epsilon$, matching the corresponding bounds in centralized settings.

Classifiers can make identical predictions yet require labels to compare their selective performance: confidence ranks weight the same errors differently. We quantify this requirement for the area under the generalized risk-coverage curve (AUGRC). A prelabel lower bound rules out insufficient budgets. With all labels known, a covering linear program bounds the minimum number of labels sufficient to fix the winner (the certificate size) within $K-1$ labels for $K$ candidates. For fixed $K$, independent uniform orders and identical predictions, the prelabel bound approaches one quarter of the pool. With iid Bernoulli errors independent of the orders, every exact acquisition policy reads almost all labels asymptotically, although a two-candidate certificate needs only half. Across 108 feature-panel comparisons on nine datasets, disagreement labels settle every accuracy choice but no AUGRC choice. A 20% budget is ruled out in 96 conditions; certificates need 56-57% on average. On ten conditions with pretrained image classifiers, confidence-score choice reads 68-91% of 10,000 labels for exact selection and 50-67% with AUGRC tolerance $5\times10^{-4}$. An exact stopping test works with any acquisition order. Together, these results link confidence ranks to label budgets and certified model comparison.

Randomization or design-based inference is becoming an increasingly popular tool for analysing data from randomized experiments: It does not require modelling assumptions on the distribution of outcomes or covariates, and hypothesis testing and estimation are respectively valid and unbiased in finite samples. Yet, confidence intervals for the sample average treatment effect (SATE) are still constructed via finite-population central limit theorems and their coverage is only asymptotic. In this work, we explore an alternative approach: We use concentration inequalities to construct confidence intervals for the SATE with non-asymptotic guarantees. We develop this approach for the most common experimental designs (Bernoulli trials and completely randomized experiments) and provide Hoeffding and Bernstein-type confidence intervals. Moreover, we extend these results to matched-pair, cluster and stratified randomized experiments. Our key technical contributions are a novel Bernstein-type concentration inequality for i.i.d. data points as well as a concentration result for Neyman's variance estimator.

Causal discovery aims to recover causal relationships from observed data. In various fields, exploring causal relationships among variables remains an important topic, but this task becomes challenging due to the existence of latent confounders. Ignoring such confounders can lead to false associations and incorrect edge directions. In this paper, we study the linear structural equation model with latent confounders. We propose an algorithm that iteratively identifies terminal (observed) nodes and reconstructs the directed acyclic graph of the observed variables. To do this, we recover the precision matrix of the observed variables as a sparse plus low-rank matrix: a sparse matrix captures the conditional dependencies among observed variables, while a low-rank matrix captures the combined influence of a few latent confounders. We establish that for $p$ observed variables, $r$ latent confounders and $s$ edges, our procedure correctly identifies the directed causal relationship among observed variables, for $n \gtrsim \max\{s\log p,\ r p\}$ samples. Experimental results validate our theoretical contributions.

Diffusion samplers reverse a learned noising process using either stochastic (DDPM) or deterministic (DDIM) updates, which represent endpoints of a single family controlled by a scalar noise-injection variance that is applied identically at every spatial location. This uniform approach neglects the geometry of natural images: high-curvature regions such as edges and textures, where the denoiser is uncertain, benefit from stochastic correction, whereas smooth regions, where the score is precise, are degraded by injected noise. This work investigates whether each pixel requires stochastic correction at a given timestep and introduces Spatially Adaptive Noise Injection (SANI), a novel sampling framework that dynamically adjusts noise application on a per-pixel basis. SANI integrates a probabilistic gating mechanism with a derived spatially adaptive variance, ensuring that noise is injected precisely where needed to refine complex features while preserving well-formed structures. Experimental results and decoupling ablations demonstrate that SANI consistently improves Fréchet Inception Distance (FID) over the vanilla DDPM and DDIM endpoint samplers across diverse sampling timesteps, while remaining competitive with variance-learning baselines, highlighting the importance of spatial adaptivity in diffusion sampling.

Benford's law specifies a marginal distribution for significant digits, whereas the usual first-digit Pearson $p$-value is calibrated under an independent multinomial sampling model. We separate these statements with four constructions that share the same one-time or pooled Benford target but have different joint structures. Under a wrapped-Gaussian circular Markov construction, the nominal 5% Pearson test rejects 21.6% of samples at a fixed persistence level even though every one-time marginal is exactly Benford; a sequence-level second-order moment-matching calibration reduces the rejection rate to 4.8%. The same strategy performs well in the randomized-rotation and random-composition designs. The effect persists across changes in sample size and sequence length, and a grid approximation to a continuous log-significand Cramér-von Mises discrepancy shows the same design dependence. In the wrapped-Gaussian design examined here, corrected procedures retain rejection probabilities that rise with the size of a smooth marginal departure. Finally, we distinguish calibration from an information boundary: significand-only methods cannot detect changes that leave the complete significand process unchanged, but they readily detect perturbations that alter it. Benford conformity is a marginal statement; a Benford $p$-value is valid only relative to a specified statistic, sampling law, and calibration procedure; an integrity claim requires substantive competing models.

Stochastic subspace methods have gained popularity as gradient descent based techniques for large scale optimisation problems, especially in distributed settings. In this paper, we introduce the technique of "persistence of memory" to greatly extend and improve the random subspace methods. To this end, we leverage a vector that is only weakly correlated with the gradient in order to provide a guiding structure to the generative process of the random subspace along which the descent is going to take place. This guidance vector may be fixed for a large number of iterations, only to be refreshed at wide intervals (on whose size we can provide guarantees in terms of problem parameters). In important machine learning settings, such as optimisation problems embodying sparsity or a minibatch structure, we show that the guidance vector can be obtained in an effective and computationally inexpensive manner by leveraging the structured properties of the problem. En route, we establish to our knowledge the first theoretical analysis of classical SSD methods for sparse functions. In a local neighbourhood of the optimum, we demonstrate an alignment phenomenon of our gradient estimates with a low-lying eigenvector of the Hessian, allowing a once-for-all computation of the guidance vector which renders the method computationally favourable even in scenarios with unstructured objectives.

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

Reliable agent evaluation is complicated by automatic harness optimization, which repeatedly uses a released benchmark $B_{\mathrm{rel}}$ to guide a Proposer that edits prompts, memory, retrieval, tools, and control code around a fixed foundation model. Task holdout is commonly used to guard against harness overfitting. It varies semantic tasks but leaves the benchmark protocol fixed, so a bad genius Proposer can produce a cheating harness whose improvement over the initial harness on $B_{\mathrm{rel}}$ depends on a benchmark-wide shortcut. We introduce Counterfactual Harness Search and Evolution (CHASE), which casts harness evolution as constraint generation over valid counterfactual benchmarks. After each Proposer update, a Challenger searches for an executable protocol transformation with large gain destruction. A validity firewall checks that task semantics are preserved, while a held-out confirmation set determines whether the counterfactual enters a finite archive. We formalize an ideal shortcut-neutralized benchmark $B_0$ and establish theoretical guarantees linking finite counterfactual archives to $B_0$ and characterizing sequential Challenger search. We evaluate CHASE on Syn-Ledger and OfficeQA, where CHASE retains strong released-benchmark gains while substantially reducing gain destruction under valid protocol transformations.