SCIENTIFIC REPORTS
The Mathematics of Data
02 Jan 2024–26 Jan 2024
Organizing Committee
Afonso S. BandeiraETH Zurich

Subhroshekhar Ghosh
National University of Singapore

Philippe Rigollet
Massachusetts Institute of Technology

Hong T.M. Chu
National University of Singapore

CONTENTS PAGE
Page

Pierre Alquier

Robust estimation and regression with MMD

4Formulation and resolution of inverse problems in signal
and image processing - From classical methods to hybridAI

7

Optimal transport and high-dimensional probability

9

Approximation Bounds for Sparse Programs

11

Gradient flows for empirical Bayes in high-dimensional
linear models

12

A mathematical perspective on Transformers

15(Skew) Gaussian surrogates for high-dimensional
posteriors: from tighter bounds to tighter approximations

18

Information-Theoretic Thresholds for Planted Dense
Cycles21

The geometry of the deep linear network

24

Deep Learning based algorithm for nonlinear PDEs infinance and gradient descent type algorithm fornonconvex stochastic optimization problems with ReLU
neural networks

27

Jonathan Niles-Weed

Optimal transport map estimation in general function
spaces38

Mark Rudelson

University of Michigan, USA

How to check when a system of real quadratic equations
has a solution

39

Mark Rudelson

Approximately Hadamard matrices and random frames

41

ESSEC Business School, Singapore

Caroline Chaux
CNRS@CREATE, Singapore

Sinho Chewi
Massachusetts Institute of Technology, USA

Yair Shenfeld
Brown University, USA

Alexandre d’Aspremont
École Normale Supérieure, France

Zhou Fan
Yale University, USA

Borjan Geshkovski
Massachusetts Institute of Technology, USA

Anya Katsevich
Massachusetts Institute of Technology, USA

Cheng Mao
Georgia Institute of Technology, USA

Govind Menon
Brown University, USA

Ariel Neufeld
Nanyang Technological University, Singapore

New York University, USA

University of Michigan, USA

Contents Page (Continued)

Jonathan Scarlett

Recent Developments in Group Testing: Fundamental
Limits and Algorithms

44

Piyush Srivastava

Sampling from convex bodies using multiscaledecompositions

46

Random perturbation of low-rank matrices

47

National University of Singapore, Singapore

Tata Institute of Fundamental Research, India

Ke Wang

The Hong Kong University of Science andTechnology, ChinaROBUST ESTIMATION AND REGRESSION WITH MMD
PIERRE ALQUIERClassification AMS 2020: 62F10, 62F35, 62J02, 68T05, 46E22.
Keywords: universal estimation; robust statistics; kernel methods; minimum distanceestimation.
Popular estimation methods in statistics, such as the maximum likelihood estimator(MLE), or the method of moments, require strong assumptions on the statistical modeland the data-generating process to converge. When these conditions are not met, theseestimators can become very unstable. We are interested in universal estimators, thatwould converge without assumptions on the model. For example, [15] proved that socalled minimum distance estimators converge under far more general assumptions thanthe MLE. More recently, the ρ-estimators defined in [4] are not only convergent, but
minimax-optimal, in a very large class of models. However, there are still a few limitating
assumptions in [15, 4], and the computation of these estimators might be difficult inpractice.
This talk will summarize a recent line of work on a variant of minimum distanceestimators based on the so-called maximum mean discrepancy (MMD). In the case ofi.i.d. data, these estimators converge without any assumption on the model nor on thedata generating process. Moreover, relatively efficient algorithms are available tocompute these estimators.
Let X1 , . . . , Xn be X -valued random variables i.i.d. from a probability distribution
P . Let M = (Pθ , θ ∈ Θ) be a statistial model. Note that we don’t assume P 0 ∈ M.
Let H be a reproducing kernel Hilbert space (RKHS) equipped with a scalar product⟨·, ·⟩, its associated norm ∥ · ∥ and a kernel k: there is a map φ : X → H with k(x, y) =⟨φ(x), φ(y)⟩. The kernel mean embedding (KME) is defined for any probability P on X byµ(P ) = EX∼P [φ(X)]. We refer the reader to [11] for more details on this construction.
In particular: if k is bounded, then the KME is indeed well-defined for any P ; if k ischaracteristic (see [11] for a definition), µ is one-to-one. Thus, if k is both bounded
and characteristic, Dk (P, P ′ ) = ∥µ(P ) − µ(P ′ )∥ defines a metric on probabilities over X .
Finally, [11] provides examples of kernels that are indeed bounded and characteristic:
for example, when X = Rd , the Gaussian kernel k(x, y) = exp(−∥x − y∥2 /γ).
0Definition 0.1. We define the MMD-minimum distance estimator, or MMD-MDE, asθ̂ = arg min Dk (Pθ , P̂n )
θ∈Θ

Pn
1

where P̂n = n

i=1 δXi is the empirical distribution of the sample X1 , . . . , Xn .

In the talk, I will prove the following result:

1

2

Theorem 0.2 (Theorem 3.1 in [6]). Assume k ≤ 1, then2
E Dk (Pθ̂ , P 0 ) ≤ inf Dk (Pθ , P 0 ) + √ .
θ∈ΘnThe fact that there are no assumptions on P 0 nor on M = (Pθ , θ ∈ Θ) in the theorem
is not due to lack of space: this is actually a nice feature of θ̂! Variants of Theorem 0.2,
including a bound that holds with large probability rather than in expectation, can befound in [5, 6].
Various topics on θ̂ will be covered in this talk:• the MMD distance can be written in terms of expectations, and thus, theminimization in Definition 0.1 can be done with a stochastic gradient method.
Details are discussed in [5, 7].
• Theorem 1 leads to convergence: if the model is well-specified, P 0 ∈ M, then√E[Dk (Pθ̂ , P 0 )] ≤ 2/ n → 0 when n → ∞. It also leads to robustness as definedby Huber: if P 0 = (1 − ε)Pθ0 + εQ for an arbitrary contamination Q and a small√ε, then E[Dk (Pθ̂ , P 0 )] ≤ 4ε + 2/ n. This is proven in [7], where we also prove
more difficult robustness results under adversarial contamination of the data.
• the asymptotic normality of θ̂ is studied in [5] (this requires assumptions).
• an extention of Theorem 0.2 to non-i.i.d. observations (time series) is provided
in [7], under a new mixing condition. We also prove that this mixing conditionis less restrictive than the standard β-mixing condition.
√• the term 2/ n in Theorem 0.2 cannot be improved in general. However, underassumptions on the variance of P 0 , it can actually be improved: variance-aware
versions of Theorem 0.2 are proven in [14].
• we can define Bayesian-flavored estimators by using the MMD to define a pseudolikelihood. We studied such estimators in [7] and proved their convergence using
tools from the PAC-Bayes theory [1]. Other Bayesian-inspired variants of θ̂ arebased on sampling [8] and Approximate Bayesian Computation or ABC [10].
• this estimation strategy was successfully implemented in a wide range ofapplications:generative artificial intelligence [9], quantisation andclustering [13], estimation of copulas [2], estimation of parameters in stochastic
PDEs [5] etc.
A large part of the talk will be dedicated to parametric regression (linear, or not). In
this case, the observations are pairs input-output: Xi = (Zi , Yi ). Theorem 0.2 guarantees
that we can estimate the joint distribution of these pairs. But this is not relevant inregression, where the objective is rather the estimation of the conditional distribution ofYi given Zi .
The problem is that, while the theory of KME looks simple and elegant, a rigorous
definition of conditional KME turns out to be far more difficult and cumbersome! Werefer the reader to [12] for a recent account and a general approach to solve the problem.
In our recent paper [3], we proved that, in standard regression models: linear
regression with Gaussian noise, logistic regression, Poisson regression, etc., theconditions for the existence of conditional KMEs are met. This can be used to define
consistent and robust estimation of the regression parameters in the spirit of

3

Definition 0.1 above. These estimators turn out to perform extremely when compared
to existing robust regression procedures.
R EFERENCES[1] Pierre Alquier. User-friendly Introduction to PAC-Bayes bounds. Foundations and Trends® in Machine
Learning, 17(2), 174–203, 2024.
[2] Pierre Alquier, Badr-Eddine Chérief-Abdellatif, Alexis Derumigny and Jean-David Fermanian.
Estimation of copulas via maximum mean discrepancy. Journal of the American Statistical Association,
118(543), 1997-2012, 2023.
[3] Pierre Alquier and Mathieu Gerber. Universal robust regression via maximum mean discrepancy.
Biometrika, 111(1), 71–92, 2024.
[4] Yannick Baraud, Lucien Birgé and Mathieu Sart. A new method for estimation and model selection:
ρ-estimation. Inventiones mathematicae, 207(2), 425–217, 2017.
[5] Francois-Xavier Briol, Alessandro Barp, Andrew B. Duncan, and Mark Girolami. Statistical inference
for generative models with maximum mean discrepancy, arXiv preprint arXiv:1906.05944, 2019.
[6] Badr-Eddine Chérief-Abdellatif and Pierre Alquier. MMD-Bayes: Robust Bayesian estimation via
maximum mean discrepancy. Symposium on Advances in Approximate Bayesian Inference, Proceedingsin Machine Learning Research 118, 1–21, 2020.
[7] Badr-Eddine Chérief-Abdellatif and Pierre Alquier. Finite sample properties of parametric MMD
estimation: robustness to misspecification and dependence. Bernoulli, 28(1), 181–213, 2022.
[8] Charita Dellaporta, Jeremias Knoblauch, Theodoros Damoulas and Francois-Xavier Briol. Robust
Bayesian inference for simulator-based models via the MMD posterior bootstrap. International
Conference on Artificial Intelligence and Statistics, Proceedings in Machine Learning Research 151,
943–970, 2022.
[9] Gintare Karolina Dzuigaite, Daniel Roy and Zoubin Ghahramani. Training generative neural networks
via maximum mean discrepancy optimization, arXiv preprint arXiv:1505.03906, 2015.
[10] Sirio Legramanti, Daniele Durante, and Pierre Alquier. Concentration and robustness of discrepancybased ABC via Rademacher complexity, arXiv preprint arXiv:2206.06991, 2022.
[11] Krikamol Muandet, Kenji Fukumizu, Bharath Sriperumbudur and Bernhard Schölkopf. Kernel Mean
Embedding of Distributions: A Review and Beyond. Foundations and Trends® in Machine Learning,10(1–2), 1–141, 2017.
[12] Junhyung Park and Krikamol Muandet A measure-theoretic approach to kernel conditional meanembeddings. Advances in neural information processing systems 33, 21247–21259, 2020.
[13] Onur Teymur, Jackson Gorham, Marina Riabiz, and Chris Oates. Optimal quantisation of probability
measures using maximum mean discrepancy. International Conference on Artificial Intelligence andStatistics, Proceedings in Machine Learning Research 130, 1027–1035, 2021.
[14] Geoffrey Wolfer and Pierre Alquier. Variance-aware estimation of kernel mean embedding, arXiv
preprint arXiv:2210.06672, 2022.
[15] Yannis Yatracos. Rates of convergence of minimum distance estimators and Kolmogorov’s entropy.
The Annals of Statistics, 13(2), 768–774, 1985.
ESSEC B USINESS S CHOOL , A SIA -PACIFIC CAMPUS , 5 N EPAL PARK , 139408 SINGAPORE
Email address: alquier@essec.edu

FORMULATION AND RESOLUTION OF INVERSE PROBLEMS IN SIGNAL AND
IMAGE PROCESSING - FROM CLASSICAL METHODS TO HYBRID AI
CAROLINE CHAUXClassification AMS 2020: 94A08, 68T07
Keywords: Unrolled algorithms, parameter estimation, maximum a posteriori, wavelets,
sparsity, low-rank, deconvolution, robust PCA.
We are interested here in inverse problems arising in signal and image processing.
Solving such problems firstly consists in formalising the direct problem byunderstanding the physics behind and secondly, solving the associated inverse problem
which results in solving an optimization problem. Two main ingredients are involved in
the optimization functional to be minimized: a data fidelity term accounting for themodel and a regularisation term (possibly several ones) accounting for priorinformation on the targeted solution. Regularisation parameters come into the play toguarantee the best trade-off between the two quantities.
Classical optimization-based approaches consist in, once the optimization problem
has been formulated, proposing iterative procedures (e.g. proximal algorithms [1])converging to a solution of the considered inverse problem. In such a case, most of thetime, the parameters (either regularisation or algorithm’s hyperparamaters) are fixed
empirically. More recently, unrolled or unfolded neural networks have been proposed[2]. They combine optimization and learning, constitute interpretable networks, andintegrate information about the direct model.
Inverse problems are encountered in many scientific domains such as biology, medical
imaging, chemistry, audio signal processing for which, different tasks must be tackled
such as deconvolution, restoration, unmixing, missing data reconstruction, etc.
We presented in this talk a physics-informed unrolled network [3] to automaticallychoose regularisation parameters in image deconvolution when the regularisation is
performed in the wavelet domain. The resulting optimisation problem is solved byusing an unrolled version of FISTA algorithm [4]. Indeed, hyper-parameter tuning, andespecially regularisation parameter estimation, is a challenging but essential task whensolving inverse problems. We proposed to perform their estimation under an unrolled
strategy together with the inverse problem solving. The resulting network is trained
while incorporating information on the model through Maximum a Posteriori
estimation which drastically decreases the amount of data needed for the training andresults in better estimation results.
We also presented a physics-informed unrolled network [5] to solve the problem ofRobust Principal Component Analysis (RPCA) which consists in decomposing a matrix
into the sum of a low rank matrix and a sparse matrix. We proposed a deep unrolled
algorithm based on an accelerated alternating projection algorithm [6] which aims tosolve RPCA in its nonconvex form and where hyperparameters are automatically learnt.
We demonstrated the unrolled algorithm’s effectiveness on synthetic datasets and also1

2

on a face modelling problem, where it leads to both better numerical and visualperformances.
These works have been done in collaboration with Vincent Tan, Emmanuel Soubiès,
Pascal Nguyen and Elisabeth Tan.
R EFERENCES[1] P. Combettes and J.-C. Pesquet, Proximal Splitting Methods in Signal Processing. Springer
Optimization and Its Applications , Springer New York, pp. 185-212, 2011.
[2] V. Monga, Y. Li, and Y. Eldar. Algorithm unrolling: Interpretable, efficient deep learning for signal
and image processing. IEEE Signal Process. Mag., vol. 38, no. 2, pp. 18–44, 2021.
[3] P. Nguyen, E. Soubies, C. Chaux. MAP-informed Unrolled Algorithms for Hyper-parameter
Estimation. Proc. Int. Conf. Image Process., Kuala Lumpur, Malaysia, Oct. 8-11, 2023.
[4] A. Beck and M. Teboulle. A fast iterative shrinkage thresholding algorithm for linear inverse problems.
SIAM J. Imaging Sci., vol. 2, no. 1, pp. 183–202, 2009.
[5] E. Z. C. Tan, C. Chaux, E. Soubies, and V. Y. F. Tan. Deep Unrolling for Nonconvex Robust Principal
Component Analysis. IEEE Int. Workshop Mach. Learn. Signal Process., Rome, Italy, Sep. 17-20, 2023.
[6] H. Cai, J.-F. Cai, and K. Wei. Accelerated alternating projections for robust principal component
analysis. J. Mach. Learn. Res., vol. 20, no. 1, pp. 685–717, Jan 2019.
IPAL AND CNRS IRL 2955, CNRS@CREATE, 1 C REATE WAY, #08-01 C REATE T OWER , S INGAPORE
138602
Email address: caroline.chaux@cnrs.fr

MINI-COURSE ON OPTIMAL TRANSPORT AND HIGH-DIMENSIONALPROBABILITYSINHO CHEWI AND YAIR SHENFELD

The mini-course focused on the foundations and applications of optimal transport and highdimensional probability. The course was composed of 8 lectures whose outline is given below.
Lecture 1. The notion of transportation of measures was introduced, followed by the Monge problem of transporting between two measures in an optimal way. Since the Monge problem is difficult
to solve, a better behaved relaxation, known as the Kantorovich problem, was introduced, which
involves the notion of couplings as substitutes for transport maps.
Lecture 2. The Kantorovich problem is a linear programming problem where duality plays animportant. This duality was proved, yielding characterizations of the primal and dual solutions,as well as Brenier’s theorem on the original problem of Monge. The important concept of cyclicalmonotonicity was introduced and used in the proof.
Lecture 3. The Kantorovich problem leads to the notion of the Wasserstein distance over the spaceof probability measures. It was shown that the Wasserstein distance is indeed a metric, which iscompatible with weak convergence. The formula for the linearization of the Wasserstein distance
was also derived.
Lecture 4. The Wasserstein space, the space of probability measures endowed with the Wasserstein
distance, can viewed as a formal infinite-dimensional Riemannian manifold. This perspective was
explained by introducing the dynamic formulation of optimal transport, known as the Benamou–
Brenier principle. Analogies to fluid dynamics were drawn. This led to the development of the Otto
calculus over the Wasserstein space.
Lecture 5. The Otto calculus over the Wasserstein space allows for the development of a theoryof gradient flows of functionals of probability measures. Formulas for gradients and gradient flows
were derived, emphasizing the roles of important functionals such as entropy, as well as the notionof displacement convexity.
Lecture 6. Functional inequalities capture the convergence to equilibrium of gradient flows, andsome important examples of these inequalities were introduced: the log-Sobolev inequality, Talagrand’s inequality, and the HWI inequality. It was shown how the notion of displacement convexity
can be used to prove functional inequalities, and how in turn, the property of log-concavity of aprobability measure is closely connected with the displacement convexity of the relative entropy.
Lecture 7. Functional inequalities can also be used to establish the concentration of measurephenomenon. The following notions were presented: sub-Gaussian variables, the entropy method
attributed to Herbst, the Bobkov–Götze theorem, tensorization, and transport inequalities.
1Lecture 8. The perspective of viewing flows of probability measures as gradient flows in Wassersteinspace ultimately leads to a powerful set of tools and insights from optimization which can be applied
to the analysis of sampling algorithms. This idea was developed in the context of providing a nonasymptotic convergence analysis for the Euler–Maruyama discretization of the Langevin diffusion.
School of Mathematics, Institute for Advanced Study, Princeton, NJ, USA
E-mail address: schewi@ias.edu
Division of Applied Mathematics, Brown University, Providence, RI, USA
E-mail address: Yair_Shenfeld@Brown.edu

2PARSITY, FEATURE SELECTION &THE SHAPLEY-FOLKMAN THEOREM.
ALEX D’ASPREMONTClassification AMS 2020: 65K05
Keywords: Duality, Shapley-Folkman
Due to its linear complexity, naive Bayes classification remains an attractive
supervised learning method, especially in very large-scale settings. We propose a sparseversion of naive Bayes, which can be used for feature selection. This leads to a
combinatorial maximum-likelihood problem, for which we provide an exact solution in
the case of binary data, or a bound in the multinomial case. We prove that our convex
relaxation bounds becomes tight as the marginal contribution of additional features
decreases, using a priori duality gap bounds dervied from the Shapley-Folkman
theorem. We show how to produce primal solutions satisfying these bounds. Bothbinary and multinomial sparse models are solvable in time almost linear in problem
size, representing a very small extra relative cost compared to the classical naive Bayes.
Numerical experiments on text data show that the naive Bayes feature selection method
is as statistically effective as state-of-the-art feature selection methods such as recursive
feature elimination, l1 -penalized logistic regression and LASSO, while being orders ofmagnitude faster 1
CNRS - ENS PARIS .
Email address: aspremon@ens.fr

1

A python implementation can be found at https://github.com/aspremon/NaiveFeatureSelection1

Gradient flows for empirical Bayes in high-dimensional linear
modelsZhou Fan∗, Leying Guan†, Yandi Shen∗, Yihong Wu∗

Introduced by Robbins [Rob51, Rob56] in the 1950s, empirical Bayes (EB) is a powerful framework for large-scale inference that learns and adapts to latent structure in data. In this work, we
consider empirical Bayes estimation in a linear model with Gaussian noise and i.i.d. prior,
iid

iid

εi ∼ N (0, σ 2 ).

θj ∼ g∗ ,y = Xθ + ε,

(0.1)

Assuming that the sample size n and dimension p are both large, we study estimation of the effectsize prior g∗ from the regression data (X, y). Many methods based upon variational inference orMonte Carlo EM have been proposed for this problem in the context of inferring genetic architectures of complex traits [ZQPC18, O’C21, ZZ21, SSAAP22, MCW+ 23], but the problem has only
recently been the subject of theoretical study [MSS23].
We propose and study a new gradient-flow procedure for nonparametric estimation of g∗ via thenonparametric maximum likelihood estimator (NPMLE). For a user-specified parameter τ 2 > 0,
reparametrizing the model by a smoothed regression vector φ = θ + N (0, τ 2 Id) and modified
residual εe ∼ N (0, Σ) where Σ = σ 2 Id −τ 2 XX⊤ , the negative log-likelihood F̄n (g) admits a Gibbs
variational representationF̄n (g) = min Fn (q, g),q

1
Fn (q, g) :=
p

Z 


p
X
1
⊤ −1
(y − Xφ) Σ (y − Xφ) −log[Nτ ∗ g](φj ) + log q(φ) q(φ)dφ + constant.
2
j=1

We propose to jointly optimize Fn (q, g) over posterior densities q on Rp and prior densities g on abounded support [−M, M ], via the coupled gradient flowsd2
qt (φ) = −p · gradW
q Fn (qt , gt )[φ]
dt
"

 !#
[Nτ ∗ gt ]′ (φj ) p= ∇ · qt (φ) X Σ (Xφ − y) −+ ∆qt (φ)[Nτ ∗ gt ](φj ) j=1



d
q̄t
FR
gt (θ) = −α · gradg Fn (qt , gt )[θ] = α gt (θ) Nτ ∗(θ) − 1dt
Nτ ∗ gt
⊤



−1

(0.2)
(0.3)

FR
2
where gradW
q denotes the Wasserstein-2 gradient in q, gradg denotes the Fisher-Rao gradient in g,
2Nτ ∗ gt is the convolution of the N (0, τ ) density with gt , and q̄t is the univariate averaged marginal∗†

Department of Statistics and Data Science, Yale University
Department of Biostatistics, Yale University
zhou.fan@yale.edu, leying.guan@yale.edu, yandi.shen@yale.edu, yihong.wu@yale.edu

1density of coordinates under qt . The q-flow (0.2) is simulated via a Langevin diffusion [JKO98] !
√[Nτ ∗ gt ]′ (φt,j ) p⊤ −1
dt + 2 dBt(0.4)
dφt = −X Σ (Xφt − y) +[Nτ ∗ gt ](φt,j ) j=1
with time-evolving drift. This is coupled with the g-flow (0.3) which is simulated by a fixed-griddiscretization of the domain [−M, M ] and an empirical estimate of q̄t from the coordinates of theLangevin sample φt . The reparametrization of the model by φ rather than θ ensures that theLangevin diffusion targets a continuous and sufficiently regular posterior density, even if the priorestimate gt approaches a discrete measure. This idea of bivariate optimization of a Gibbs variational
representation for maximum likelihood inference can be attributed to [NH98], and has also been
proposed and studied recently in [KLJ23, ACG+ 23] in different contexts.
Our results consist of (1) a statistical consistency guarantee for the NPMLE in this problem,
extending recent work of [MSS23], (2) a new log-Sobolev inequality for mixing of Langevin dynamics
in the linear model, deriving from an insight of [BB19], and (3) a local convergence guarantee forthe above gradient flow algorithm.
For consistency, we show the following result.
Theorem 0.1. Suppose g∗ is supported on a known interval [−M, M ]. As n, p → ∞, under adeterministic condition for the design matrix X, which in particular assumes the existence of test
vectors zj ∈ Rn for sufficiently many columns xj of X such thatz⊤j xj → 1,

2+ε
|z⊤
→ 0,
j xk |

any approximate NPMLE gb over [−M, M ] satisfies gb ⇒ g∗ almost surely in the sense of weak
convergence.
Proposition 0.2. Suppose n/p ≥ γ for any constant γ > 0, and X ∈ Rn×p has i.i.d. sub-Gaussianrows { √1n x(i) }ni=1 where x(i) has mean 0, covariance ΣX ∈ Rp×p satisfyingc ≤ λmin (ΣX ) ≤ λmax (ΣX ) ≤ C,
and bounded sub-Gaussian norm. Then the assumptions for X needed in the above theorem holdalmost surely as n, p → ∞.
For mixing of the Langevin diffusion, let ν[g] be the posterior density of φ given (X, y) under
the prior g for θ. We prove the following uniform log-Sobolev inequality, which implies exponential
contraction in KL-divergence for the Langevin diffusion with ν[g] as its stationary law.
Theorem 0.3. Let g be supported on [−M, M ] and fix a constant δ > 0. Suppose that σ 2 /∥X∥2op >M 2 + δ. Then for any τ 2 ∈ (0, σ 2 /∥X∥2op − δ), ν[g] satisfies a LSI with constant C = C(M, τ, δ) > 0
that is uniform over all priors g on [−M, M ].
Finally, for local convergence of the bivariate gradient flow, we show the following.
Theorem 0.4. Suppose that ν[g] satisfies a LSI with some constant C = C(M, τ ) > 0, uniformly
over priors g supported on [−M, M ]. Let {(qt , gt )}t≥0 be the solution to the bivariate gradient flow
(0.2–0.3).
If F̄n is convex over the sub-level set {g : F̄n (g) ≤ Fn (q0 , g0 )} and additional (mild) technicalconditions hold for the initialization (g0 , q0 ), then for any ε > 0, some constant C(ε) > 0, and sometime t ≤ C(ε)(n + p), it holds that F̄n (gt ) − F̄n (g∗ ) ≤ ε.
We leave for future work several interesting open questions, including a better theoretical characterization of the optimization landscape of F̄n (g) and conditions for global convergence, a theoreticalanalysis of time discretization and consistency of estimating q̄t within the empirical implementationof the g-flow equation (0.3), and extensions to more complex models and sampling algorithms.
2

References[ACG+ 23] Ö Deniz Akyildiz, Francesca Romana Crucinio, Mark Girolami, Tim Johnston, andSotirios Sabanis, Interacting particle Langevin algorithm for maximum marginal likelihood estimation, arXiv preprint arXiv:2303.13429 (2023).
[BB19]

Roland Bauerschmidt and Thierry Bodineau, A very simple proof of the LSI for hightemperature spin systems, J. Funct. Anal. 276 (2019), no. 8, 2582–2588.

[JKO98]Richard Jordan, David Kinderlehrer, and Felix Otto, The variational formulation of
the Fokker-Planck equation, SIAM J. Math. Anal. 29 (1998), no. 1, 1–17.

[KLJ23]

Juan Kuntz, Jen Ning Lim, and Adam M Johansen, Particle algorithms for maximum
likelihood training of latent variable models, International Conference on Artificial Intelligence and Statistics, PMLR, 2023, pp. 5134–5180.

[MCW+ 23] Fabio Morgante, Peter Carbonetto, Gao Wang, Yuxin Zou, Abhishek Sarkar, andMatthew Stephens, A flexible empirical Bayes approach to multivariate multiple regression, and its improved accuracy in predicting multi-tissue gene expression fromgenotypes, PLoS Genetics 19 (2023), no. 7, e1010539.
[MSS23]

Sumit Mukherjee, Bodhisattva Sen, and Subhabrata Sen, A mean field approachto empirical bayes estimation in high-dimensional linear regression, arXiv preprint
arXiv:2309.16843 (2023).

[NH98]

Radford M Neal and Geoffrey E Hinton, A view of the EM algorithm that justifies
incremental, sparse, and other variants, Learning in graphical models, Springer, 1998,
pp. 355–368.

[O’C21]

Luke J O’Connor, The distribution of common-variant effect sizes, Nature Genetics 53(2021), no. 8, 1243–1249.

[Rob51]

Herbert Robbins, Asymptotically subminimax solutions of compound statistical decision
problems, Proceedings of the Second Berkeley Symposium on Mathematical Statistics
and Probability, The Regents of the University of California, 1951.

[Rob56], An empirical Bayes approach to statistics, Proceedings of the Third Berkeley
Symposium on Mathematical Statistics and Probability, Volume 1: Contributions to
the Theory of Statistics, The Regents of the University of California, 1956.

[SSAAP22] Jeffrey P Spence, Nasa Sinnott-Armstrong, Themistocles L Assimes, and Jonathan K
Pritchard, A flexible modeling and inference framework for estimating variant effectsizes from GWAS summary statistics, BioRxiv (2022), 2022–04.
[ZQPC18]

Yan Zhang, Guanghao Qi, Ju-Hyun Park, and Nilanjan Chatterjee, Estimation ofcomplex effect-size distributions using summary-level statistics from genome-wide association studies across 32 complex traits, Nature Genetics 50 (2018), no. 9, 1318–1326.

[ZZ21]

Geyu Zhou and Hongyu Zhao, A fast and robust Bayesian nonparametric method forprediction of complex traits using summary statistics, PLoS genetics 17 (2021), no. 7,
e1009697.

3

ON THE EMERGENCE OF CLUSTERS IN SELF-ATTENTION DYNAMICSBORJAN GESHKOVSKIClassification AMS 2020: 34D05, 34D06, 35Q83
Keywords: Transformers, self-attention, interacting particle systems, clustering, gradient
flows.
1. I NTRODUCTION
The introduction of Transformers in 2017 marked a milestone in the development ofneural network architectures. Central to this is self-attention, a novel mechanism which
distinguishes Transformers from traditional architectures, and which plays a substantial
role in their superior practical performance. We present a mathematical framework foranalyzing Transformers based on their interpretation as interacting particle systems, inwhich time plays the role of layers. Our analysis reveals that clusters emerge in long time,
confirming previous empirical findings, and shedding light on the role of the attention
mechanism.
2. M AIN RESULTAs first done in [3, 8], we define an idealized model of the Transformer architecturethat consists in viewing the discrete layer indices as a continuous time variable, andwhich focuses exclusively on two key components of the Transformers architecture: selfattention and layer normalization1. This results in the dynamics!
nX1eβ⟨xi (t),xj (t)⟩ xj (t)(SA)
ẋi (t) = P⊥
xi (t)
Zβ,i (t) j=1
for i ∈ [n] and t ≥ 0, where(2.1)

Zβ,i (t) =nXeβ⟨xi (t),xk (t)⟩

k=1
⊤
d−1
and P⊥. We prove the following resultx = Id − xx is the orthogonal projector to Tx S

(see [1, 2]).
Theorem 2.1. Let d, n ≥ 2 and β ≥ 0, and suppose that either d ≥ n, or β ≳d n2 , orβ ≲ n−1 . Consider the unique solution (xi (·))i∈[n] ∈ C 0 (R≥0 ; (Sd−1 )n ) to the Cauchy problem
for (SA), corresponding to an initial sequence of points (xi (0))i∈[n] ∈ (Sd−1 )n distributed
uniformly at random. Then almost surely there exists x∗ ∈ Sd−1 such thatlim xi (t) = x∗t→+∞
1

Strictly speaking, we are using root-mean-square (RMS) norm. To the best of our knowledge, theoriginal layer normalization, which consists in standardizing every component of every token at every
layer has not yet been addressed in the continuous-time formulation of Transformers.

1

2

for all i ∈ [n].
t = 0.0

t = 15.0

F IGURE 1. The evolution of trajectories can be fully described in very large
dimension (d ≫ poly n), beyond solely long time asymptotics with rates,
as seen in the phase diagram above. See [2].
Note that when d = 2 and β = 0, we recover the celebrated Kuramoto model ofoscillators [7] (see [6] for mathematical details in the case d ≥ 2).
Analyzing the self-attention model with a view toward understanding fine-grained andempirically observed phenomena in Transformers has borne fruit to several interesting
results in the past year. Stability of the asymptotic results presented in [1] with respect
to perturbations of the parameters has been addressed in [5]. A mean-field version of
the masked attention mechanism, as well as a thorough study of the Lipschitz constantof attention is presented in [4]. Steering ensembles of empirical measures to ensemblesof empirical measures via parametrized self-attention dynamics has been done in [9].
R EFERENCES[1] B. Geshkovski, C. Letrouit, Y. Polyanskiy, and P. Rigollet, “The emergence of clusters in self-attention
dynamics,” Advances in Neural Information Processing Systems, vol. 36, 2024.
[2] B. Geshkovski, C. Letrouit, Y. Polyanskiy, and P. Rigollet, “A mathematical perspective on
Transformers,” arXiv preprint arXiv:2312.10794, 2023.
[3] M. E. Sander, P. Ablin, M. Blondel, and G. Peyré, “Sinkformers: Transformers with doubly stochastic
attention,” in International Conference on Artificial Intelligence and Statistics, 2022, pp. 3515–3530,
PMLR.
[4] V. Castin, P. Ablin, and G. Peyré, “Understanding the Regularity of Self-Attention with OptimalTransport,” arXiv preprint arXiv:2312.14820, 2023.
[5] H. Koubbi, M. Boussard, and L. Hernandez, “The Impact of LoRA on the Emergence of Clusters inTransformers,” arXiv preprint arXiv:2402.15415, 2024.
[6] A. Frouvelle and J.-G. Liu, “Long-time dynamics for a simple aggregation equation on the sphere,” inStochastic Dynamics Out of Equilibrium: Institut Henri Poincaré, Paris, France, 2017, Springer, 2019,
pp. 457–479.
[7] Y. Kuramoto, “Self-entrainment of a population of coupled non-linear oscillators,” in International
Symposium on Mathematical Problems in Theoretical Physics: January 23–29, 1975, Kyoto University,
Kyoto/Japan, Springer, 1975, pp. 420–422.
[8] Yiping Lu, Zhuohan Li, Di He, Zhiqing Sun, Bin Dong, Tao Qin, Liwei Wang, and Tie-Yan Liu.
Understanding and improving transformer from a multi-particle dynamic system point of view. arXiv
preprint arXiv:1906.02762, 2019.
[9] Agrachev, Andrei, and Cyril Letrouit. Generic controllability of equivariant systems and applications toparticle systems and neural networks. arXiv preprint arXiv:2404.08289 (2024).

3

I NRIA & L ABORATOIRE JACQUES -L OUIS L IONS , S ORBONNE U NIVERSIT É , 4 P LACE J USSIEU, 75005 PARIS
Email address: borjan@mit.edu

THE LAPLACE APPROXIMATION IN HIGH-DIMENSIONAL BAYESIAN INFERENCE
ANYA KATSEVICHClassification AMS 2020: 62E17, 62F15
Keywords: Laplace approximation, high dimensional Bayesian inference
Developing cheap and accurate computational techniques for Bayesian inference is an
important goal, as Bayesian inference tasks can be very computationally intensive.
These tasks include constructing posterior credible sets, computing the posterior meanand covariance, and computing the posterior predictive distribution of a new data point.
Computing all of
R these quantities involves either sampling from the posterior π, ortaking integrals gdπ against the posterior. When the dimensionality of the parameterof interest is large, these tasks can be very expensive.
The most common approach in Bayesian inference to both generate samples from,and integrate against π is Markov Chain Monte Carlo (MCMC) [3]. In principle, MCMC
schemes can be made arbitrarily accurate by tuning parameters, e.g. by extending thesimulation time of the Markov chain or by decreasing step sizes. However, MCMC iscomputationally intensive in high dimensions and it also has other disadvantages. For
example, it can be difficult to identify clear-cut stopping criteria for the algorithm [5].
Another popular approach is to find a simple distribution γ̂ which approximates π,and to use this distribution as a proxy for π to do all of one’s inference tasks. In theideal scenario, many integrals against γ̂ are computable in closed form, and it is cheap
to sample from γ̂. The idea of using an approximation γ̂ to π is at the heart ofapproximate Bayesian inference methods such as variational inference [2, 16],expectation propagation [12], and the Laplace approximation, a Gaussian
approximation first introduced in the Bayesian inference context by [14].
The Laplace approximation (LA) exploits large sample properties of the posterior
which have been established in the Bernstein von-Mises theorem (BvM). The BvM is afundamental result which in its classical form states that if the model is well-specified,then the posterior contracts around the ground truth parameter and becomes
asymptotically normal in the large sample limit [15, Section 10.2]. Despite thetheoretical and philosophical importance of the BvM (it has been used as a justificationof Bayesian procedures from the perspective of frequentist inference), this result does
not give an implementable Gaussian approximation to the posterior. This is where theLA comes in.
To explain the LA construction, consider a posterior π whose mass concentrates in asmall neighborhood of the mode, which we call x̂. When the conditions of the BvM aresatisfied, this highest mode should be unique; otherwise, concentration cannot occur.
Since most of the mass of π is near x̂, we should incur only a small error by replacing thelog posterior with its second order Taylor expansion about x̂. This gives rise to the LA,
Date: May 17, 2024.

1

2

the Gaussian density γ̂ which is given by(0.1)

γ̂ = N x̂, ∇2 V (x̂)−1 ,x̂ = arg min V (x)x∈Θwhere π ∝ e−V is a density on Θ ⊆ Rd . The LA has proved to be an invaluable tool for
Bayesian inference in applications ranging from deep learning [6] to inverse
problems [4] to variable selection in high-dimensional regression [1].
Unlike methods such as MCMC, which can be made arbitrarily accurate, there are no
parameters to tune in the LA — it incurs some fixed approximation error. Quantifying
the LA’s error as a function of dimension d, sample size n, and model parameters, is aworthy task given its widespread use. It is also a challenging theoretical endeavor whendimension d is large, and currently a very active research area. Major contributions havebeen made e.g. by [9, 13, 8]. We also study the LA error in this work.
But arguably, it is even more important to go beyond the LA to develop new, moreaccurate approximations which better capture the complexity of the posterior π. For
example, a known downside of the LA γ̂ is that it is symmetric about the mode andtherefore cannot capture skewness of π. Instead of constructing an entirely new kind of
approximation, a natural idea is to correct the LA in some way to get a higher-orderaccuracy approximation.
Non-rigorous skew normal approximations have beendeveloped in [17], but we are aware of only a single work [7] that rigorously derives ahigher-order accurate LA. However, the corrected LA obtained by [7] is only shown tobe accurate in constant dimension d.
In modern applications involving very
high-dimensional parameters, d cannot be considered constant relative to sample size n.
So far, no prior work has obtained a higher-accuracy LA which is rigorously justified inhigh dimensions.
In this work, we develop a powerful technique to analyze the Laplace approximation
more precisely than was possible before. This technique leads us to derive the first evercorrection to the LA which provably improves its accuracy by an order of magnitude, in highdimensions. Our approach allows us to prove error bounds on this corrected LA in termsof a variety of error metrics discussed below. It also enables us to prove both tighter
upper bounds and the first ever lower bounds on the standard LA in high dimensions. In
particular, we prove that d2  n is in general necessary for accuracy of the LA. Finally, weapply our theory in two example models: a Dirichlet posterior arising from a multinomial
observation, and logistic regression with Gaussian design. In the latter setting, we provehigh probability bounds on the accuracy of the LA and skew-corrected LA in powers of√d/ n alone.
R EFERENCES[1] Rina Foygel Barber, Mathias Drton, and Kean Ming Tan. Laplace approximation in high-dimensional
Bayesian regression. In Statistical Analysis for High-Dimensional Data: The Abel Symposium 2014,
pages 15–36. Springer, 2016.
[2] David M Blei, Alp Kucukelbir, and Jon D McAuliffe. Variational inference: A review for statisticians.
Journal of the American statistical Association, 112(518):859–877, 2017.
[3] William M Bolstad. Understanding computational Bayesian statistics, volume 644. John Wiley & Sons,
2009.
[4] Tan Bui-Thanh, Carsten Burstedde, Omar Ghattas, James Martin, Georg Stadler, and Lucas C. Wilcox.
Extreme-scale uq for Bayesian inverse problems governed by PDEs. In SC ’12: Proceedings of the

3International Conference on High Performance Computing, Networking, Storage and Analysis, pages
1–11, 2012.
[5] Mary Kathryn Cowles and Bradley P. Carlin. Markov chain Monte Carlo convergence diagnostics: A
comparative review. Journal of the American Statistical Association, 91(434):883–904, 1996.
[6] Erik Daxberger, Agustinus Kristiadi, Alexander Immer, Runa Eschenhagen, Matthias Bauer, andPhilipp Hennig. Laplace redux-effortless Bayesian deep learning. Advances in Neural Information
Processing Systems, 34:20089–20103, 2021.
[7] Daniele Durante, Francesco Pozza, and Botond Szabo. Skewed Bernstein-von Mises theorem andskew-modal approximations. arXiv preprint arXiv:2301.03038, 2023.
[8] Tapio Helin and Remo Kretschmann. Non-asymptotic error estimates for the Laplace approximationin Bayesian inverse problems. Numerische Mathematik, 150(2):521–549, 2022.
[9] Mikolaj J Kasprzak, Ryan Giordano, and Tamara Broderick. How good is your Gaussian
approximation of the posterior? Finite-sample computable error bounds for a variety of useful
divergences. arXiv preprint arXiv:2209.14992, 2022.
[10] Bas JK Kleijn and Aad W van der Vaart. The Bernstein-von-Mises theorem under misspecification.
Electronic Journal of Statistics, 6:354–381, 2012.
[11] Jeffrey W Miller. Asymptotic normality, concentration, and coverage of generalized posteriors.
Journal of Machine Learning Research, 22(168):1–53, 2021.
[12] Thomas P. Minka. Expectation propagation for approximate Bayesian inference. In Proceedings of the
17th Conference in Uncertainty in Artificial Intelligence, UAI ’01, pages 362–369, San Francisco, CA,
USA, 2001. Morgan Kaufmann Publishers Inc.
[13] Vladimir Spokoiny. Dimension free nonasymptotic bounds on the accuracy of high-dimensional
laplace approximation. SIAM/ASA Journal on Uncertainty Quantification, 11(3):1044–1068, 2023.
[14] Luke Tierney and Joseph B Kadane. Accurate approximations for posterior moments and marginaldensities. Journal of the American Statistical Association, 81(393):82–86, 1986.
[15] A. W. van der Vaart. Asymptotic Statistics. Cambridge Series in Statistical and Probabilistic
Mathematics. Cambridge University Press, 1998.
[16] Martin J Wainwright, Michael I Jordan, et al. Graphical models, exponential families, and variational
inference. Foundations and Trends® in Machine Learning, 1(1–2):1–305, 2008.
[17] Jackson Zhou, Clara Grazian, and John T Ormerod. Tractable skew-normal approximations via
matching. Journal of Statistical Computation and Simulation, 94(5):1016–1034, 2024.
77 M ASSACHUSETTS AVENUE , B UILDING 2, R OOM 350B, C AMBRIDGE MA 02139
Email address: akasevi@mit.edu

INFORMATION-THEORETIC THRESHOLDS FOR PLANTED DENSE CYCLESCHENG MAO, ALEXANDER S. WEIN, AND SHENDUO ZHANGClassification AMS 2020: 62B10, 05C82
Keywords: Planted dense cycle, small-world network, information-theoretic thresholds
The Watts–Strogatz small-world model has been an influential random graph model
since its proposal in 1998 due to the ubiquity of the small-world phenomenon in complex
networks [12, 11]. In this model, there are n vertices with latent positions on a circle,
and the vertices are more likely to be connected to their k-nearest geometric neighbors
than to more distant vertices. In other words, a denser cycle of length n and width k is“planted” in the sparser ambient random graph on n vertices.
While there has been extensive literature on small-world networks and geometricgraphs in general, the associated statistical problems, such as detection and recovery of
the latent geometry from the observed random graph, have only gained attention more
recently.
The information-theoretic thresholds and efficient algorithms for thesmall-world model are studied in [2], but there remain several gaps between upper andlower bounds that are unknown to be inherent or not. Much sharper characterizationsof the recovery thresholds are given in [1, 3] but only for a small bandwidth parameterk = no(1) . From the perspective of graphon estimation, there have also been algorithms
and statistical analyses introduced for related models recently [6, 10, 5]. Moreover, thestudy of random geometric graphs has received broad interests in recent years; see thesurvey [4] and the works [7, 8] for high-dimensional random geometric graphs withedge noise similar to what we consider.
Since the model of interest consists of a hidden dense cycle planted in a sparser randomgraph, the recent work [9] refers to it as the planted dense cycle problem, following
the etymology of planted clique and planted dense subgraph problems. The work [9]studies the problem in the framework of low-degree polynomial algorithms, a frameworkthat has proved to be successful at predicting computational thresholds. However, themore fundamental information-theoretic (or statistical) thresholds—where no constraints
are placed on computation time—remain open, and we aim to address them in thiswork. More specifically, suppose that in an ambient random graph on n vertices withedge density q, there is a hidden cycle with expected width k = nτ and edge density
p. We find the information-theoretic thresholds for detecting the presence of the cycleand for recovering the location of the cycle, in terms of the parameters n, τ, p, q. In
particular, the information-theoretic thresholds for detection and recovery both differ
from the computational thresholds given in [9], justifying the existence of statistical-tocomputational gaps for this problem.
Problem setup. The planted dense cycle model can be described as follows. For any
a, b ∈ [0, 1], define d(a, b) := min{|a − b|, 1 − |a − b|}. In other words, d(a, b) is the distance
between a and b on a circle of circumference 1. Throughout the paper, we consider the1

2

setting where the number of vertices n grows, and other parameters p, q, r, τ ∈ [0, 1] may
depend on n. In the models to be defined, p and q will be the average edge densities on
and off the planted cycle respectively, r is the average edge density of the entire graph,and nτ is the bandwidth of the cycle, satisfying(1)n → ∞,

0 < q < r < p ≤ 1,

0 < τ < 1/2,

r = τ p + (1 − τ )q.

Definition 1 (Model P, planted dense cycle). Let z ∈ [0, 1]n be a latent random vector
whose entries z1 , . . . , zn are i.i.d. Unif([0, 1]) variables. Let Xij := 1d(zi ,zj )≤τ /2 for all (i, j) ∈[n]. For concreteness, we let Xii = 0 and Xji = Xij for i ̸= j, so that X ∈ Rn×n is the2adjacency matrix of the underlying cycle. We observe an undirected graph with adjacency
matrix A ∈ Rn×n whose edges, conditional on z1 , . . . , zn , are independently sampled asfollows: Aij ∼ Bern(p) if Xij = 1 and Aij ∼ Bern(q) if Xij = 0 for (i, j) ∈ [n]. We write2
A ∼ PA and (A, X) ∼ P (or, equivalently, (A, z) ∼ P).
Definition 2 (Model Q, Erdős–Rényi graph). We observe a G(n, r) Erdős–Rényi graph withadjacency matrix A ∈ Rn×n . We write A ∼ Q.
We now formulate the detection and recovery problems of interest.
Definition 3 (Detection). Let P and Q be the models from Definitions 1 and 2 respectively,
with parameters in (1). Observing the adjacency matrix A ∈ Rn×n of a graph, we test
H1 : A ∼ PA against H0 : A ∼ Q. We say that a test Φ, a {0, 1}-valued measurable functionof the observation A, achieves• strong detection, if limn→∞ [P{Φ(A) = 0} + Q{Φ(A) = 1}] = 0;
• weak detection, if lim supn→∞ [P{Φ(A) = 0} + Q{Φ(A) = 1}] < 1.
Definition 4 (Recovery). Let P be the model from Definition 1 with parameters in (1).
For (A, X) ∼ P, observing A, we aim to estimate X with an estimator X̂ ∈ Rn×n that ismeasurable with respect to A.
Consider the mean squared errorP
2
R(X̂, X) := 1≤i<j≤n E[(X̂ij − Xij ) ], where the expectation is with respect to (A, X) ∼ P.
We say that an estimator X̂ achieves• strong recovery, if limn→∞ nR(τX̂,X)= 0;( 2 ) (1−τ )• weak recovery, if lim supn→∞ nR(τX̂,X)< 1.
( 2 ) (1−τ )Note that each Xij is marginally a Bern(τ ) random variable, so estimating Xij by its mean
τ for all (i, j) ∈ [n]yields a trivial mean squared error n2 τ (1 − τ ), which justifies the2above definition.
Main results. Our main results are summarized in the following theorem.
Theorem 5 (Information-theoretic thresholds). Consider the detection and recovery
problems in Definitions 3 and 4 respectively, with parameters n, τ, p, q, r in (1).
2Furthermore, suppose that (log n)3 ≤ nτ ≤ (lognn)2 and logn n ≤ r ≤ 21 . Define λ := (p−q)
r(1−r)
which can be seen as the signal-to-noise ratio of the problem. Then we have:• If nτ λ → 0 as n → ∞, then no test achieves weak detection, and no estimatorachieves weak recovery.
(p−r)
nτ λ
• If log→ ∞ and nτlog→ ∞ as n → ∞, then there is a test that achieves strongnnrecovery, and there is an estimator that achieves strong recovery.

3

Our conditions for positive and negative results match up to a logarithmic factor inmost cases. The threshold for both detection and recovery is nτ λ = Θ̃(1), where nτ isthe bandwidth of the planted dense cycle, and λ is the edge-wise signal-to-noise ratio.
Remark 6 (Statistical-to-computational gap). Our information-theoretic results for thedetection and recovery of a planted dense cycle complement the work [9] which studies lowdegree polynomial algorithms for the same model. Suppose that the parameters in (1) satisfy
Cq ≤ p ≤ C ′ q for some constants C ′ > C > 1. As shown in [9], the detection threshold forthe class of low-degree polynomial algorithms is n3 p3 τ 4 = no(1) , and the recovery thresholdfor the class of low-degree polynomial algorithms is npτ 2 = no(1) . In this regime, the signal2
to-noise ratio λ = (p−q)
has the same order as p, so the information-theoretic thresholdr(1−r)
from Theorem 5 can be expressed as npτ = Θ̃(1). Therefore, Theorem 5 and [9] together
suggest that there are statistical-to-computational gaps for both detection and recovery of aplanted dense cycle.
R EFERENCES[1] Vivek Bagaria, Jian Ding, David Tse, Yihong Wu, and Jiaming Xu. Hidden hamiltonian cycle recovery
via linear programming. Operations research, 68(1):53–70, 2020.
[2] Tony Cai, Tengyuan Liang, and Alexander Rakhlin. On detection and structural reconstruction ofsmall-world random networks. IEEE Transactions on Network Science and Engineering, 4(3):165–176,
2017.
[3] Jian Ding, Yihong Wu, Jiaming Xu, and Dana Yang. Consistent recovery threshold of hidden nearest
neighbor graphs. In Conference on Learning Theory, pages 1540–1553. PMLR, 2020.
[4] Quentin Duchemin and Yohann De Castro. Random geometric graph: Some recent developments
and perspectives. arXiv preprint arXiv:2203.15351, 2022.
[5] Christophe Giraud, Yann Issartel, and Nicolas Verzelen. Localization in 1d non-parametric latent
space models from pairwise affinities. Electronic Journal of Statistics, 17(1):1587–1662, 2023.
[6] Jeannette Janssen and Aaron Smith. Reconstruction of line-embeddings of graphons. Electronic
Journal of Statistics, 16(1):331–407, 2022.
[7] Suqi Liu and Miklós Z Rácz. Phase transition in noisy high-dimensional random geometric graphs.
Electronic Journal of Statistics, 17(2):3512–3574, 2023.
[8] Suqi Liu and Miklós Z Rácz. A probabilistic view of latent space graphs and phase transitions.
Bernoulli, 29(3):2417–2441, 2023.
[9] Cheng Mao, Alexander S. Wein, and Shenduo Zhang. Detection-recovery gap for planted dense
cycles. In Proceedings of Thirty Sixth Conference on Learning Theory, volume 195 of Proceedings of
Machine Learning Research, pages 2440–2481, 12–15 Jul 2023.
[10] Amine Natik and Aaron Smith. Consistency of spectral seriation. arXiv preprint arXiv:2112.04408,
2021.
[11] Duncan J Watts. Small worlds: the dynamics of networks between order and randomness, volume 36.
Princeton university press, 2004.
[12] Duncan J Watts and Steven H Strogatz. Collective dynamics of ‘small-world’ networks. Nature,393(6684):440–442, 1998.
S CHOOL OF M ATHEMATICS , G EORGIA I NSTITUTE OF T ECHNOLOGY, ATLANTA , GA, USA
Email address: cheng.mao@math.gatech.edu
D EPARTMENT OF M ATHEMATICS , U NIVERSITY OF C ALIFORNIA , D AVIS , CA, USA
Email address: aswein@ucdavis.edu
S CHOOL OF M ATHEMATICS , G EORGIA I NSTITUTE OF T ECHNOLOGY, ATLANTA , GA, USA
Email address: szhang705@gatech.edu

GEOMETRY AND DYNAMICS OF THE DEEP LINEAR NETWORKGOVIND MENONClassification AMS 2020: 68Q32, 68T05, 93E35
Keywords: Deep learning, gradient flows, entropy of overparametrization.

1. G RADIENT DYNAMICS IN DEEP LEARNING
These lectures study the training dynamics of the deep linear network (DLN) using
the geometric theory of dynamical systems. The DLN is deep learning restricted tolinear functions.
However, it retains two essential features of deep learning:
overparametrization and degenerate loss functions. Overparametrization provides afoliation of phase space by invariant manifolds. Of these, there is a fundamental
invariant manifold, the balanced manifold, which is itself foliated by group orbits. Thisgeometric structure allows us to define a natural Boltzmann entropy (the logarithm of
the volume of a group orbit) that may be computed explicitly. Our approach unifies thework of several authors into a thermodynamic framework.
2. T HE MODEL
We fix two positive integer d and N referred to as the width and depth of the network.
The state space for the DLN is MN
d , where Md denotes the space of d × d real matrices.
N
Each point W ∈ Md is denoted by W = (WN , WN −1 , . . . , W1 ). We equip Md with
 the
PN
T
N
2
Frobenius norm so that Md is Euclidean with the norm ∥W∥ = p=1 Tr Wp Wp . We
define the projection ϕ : MN
d → Md and end-to-end matrix X by
ϕ(W) = WN WN −1 · · · W1 =: W.

(2.1)

We assume given an energy E : Md → R. The training dynamics are described by thegradient flow in MN
d of the ‘lifted’ loss function L = E ◦ ϕ
Ẇ = −∇W L(W).

(2.2)

A computation then reduces equation (2.2) to a collection of N equations in Md(2.3)

Ẇp = −(WN · · · Wp+1 )T E ′ (X)(Wp−1 · · · W1 )T ,Here E ′ (X) denotes the d × d matrix with entries(2.4)

E ′ (X)ij =

∂E
,
∂Xjk

Date: May 20, 2024.

1

1 ≤ j, k ≤ d.

p = 1, . . . , N.

2

3. O VERPARAMETRIZATION AND DEGENERATE LOSS FUNCTIONS
Let fp : Rd → Rd define the linear function fp (x) = Wp x for 1 ≤ p ≤ N . Then the linear
function f : Rd → Rd corresponding to the matrix f (x) = Xx is f = fN ◦ fN −1 . . . ◦ f1 .
The output function f depends only on the end-to-end matrix X. However, the samefunction f may be represented by N d2 choices of training parameters W. Thus, despite
the absence of the nonlinear activation element and shifts, the choice of variables in the
DLN models overparametrization in deep learning.
Natural learning tasks, such as matrix completion, give rise to degenerate loss functions.
Assume given a subset S ⊂ {(i, j)}1≤i,j≤d and assume given the values of Xij , for (i, j) ∈S, say Xij = aij . The task in matrix completion task is to obtain a principled answer tothe question: how do we reconstruct X from the partial observations Xij for (i, j) ∈ S?
The DLN is used to study this question as follows. Introduce the quadratic loss function
1 X
(3.1)
E(X) =|Xij − aij |2 ,2(i,j)∈S

and seek the limit of X(t) = ϕ(W(t)) as t → ∞ when W(t) solves the gradient flow (2.3).
The loss function E(X) is degenerate because it does not depend on Xij when (i, j) ̸∈ S.
Thus, E is minimized on the affine subspace
S = {X ∈ Md : Xij = aij ,(i, j) ∈ S}.

In particular, E does not have compact sublevel sets, and we cannot apply La Salle’s
invariance principle. However, a new mathematical structure emerges.
4. G EOMETRIC FEATURES OF TRAINING DYNAMICS(1) Invariant varieties. Equation (2.3) shows that each Ẇp is obtained by a linear
transformation of E ′ (X). Thus, the space of gradients of L = E ◦ ϕ at W has
2only d2 dimensions, whereas the dimension of the tangent space TW MN
d is N d .
A simple calculation then yields (N − 1)d(d + 1)/2 conserved quantities [1, 2, 3](4.1)

T
Gp = Wp+1
Wp+1 − Wp WpT ,1 ≤ p ≤ N − 1.

The solution set to these equations is a conic section in MN
d parametrized by N −1
symmetric matrices. We call these the G-balanced varieties, denoted MG , whereG = (GN −1 , . . . , G1 ). Each of these varieties is invariant under the flow (2.2).
(2) Riemannian submersion of M When G = 0 and X has full rank, we obtain afundamental invariant manifold termed the balanced manifold, M. The dynamicsof X(t) ∈ Md are slaved to the dynamics of W(t) ∈ MG . However, when W(t)lies on M, then X(t) satisfies the Riemannian gradient flow [4](4.2)

Ẋ = −gradgN E(X),

where the metric g N is obtained by Riemannian submersion of the metric on MNdrestricted to M.
(3) Group orbits and an entropy formula Given X ∈ Md , we may lift it to an O(d)N −1
orbit OX ∈ MN
d , such that ϕ(W) = X for each W ∈ OX . We may then quantify
overparametrization
withaBoltzmannentropyoftheform
S(X) = log vol(OX ) [7], improving the main result in [5].

3(4) From gradient descent to the Riemannian Langevin equation (RLE). We augment
the gradient dynamics with natural stochastic dynamics, using the geometry of
overparametrization to add noise in the ‘null directions’. This extends the idea of
optimization by the DLN to that of Gibbs sampling, providing a thermodynamic
framework in which the energy E(X) is replaced by a Helmholtz free energy
F (X) = E(X) − β −1 S(X) at inverse temperature β > 0 [6, 7].
The volume formula is as follows. Given X ∈ Md with full rank, let its SVD be X =QN ΣQT0 , where QN , Q0 ∈ O(d), Σ = diag(σ1 , ..., σd ). Then
v
s
u 2
2
Y
u σi − σj2van(Σ )
N −1
N −1
t 2
(4.3)
vol(OX ) = cd
= cd2
2van(Σ N )σiN − σjN
1≤i<j≤dr
Q
1
where cd = vol(O(d)) = 2 2 d(d+3) dr=1 Γ(π 2r ) is the volume of O(d)van(Λ) denotes the2
Vandermonde determinant associated to a diagonal matrix Λ.
The formulation of the RLE, especially the definition of the noise in the null directions,requires greater detail and is contained in the forthcoming work [6, 7].
R EFERENCES[1] Arora, Sanjeev and Cohen, Nadav and Hazan, Elad. On the optimization of deep networks: Implicit
acceleration by overparameterization. International conference on machine learning, 244-253, 2018.
[2] Arora, Sanjeev and Cohen, Nadav and Hu, Wei and Luo, Yuping. Implicit regularization in deep
matrix factorization. Advances in Neural Information Processing Systems, 32, 2019.
[3] Arora, Sanjeev and Cohen, Nadav and Hu, Wei and Luo, Yuping. A convergence analysis of gradientdescent for deep linear neural networks. arXiv:1810.02281, 2018.
[4] Bah, Bubacarr and Rauhut, Holger and Terstiege, Ulrich and Westdickenberg, Michael. Learning
deep linear neural networks: Riemannian gradient flows and convergence to global minimizers.
Information and Inference: A Journal of the IMA, 11, 307–353, 2022.
[5] Cohen, Nadav and Menon, Govind and Veraszto, Zsolt. Deep Linear Networks for Matrix
Completion—an Infinite Depth Limit. SIAM Journal on Applied Dynamical Systems, 22, 3208–3232,
2023.
[6] Menon, Govind. The geometry of the deep linear network. Preprint, 2024.
[7] Menon, Govind and Yu, Tianmin. An entropy formula for the deep linear network. Preprint, 2024.
B ROWN U NIVERSITY, 182 G EORGE S T., P ROVIDENCE , RI 02912, USA
Email address: govind menon@brown.eduOPTIMAL TRANSPORT MAP ESTIMATION IN GENERAL FUNCTION SPACESJONATHAN NILES-WEED

Classification AMS 2020: 62G05
Keywords: Optimal transportation, Brenier maps, empirical processes
We study the problem of estimating a function T given independent samples from a
distribution P and from the pushforward distribution T♯ P . This setting is motivated byapplications in the sciences, where T represents the evolution of a physical system overtime, and in machine learning, where, for example, T may represent a transformation
learned by a deep neural network trained for a generative modeling task.
To ensure identifiability, we assume that T = ∇ϕ0 is the gradient of a convex
function, in which case T is known as an optimal transport map. The estimation of suchmaps was inaugurated by Hütter and Rigollet [2], and has been subsequently studied ina number of works [2, 1, 3, 4, 5]. These works all study estimation of T under theassumption that it lies in a Hölder class, but general theory is lacking. We present a
unified methodology for obtaining rates of estimation of optimal transport maps ingeneral function spaces. Our assumptions are significantly weaker than those appearing
in the literature: we require only that the source measure P satisfy a Poincaré inequality
and that the optimal map be the gradient of a smooth convex function that lies in aspace whose metric entropy can be controlled. As a special case, we recover known
estimation rates for Hölder transport maps, but also obtain nearly sharp results in manysettings not covered by prior work. For example, we provide the first statistical rates of
estimation when P is the normal distribution, between log-smooth and strongly
log-concave distributions, and when the transport map is given by an infinite-width
shallow neural network.
R EFERENCES
[1]
[2]
[3]
[4]
[5]

N. Deb, P. Ghosal, and B. Sen. Rates of estimation of optimal transport maps using plug-in estimators
via barycentric projections. arXiv preprint arXiv:2107.01718, 2021.
J.-C. Hütter and P. Rigollet. Minimax estimation of smooth optimal transport maps. The Annals of
Statistics, 49(2):1166–1194, 2021.
T. Manole, S. Balakrishnan, J. Niles-Weed, and L. Wasserman. Plugin estimation of smooth optimaltransport maps. arXiv preprint arXiv:2107.12364, 2021.
B. Muzellec, A. Vacher, F. Bach, F.-X. Vialard, and A. Rudi. Near-optimal estimation of smoothtransport maps with kernel sums-of-squares. arXiv preprint arXiv:2112.01907, 2021.
A.-A. Pooladian and J. Niles-Weed. Entropic estimation of optimal transport maps. arXiv preprint
arXiv:2109.12004, 2021.

C OURANT I NSTITUTE FOR M ATHEMATICAL S CIENCES , N EW Y ORK U NIVERSITY, N EW Y ORK , NY
Email address: jnw@cims.nyu.edu

1WHEN A SYSTEM OF REAL QUADRATIC EQUATIONS HAS A SOLUTIONMARK RUDELSONClassification AMS 2020: Primary: 14P05. Secondary: 14Q30, 90C22.
Keywords: positive semidefinite relaxation, quadratic equations, algorithms
This talk is based on a joint paper with Alexander Barvinok [1]. It discussed a problem
of providing a computationally efficient certificate for the existence of a real solution of asystem of m real quadratic equations with n variables. The question of feasibility of thistype of a system of equations systems appears in various computer science contexts, see,for example, [2], [3].
If m and n are both allowed to grow, the problem becomes computationally hard.
Unless the computational complexity hierarchy collapses, no polynomial time algorithm
provides a necessary and sufficient condition for the existence of solutions of such asystem. In fact, testing the feasibility of an arbitrary system of real polynomial equations
can be easily reduced to testing the feasibility of a system quadratic ones. First, we
gradually reduce the degree of polynomials by repeatedly introducing new variables andequations of the type ξij − ξi ξj = 0, which allows us to replace the product ξi ξj of oldvariables by a single new variable ξij , and hence eventually reduce a given polynomial
system to a systemqi (x) = 0 for i = 1, . . . , m,where qi are quadratic, not necessarily homogeneous, polynomials. Then we introduce
another, more delicate change of variables which is based on a semi-definite relaxation.
This semi-definite relaxation problem can be efficiently solved in polynomial time. It
allows to further reduce the system of general quadratic equations to a systemqi (x) = ⟨Qi x, x⟩ = tr(Qi ) for i = 1, . . . , m,

(0.1)

where Qi are n × n symmetric matrices and⟨x, y⟩ =

n
X

ξi ηi

for x = (ξ1 , . . . , ξn )and y = (η1 , . . . , ηn )i=1

is the standard inner product in Rn and tr(A) denotes the trace of a matrix A.
We present computationally simple sufficient criterion for (0.1), to have a solution.
First, note that any linear combination of the equations of the form (0.1) has the sameform. Using this observation, we can choose the simplest set of matrices Q1 , . . . , Qm
which generates the same system of equations. To this end, we introduce the standard
inner product of matrices⟨A, B⟩ = tr(A⊤ B).
Then without loss of generality, we can assume that the matrices Q1 , . . . , Qm form an
orthonormal system with respect to this product. Computing such an orthonormal system
is also efficient and can be done, for example, by the Gram-Schmidt procedure.

1

2

For an n × n real symmetric matrix Q, we denote by ∥Q∥op the operator norm of Q,
that is, the largest absolute value of an eigenvalue of Q.
Now, we can formulate the main result.
Theorem 0.1. There is an absolute constant η > 0 such that the following holds. Let
Q1 , . . . , Qm , m ≥ 3, be linearly independent n × n symmetric matrices. Suppose that
m
X
i=1

A2i

η
m

≤
op

for some (equivalently, for any) orthonormal basis A1 , . . . , Am of the linear subspacespanned by Q1 , . . . , Qm . Then the system of quadratic equations⟨Qi x, x⟩ = tr(Qi ) fori = 1, . . . , m

has a solution x ∈ Rn .
We also show that while the choice of the basis
Pm A21 , . . . , Am depends on theorthogonalization procedure, the key quantity ∥ i=1 Ai ∥op is independent of it.
Moreover, since this quantity is the largest eigenvalue of a positive definite matrix, it
can be computed efficiently.
To illustrate the applicability of Theorem 0.1, we show that if Q1 , . . . , Qm are symmetricmatrices with independent identically distributed entries, then with high probability, the√system is feasible provided that m ≤ c n for some absolute constant c > 0.
While the condition we obtain is of an algebraic nature, the proof relies on analytic
tools including Fourier analysis and measure concentration.
R EFERENCES[1] A. Barvinok, M. Rudelson. When a system of real quadratic equations has a solution. Adv. Math. , 403(2022), Paper No. 108391, 38 pp.
[2] D. Bienstock. A note on polynomial solvability of the CDT problem. SIAM Journal on Optimization,
26 (2016), no. 1, 488–498.
[3] L. Liberti, C. Lavor, N. Maculan and A. Mucherino. Euclidean distance geometry and applications.
SIAM Review , 56 (2014), no. 1, 3–69.
D EPARTMENT OF M ATHEMATICS , U NIVERSITY OF M ICHIGAN , A NN A RBOR , MI 48109-1043, USA
Email address: rudelson@umich.edu

APPROXIMATELY HADAMARD MATRICESAND RIESZ BASES IN RANDOM FRAMESMARK RUDELSONClassification AMS 2020: 60B20.
Keywords: random frames, random matrices, Hadamard matrices.
This talk was based on a joint paper with Xiaoyu Dong [4]. We considered a problem
which originated in signal processing. This problem was quickly reduced to a question
about sub-matrices of a random matrix. We showed that this probabilistic question canbe further reduced to a completely deterministic problem about the existence ofapproximately Hadamard matrices. The solution of the last problem was achieved by a
combination of number theoretic and probabilistic tools, thus going back to the realm of
probability.
Let n < N be natural numbers. A set of vectors X1 , . . . , XN ∈ Rn is called a frame if(0.1)

K(n, N ) ∥x∥22 ≤N
X⟨x, Xj ⟩2 ≤ RK(n, N ) ∥x∥22

j=1

for all x ∈ Rn . Here R ≥ 1 is called the frame constant, and K(n, N ) > 0 is some functionof n and N . The frame is considered to be good if its frame constant is relatively small.
The notation ∥x∥2 stands for the Euclidean norm of the vector x = (x1 , . . . , xn ).
In the last 40 years, frame theory became a well-developed area of appliedmathematics, see [1], [2], [3], and the references therein. A frame can intuitively be
regarded as overcomplete basis in Rn . Because of this property, frames became avaluable tool in signal transmission. A signal which is viewed as an n-dimensional
vector can be encoded by the sequence of its inner products with the frame vectors. If
this sequence is transmitted over a communication line, then the original signal can bereconstructed even if part of the coefficients is lost or corrupted in the process oftransmission. Moreover, this encoding is robust, which means that if the inner products
are evaluated with some noise, then the reconstructed version will be close to theoriginal one with the error depending on the noise magnitude.
One of the most popular classes of frames in algorithmic applications is the set of
random frames. Such frames became also the method of choice in compressed sensing
where one needs to reconstruct a low complexity signal from a small number of linear
measurements, see, e.g., [5]. For example, if complexity is measured as the size of the
support, and the support itself is unknown, the random frames provide robust recovery
with optimal or almost optimal theoretical guarantees.
We consider a problem when a random frame contains a copy, or many copies of a“nice” basis. This problem can be conveniently translated to the language of randommatrices. Define the condition number of a matrix A asmax∥x∥2 =1 ∥Ax∥2
κ(A) =.
min∥x∥2 =1 ∥Ax∥2

1

2

With this notation, the frame property (0.1) can be rewritten as κ(An,N ) ≤ C where
An,N is the n × N matrix with columns X1 , . . . , XN . Thus, the problem of existence ofa ‘nice” basis in a random frame can be recast as the question of existence of one or
many well-conditioned square n × n sub-matrices of an n × N random matrix An,N with
i.i.d. entries. Our main result shows that the probability of finding such a sub-matrix
undergoes a phase transition when N is exponential in terms of n. Since the upper
and the lower bound hold under somewhat different assumptions, we formulate themseparately.
Denote by [N ] the set {1, . . . , N }. Let A be an n × N matrix. If I ⊂ [N ], denote by AI
the sub-matrix of A whose columns belong to I. The following theorem shows that if N
is exponential in n, then with high probability, the n×N random matrix has many squaresubmatrices with uniformly bounded condition numbers. In the language of frames, itmeans that a random frame with exponentially many vectors contains a large number ofbases whose frame constants are uniformly bounded.
Theorem 0.1. Let A be an n × N matrix with i.i.d. symmetric non-degenerate entries.
Then there exist constants c, C, α, β > 0 depending on the distribution of entries of A withthe following property.
Assume that N ≥ exp(Cn). Then there exists L ≥ exp(cn) such that
P exist disjoint subsets I1 , . . . , IL of [N ] with |Ij | = n and κ(AIj ) < α for all j ∈ [L]≥ 1 − exp (− exp(βn)) .
The strategy of proving Theorem 0.1 relies on finding columns of A which are close
to the columns of a certain deterministic n × n matrix V having a bounded condition
number. The key to this strategy is a successful choice of the pattern matrix V . The
requirement that a column of A can be close to a column of V with a non-negligible
probability forces us to look for a matrix V which is a scaled copy of a matrix with ±1
entries. Such matrices are known in some cases. For instance, the condition number ofany Hadamard matrix is one. An n × n matrix H is called Hadamard if n−1/2 H is an
isometry. Hadamard matrices is a well-studied subject, and a number of constructionsof such matrices are available. Yet, the dimensions in which Hadamard matrices wereconstructed are rare, so we have to use approximately Hadamard matrices instead.
The existence of an approximately Hadamard matrix in any dimension is establishedin the next theorem.
Theorem 0.2. There exists a constant C ≥ 1 such that for any n ∈ N, one can find an n × nmatrix V with ±1 entries satisfying
κ(V ) ≤ C.
The proof of Theorem 0.2 relies on Vinogradov’s theorem from analytic number theory
and combines number-theoretic and probabilistic ideas.
The conclusion of Theorem 0.1 holds under minimal assumptions on the distributionof entries. If we assume that the entries of the matrix are sub-gaussian, then the boundof Theorem 0.1 becomes sharp. Recall
 that a random variable X is called subgaussianX2
if there is a > 0 such that E exp a2 ≤ 2. Subgaussian random variables form a large
family containing many naturally arising ones, see, e.g. [5].

3The next theorem shows that finding a submatrix with a bounded condition number
requires an exponential number of columns for matrices with subgaussian entries.
Theorem 0.3. Let X be a centered subgaussian random variable.
C, c, c̃, t0 > 0 with the following property. Let t > t0 , and assume thatc̃N ≤ exp 4 n .
tThen there exist

Let A be an n × N matrix whose entries are independent copies of X. Thenn2
P (∃I ⊂ [N ] |I| = n and κ(AI ) < t) ≤ exp −c 4 .
t
R EFERENCES[1] P. Casazza, G. Kutyniok, F.Philipp. Introduction to finite frame theory. Finite frames, 1–53, Appl.
Numer. Harmon. Anal., Birkhäuser/Springer, New York, 2013.
[2] O. Christensen. Frames and bases. An introductory course. Applied and Numerical Harmonic
Analysis. Birkhäuser Boston, Inc., Boston, MA, 2008.
[3] O. Christensen. An introduction to frames and Riesz bases. Second edition. Applied and Numerical
Harmonic Analysis. Birkhäuser/Springer, 2016.
[4] X. Dong, M. Rudelson. Approximately Hadamard matrices and Riesz bases in random frames. Int.
Math. Res. Not., (2024), no. 3, 2044–2065.
[5] R. Vershynin. High-dimensional probability. An introduction with applications in data science. Witha foreword by Sara van de Geer. Cambridge Series in Statistical and Probabilistic Mathematics, 47.
Cambridge University Press, Cambridge, 2018.
D EPARTMENT OF M ATHEMATICS , U NIVERSITY OF M ICHIGAN , A NN A RBOR , MI 48109-1043, USA
Email address: rudelson@umich.edu

RECENT DEVELOPMENTS IN GROUP TESTING:
FUNDAMENTAL LIMITS AND ALGORITHMS
JONATHAN SCARLETTClassification AMS 2020:
Keywords: Group testing, sparsity, information-theoretic limits, discrete algorithms
The group testing problem concerns discovering a small number of defective items
within a large population by performing tests on pools of items. A test is positive if
the pool contains at least one defective, and negative if it contains no defectives. Thisis a sparse inference problem with a combinatorial flavor, with applications in medical
testing, biology, multi-access communication, database systems, and more.
In this talk, I reviewed recent advances in the mathematics of group testing, including
both information-theoretic limits and performance bounds for practical algorithms, with
an emphasis on the following defining features:• Non-adaptive testing (all tests must be designed in advance) vs. adaptive testing(tests are designed sequentially based on previous outcomes)• Noiseless testing (tests are perfectly reliable) vs. noisy tests (some test outcomesare corrupted)
These two features lead to 4 distinct settings with varying degrees of difficulty. The
over-arching goal is to determine the smallest possible number of tests while
maintaining reliable recovery of the defective set. This question has been studied from awide variety of perspectives, including high-dimensional statistics, information theory,
discrete algorithms, combinatorics, error-correcting codes, and others.
The results that I surveyed in this talk are summarized as follows:• The noisy adaptive setting has long been very well-understood, with a prominent
approach being Hwang’s adaptive binary splitting algorithm [9].
• The noiseless non-adaptive setting was extensively studied over the last decade orso (by myself and others), eventually leading to exact information-theoretic limits
and a method for matching them in with an efficient algorithm [1, 12, 2, 10, 6, 7].
• The noisy non-adaptive setting has had analogous developments and generally
lagged behind [4, 13, 14], but very recent works (including ours) works havesubstantially closed these gaps [5, 8].
• The noisy adaptive setting was seemingly overlooked for a long time, and myworks showed that it can offer significant reductions in the number of tests
compared to non-adaptive methods [11, 15].
A detailed survey covering all of these settings can be found in [3].

1

2

R EFERENCES[1] M. Aldridge, L. Baldassini, and O. Johnson, “Group testing algorithms: Bounds and simulations,”
IEEE Trans. Inf. Theory, vol. 60, no. 6, pp. 3671–3687, June 2014.
[2] M. Aldridge, “The capacity of Bernoulli nonadaptive group testing,” IEEE Trans. Inf. Theory, vol. 63,
no. 11, pp. 7142–7148, 2017.
[3] M. Aldridge, O. Johnson, and J. Scarlett, “Group testing: An information theory perspective,” Found.
Trend. Comms. Inf. Theory, vol. 15, no. 3–4, pp. 196–392, 2019.
[4] C. L. Chan, P. H. Che, S. Jaggi, and V. Saligrama, “Non-adaptive probabilistic group testing withnoisy measurements: Near-optimal bounds with efficient algorithms,” in Allerton Conf. Comm., Ctrl.,
Comp., Sep. 2011, pp. 1832–1839.
[5] J. Chen and J. Scarlett, “Exact thresholds for noisy non-adaptive group testing,” arXiv preprint
arXiv:2401.04884, 2024.
[6] A. Coja-Oghlan, O. Gebhard, M. Hahn-Klimroth, and P. Loick, “Information-theoretic and algorithmic
thresholds for group testing,” in Int. Colloq. Aut., Lang. and Prog. (ICALP), 2019.
[7] A. Coja-Oghlan, O. Gebhard, M. Hahn-Klimroth, and P. Loick, “Optimal group testing,” in Conf. Learn.
Theory (COLT), 2020.
[8] A. Coja-Oghlan, M. Hahn-Klimroth, L. Hintze, D. Kaaser, L. Krieg, M. Rolvien, and O. Scheftelowitsch,“Noisy group testing via spatial coupling,” arXiv preprint arXiv:2402.02895, 2024.
[9] F. Hwang, “A method for detecting all defective members in a population by group testing,” J. Amer.
Stats. Assoc., vol. 67, no. 339, pp. 605–608, 1972.
[10] O. Johnson, M. Aldridge, and J. Scarlett, “Performance of group testing algorithms with near-constanttests-per-item,” IEEE Trans. Inf. Theory, vol. 65, no. 2, pp. 707–723, Feb. 2019.
[11] J. Scarlett, “Noisy adaptive group testing: Bounds and algorithms,” IEEE Trans. Inf. Theory, vol. 65,
no. 6, pp. 3646–3661, June 2019.
[12] J. Scarlett and V. Cevher, “Phase transitions in group testing,” in ACM-SIAM Symp. Disc. Alg. (SODA),
2016.
[13] J. Scarlett and V. Cevher, “Near-optimal noisy group testing via separate decoding of items,” IEEE J.
Sel. Topics Sig. Proc., vol. 2, no. 4, pp. 625–638, 2018.
[14] J. Scarlett and O. Johnson, “Noisy non-adaptive group testing: A (near-) definite defectives
approach,” IEEE Trans. Inf. Theory, vol. 66, no. 6, pp. 3775–3797, 2020.
[15] B. Teo and J. Scarlett, “Noisy adaptive group testing via noisy binary search,” IEEE Trans. Inf. Theory,
vol. 68, no. 5, pp. 3340–3353, 2022.
NUS S CHOOL OF C OMPUTING (O FFICE COM3-02-57), 13 C OMPUTING D RIVE , 117417
Email address: scarlett@comp.nus.edu.sg

SAMPLING FROM CONVEX SETS WITH A COLD START
USING MULTISCALE DECOMPOSITIONSHARIHARAN NARAYANAN, AMIT RAJARAMAN, AND PIYUSH SRIVASTAVAA BSTRACT. A standard approach for sampling approximately uniformly from a convex
body K ⊆ Rn is to run a random walk within K. The requirement is that starting froma suitable initial distribution, the random walk should “mix rapidly”, i.e., after a numberof steps that is polynomial in n and the aspect ratio R/r (here, K is assumed to containa ball of radius r and to be contained within a ball of radius R), the distribution of therandom walk should come close to the uniform distribution πK on K. Different random
walks differ in aspects such as the ease of implementation of each step, or suitability fora specific class of convex bodies. Therefore, the rapid mixing of a wide variety of random
walks on convex bodies has been studied.
Many proofs of rapid mixing of such random walks however require that the initial
distribution of the random walk is not too different from the target distribution πK . In
particular, they require that the probability density function of the initial distribution
with respect to the uniform distribution πK on K must be bounded above by poly(n):this is called a warm start. Achieving such a warm start often requires a non-trivial
pre-processing step before the random walk can be started. This motivates the problem
of proving rapid mixing from “cold starts”, i.e., when the density of the initial distribution
with respect to πK can be as high as exp(poly(n)). In contrast to warm starts, a cold startis usually trivial to achieve. However, rapid mixing from a cold start may not hold forevery random walk, e.g., the well-known “ball walk” does not have rapid mixing from anarbitrary cold start. On the other hand, for the “hit-and-run” random walk, Lovász andVempala proved rapid mixing from a cold start. For the related coordinate hit-and-run(CHR) random walk, which has been found to be promising in computational
experiments, a rapid mixing result starting from a warm start was proven only recently,
while the question of whether CHR mixes rapidly from a cold start remained open.
In this paper, we construct a family of Markov chains inspired by classical multiscaledecompositions of subsets of Rn into countably many axis-aligned cubes. We show thateven with a cold start, the mixing times of these chains are bounded by a polynomial
in n and the aspect ratio of the body. Our main technical ingredient is an isoperimetricinequality for K for a metric that magnifies distances between points that are close to theboundary of K. As a byproduct of the analysis of this new family of chains, we show thatthe coordinate hit-and-run (CHR) random walk also mixes rapidly from a cold start, andalso from any point that is not too close to the boundary of the body.

Classification AMS 2020: 60J22, 68Q25
Keywords: Convex bodies, Markov chains, Isoperimetric inequalities.
A full version of the paper is available at arXiv:2211.04439. An extended abstract
appears in Proceedings of the 55th Annual ACM Symposium on Theory of Computing(STOC), June 2023, pp. 117–130.

1

RANDOM PERTURBATION OF LOW-RANK MATRICESKE WANGClassification AMS 2020: 60B20, 62H25, 62H30.
Keywords: Singular vector perturbation, singular subspace perturbation, low-rank
structures, random matrices.
We consider an unknown N × n data matrix A. Suppose we cannot observe A directly
e given by
but instead have access to a corrupted version A
e := A + E,
A
where E represents the noise matrix. A classical question is to estimate the leading
singular values and singular vectors (alternatively, eigenvalues and eigenvectors) of Ae This problem is crucial across various fields including engineering, statistics,
using A.
machine learning, computer science, and mathematics. Generally, it is assumed that A is
not arbitrary but exhibits specific structural characteristics, such as low-rank, which is acommon assumption in numerous applications.
Assume that A has rank r ≥ 1. The singular value decomposition (SVD) of A takes theform A = U ΣV T , where Σ = diag(σ1 , . . . , σr ) is a diagonal matrix containing the non-zerosingular values σ1 ≥ σ2 ≥ · · · ≥ σr > 0 of A; the columns of the matrices U = (u1 , . . . , ur )and V = (v1 , . . . , vr ) are the orthonormal left and right singular vectors of A, respectively.
In other words, ui and vi are the left and right singular vectors corresponding to σi . It
follows that U T U = V T V = Ir , where Ir is the r × r identity matrix. For conveniencee similarly by A
e=UeΣe Ve T , where
we will take σr+i = 0 for all i ≥ 1. Denote the SVD of A
e are the singular values σ
the diagonal entries of Σe1 ≥ σe2 ≥ · · · ≥ σemin{N,n} ≥ 0, and thee and Ve are the orthonormal left and right singular vectors, denoted by ucolumns of Ueiand vei , respectively.
The famous Davis–Kahan bound addresses this question for the eigenvectors ofdeterministic real symmetric matrices. An analogous version for the singular vectors of
non-square matrices was established by Wedin. The key parameters in this bound are
the gap (or separation) δ1 between the largest singular values of A given byδ1 := σ1 − σ2 and the spectral norm of E defined by kEk := maxkuk=1 kEuk, where kuk
denotes the Euclidean norm of the vector u. The classical Wedin’s bound givessin ∠(u1 , ue1 ) ≤ C

kEk,
δ1

where ∠(u1 , ue1 ) is the acute angle between u1 and u
e1 taken in [0, π/2] and C > 0 is anabsolute constant. The same bound holds for sin ∠(v1 , ve1 ).
For the case when E is a random matrix, an earlier paper by O’Rourke, Vu and myself[1] proved a version of bounds of the form
(0.1)

sin ∠(u1 , ue1 ) .

r1/α kEk kEk2
+
+
,
δ1
σ1
σ 1 δ11

2which holds with high probability. Here, α > 0 is a parameter that depends on the
distribution of the random matrix E. When the rank r of A is sufficiently small comparedto the dimensions and σ1 , δ1 are sufficiently large, this bound improves upon the boundin Wedin’s theorem. The first term on the right-hand side of (0.1) was conjectured asa consequence of the true dimension being actually r. The second term represents thesignal-to-noise ratio. However, the third term on the right-hand side of (0.1) appears
unnatural and unnecessary. By a completely different method, which uses tools fromrandom matrix theory, in a recent joint work with O’Rourke and Vu [2], we removed thethird term from (0.1) when the entries of E are independent and identically distributed(i.i.d.) copies of a standard normal random variable. In particular, we obtained that withhigh probabilitypr log(N + n) kEk
sin ∠(u1 , ue1 ) .
+.
δ1
σ1
My talk is concerned with the most recent progress in this direction of research. In [3],we have enhanced the r-dependence in the bounds obtained from our previous work [2]and have eased some technical assumptions. For instance, for the largest singular vector,
we have obtained thatpr + log(N + n) kEk
sin ∠(u1 , ue1 ) .
+.
δ1
σ1
Furthermore, we have presented a comprehensive extension of the Davis-Kahan-Wedinsin Θ theorem. This extension applies to any unitarily invariant norm and operates
under the assumption that E contains i.i.d. standard normal entries. Moreover, we havederived precise `∞ bounds for the perturbed singular vectors and the `2,∞ bounds forthe perturbed singular subspaces of A + E. Beyond these specific bounds, we have alsoestablished results pertaining to the generalized components - also known as linear andbilinear forms - of the perturbed singular vectors and singular subspaces. We further
investigate the `2,∞ bounds on the perturbed singular vectors, taking into account theweighting by their respective singular values. These fine-grained analysis is motivatedby the substantial impact and wide-ranging applications these analyses offer in statistics
and machine learning.
R EFERENCES[1] Sean O’Rourke, Van Vu and Ke Wang. Random perturbation of low rank matrices: improving classicalbounds. Linear Algebra Appl., 540:26–59, 2018.
[2] Sean O’Rourke, Van Vu and Ke Wang. Matrices with Gaussian noise: optimal estimates for singularsubspace perturbation. IEEE Transactions on Information Theory, 70(3):1978–2002, 2024.
[3] Ke Wang. Analysis of singular subspaces under random perturbations. arXiv:2403.09170.
D EPARTMENT OF M ATHEMATICS , H ONG KONG U NIVERSITY OF S CIENCE AND T ECHNOLOGY, H ONG KONG
E-mail address: kewang@ust.hk