Papers

We develop an exact posterior-sampling framework for Bayesian inverse problems when selected bounded forward observables can be accessed through Bernoulli events. A Bernstein–Poisson construction converts these forward coins into scaled Gaussian likelihood coins, and thinning an inflated prior Poisson point process yields posterior atoms that are iid conditional on their number. For independent Gaussian observations, we derive an exact mean-work identity for the implemented early-stopped factory and sharp small-noise complexity laws governed by local prior-predictive mass near the exact-fit set; a factorial-moment construction extends the likelihood factory to correlated Gaussian errors. The posterior algorithm is model-agnostic once Bernoulli access is available. As one continuum realization, we use Feynman–Kac sampling, Poisson killing, and lazy random-series evaluation for a bounded elliptic resolvent problem with a function-valued coefficient. Numerical experiments validate the forward and likelihood coins, support the predicted work regimes, and demonstrate posterior sampling without deterministic spatial discretization or fixed parameter truncation in the target. A matched-accuracy benchmark against finite-difference prior rejection illustrates how deterministic discretization bias changes the posterior accuracy–cost balance.

Isotonic regression and mean estimation over a simplex are among the most classical denoising problems under geometric constraints. Chatterjee and Lafferty generalized both to the recovery of a flow on a rooted tree from Gaussian noise, with isotonic regression on a path and the simplex on a star. Their results reach only special trees; on a general tree the minimax risk remained unknown for a decade. We settle the question in full. On every tree with $n$ vertices, at every budget $V$ and noise level $\sigma$, the minimax risk is of order $\min\{V^2H_T,\ \sigma^2k_0\}$, with $H_T$ the diameter and $k_0$ the crossing index of a truncated ancestor profile. A deterministic estimator, an exponentially weighted aggregate over integer states, attains it in $O(n^2\log(n+1))$ arithmetic operations. Least squares at a known budget, the natural convex program for the model, is suboptimal in the worst case by a factor of order $(\log n)^{2/5}$. One functional, the truncated ancestor profile, carries the rate, the estimator, and the least squares comparison. It serves as covering radius, as packing entropy in the positive cone, and as the range of the one-dimensional state behind the exact computation. Two classical theories become one on trees, sharp and computable throughout.

Consider Gaussian denoising with a fixed prior and varying noise variance. When the prior has a density, maximum a posteriori estimation uses its negative logarithm as a fixed regularizer; only the weight of the quadratic data term changes with the noise variance. We ask whether a fixed regularizer can also recover the posterior mean. Although a suitable regularizer exists at each variance, we prove that no single regularizer works at two distinct positive variances for any compactly supported prior other than a point mass. This holds even for nonconvex and nonsmooth regularizers. A residual computed from two denoiser evaluations gives positive lower bounds on the uniform approximation error and excess mean squared error of estimators using a common regularizer. Numerical examples illustrate these bounds and the cost of reusing a regularizer across noise levels.

We investigate the nonparametric estimation of the diffusion matrix in stochastic differential equations featuring multidimensional, strong mixing covariate processes. We propose a flexible statistical framework based on general function classes that does not require a linear basis representation, rendering our results directly applicable to deep neural network estimators. Our approach employs a two-step estimation procedure: constructing a preliminary nonparametric quasi-likelihood estimator and subsequently regularizing it via a $\beta$-Hölder class approximation. We establish general risk comparison theorems between empirical and generalization risks for arbitrary estimators without relying on a specific probabilistic structure of the underlying process. In diffusion matrix learning based on $n + 1$ observations over the time interval $[0,T]$, the derived upper bounds capture the intrinsic interplay between the $T$-rate, associated with the mixing behavior of the covariate process, and the intrinsic $n$-rate governing the volatility estimation.

Developing an asymptotic theory for restricted maximum likelihood (REML) estimators in generalized linear mixed models (GLMMs) has remained an open problem for decades. In this paper, we establish the asymptotic distributions of the REML and maximum likelihood (ML) estimators for GLMMs when both the number of clusters and the cluster sizes tend to infinity. Under very mild conditions, requiring only finite-moment assumptions on the random effects rather than normality and imposing no restriction on the relative rates at which the number of clusters and the cluster sizes diverge, both estimators are proved to be asymptotically normal with an explicit block-diagonal covariance structure, and the bias-correction property of REML is established. Simulation studies support the asymptotic results, and a genetic data analysis illustrates the methodology.

The Swing filter extends a feasible line through its previous recording and restarts from a fitted value at the last accepted sample. We show that this inherited endpoint causes a singular loss of compression on rough input. For uniformly sampled Brownian motion at fixed tolerance, every fixed number of recording times converges to the first continuous feasibility boundary, with alternating boundary errors. The longest post-first segment tends to zero in probability. The segment-count ratio against a continuous approximation through true sampled knots diverges almost surely whenever that boundary precedes the observation horizon. For unit Gaussian increments, constrained least squares and no maximum lag, the first mean segment length has quadratic order in the tolerance $E$, the second has linear order, and every fixed mean from the third onward has order $E\log E$. We explain the logarithm through the conditional survival tail and prove a geometric bound on the subsequent approach to the stationary coefficient in the large-tolerance limit. The stationary mean satisfies $C_E\sim\gamma_{\mathrm{Swing}}E\log E$. An explicit one-dimensional kernel specifies the single coefficient $\gamma_{\mathrm{Swing}}\approx1.8120703$ and admits a $3/4$ Wasserstein contraction and convergent evaluation brackets. The decimal is not a certified error interval. A half-budget slope corridor, true endpoint chords and one-sample bridges repair the degeneration while preserving the original error bound, continuous output and constant working memory. Including bridge costs, the long-stream mean span per output segment is asymptotic to $21\zeta(3)E^2/(8\pi^2)$. Its explicit coefficient follows from the known Brownian anchored lifetime. Thus a change to the recording rule restores the quadratic scale.

In this paper, we introduce a variance-based inaccuracy measure for order statistics based on the concept of varextropy. The proposed measure provides a framework for quantifying the variability of the inaccuracy function between an order statistic and its parent distribution. We establish a transformation property of the proposed measure under strictly increasing differentiable transformations. We establish theoretical upper and lower bounds for this measure. We further develop a dynamic residual version of the proposed inaccuracy measure and establish characterization results for the underlying distribution. Furthermore, we derive bounds for the dynamic residual inaccuracy measures based on extropy and varextropy. Finally, a nonparametric estimator of the proposed measure is developed and applied to real data to assess the goodness-of-fit of candidate distributions and facilitate model selection. 1.

A switchback experiment alternates an entire system between treatment and control over time. It is especially useful when interactions between units can undermine standard unit-level experiments. Switchback experiments are commonly implemented on fixed temporal grids, which impose highly structured restrictions on when treatment can switch. We study a broader class of random-duration switchbacks, in which treatment and control alternate across runs whose durations are drawn from a common distribution. Under finite carryover, we show that the mean squared error of the Horvitz-Thompson estimator has a simple frequency-domain representation governed by how temporal outcome patterns align with design-induced imbalance and contamination patterns. This representation yields a tractable worst-case design criterion that can be computed directly from the duration distribution. Optimizing over even a simple two-rate family produces a run-age-dependent switching rule that reduces the asymptotic worst-case mean squared error by at least 32% relative to the standard independently randomized fixed-block switchback.

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.

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.

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.

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.

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.

In this work, we provide a new perspective on semiparametric efficiency theory. In particular, we reformulate the questions of optimal efficiency and its attainment as a variational problem: minimize variance over estimating functions subject to robust unbiasedness constraints. We develop this theory for marginal-integral functionals without requiring pathwise differentiability. We study estimating functions whose expectations remain equal to the target when any one specified nuisance component is misspecified and the others are correct. We characterize the infimum of their variances through conditional variance minimization. For maxima of affine functions of treatment-specific conditional means, we obtain explicit optimal weights and construct cross-fitted estimators that attain the bound under conditions on nuisance estimation and the probability of near ties. For analogous maxima based on jointly observed quantities, the optimal weights use the full conditional covariance matrix. We also identify conditions under which the variance bound agrees with a classical convolution bound for parametric perturbations that preserve ties to first order. Examples include optimal policy values, $L^1$ calibration error, Balke–Pearl bounds, and mediation parameters. In an application to the National Longitudinal Survey of Young Men, covariance weighting reduces the median estimated variance of the cross-fitted estimating function relative to equal weighting by 12.9% for the lower Balke–Pearl endpoint and 16.5% for the upper, with improvements in all 20 repeated cross-fitting splits.

We estimate the scale of an exponential lifetime from observations censored randomly from both sides, under an informative model in which each censoring law is a power of the lifetime survival function – a two-sided generalization of the Koziol–Green proportional hazards model. Substituting the empirical distribution function into the likelihood equation yields a closed-form pseudo-maximum-likelihood estimator of the scale that is computed in one pass over the order statistics and needs neither iteration nor numerical optimization. A digamma identity shows the estimating equation to be exactly, not merely asymptotically, Fisher-consistent. Strong consistency and asymptotic normality are established with no restriction on the censoring depth; the key is an identity that rewrites the underlying $L$-statistic without its unbounded score function and so removes the condition that classical limit theorems would impose. The influence function and the asymptotic variance are explicit in polygamma functions, so confidence intervals require no numerical integration. Against the information bound of the observed-data model the efficiency exceeds $98.6\%$ in the designs considered and stays above $96.2\%$ over a wide sweep. Only the entry half of the model is needed: the estimator remains exactly Fisher-consistent and strongly consistent whatever the right-censoring law, and exponentiality is shown to be the precise price of that freedom. The construction extends verbatim to a proportional hazards class with known baseline. Three data sets illustrate the procedure, including a cohort in which left censoring is genuine.

We establish the stochastic-gradient query complexity of sampling smooth strongly log-concave distributions in every fixed dimension $d\geq1$. For $\mu$-strongly convex, $L$-smooth potentials with unknown minimizers $x_f^\star$ in the ball $\mathbb{B}(0,\mu^{-1/2})$, and unbiased gradient oracles with variance at most $\sigma^2$, the minimax worst-case expected query complexity is $\Theta\left(\log(1+\kappa)+\frac{\sigma^2}{\mu\varepsilon}\right)$, jointly optimal for the condition number $\kappa:=L/\mu$, variance $\sigma^2\ge 0$, and TV accuracy $0<\varepsilon\leq 1/10$. In the noiseless setting $\sigma=0$, this tight complexity $\Theta\left(\log(1+\kappa)\right)$ is independent of $\varepsilon$. Moreover, a gradient-only sampler generates an exact sample with ${O}(\log(1+\kappa))$ worst-case expected queries. Without previous ball $\mathbb{B}(0,\mu^{-1/2})$, exact sampling from any initial point $x_0$ can be implemented with expected cost ${O}\left(\log(1+\kappa)+\log(1+\sqrt{\mu}\|x_0-x_f^\star\|)\right)$, without knowing the initial distance. We also show that dependence on initial distance is generally unavoidable.

Community detection in bipartite networks is a fundamental problem in modern data analysis, with applications in recommendation systems, biological networks, and social network analysis. Unlike conventional unipartite graphs, bipartite networks consist of two distinct types of nodes with edges only connecting across types, so recovering latent communities requires estimating labels on the two node types. The stochastic co-blockmodel is a classical probabilistic framework for such networks, yet theoretical guarantees for exact community recovery in this setting remain limited, especially when the number of communities grows, the community sizes are unbalanced, or the degrees are heterogeneous. In this work, we prove that a simple spectral clustering algorithm based on the diagonal-deleted Gram matrix achieves exact recovery with high probability under mild conditions on sparsity, community balance, and the number of clusters. We further extend the result to the degree-corrected stochastic co-blockmodel, where each node carries its own degree heterogeneity parameter, and show that a row-normalized version of the same algorithm maintains the exact recovery guarantee. Extensive experiments validate our theoretical findings.

We study an online sketched Newton method that approximates the Newton direction at each step via a state-of-the-art sketching solver, called the generalized accelerated sketch-and-project solver (GAS), thereby mitigating the computational bottleneck of classical second-order methods. The GAS solver improves upon vanilla, unaccelerated sketch-and-project solvers by achieving accelerated convergence through Nesterov momentum updates, and accommodates a flexible projection metric whose proper choice further reduces computational cost. Building on this design, we establish asymptotic normality of the averaged sketched Newton iterates and characterize their limiting covariance matrix. The resulting covariance recovers that of the unaccelerated sketched Newton method under a specific choice of acceleration parameters, converges more rapidly (in the number of sketching steps) to the minimax-optimal covariance in general, and is smaller than that of the last iterate produced by the accelerated method. Finally, we strengthen these results by establishing a functional central limit theorem for the Newton iterates, which allows us to bypass explicit covariance estimation and develop an online inference procedure based on random scaling. Specifically, we construct a pivotal test statistic by appropriately rescaling the averaged iterates, so that its limiting distribution is free of any unknown parameters, enabling asymptotically valid online inference. Numerical experiments demonstrate superior performance of the proposed inference procedure.

We study the problem of estimating unknown parameters in PDE-constrained inverse problems from noisy observations, where the PDE solution is approximated using Physics-Informed Neural Networks (PINNs). While PINNs have demonstrated remarkable empirical success, existing estimators often inherit the slow nonparametric convergence rate of the neural-network solution, leading to biased and statistically inefficient inference for the finite-dimensional parameters of interest. To address this, we propose a two-step debiased estimation procedure that combines neural-network-based nonparametric estimation with an influence-function-based bias correction. By eliminating the first-order sensitivity of the estimator to errors in the nuisance function, our procedure yields a $\sqrt{n}$-consistent and asymptotically normal estimator without requiring undersmoothing of the neural network component. We further extend this framework to Bayesian inference by replacing the original likelihood with a debiased quasi-likelihood and establish a Bernstein-von Mises theorem showing that the resulting posterior contracts at the $\sqrt{n}$-rate with an asymptotic covariance matching that of the frequentist estimator. As a by-product of our analysis, we establish near-minimax optimal convergence rates for estimating a nonparametric regression function and its derivatives in Sobolev spaces using neural networks. Extensive numerical experiments corroborate our theoretical findings and demonstrate the necessity of the proposed debiasing procedure for valid statistical inference in PDE-constrained inverse problems.

Sampling error affects the finite-sample values of residual-based fit indices. We show that the unweighted comparative fit index ($CFI_u$) rests on a quadratic-kernel maximum mean discrepancy and that a Gaussian-kernel comparative fit index (KCFI) converges to $CFI_u$ as bandwidth increases. Because kernel discrepancies lack an analytic sampling correction, we extend expected-discrepancy bootstrap methods to a fitted-correlation reference state: the sample is transformed to reproduce the fitted correlation matrix, then resampled and refitted to estimate the exact-fit sampling component. Across 18 confirmatory factor analysis conditions, this correction reduced the root mean square error (RMSE) of both indices by 62%-87%. Under normality it matched an analytic residual-covariance estimator of $CFI_u$ within .0001, and in small samples it outperformed ordinary bootstrap correction in accuracy and availability. Corrected KCFI had lower RMSE than identically corrected $CFI_u$ in all 29 nonnormal conditions, even against a common target, but in six heavy-tailed small-sample conditions the analytic estimator was more accurate than both. The results validate the reference-state correction where it is needed and delimit where the analytic estimator is preferable.

We study estimation of a linear functional in random-design linear models when $\gamma=p/n \in (1, \infty)$ under the dense regime where the design covariance is unknown, and neither the regression vector nor the precision matrix is assumed sparse. Suppose the spectrum of the design covariance belongs to $[\ell,L]$, we establish the minimax mean squared rate for the linear functional has the order $n^{-\alpha(\kappa,\gamma)}$, where the exponent $\alpha(\kappa,\gamma)$ has an explicit elliptic-integral representation and $\alpha(\kappa,\gamma) < 1$ when both $\kappa > 1$ and $\gamma > 1$. We propose the Bocca estimator, which combines empirically calibrated resolvents through a complex contour integral to attain the optimal rate.

One method of generalizing a distribution is to add new parameter to the baseline distribution. In this paper the author discusses this methodology to create a new six parameters generalized Weibull-Kumaraswamy-{Lomax) Distribution. The author utilizes the approach of the T-X{Y} family to add new parameters the baseline distribution Kumaraswamy distribution defined on the unit interval with two shape parameters a & b>0. The author used the Lomax distribution with two parameters {\alpha} & {\lambda}>0 , scale and shape parameters respectively, to link the baseline distribution with the generator distribution. The generator distribution is the Weibull distribution with the two parameters, c & {\gamma}>0 , shape and scale parameters respectively as will be explained in the paper. The result is a new distribution with six parameters. The author derived the PDF, CDF, and quantile functions. The author used a new technique to expand the pdf in a linear series expansion that is simple and does not require special function. This makes the computational statistics for the raw moments easier and simpler. The author discusses the Shannon entropy and the paper contains figures for the PDF, CDF and the quantile functions.

We establish a sharp pathwise signal-to-noise criterion for quasi-maximum-likelihood (QML) estimation of a dominant breakpoint in the second-moment structure of a multivariate time series: the QML estimator is consistent whenever the between-regime contrast exceeds the within-regime fluctuation by an explicit factor, and below this threshold global recovery can fail. Two innovations drive the result. First, the framework is pathwise: no stochastic model is imposed on the data. Second, the log-determinant objective carries an additive ridge regularization: it removes the endpoint boundary layers, so no trimming of the candidate set is required, and tuning the ridge weakens the consistency condition. As an application of the main theorem, we establish consistency of QML breakpoint estimation in pervasive factor models whose dimension diverges with the sample size, for completely general error terms – not necessarily independent, idiosyncratic, or even random.

In this work we consider numerical models with random input parameters, where the model output is a spatially distributed random field. We carry out a statistical sensitivity analysis of spatial sets which arise for instance when the model output exceeds a critical threshold. We consider two approaches: (i) a kernel-based sensitivity analysis working with the Hilbert–Schmidt Independence Criterion, and (ii) a function-valued sensitivity analysis working with generalized Sobol' indices, where the model output is the indicator function of the spatial set of interest. We develop efficient estimators for the sensitivity indices, especially for finite-element-based numerical models, where we approximate volume integrals in the sensitivity measures by finite element quadrature, thereby avoiding Monte Carlo sampling over the spatial domain. This improves the state-of-the-art for the kernel-based indices and makes these methods practically accessible for expensive numerical models. For the generalized Sobol' indices we utilize the same finite element quadrature for efficient computation. We present numerical experiments for a hydrogen combustion process, modeled by a coupled system of convection-diffusion-reaction equations in a two-dimensional spatial domain, and ask whether the temperature remains below a critical value in selected regions of the combustion domain. Our results show that in this example the kernel-based and the generalized Sobol' indices give qualitatively similar input importance rankings. Moreover, the ranking depends on the chosen observation window and can flip between different windows. The algorithms for implementing the proposed methods are provided in the supplementary materials and in an open source code repository.

The Lehmer transform of a positive signal is the ratio of its power sums at consecutive orders, read as a function of the order, and it is a strictly increasing analytic curve from the smallest observation to the largest. Composing this curve with a map whose inverse is the Lambert function turns it into a probability law on the real line, the Lehmer-Lambert distribution, which we introduce and study. Its distribution, density and quantile functions are closed-form, its random variates are exact, and the law is invariant under changes of gain and of power-law calibration of the signal. Its modes lie at orders at which the signal changes scale, and for a sample its tails are exponential with rates equal to the two extreme logarithmic gaps of the sample. We derive moments, generating functions and limit theorems, and we investigate applications of the Lehmer-Lambert distribution in the analysis of electroencephalograms of patients with major depressive disorder.

Generalized latent factor models provide a flexible framework for analyzing high-dimensional non-Gaussian data, but principled estimation and uncertainty quantification under missingness remain substantially less developed. We develop a theory that connects a computationally tractable nonconvex procedure directly to statistical inference for nonlinear latent factor models with exponential-family links and partially observed entries. Our procedure combines a link-aware double-SVD initialization, a unilateral refinement that achieves rowwise consistency, and vanilla gradient descent. We show that the refined initializer enters a region of incoherence and contraction and that gradient descent remains in this region through implicit regularization, contracting rapidly down to the statistical estimation error without explicit incoherence or balancing regularization. Our central result is a uniform rowwise linear approximation for the actual output of gradient descent that isolates the leading score fluctuations from higher-order estimation and optimization errors. These expansions yield asymptotically valid individual and Gaussian multiplier-bootstrap simultaneous inference for latent factors, together with simultaneous confidence bands for missing-entry means, without requiring an additional debiasing step. The resulting estimation rate matches a restricted-class minimax lower bound up to logarithmic factors, while the theory accommodates severe missingness, weak low-rank signals, and diminishing local curvature. Simulations support the theoretical findings, and an application to large language model evaluation illustrates uncertainty-aware estimation and ranking of latent model capabilities.

We introduce statistical theory for the matching of finitely many objects, represented as metric measure spaces (mm-spaces). The approach is based on the second lower bound (SLB) of the Gromov-Wasserstein distance and thus is able to identify deviations in the distributions of the (pairwise) distances within each mm-space. We introduce a surrogate of the SLB barycenter which can be easily computed and expressed explicitly in terms of the distance distributions of each object. When comparing $m$ mm-spaces for $n$ randomly drawn samples in each space, the resulting statistic then can be calculated efficiently in $O(m \cdot n^2 \log(n))$ basic operations. We derive the asymptotic distribution and finite-sample bounds of the proposed test statistic, which serves as a basis for a variety of tools for statistical inference, specifically an asymptotic test for pose-invariant object discrimination and a classification method (based on the SLB barycenter) with controlled error rates. These methods are investigated in simulations and applied to the structural comparison of protein domains.

Random geometric graph detection asks whether a graph generated from $n$ independent uniform points on $S^{d-1}$, with edge probabilities depending on inner products, can be distinguished from an Erdos-Renyi graph of the same edge density using only its adjacency matrix. Mao, Wu, and Xu conjectured that the cubic trace of the standardized kernel determines the detection boundary: $n^3(\mathrm{tr}\,K^3)^2$ tending to zero or infinity should imply impossibility or strong detection, respectively. We construct counterexamples and prove a uniform result for spherical polynomial kernels of uniformly bounded degree with mean edge density bounded away from zero and one. For the integral operator $K$ of the standardized centered kernel, define \[ n_* = \min\left\{ \frac{1}{|\mathrm{tr}\,K^3|^{2/3}},\quad \frac{1}{\sqrt{\mathrm{tr}\,K^4}},\quad \frac{d}{\mathrm{tr}\,K^2} \right\}, \] with zero denominators interpreted as infinity. The total variation distance between the geometric and independent-edge models tends to zero if $n/n_* \to 0$ and to one if $n/n_* \to \infty$. The three scales correspond to triangles, four-cycles, and global geometry. The pure degree-four spherical harmonic kernel has only positive nonzero eigenvalues, yet its detection scale is $n \asymp d^5$. Throughout $d^5 \ll n \ll d^{16/3}$, the triangle and four-cycle signals and every fixed-degree standardized polynomial mean gap vanish, while the full graph remains strongly detectable. Thus the cubic-trace conjecture fails even without spectral sign cancellation. The impossibility proof combines higher-order spherical integration by parts, a graph expansion of cumulants, and a finite-order relative entropy comparison.

A novel, non-parametric clustering algorithm is developed, based on flows in the space of soft assignments $P_k^i$, representing the probability that point $x^i$ belongs to class $k$. These flows have two drivers: a nonlinear reaction component that relaxes $P$ to an evolving Bayesian-like posterior, and a diffusive component that relaxes $P$ to its local average, with a diffusivity $\nu$ that adapts dynamically so as to balance the two components. The methodology extends to semi-supervised classification, where a subset of the labels are known. In the limit of infinitely many observations, it yields a set of non-standard reaction-diffusion equations, which produces sharp boundaries between species that separate spontaneously into well-balanced domains. Including more than one diffusive network opens the way to a broader class of applications, including the detection of regime changes in time series and a clustering procedure based on numerous features that defeats the curse of dimensionality. In the continuous limit, this extension results in a novel, non-local class of diffusive operators.

Identifiability of tensor decompositions is often established through linear-algebraic conditions on the factor families. For nonnegative decompositions, however, positivity provides additional information that is not captured by dimension and independence alone: nonnegative terms cannot cancel, and their supports constrain competing decompositions. We introduce a positive scattering term that quantifies this additional source of identifiability and combine it with the dimension budget underlying the Lovitz–Petrov generalization of Kruskal's theorem. For every subset of components, we obtain two sufficient conditions: a threshold of $2|S|-2$ guarantees minimality and nonnegative rank, while the stronger threshold $2|S|-1$ guarantees uniqueness among nonnegative decompositions of the same length. The key result is a positive splitting inequality for irreducible exchanges of nonnegative rank-one tensors, which combines the dimension constraint with support-induced geometric rigidity. Although the scattering term is defined through an optimization over intermediate factor spaces, we show that its mode costs are exactly $0$, $1$, or $+\infty$, yielding an exact activation characterization in terms of graph connectivity. The resulting criterion can strictly certify sparse nonnegative tensor decompositions beyond the reach of Kruskal and Lovitz–Petrov conditions, including examples for which those conditions fail even after reshaping. In the matrix case, the two criteria reduce respectively to full-rank factorization and two-sided separability.