Papers

We study the estimation of the interaction parameter in singular interacting particle systems from continuous and discrete observations over a fixed time horizon. The singularity of the interaction drift prevents the direct application of standard likelihood and contrast-based methods, which rely on regularity assumptions that fail in this setting. Our approach exploits an exact reduction of the particle dynamics to a one-dimensional stochastic differential equation for a quadratic empirical moment. By the symmetry of the interaction kernel, the singular terms cancel, yielding an explicit estimator in which the unknown parameter enters linearly. We establish consistency and asymptotic normality as the number of particles, \(N\), tends to infinity. Letting \(\varepsilon_N\) denote the microscopic noise scaling of the system, we obtain asymptotic normality at the rate \(\sqrt{N/\varepsilon_N}\) for strictly positive interaction kernels, under an additional condition on the limiting empirical information. We further extend these results to degenerate interaction kernels, achieving consistency and asymptotic normality at rate \(N\) under suitable structural assumptions and the canonical scaling \(\varepsilon_N=N^{-1}\). In both regimes, we provide conditions on the observation mesh ensuring that the statistical limits remain valid under discrete sampling, in high frequency. Finally, we illustrate the applicability of our methodology through Dyson Brownian motion, the Dyson Ornstein–Uhlenbeck process, and Wishart–Laguerre particle systems, alongside a numerical investigation of the estimators' finite-sample performance.

It is folklore that passing "convergence" diagnostics does not guarantee that a Markov chain has converged. We prove this as an impossibility for any finite checklist of diagnostics. Running more chains, longer chains, and smarter diagnostics does not avoid this impossibility in certifying convergence. Under regularity assumptions, we determine how the worst-case undetected non-convergence decays as more diagnostics are added, and how to choose diagnostics optimally. Folklore suggests that the sampler's slowest mixing modes are the "dangerous" ones that must be monitored; we prove that the slow modes yield nearly optimal diagnostics as sampling proceeds. Given a source class, we also give an anytime-valid bound on the remaining discrepancy and show how undetected bias caps the effective sample size that a finite checklist can certify.

We study the minimax estimation of the pseudo-spectral gap $\gamma_{ps}$ of a Markov chain from the observation of a single trajectory. Extending earlier works on finite-state space chains by Wolfer and Kontorovich (2024), we provide an estimator in the case of a compact, continuous state-space. Our estimator appears to be minimax-optimal on the class of chains with a transition operator with Lipschitz density, up to logarithmic terms. Our estimator, combined with Paulin's version of Bernstein inequality, allows to derive empirical concentration inequalities for Markov chains. This unlocks many applications, from confidence intervals, to active learning and generalization bounds. In particular, for rested bandits with Markov rewards, we provide a version of UCB with a regret bound similar to the one in the iid case, without assuming any knowledge on the pseudo-spectral gap of the chains.

This paper introduces a new cumulative residual Rényi entropy (CRRE), which addresses limitations of the existing cumulative residual Rényi entropy. The conventional measure suffers from two primary drawbacks: its failure to converge to the cumulative residual entropy (CRE) as the order parameter approaches unity, and its lack of guaranteed positivity. Our proposed CRRE resolves these issues, offering a theoretically consistent and always nonnegative measure of uncertainty in residual lifetimes. We explore its mathematical properties, including its relationship with existing information measures. Furthermore, the concept is extended to its dynamic version, and its time-dependent characteristics are investigated. A non-parametric estimator of the proposed CRRE is developed, its almost sure convergence is established, and is validated through Monte Carlo simulation. The effectiveness of CRRE as a tool for quantifying uncertainty is confirmed through its application in chaotic theory and financial data.

The maximum likelihood degree of a Gaussian graphical model measures the algebraic complexity of its likelihood equations. This paper studies maximum likelihood degrees of Gaussian graphical models through chordal completion. When a graph becomes chordal after one missing edge is added, its maximum likelihood degree is determined by the maximal cliques and clique-tree separators containing that edge. This recovers Uhler's formula for complete bipartite graphs with one part of size two, and extends it to deleted-spine clique books, whose degree is unchanged when a page is thickened to a clique. For several fill edges, the same method gives a product formula when no clique or separator bag contains two fill edges. We also determine four classes of vertices at which clique blow-up preserves the degree. Two follow from known results, and two are established here.

Individualized survival predictions are often summarized by scalar risk scores or event probabilities at selected time horizons, although clinically relevant risk may also differ in its temporal organization. We consider finite-dimensional representations of individual hazard trajectories through interpretable hazard functionals, such as peak timing, cumulative hazard over prespecified intervals, and persistence of elevated risk. A central issue is that the richness of such representations is constrained by the hazard geometry supported by the underlying survival model. We formalize this through model-induced tangent spaces and functional ranks, derive an exact decomposition of model-supported variation into directions retained and lost by a chosen representation, and introduce a model-relative notion of local completeness. We further distinguish intrinsic model capacity from the variation actually activated by covariates. These results motivate a structured extension of AFT time scaling through a positive-affine deformation on the log-time scale, yielding a power-warp mechanism with an additional temporal deformation direction. Combined with a log-logistic reference hazard, this gives a parsimonious likelihood-based working specification for which peak timing and peak height provide complete coordinates over the unimodal regime. We also establish plug-in recovery and asymptotic inference for hazard functionals. Simulations and a lung-cancer application illustrate how model-induced restrictions propagate to hazard-functional representations.

Differential Privacy (DP) provides a rigorous framework for quantifying privacy guarantees. While most methods achieve DP by injecting calibrated random noise, the intrinsic randomness of certain procedures such as sampling can yield DP guarantees “for free”. In this work, we develop a unified Rényi divergence framework to characterize the inherent DP guarantees of sampling from a Bayesian posterior distribution. We show that the privacy loss incurred by One-Posterior-Sampling (OPS) is governed by posterior exponential moments of the single-record log-likelihood ratio, and establish explicit DP guarantees under uniformly bounded, sub-Gaussian, or sub-exponential tail regimes. We apply the framework to representative Bayesian models, including categorical likelihoods with arbitrary priors, Gaussian-Gaussian conjugate pair, generalized linear models from the exponential family and linear regression with Gaussian likelihood, both with Gaussian priors on the regression coefficients. Our results elucidate sample size, model structure, neighboring relations in DP (bounded or unbounded), and how prior if applicable, jointly determine the inherent DP guarantees of OPS. Furthermore, our results recover the existing DP guarantees for OPS as special cases – often with tighter privacy loss bounds – and substantially broaden the class of Bayesian models for which the inherent DP guarantees of OPS can be rigorously characterized.

Softplus affine\((1,1)\) models keep the conditional mean positive while allowing coefficients of either sign. We study them under Poisson, negative binomial, gamma and related conditional distributions. A stationary solution with moments of all positive orders exists whenever the coefficient sum is less than one, strictly enlarging previously known sufficient regions. At the unit root, the coefficient sum equals one. Under an invertibility condition, a stationary solution exists in this case if and only if the intercept is negative. At a stationary unit root, the conditional mean has a power-law tail with an explicit index, and the observations can exhibit long memory. For a fixed softplus scale, the Poisson quasi-maximum likelihood estimator is consistent and asymptotically normal at the usual \(\sqrt n\) rate below and at the unit root, with the coefficient sum left unrestricted in estimation. At the unit root, these results require only a finite mean and yield standard Wald inference even when the observations have infinite variance or long memory. Simulations and two applications examine finite-sample inference and sensitivity to the softplus scale.

In this paper we study the behavior of transformer networks generating text over long periods of time. We show that even a three-token, single-layer attention-only transformer can fail to converge. We then study a rich class of transformers with arbitrary finite depth, layer normalization, MLPs, RoPE encodings, and an arbitrary mix of local and global layers. Every transformer in this class has a softmax temperature threshold above which its empirical cache distribution converges to a unique limit, independent of initial conditions. We show further that, while each transformer has a temperature threshold, the supremum of this set is infinite, so there is no global minimal temperature that forces all transformers to converge.

We study inference for continuous linear functionals of solutions to ill-posed inverse problems, including Euclidean parameters arising in nonparametric instrumental variables, proximal causal inference, and shadow-variable models. Although pointwise root-$n$ inference is possible under standard conditions, we show that this need not translate into locally reliable inference. Under a model imposing only minimal regularity conditions, the target is Hellinger discontinuous when both inverse problems are ill-posed: arbitrarily small changes in the data-generating law can produce arbitrarily large changes in the target, ruling out informative locally uniform inference. Imposing a uniform bound on a minimum-norm solution restores Hellinger continuity under regularity conditions, but does not by itself ensure stable root-n inference. If, in addition, both inverse problems are ill-posed and the solution whose norm is uniformly bounded is nonunique at a law $P_0$, then under regularity conditions the target fails to be pathwise differentiable; consequently, asymptotically linear estimators exhibit nonvanishing root-$n$ bias under some local alternatives, causing Wald confidence intervals to undercover asymptotically. Under regularity conditions, pathwise differentiability is restored if the uniformly bounded solution is unique at $P_0$. It is also restored, regardless of uniqueness, under the stronger requirement that the corresponding solution satisfy a uniformly bounded $\beta$-source condition for some $\beta$>0. We relate these regimes to very weak, weak, and strong identification in instrumental variables models.

We consider nonparametric regression with random design and non-Gaussian noise, where the regression function has the multi-index form $F(x)=f(\Theta x)$ with an unknown link $f:\RR^m\to\RR$ in an anisotropic Hölder class and an unknown matrix $\Theta$ with $m\le d$ orthonormal rows. For estimation of $F$ at a fixed point we prove that the minimax risk is bounded from below by a multiple of $L^{1/(2\gamma+1)}(n^{-1}\ln n)^{\gamma/(2\gamma+1)}$, where $\gamma$ is the harmonic mean of the smoothness indices of the link divided by $m$. For $m<d$ this extends a Gaussian white noise result to regression. The main new case is $m=d$, where no dimension reduction takes place: the logarithmic factor still appears as soon as the link is anisotropic and the orientation of its axes is unknown. In contrast, when the link is isotropic with smoothness index at most one, the class is rotation invariant up to a constant and the classical rate without logarithm is attained. The logarithmic loss is therefore caused by the unknown orientation of the directions in which the function is least smooth

Selecting a budget-constrained set of documentary fragments, for instance the knowledge base of a retrieval-augmented generation system, raises a statistical question: how much expected query coverage is lost relative to the best feasible selection when the query distribution is unknown? We construct finite-sample, distribution-free upper confidence bounds for this coverage regret that remain valid after data-dependent selection. A pilot stage may select the document set, screen competitors, build comparison networks and learn an annotation policy. Then an independent audit acquires relevance labels with probabilities that are constant within observable strata. The main result is a simultaneous bound for all retained coverage contrasts. The bound is made observable by empirical-Bernstein variance calibration or by a one-stage construction without calibration, and it accounts for estimated stratum proportions. For maximum-type coverage we prove variance contractions for swaps of redundant fragments, give a condition under which multiscale radii improve on a direct union bound, and show that such redundancy cannot be certified from samples alone. A simulation study confirms validity and quantifies the gains from variance calibration and certified redundancy.

We ask when the law of an infinite exchangeable sequence is uniquely determined by a finite-dimensional marginal. For mixtures of a probability kernel, a weak localization condition makes three coordinates determine the directing law. Completeness alone makes four coordinates sufficient. Neither order can be lowered in general. For exponential kernels with a nonnegative sufficient statistic, a Müntz-Szász support condition implies three coordinates suffice. In the scalar case, this result is sharp. We also characterize when completeness holds simultaneously for all scalar mixing laws by a weighted Müntz density condition and show that statistic support alone need not decide completeness. These methods yield sharp finite-dimensional refinements of Freedman's theorem on rotatable sequences, namely that spherical symmetry of three coordinates makes an infinite exchangeable sequence rotatable while symmetry of the first two coordinates is in general insufficient. Finally, for energy-radial densities we give a complete-monotonicity criterion for infinite extendibility and determine the exact uniqueness thresholds for $\ell_q$-norm exponential laws.

We present an optimal tester for support size that determines whether an unknown probability distribution is supported on $n$ elements or $\varepsilon$-far away (in total variation distance) with sample complexity $\Theta(\frac{n}{\varepsilon\log(n)})$, improving the prior best $O(\frac{n}{\varepsilon\log(n)} \log(\min\{\frac{1}{\varepsilon}, n\}))$ due to Ferreira Pinto Jr. and Harms (STOC 2025). As a generalization, we provide a pretty good lower bound on support size.

We consider the problem of estimating the expectation of a function (observable) with respect to a high-dimensional probability distribution, under the assumption that a low-dimensional collective variable map captures all the directions of slow mixing for the distribution. More precisely, under the assumption that the fibers of the collective variable map satisfy a log-Sobolev inequality, we establish non-asymptotic relative-error guarantees for this problem using the eigenvector method for umbrella sampling (EMUS). Our analysis applies to both log-Lipschitz observables and indicator of convex sets, and gives complexity bounds that are polynomial in regularity parameters and inverse relative accuracy, and exponential in the dimension of the collective variable space. As a relative error guarantee, our bounds apply to rare event estimation. We consider Gaussian biasing windows defined on a low-dimensional collective variable space and analyze a modified sparse matrix version of EMUS. The bounds quantify the effect of the map through mixing of the conditional distribution and fluctuation of the observable along the fibers.

We construct confidence sets for the entire Fréchet mean set by inverting studentized pairwise risk inequalities. Gaussian multiplier calibration localized around empirical near minimizers gives asymptotic coverage at least the nominal level, uniformly over distribution classes satisfying empirical-process conditions. The construction requires neither unique ness nor a manifold structure. We characterize the limiting binding directions and give conditions for exact asymptotic coverage. On the intrinsic circle, exactness holds for every fixed non-point-mass law and uniformly over a family in which two means merge at a quar tic minimum. With ideal multiplier quantiles, the construction reproduces the Euclidean Wald set under both nonsingular and singular covariance. Under regular local conditions, it is first-order equivalent to intrinsic Wald inference.

We study the role of design smoothness in the minimax estimation of conditional average treatment effects at a fixed covariate point. For design densities bounded above and away from zero, without further regularity, we derive upper and lower risk bounds with an arbitrarily small gap in the convergence-rate exponent. Under uniform design, we obtain matching upper and lower bounds at interior points in the regime of low nuisance smoothness. These results show that rough covariate designs can lead to slower rates than uniform designs, even under the same smoothness assumptions on the treatment effect and nuisance functions, highlighting the role of the covariate distribution in pointwise causal estimation.

We characterize several properties of higher-arity function classes, specifically finite $\VC_k$-dimension, slicewise finite $\VC$-dimension, and finite slicewise Littlestone dimension, in terms of generalizations of Rademacher complexity. These characterizations sets up simple new proofs of several results on finite $\VC_k$-dimension and slicewise finite $\VC$-dimension. These include the preservation of finite $\VC_k$-dimension under taking expectations and uniform law of large numbers statements for both properties. All of these results apply naturally to real-valued function classes as well as boolean-valued classes. We also prove a new adversarial uniform law of large numbers for classes with finite slicewise Littlestone dimension.

In the past two decades, the Bayesian approach to inverse problems in infinite dimensional spaces has seen significant developments in theory, Markov chain Monte Carlo (MCMC) methods, and in applications. The Bayesian approach in finite dimension produces posterior distributions of finite dimensional random variables. In infinite dimensions, the desired objects are random fields and their laws as measures on Hilbert or Banach spaces. There are many textbooks and tutorials available that explain Bayesian methods and MCMC in finite dimensions. To enter the field of Bayesian inverse problems in infinite dimensions, a newcomer has to engage with a quite diverse and extensive literature. The purpose of this article is to provide a basic path to the subject that is easily accessible to a reader with general knowledge in probability and mathematical analysis. We start with a swift introduction to the finite dimensional case. The infinite dimensional case is introduced for Gaussian priors on Hilbert spaces and Gaussian-type likelihood functions. In this situation, mere continuity of the forward map suffices to prove that the posterior measure is well-defined. We also give some new and simple stability results. Two application examples are elaborated, one finite dimensional, one infinite dimensional.

We propose a multilevel main effects matrix factor model (MMEFM) to extract meaningful row and column effects from groups of matrix data. As a generalization of a multilevel matrix factor model, we provide rigorous definitions of MMEFM and identifications, together with iterative algorithms for the estimators of all global and local components in the main effects and common component, including global/local main effects and core ranks, with rates of convergence spelt out. An important feature of the model is that all main effects can be nonstationary and still be consistently estimated at each time point. At the same time, our estimated global common factors are alignment independent, and hence are truly global, and can be presented and given interpretations more easily. We propose a Bootstrap procedure for detecting the existence and identifying the groups of data that are driven by the global common factors. The power of the test, together with the accuracy of various estimators, are demonstrated in extensive simulation settings. Estimation performance is also compared to other state-of-the-art methods for matrix time series by simulations and a set of U.S. equity data. All numerical works are done using the R package MMEFM, which is available at CRAN.

In many applications, covariates are naturally organized into groups, yet conventional group-selection methods often retain all variables within selected groups, failing to achieve within-group sparsity. To address this limitation, we propose an adaptive two-stage sparse-group selection procedure for high-dimensional quantile regression. The proposed method integrates convolution-smoothed quantile loss with $\ell_{2,0}$ and $\ell_0$ penalties to identify important groups and relevant variables within them, achieving sparsity at both levels. We establish theoretical guarantees for sure screening and exact support recovery at both stages, as well as consistent model-size selection using information criteria, with explicit bounds on failure probabilities. The proposed procedure is computationally efficient by restricting optimization to selected groups and active variables. Simulation studies and a breast cancer gene-expression analysis demonstrate its competitive performance in variable selection, estimation accuracy, and computational efficiency, while identifying biologically relevant genes and delivering accurate quantile predictions.

In the literature on multiple testing methods that control the familywise error rate, a fundamental concept is "coherence". A procedure is coherent if it is "logically consistent" (Hochberg and Tamhane, 1987) in the sense that if a hypothesis is rejected, then all hypotheses that imply it are also rejected. If a method is not coherent, then it can be made coherent without increasing the familywise error rate. Recently, the concept of coherence was generalized to a much broader class of multiple testing procedures, namely all methods that provide simultaneous confidence bounds for false discovery proportions. Such procedures were defined to be coherent if they cannot be further improved through a process called "interpolation". However, we show that such a coherent procedure can sometimes be further improved by exploiting remaining logical inconsistencies. This leads to what we call a "perfectly coherent" method. Thus, interpolation does not necessarily lead to perfect coherence. However, we provide conditions under which it does. Fortunately, many methods based on interpolation turn out to be perfectly coherent.

In this paper, we study the estimation of the signal subspace dimension $k$ in a $p$-variate signal-plus-noise complex Gaussian model, when the noise covariance matrix is assumed to be diagonal or Toeplitz. We consider the high-dimensional asymptotic regime in which both the data dimension $p$ and the available sample size $n$ diverge to infinity at the same rate, while we also allow $k$ to diverge at a slower speed. In this context, the classical estimates of the literature, consistent in the high-dimensional regime, often rely on the a priori knowledge of an upper bound $k_*$ of $k$, which can be misspecified in practice. Our proposed approach relies on building estimates of the tracial moments of the noise covariance matrix using the $p-k_*$ smallest eigenvalues of the sample covariance matrix (SCM), which are in turn used to build a threshold above which the $k$ largest eigenvalues are isolated. The resulting estimate of $k$ is proved to be consistent in the high-dimensional regime, even for the misspecified situations where $k_* < k$. Numerical experiments in different challenging scenarios illustrate the improvement of the proposed approach over the aforementioned methods.

Sampling Gaussian measures restricted to nonlinear manifolds in high dimensions requires a proposal kernel that respects both geometric constraints and dimensional scaling requirements. We introduce tangent preconditioned Crank-Nicolson (pCN), which combines a manifold projection with a cancellation of the ambient Gaussian contribution, leaving only a geometric residual. Considering a shifted sphere with $\sqrt d$ scaling of radius and centre as a tractable benchmark embedded manifold, we prove that uncorrected tangent pCN results in a proposal noise scale of order $d^{-1/6}$ (compared to $d^{-1/2}$ for a tangent random-walk). We establish asymptotic optimal tuning for the geodesic expected squared jump distance which is maximised at an acceptance probability of 0.574 for the family of drift-corrected tangent proposals containing tangent pCN and tangent MALA. We also prove a joint diffusion limit for the parallel and orthogonal coordinates to the sphere's centre vector of tangent pCN. The orthogonal diffusion speed matches the limiting geodesic expected squared jump distance and is maximised at acceptance probability 0.574, whereas the parallel speed is maximised at 0.712, demonstrating that optimal tuning depends on the observable. We then identify an intrinsic precision correction characterised by potential Hessian and Ricci curvature that yields a proposal kernel with dimension-independent fixed step size limit. We discuss a generalisation to arbitrary hypersurfaces and present nondegenerate scaling of tangent pCN applied to conditioned diffusions under resolution refinement.

We study the geometry and resulting properties of Fréchet functions and barycenters in metric spaces with curvature bounded above by a positive constant in the sense of Alexandrov. Existing results typically rely on the assumption that the support is small enough that the squared distance, and thus the Fréchet function, is strongly geodesically convex. We work in the sharper regime where existence and uniqueness of the barycenter are guaranteed, but where convexity typically fails. We develop new geometric properties of these spaces, including a quantitative four-point inequality and a geodesic stability result. Using these tools, we establish new finite sample error bounds for two statistical estimators: the empirical barycenter and the iterated barycenter, as well as new results on sparse representations and approximations of barycenters.

Watermarking embeds an imperceptible signal in text from large language models that is statistically detectable with a secret key. The recent textGrain watermark (Li et. al. (2026b)) calibrates this signal against the model's sampling entropy through an explicit entropy budget, combining the stochasticity of Red-Green and the unbiasedness of Gumbel-max watermarks, and its accompanying report proposes a detector that sums pivotal statistics. We study when, and how efficiently, textGrain can be detected after human edits, modeled as a sparse mixture in which a fraction $\varepsilon_n\asymp n^{-p}$ of positions retain the watermark and token uncertainty decays as $n^{-q}$. First, entropy calibration comes at a cost: the minimax detection boundary of textGrain is $p+q=1/2$, so its detectable region is strictly smaller than that of Gumbel-max and maximal-coupling watermarks, whose boundary is $2p+q=1$. This cost stems from textGrain's construction rather than from entropy calibration itself: an unbiased scheme with the same entropy loss remains detectable where textGrain is not. Second, perhaps surprisingly, this loss makes the sum detector optimal at the boundary when the number of transport columns is fixed, unlike for Gumbel-max. When the number of columns grows and the entropy budget does not vanish, the detectable region expands while that of the sum detector does not. Third, for a fixed edit fraction, we characterize the optimal exponential rate of type II error decay and show that the sum detector falls strictly short of it at every fixed alternative. Finally, we propose a weighted multiscale binomial scan, valid in finite samples and needing no knowledge of the next-token distributions, entropy budget, or edit fraction; it attains every detection boundary above and the optimal exponent.

Many systems in engineering applications may experience complex multi-dimensional degradation processes involving random jumps. Analyzing the reliability of such degradation processes presents significant challenges. Previously developed models have utilized multi-dimensional Lévy subordinators and Lévy copulas to describe these processes, with reliability functions expressed through characteristic functions and Laplace transforms (Shu et al., 2015, 2016, Shi et al., 2020, Shi and Feng, 2024). However, these models are restricted to stationary linear paths and homogeneous Lévy subordinators. In real-world scenarios, increments in a degradation process are often time-dependent and nonhomogeneous. To address this limitation of previous studies, this work extends existing degradation models by introducing nonhomogeneous Lévy subordinators, which are designed to capture time-varying complexities inherent in degradation processes. Unlike conventional degradation models that often assume a constant rate of jumps, nonhomogeneous Lévy subordinators aim to provide a more flexible representation of degradation behaviors that evolve under varying external conditions and internal system dynamics. Additionally, the joint behaviors of multiple degradation processes are characterized using Lévy copulas, enabling a more comprehensive analysis of system reliability. Numerical examples are provided to illustrate the applicability of these models in lifetime analysis.

Bridge sampling estimates of marginal likelihoods are routinely reported with a delta-method Monte Carlo standard error (MCSE). Existing theory is asymptotic and cannot say whether a given run is long enough to trust the estimate or its MCSE. We give finite-sample answers. The optimal bridge estimator fits a logistic regression that separates posterior from proposal draws, and its accuracy is governed by their overlap, within a factor of two of the classification's Bayes error. Because the estimating function is a sum of bounded terms, no moment conditions are needed, and all bounds hold uniformly over targets, proposals and dimensions. A number of draws inversely proportional to the overlap suffices and is necessary in the worst case, paralleling the exponential-in-Kullback–Leibler requirement of importance sampling. A Berry–Esseen bound with an explicit constant holds at the optimal rate, and the MCSE is accurate once the effective overlap is large; it is estimated at no cost. For Markov chain draws, the required sample size grows in proportion to the sampler's relaxation time, unavoidably in the worst case. The estimated bridge, the fixed-point iteration and the fitted proposal matter only at second order, and we give free diagnostics for sampler failure and misestimated autocorrelation.

For the counts $(X_1,\dots,X_n)\sim\mathrm{Mult}(N;1/n,\dots,1/n)$ of $N$ observations in $n$ ordered bins we study the variance of the discrete total variation $d=\sum_{i=1}^{n-1}|X_{i+1}-X_i|$. An exact decomposition splits $\mathrm{Var}(d)$ into per-edge, adjacent-pair and disjoint-pair terms and yields, for fixed $n$ as $N\to\infty$, $\mathrm{Var}(d)/N=\alpha(n)-\beta(n)/N+\gamma(n)/N^2+o(N^{-2})$. For every $n\ge 2$, $N\alpha(n)$ is the classical normal-theory variance (Hoel, Kamat) of the sum of absolute successive differences of $n$ independent $N(0,N/n)$ variables; the correction $\beta(n)$ is obtained in closed form for $n\ge 3$ (the case $n=2$ is treated separately). Two disjoint edge differences $D_1=X_2-X_1$ and $D_3=X_4-X_3$ are exactly uncorrelated and asymptotically independent, yet $\mathrm{Cov}(|D_1|,|D_3|)=-1/(\pi n)+O(1/N)$: a correlation of relative order $1/N$, invisible in the Gaussian limit and absent under Poissonized sampling, whose constant is proportional to the exact joint fourth cumulant $\kappa_{2,2}(D_1,D_3)=-4N/n^2$ and is proved here on the lattice. The adjacent pair contributes the lattice-Edgeworth constant $a_*=-2\sqrt{3}/(9\pi)$, of which one quarter is a smooth Edgeworth term and three quarters a lattice-corner term; the disjoint skeleton $h(m)$ is closed through a Bessel generating function; and $\gamma(n)$ ($n\ge 5$) is expressed through a single triangular-lattice constant whose closed form we conjecture. Together with the mean of Part 1, the variance supplies the two moments used by moment-matched approximations to the null law of $d$. This is Part 2 of a series; the cyclic statistic $d^{c}$ will be treated in Part 5.

We study whether an unnormalized graph Laplace operator can determine both the geometry of a manifold and its sampling density. On a fixed compact smooth manifold without boundary, we use a known positive, continuous, strictly decreasing radial kernel at a fixed positive bandwidth. For positive sampling densities, the operator built from geodesic distances determines the Riemannian metric and identifies the density up to equality almost everywhere. When the metric is induced by a Euclidean embedding, the operator built from ambient distances also determines the embedding up to rotations, reflections and translations, and the density. On compact parameter classes with positive continuous densities, minimizing the distance between empirical and popula- tion operators gives estimates that converge to the true parameters under independent, identically distributed sampling. This recovery result assumes that the empirical operator can be evaluated at arbitrary points of the manifold.