Optimal rates of estimation for
multi-reference alignment
Afonso S. Bandeira∗ , Philippe Rigollet† , and Jonathan Weed‡

arXiv:1702.08546v1 [] 27 Feb 2017

Courant Institute of Mathematical Sciences, New York University
Massachusetts Institute of Technology
Massachusetts Institute of Technology

Abstract. This paper describes optimal rates of adaptive estimation of
a vector in the multi-reference alignment model, a problem with important applications in fields such as signal processing, image processing,
and computer vision, among others. We describe how this model can
be viewed as a multivariate Gaussian mixture model under the constraint that the centers belong to the orbit of a group. This enables us
to derive matching upper and lower bounds that feature an interesting
dependence on the signal-to-noise ratio of the model. Both upper and
lower bounds are articulated around a tight local control of KullbackLeibler divergences that showcases the central role of moment tensors
in this problem.
AMS 2000 subject classifications: Primary Statistics; secondary Invariant Theory, Signal Processing.
Key words and phrases: Multi-reference alignment, Orbit retrieval, Mixtures of Gaussians.
1. INTRODUCTION
The multi-reference alignment problem and its variants arise in various scientific and engineering applications such as structural biology [SVN+ 05, TS12,
Sad89], image recognition [Bro92], and signal processing [ZvdHGG03]. A striking
feature of this class of problems is that each observation is not only observed in a
noisy setting but is also altered by an latent transformation that reflects underlying heterogeneity of the data. The precise nature of this transformation depends
on the specific application, but it can often be characterized as the action of the
unknown element of a known group.
1.1 The multi-reference alignment problem
Consider the multi-reference alignment problem on IRd under cyclic shifts [BCSZ14].
The goal in this problem is to produce an estimate of a vector θ ∈ IRd (of∗

Part of this work was done while A. S. Bandeira was with the Mathematics Department at
MIT and supported by NSF Grant DMS-1317308.
†
This work was supported in part by NSF CAREER DMS-1541099, NSF DMS-1541100,
DARPA W911NF-16-1-0551, ONR N00014-17-1-2147 and a grant from the MIT NEC Corporation.
‡
This work was supported in part by NSF Graduate Research Fellowship DGE-1122374.

1

2

BANDEIRA, RIGOLLET AND WEED

ten thought of as a signal) from noisy cyclic shifted observations. More concretely, one observes n independent random vectors Y1 , . . . , Yn ∈ IRd given by
Yi = R`i θ + σξi , where `i is an unknown parameter (shift) in [d] := {1, . . . , d};
R`i is a latent cyclic shift by `i coordinates: the jth coordinate of R`i θ ∈ IRd is
given by R`i θ)j = θj+`i (mod d) ; and ξi ∼ N (0, Id ), i.i.d.
Note that R`i is a linear operator from IRd to IRd that can be represented by
a matrix. For example when d = 3, we have






0 1 0
0 0 1
1 0 0
R1 =  0 0 1  , R2 =  1 0 0  , and R3 =  0 1 0 
1 0 0
0 1 0
0 0 1
However, since cyclic shifts in IRd form a group isomorphic to quotient ring Zd :=
Z/dZ of integers modulo d, we write R`i ∈ Zd for simplicity.
The action of Zd on IRd also has a simple representation in the Fourier domain,
where the group acts on the phases of the Fourier coefficients. Denoting by θ̂ the
discrete Fourier transform of the vector θ, we obtain that the rotated vector R` θ
satisfies
(1.1)

2πj`i/d
d
R
θ̂j
` θj = e

for −bd/2c ≤ j ≤ dd/2e.

In practical applications, the signal θ arises as the discretization of some underlying continuous signal. Therefore, instead of focusing on the group Zd , we consider
instead the action of the circle group U (1) of unit-norm complex numbers, where
given z ∈ U (1) we define the operator Rz on IRd by its action on the Fourier
transform θ̂:
(1.2)

j
d
R
z θ j = z θ̂j

for −bd/2c ≤ j ≤ bd/2c.

We define the group of such operators by F and call Rz a fractional cyclic shift.
We assume throughout that d is odd, since for real signals not all fractional cyclic
shifts are well defined when j = d/2. This group is both slightly easier to analyze
and better corresponds to the situation in practice. Moreover, a comparison of
Equations 1.1 and 1.2 shows when d is large the groups Zd and F are essentially
equivalent. For the sake of exposition, in the sequel, we will focus on F and omit
the adjective “fractional” when referring to shifts.
Multi-reference alignment is directly used in structural biology [Dia92, TS12];
radar [ZvdHGG03]; crystalline simulations [SSK13]; and image registration in a
number of important contexts, such as in geology, medicine, and paleontology
[SSK13, FZB02]. As a result, variants of this problem for groups other than F
have received some attention, but rarely in statistics. We note parallel research
efforts that have investigated a Boolean version of this problem [APS17].
Another important related problem is that of molecule reconstruction in CryoElectron Microscopy (Cryo-EM). Cryo-EM is an important technique used to
determine three-dimensional structures of biological macromolecules. (It was considered the Nature method of the year in 2015 [Edi15, Nog15]). As in the multireference alignment problem described above, one of the main difficulties with
this imaging technique is that these molecules are imaged at different unknown
orientations and each molecule can only be imaged once due to the destructive nature of the process. More precisely, each measurement consists of a tomographic projection of a rotated (by an unknown rotation in SO(3)) copy of

OPTIMAL RATES FOR MRA

3

the molecule. The task is then to reconstruct the molecule density from many
such measurements. This reconstruction problem has also received significant attention, primarily from computational perspectives, but its statistical properties
remain largely unexplored.
We present in this paper the first statistical analysis for Euclidean MRA and
describe natural classes of signals whose complexity spans the spectrum of statistical rates achievable in this context. Although we focus on F or its subgroups
rather than SO(3), the two frameworks share many features. In particular, our
results exhibit non-trivial minimax rates that are inherent to this class of problems.
1.2 The Synchronization Approach
The difficulty of the multi-reference alignment problem resides in the fact that
both the signal θ ∈ IRd and the shifts z1 , . . . , zn ∈ F are unknown and the
latter are therefore latent variables. If the shifts were known, one could easily
estimate θ by taking the average of Rz−1
Yi , i = 1, . . . , n. In fact, this simple obi
servation is the basis of the leading current approach called the “synchronization
approach” [BCSZ14, BCS15] to this problem. Specifically, synchronization aims
at recovering the latent variables Rzi by solving a problem of the form
X
2
−1
(1.3)
min
Rz−1
.
0 Yi − Rz 0 Yj
0
0 ∈U (1)
z1 ,...,zn

1≤i,j≤n

i

j

Denoting by R̂zi the solutions of (1.3), one can then estimate θ by the average
of R̂z−1
Yi , i = 1, . . . , n.
i
Synchronization problems can be formulated as estimation problems on a
graph. More precisely, one can associate each observation Yi to a graph node,
each of which has a hidden label gi ∈ G for some group G of transformations.
(In our case, G ∼
= F.) The pairwise data, which we identify with edges of the
graph, reveals information about ratios gi (gj )−1 . In the the context of (1.3), this
−1
2
information is simply kR`−1
0 Yi −R`0 Yj k . Despite synchronization problems being
i

j

computationally hard in general, certain theoretical guarantees have been derived
under specific noise models that are unfortunately not realistic for the problems of
interest in this paper. For example, it is often assumed that the edge observations
are independent instead of the more relevant independence of vertices. Among
the most prominent methods are spectral methods [Sin11, BSS11], semidefinite
relaxations [BCSZ14, BCS15, ABBS14, BBS16, JMRT16, BBV16], and methods
based on Approximate Message Passing [PWBM16] and other modified power
methods [Bou16, CC16]. Synchronization also enjoyed many interesting connections with geometry (see [GBM16] for an example).
When the noise is smaller entrywise than the signal θ, synchronization approach yields acceptable results because macroscopic features underlying signal
are still visible. However, in applications, the noise level is often significantly
larger than signal, and it is on this regime that we focus our attention. An illustration of the difference between these regimes appears in Figure 1.
Unfortunately, for the high noise regimes we are interested here the shifts are
impossible to reliably estimate, regardless of the number of samples [ADBS16,
WW84]. In fact, it is not difficult to see that even if we had access to θ, one would
still not be able to reliably estimate the shifts in a high noise regime.

4

BANDEIRA, RIGOLLET AND WEED

2
1
0
1
2
2
1
0
1
2
2
1
0
1
2

15
10
5
0
5
10
15
15
10
5
0
5
10
15
15
10
5
0
5
10
15

Figure 1: Instances of the multi-reference alignment problem, at low (left column)
and high (right column) noise levels. The true underlying signal appears in gray,
and the noised version appears in red. When the noise level is low, large features
of the signal are still visible despite the noise; in the presence of large noise,
however, the signals cannot reliably be synchronized.

1.3 Notation
Denote by W ∈ C d×d the discrete Fourier transform matrix, with entries given
by
1
Wjk = √ exp(2πijk/d) ,
1 ≤ j, k ≤ d .
d
√
The normalization factor 1/ d is chosen so that the resulting matrix W is unitary.
Given θ ∈ IRd , let θ̂ = W θ be its Fourier transform. We will index θ̂ from
−bd/2c to bd/2c.
The symbol k · k denotes the `2 norm on IRd . For any positive integer d, we
write [d] = {1, . . . , d}.
Given a vector t, let t⊗k denote the order-k tensor formed by taking the kfold tensor product of t with itself. Denote by kAk the Hilbert-Schmidt norm of
a tensor A, defined by kAk2 = hA, Ai, where h·, ·i denotes the entrywise inner
product.
A tensor A is symmetric if Ai1 ...ik = Aiπ(1) ...iπ(k) for any permutation π of [k].
For such tensors, the value Ai1 ...ik depends only on the multiset {i1 , . . . , ik }.
We define a cyclic shift of θ ∈ IRd by a unit-norm complex number z by

OPTIMAL RATES FOR MRA

5

specifying its action on the Fourier transform of θ:
k
\
R
z (θ)k = z θ̂k for k = −bd/2c, . . . , bd/2c.

The group of cyclic shifts is known as the circle group and we denote it by F.
Note that it is isomorphic to several well known groups such as SO(2) and IR/Z.
Recall that the Kullback-Leibler (KL) divergence between two distributions P
and Q such that P  Q is given by
Z
dP 
D(P k Q) = log
dP .
dQ
It is well known that D(P k Q) ≥ 0, with equality holding iff P = Q almost surely.
1.4 Organization of the paper
In Section 2, we propose a new approach to MRA using Gaussian mixture models. Then, in Section 3 we present our main result, Theorem 1, providing minimax
rates for the multi-reference alignment problem under cyclic shifts. In Section 4
we draw precise connections between the Kullback-Leibler (KL) divergence of
sampled distributions from two different signals and the distance between their
invariant moment tensors. This relation does not depend on the precise nature
of the cyclic shifts and holds for any compact subgroup of the orthogonal group.
In Section 5 we give guarantees for the maximum likelihood estimator (MLE):
Theorem 4 gives a general guarantee for the MLE under conditions on the KL
divergence that were obtained in Section 4 and in Section 5.2 a modified MLE
is developed for the specific problems of shifts. This modified MLE provides the
upper bounds in Theorem 1. Section 6 concludes by establishing the lower bound
in Theorem 1; the proof involves finding pairs of different signals with several
matching invariant moment tensors.
2. GAUSSIAN MIXTURE MODEL
We propose an alternative to the synchronization approach discussed above
that completely bypasses the estimation of the shifts `1 , . . . , `n in favor of estimating θ directly. To do so, we recast our model as a continuous mixture of
Gaussians whose centers are algebraically constrained. Since the Gaussian distribution is invariant under cyclic shift, we can without loss of generality assume
that the shifts R`i are independent and uniformly distributed over F. Indeed, to
achieve this setup, we transform the observations Yi into RUi Yi , where Ui is uniformly distributed over U (1) and independent of all other random variables. Since
Ui + `i is also uniformly distributed over U (1), these new observations are drawn
from a mixture of Gaussians with uniform mixing weights and centers given by
R` θ, ` ∈ U (1). In particular, these centers are linked together by a rigid algebraic
structure: they are the orbit of θ ∈ IRd under the action of F. We summarize this
observation into our basic model.
Let Y1 , . . . , Yn ∈ IRd be n independent copies of Y ∈ IRd where
(2.1)

Y = Rθ + σξ,

where θ ∈ IRd is the unknown parameter of interest, R is drawn uniformly (i.e.,
according to the Haar probability measure) from F and ξ ∼ N (0, Id ) is independent Gaussian noise. Denote by Pθ the distribution of a random variable Y that
satisfies (2.1).

6

BANDEIRA, RIGOLLET AND WEED

Define the maximum likelihood estimator (MLE) θ̃n by
(2.2)

θ̃n = argmax
θ∈IRd

n
X
i=1



1
log IER exp − 2 kYi − Rθk2 .
2σ

We postpone the discussion of computational efficiency to a companion paper [PWB+ 17] and focus here on the statistical properties of this estimator. In
particular, we show that a small modification of it optimally solves the multireference alignment problem in a certain sense.
Although the goal of this problem is to recover the unknown parameter θ, it
is not hard to see that we can only hope to identify θ up to a global cyclic shift.
We therefore define a (pseudo-)metric on IRd that deems two vectors θ, φ ∈ IRd
to be close if they are close up to some element of F. To this end, define
ρ(θ, φ) = min kθ − Rφk .
R∈F

We assume throughout that the noise variance σ 2 is known. This assumption is
realistic in many applications such as imaging or signal processing. In other circumstances, it may be calibrated using cross-validation, for example. Throughout
this paper, we assume that c ≤ kθk ≤ 1, where c > 0 is a universal constant so
that σ captures entirely the (inverse) signal-to-noise ratio. We note that, comparatively, IEkσξi k2 = σ 2 d.
Gaussian mixture models have been extensively studied in the statistical literature since their introduction by Pearson [Pea94] in the nineteenth century
(see e.g. [MP00] for an overview). As illustrated by the extant literature, mixture
models are quite rich and broadly applicable to a variety of statistical problems
ranging from clustering to density estimation. As a result, various statistical perspectives may apply when studying the performance of estimators. In the context
of multi-reference alignment, we are naturally interested in estimating θ or, equivalently, the centers. It has been established that in general, the rate of estimation
of the centers may scale like n−C/d for d centers that can be arbitrarily close
(see for example [MV10] and more recently [HK15] for an interesting explanation
from the point of view of model misspecification). This curse of dimensionality
arises from the minimax point of view where centers may be arbitrarily close,
with a distance that may depend on the number of observations, thus leading to
nonparametric rates of estimation. Instead, when the centers are well separated,
the general conditions of [Wal49] are satisfied so that the MLE converges at the
√
parametric rate 1/ n . In our case, when the mixture is continuous, Assumption 1, below, plays an analogous role in ensuring convergence at the parametric
rate. When such conditions are met, the Fisher information determines how the
statistical performance of the MLE scales with the noise level σ. This question
is central to signal processing problems such as multi-reference alignment, where
σ is quite large, since it determines the order of magnitude of the sample size n
required to achieve a certain accuracy.
Exact computation of the Fisher information matrix in this model is out of
√
reach and instead, we focus on the scaling of the quantity nρ(θ̃n , θ) with the
signal-to-noise ratio of the problem, where θ̃n is the MLE (2.2).

OPTIMAL RATES FOR MRA

7

3. MAIN RESULTS
√

As mentioned above, nρ(θ̃, θ) depends asymptotically on the Fisher information of the model, which can be related to the curvature of the KullbackLeibler divergence around its minimum. Conversely, (lack of) curvature of the
Kullback-Leibler divergence around its minimum is what controls minimax lower
bounds. To address both problems at once, we follow an idea originally introduced in [LNS99] in the context of functional estimation and further developed
by [CL11, WY16]. In the multi-reference alignment model, this approach allows
us to relate Kullback-Leibler divergence to moment tensors, which can in turn
be controlled using Fourier-theoretic arguments. It will follow from this analysis
that the difficulty of estimating a particular signal θ depends on the support of
the Fourier transform θ̂ of θ. In particular, define the positive support psupp(θ̂)
of θ̂ by
psupp(θ̂) = {j | j ∈ {1, . . . , d/2}, θ̂j 6= 0} .
We make the following assumption, which guarantees that the MLE converges
at a parametric rate.
Assumption 1. There exists an absolute constant c > 1, not depending on
n, such that c−1 ≤ |θ̂j | ≤ c for all j ∈ psupp(θ̂). We denote by T the set of such
vectors.
We emphasize that this is the situation of most interest to practitioners: the
existence of very small, but non-zero, coordinates whose values approaches 0
with n should rightly be considered pathological. In a way, Assumption 1 rules
out certain artificial difficult situations analogous to classical difficulties arising in
estimating mixtures of Gaussians, such as distinguishing the mixture .5N (+ε, 1)+
.5N (−ε, 1) from the single Gaussian N (0, 1) for very small ε. Moreover, it enables
us to restrict our attention to a compact set of parameters and bypass some
technical complications.
The following theorem reveals a surprising phenomenon: even under Assumption 1, the multi-reference alignment problem suffers from the curse of dimensionality.
Theorem 1. Let 2 ≤ s ≤ bd/2c. Let Ts be the set of vectors θ ∈ T satisfying
Assumption 1 and psupp(θ̂) ⊂ [s]. Then,
(3.1)

σ 2s−1
inf sup IEθ [ρ(Tn , θ)]  √ (1 + on (1)) ,
Tn θ∈Ts
n

where the infimum is taken over all estimators Tn of θ and where the symbol 
hides constants depending on d but on no other parameter. A modified MLE θ̌n
defined in (5.11) achieves this rate.
A few remarks are in order. Note first that the curse of dimensionality is inherent to the minimax paradigm. Indeed, our proof of the lower bound describes
a class of signals that satisfy Assumption 1 but have a very specific Fourier spec√
trum. Such signals drive the worst case bound of order σ 2s−1 / n. This limitation
is overcome in a companion paper [PWB+ 17], where we show that even Fourier
dense signals θ may be estimated as the same rate as signals in T2 as long as they

8

BANDEIRA, RIGOLLET AND WEED

are generic enough. Second, our proof techniques do not allow us to remove the σ
dependence of the term on (1). In particular, for small values of n, this term may
actually dominate the upper bound. We conjecture that this issue is an artifact
of the proof technique and note that preliminary numerical results in [PWB+ 17]
support this claim.
The rest of this paper is devoted to the proof of the main results in Theorem 1.
4. INFORMATION GEOMETRY
In this section, we develop several new tools to obtain precise bounds on the
divergence D(Pθ k Pφ ) for pairs of signals θ and φ. On the one hand, it is well
known that upper bounds on the divergence translate into minimax lower bounds
using LeCam’s method [LeC73]. On the other hand, we also show how to transform lower bounds on D(Pθ k Pφ ) into uniform upper bounds on the performance
of the MLE. Note that this analysis departs from the classical pointwise rate of
convergence for MLE that guarantees a rate of convergence n−1/2 for each fixed
choice of parameter as n → ∞. Our tools strengthen this result considerably.
Indeed, we show that for reasonable choices of θ, the MLE achieves a rate of
n1/2 uniformly over all choices of θ. We refer the reader to [HK15] for examples
of other Gaussian mixture problems where the pointwise and uniform rates of
estimation differ.
For convenience, we abbreviate D(Pθ k Pφ ) by D(θ k φ). The following Lemma
collects several useful facts about the function D(θ k φ). A proof appears in Appendix A.
Lemma 2. Let R have uniform distribution over any subgroup G of the orthogonal group in d dimensions. Fix θ, φ ∈ IRd . The following holds
(i) If ϑ = θ − IERθ and ϕ = φ − IERφ, then
D(θ k φ) = D(ϑ k ϕ) +

1
kIE[Rθ] − IE[Rφ]k2 .
2σ 2

(ii) If ρ(θ, φ) = ε and kθk, kφk ≤ 1, then
D(θ k φ) =

1
1
O(ε2 )
2
⊗2
⊗2 2
kIE[Rθ]
−
IE[Rφ]k
+
kIE[(Rθ)
]
−
IE[(Rφ)
]k
+
,
2σ 2
4σ 4
σ6

where O(ε2 ) hides a constant, which may depend on d but is otherwise independent
of θ, φ, and σ.
The following Theorem provides a tight bound on the quantity D(θ k φ) in
terms of Hilbert-Schmidt distance between the moment tensors IE[(Rθ)⊗k ] and
IE[(Rφ)⊗k ]. Specifically, we show that the divergence D(θ k φ) is of order σ −2k ε2
where k is the smallest natural number such that the moment tensors IE[(Rθ)⊗k ]
and IE[(Rφ)⊗k ] differ significantly. This theorem is not specific to cyclic shifts
and holds for any compact subgroup of the orthogonal group.
Theorem 3. Let θ be a fixed vector in IRd such that kθk ≤ 1. Let φ ∈ IRd
be such that ρ(θ, φ) = ε ≤ kθk. Let R be a random element drawn according to
the Haar probability measure on any compact subgroup G of the orthogonal group

9

OPTIMAL RATES FOR MRA

in d dimensions. For all m ≥ 1, let ∆m = IE[(Rθ)⊗m − (Rφ)⊗m ]. If there exists
k ≥ 1 such that, as ε → 0,
k∆m k = o(ε) for m = 1, . . . , k − 1,

k∆k k = Ω(ε) ,

and

then k∆k k = Θ(ε). Moreover, for σ ≥ 1 there exist universal constants c and C̄
and constant C d that depends only on d, all positive and such that
ck
σ 2k k!

k∆k k2 − C d

ε2
σ 2k+2

≤ D(θ k φ) ≤

2
σ 2k k!

k∆k k2 + C̄

ε2
σ 2k+2

.

In particular, there exists positive σ0 , ε0 that depend on d such that for all
σ ≥ σ0 , and θ, φ such that kθk ≤ 1, ρ(θ, φ) ≤ ε0 , it holds
D(θ k φ)  σ −2k ρ2 (θ, φ) ,
where the symbol  hides constants depending on d but on no other parameters.
Proof. The divergence D(θ k φ) and the tensors IE[(Rφ)⊗m ] are unaffected
if we replace φ by Rφ for any R ∈ G. Hence without loss of generality, we can
assume that kθ − φk = ρ(θ, φ) = ε.
We first establish that for any m ≥ 1, the quantities k∆m k are indeed of order
at most ε. Specifically, we prove that
(4.1)

kIE[(Rθ)⊗m − (Rφ)⊗m ]k2 ≤ m2 2m−1 ε2

for all m ≥ 1 and ε ∈ (0, 1). This implies that k∆k k = Θ(ε).
To prove (4.1), note that by Jensen’s inequality,
kIE[(Rθ)⊗m − (Rφ)⊗m ]k2 ≤ IEk(Rθ)⊗m − (Rφ)⊗m k2 = kθ⊗m − φ⊗m k2 .
Expanding the norm yields
kθ⊗m − φ⊗m k2 = kθk2m − 2hθ, φim + kφk2m
= kθk2m (1 − 2(1 + γ)m + (1 + 2γ + δ 2 )m ) ,
where δ = ε/kθk and γ = hθ, φ − θi/kθk2 is such that |γ| ≤ δ by Cauchy-Schwarz.
Using Lemma B.1 in the Supplementary Materials, we conclude that
kθk2m − 2hθ, φim + kφk2m ≤ kθk2m m2 2m−1 δ 2 ≤ m2 2m−1 ε2 ,
as desired.
If k < 3, then the claim follows directly from Lemma 2 (ii).
We therefore assume k ≥ 3. Then kIE[Rθ−Rφ]k = o(ε), so Lemma 2 (i) implies
D(θ k φ) = o(ε2 ) + D(ϑ k ϕ) ,
where ϑ = θ − IERθ and ϕ = φ − IERφ. Hence we can replace θ and φ by ϑ and ϕ
without affecting D(θ k φ) by more than o(ε2 ). We therefore assume without loss
of generality that IERθ = IERφ = 0.

10

BANDEIRA, RIGOLLET AND WEED

We first show the upper bound. Denote by g the density of a standard ddimensional Gaussian random variable. For all ζ ∈ IRd , let fζ denote the density
of Pζ . Then
kζk2
y > Rζ
1
1
fζ (y) = IER g(σ −1 (y − Rζ)) = g(σ −1 y)e− 2σ2 IER e σ2 .
σ
σ

Let χ2 (θ, φ) denote the χ2 -divergence between Pθ and Pφ , defined by
Z

2

χ (θ, φ) =

(fθ (y) − fφ (y))2
dy .
fθ (y)

Since IERθ = 0 by assumption, Jensen’s inequality implies
fθ (y) ≥

kθk2
kθk2
y > Rθ
1
1
g(σ −1 y)e− 2σ2 eIER σ2 = g(σ −1 y)e− 2σ2 .
σ
σ

Hence
kθk2
kθk2
kφk2
y > Rθ
y > Rφ 2 1

(fθ (y) − fφ (y))2
≤ e 2σ2 e− 2σ2 IER e σ2 − e− 2σ2 IER e σ2
g(σ −1 y) .
fθ (y)
σ

Integrating this quantity with respect to y yields a bound on the χ2 divergence.
Let ξ ∼ N (0, Id ) and observe that
i
h kθk2
kθk2
kφk2
σξ> Rθ
σξ> Rφ 2
χ2 (θ, φ) ≤ IEξ e 2σ2 e− 2σ2 IER e σ2 − e− 2σ2 IER e σ2
h kθk2
kθk2 ξ> (R+R0 )θ
kθk2 +kφk2 ξ> (Rθ+R0 φ)
σ
σ
= IEξ,R,R0 e 2σ2 e− σ2 e
− 2e− 2σ2 e
kφk2 ξ> (R+R0 )φ i
σ
+ e− σ2 e
,
where R and R0 are independent elements selected uniformly from G.
Taking expectations with respect to ξ yields
θ > Rφ
φ> Rφ 
 θ> Rθ
χ2 (θ, φ) ≤ 2IER e σ2 − 2e σ2 + e σ2 ,
2

2

Where we used that ekθk /(2σ ) ≤ 2 for σ ≥ 1 and kθk2 ≤ 1.
The random variables Rθ and Rφ have moment generating functions that
converge in a neighborhood of the origin, hence
h
i
X
2
>
m
>
m
>
m
χ2 (θ, φ) ≤
IE
(θ
Rθ)
−
2(θ
Rφ)
+
(φ
Rφ)
R
σ 2m m!
m≥0
X
2
k∆m k2
=
σ 2m m!
m≥0

≤
≤

2
σ 2k k!

k∆k k2 + ε2

X
m≥k+1
ε2

2m2
+ o(ε2 )
σ 2m m!

2
k∆k k2 + C̄ 2k+2 .
σ 2k k!
σ

The upper bound then follows from the inequality D(θ k φ) ≤ χ2 (θ, φ) [Tsy09].

OPTIMAL RATES FOR MRA

11

We now turn to the lower bound. As before, we assume that k ≥ 3 and define
Z
d1 (θ, φ) = |fθ (y) − fφ (y)| dy .
Recall that by Pinsker’s Inequality, we have D(θ k φ) ≥ 12 d21 (θ, φ). Moreover,
Z
ky−Rθk2
ky−Rφk2
1
−
−
2σ 2
2σ 2
d1 (θ, φ) =
IEe
−
IEe
dy
(2π)d/2 σ
 2σξ> Rφ−kφk2 
 2σξ> Rθ−kθk2 
2σ 2
2σ 2
ξ ∼ N (0, Id )
ξ − IER e
ξ ,
= IEξ IER e
2
2
2
kφk
 ξ> Rθ  kφk −kθk
 ξ> Rφ 
= e− 2σ2 IEξ IER e σ ξ e 2σ2
− IER e σ ξ .
We now show that |kθk2 − kφk2 | = o(ε). Indeed
kθk2 − kφk2 = IE

d
X


(Rθ)2i − (Rφ)2i



i=1

≤

≤

d
X
i=1
d
X

IE[(Rθ)2i ] − IE[(Rφ)2i ]
IE[(Rθ)i (Rθ)j ] − IE[(Rφ)i (Rφ)j ]

i,j=1

≤ dkIE(Rθ)⊗2 − (Rφ)⊗2 k = dk∆2 k = o(ε) ,
where we applied Cauchy-Schwarz to get the second inequality and last equality
follows from the assumption of the theorem when k ≥ 3. Since kθk ≤ 1 and σ ≥ 1,
we have
kφk2 −kθk2
kφk2
= 1 + o(ε) and
≤1
e 2σ2
2σ 2
for ε small enough. It yields

 ξ> Rθ 
 ξ> Rφ 
d1 (θ, φ) ≥ e−1 IEξ IER e σ ξ 1 + o(ε) − IER e σ ξ
kθk2
 ξ> Rθ 
 ξ> Rφ 
≥ e−1 IEξ IER e σ ξ − IER e σ ξ − e 2σ2 −1 o(ε)
X 1
= e−1 IEξ
h∆m , ξ ⊗m i − o(ε)
σ m m!
m≥1
X
1
1
≥ e−1 IEξ k h∆k , ξ ⊗k i − e−1
IEξ m h∆m , ξ ⊗m i − o(ε)
σ m!
σ k!
m6=k
p
X (d + m)m
ck
≥ √ k∆k k − e−1
k∆m k − o(ε)
σ m m!
σ k k!
m6=k

≥

ck
σk

ε
√ k∆k k − Cd k+1
σ
k!

where in the penultimate inequality, we employ a Khinchine-type inequality due
to [Bob00], details of which appear as Lemma B.2 in the Supplementary Materials.
Together with Pinsker’s inequality and the fact that k∆k k = Θ(ε) where ε =
ρ(θ, φ), this completes the proof of the lower bound.

12

BANDEIRA, RIGOLLET AND WEED

5. MAXIMUM LIKELIHOOD ESTIMATION
Let Y1 , . . . , Yn be i.i.d observations from the model (2.1) and consider the MLE
θ̃n that was defined in (2.2). In this section, we prove our main statistical result,
that is a uniform upper bound on the rate of convergence of the MLE.
Our proof technique extends beyond the framework of the MRA model and
can be broadly applied to derive uniform rates of convergence for the MLE from
tight bounds on the KL divergence like Theorem 3. While similar ideas are often
employed to obtain pointwise rates of convergence, extension to uniform rates
requires novel elements. From here on, positive constants may depend on d unless
noted otherwise.
We first establish an upper bound for the MLE under a general lower bound for
the KL divergence. Then we specialize this result to obtain the minimax upper
bounds over Ts that are presented in Theorem 1.
5.1 A general upper bound
Theorem 4. Let R be a random element drawn according to the Haar probability measure on any compact subgroup G of the orthogonal group in d dimensions. Assume that there exist k ≥ 1 and positive σ0 , ε0 , c0 that depend on d such
that the following hold: for all σ ≥ σ0 , and θ, φ such that θ ∈ T ,
(5.1)

D(θ k φ) ≥ Cσ −2k ρ2 (θ, φ)

∀ φ ∈ IRd ρ(θ, φ) ∈ [0, ε0 ) ,

(5.2)

D(θ k φ) ≥ Cσ −`

∀ φ ∈ IRd ρ(θ, φ) ∈ [ε0 , c0 σ) ,

(5.3)

D(θ k φ) ≥ Cσ −2k ρ2 (θ, φ)

∀ φ ∈ IRd ρ(θ, φ) ∈ [c0 σ, ∞) .

Then the MLE θ̃n satisfies uniformly over such θ,
log n
σk
,
IEθ [ρ(θ̃n , θ)] ≤ C √ + Cσ
n
n

(5.4)
where Cσ ≤ Cσ 2`+10k .

Proof. The symbols c and C denote constants whose value may change from
line to line. In the rest of this proof, we write θ̃ = θ̃n to denote the MLE.
The main goal of the proof is to combine control of the curvature of the function
D with control of the deviations of the log-likelihood function.
Define the event E = {ρ(θ̃, θ) ≤ ε} where ε ≤ ε0 . Since D is invariant under
the action of G, we can assume without loss of generality that ρ(θ̃, θ) = kθ̃ − θk.
We first establish that on this event kθ̃ − θk can be controlled in terms of the
metric induced by the Hessian of D at θ.
Fix θ ∈ IRd and denote by H the Hessian of √
the function φ 7→ D(θ k φ) evalud
ated at φ = θ. For any u ∈ IR , define kukH = u> Hu.
On E, D(θ k θ̃) ≥ Cσ −2k ρ2 (θ, θ̃) and ρ(θ, θ̃) = kθ̃ − θk so using a third order
Taylor expansion, we get
1
kθ̃ − θk3
D(θ k θ̃) − kθ̃ − θk2H ≤ C
≤ Cεσ 2k−3 kθ̃ − θk2H .
2
σ3
Therefore, there exists a C such that if ε < Cσ 3−2k , we get
(5.5)

1
D(θ k θ̃) ≥ kθ̃ − θk2H ≥ Cσ −2k kθ̃ − θk2 .
4

OPTIMAL RATES FOR MRA

13

We now control the geometry of the log-likelihood function near θ. Define
n

Dn (θ k φ) =

fθ
1X
log (Yi ) ,
n
fφ
i=1

where Yi are i.i.d from Pθ and we recall that fζ is the density of Pζ , ζ ∈ IRd . Note
that Dn (θ k θ) = 0, so Dn (θ k θ̃) ≤ 0.
Since θ is held fixed throughout the proof, we abbreviate D(θ k φ) and Dn (θ k φ)
as D(φ) and Dn (φ), respectively.
Using a second order Taylor expansion, we get
1
D(θ̃) − Dn (θ̃) = −∇Dn (θ)> h + h> ∇2 (D − Dn )(θ̄)h ,
2
where h = θ̃ − θ and θ̄ lies on a segment between θ̃ and θ.
For all ζ ∈ T , write Hn (ζ) for the Hessian of Dn (φ) evaluated at φ = ζ, and
similarly let H(ζ) be the Hessian of D(φ) evaluated at φ = ζ. We obtain
(5.6)

1
1
khk2H ≤ D(θ̃) − Dn (θ̃) ≤ −∇Dn (θ)> h + h> (H(θ̄) − Hn (θ̄))h .
4
2

For the first term, we note that if H were invertible (and hence k · kH a genuine
metric), then it is well known (see, e.g. [HUL01]) that
sup

−∇Dn (θ)> u = k∇Dn (θ)k∗H = k∇Dn (θ)kH −1 ,

u:kukH =1

where k · k∗H , the dual norm to k · kH , is such that k · k∗H = k · kH −1 . In general,
H is not invertible, but we still have k · k∗H = k · kH † , where H † denotes the
Moore-Penrose pseudo-inverse of the matrix H.
To control the second, note first that by (5.5), it holds khk ≤ Cσ k khkH . Therefore,
h> (H(θ̄) − Hn (θ̄))h ≤ Cσ k khkH khk sup kH(φ) − Hn (φ)kop ,
φ∈Bε

IRd

where Bε := {φ ∈
: ρ(φ, θ) ≤ ε}.
Combining the above three bounds and dividing by khkH , we get that on E,
1
khkH ≤ kDn (θ)kH † + Cσ 2k sup kH(φ) − Hn (φ)k2op + khk2 ,
4
φ∈Bε
where we applied Young’s inequality.
On the other hand, using (5.5), we have khkH ≥ Cσ −k ρ(θ̃, θ). Therefore, applying the Cauchy-Schwarz inequality and Chebyshev’s inequality, we get
IE[ρ(θ̃, θ)] = IE[ρ(θ̃, θ)1IE ] + IE[ρ(θ̃, θ)1IE c ]
1/2
1/2
≤ Cσ k IE[khkH 1IE ] + IE[ρ2 (θ̃, θ)]
IP(E c )
≤ Cσ k IE[khkH 1IE ] + ε−1 IE[ρ2 (θ̃, θ)] .
Thus,
(5.7)


IE[ρ(θ̃, θ)] ≤ C σ k IEkDn (θ)kH † +σ 3k IE sup kH(φ)−Hn (φ)k2op +(σ k +ε−1 )IE[ρ2 (θ̃, θ)]
φ∈Bε

14

BANDEIRA, RIGOLLET AND WEED

It suffices to control the right side of the above inequality. The main term is
the first one: Jensen’s inequality and the second Bartlett identity imply
r
r
q
1
d
†
>
†
tr(H H) ≤
.
(5.8)
IEkDn (θ)kH † ≤ tr(H IE[Dn (θ)Dn (θ) ]) =
n
n
Standard matrix concentration bounds can be applied to show
IE sup kH(φ) − Hn (φ)k2op ≤ C

(5.9)

φ∈Bε

log n
.
nσ 4

A proof of (5.9) appears as Lemma B.4 in the Supplementary Materials.
Likewise, a standard slicing argument, Lemma B.5 in the Supplementary Materials, implies
IE[ρ2 (θ̃, θ)] ≤ C

(5.10)

σ (8k−4)∧(2`+2)
.
n

Plugging (5.8), (5.9), and (5.10) into (5.7), we get
 σk
σ 5k−4 log n σ (9k−4)∧(2`+k+2) 
IE[ρ(θ̃, θ)] ≤ C √ +
+
,
n
nε
n
as desired.
5.2 A modified MLE
The MLE itself may not achieve the optimal rate of convergence because the
lower bound (5.1) may not be satisfied with the optimal choice of k when φ
is a specific perturbation of θ. Namely, in the specific case of cyclic shifts, the
divergence D(φ) is not curved enough in directions that perturb a null Fourier
coefficient of θ. To overcome this limitation, we split the sample Y1 , . . . , Yn into
two parts: with the first part we estimate the support of θ̂ under Assumption 1 and
with the second part, we compute a maximum likelihood estimator constrained
to have the estimated support.
Specifically, assume for simplicity that we have a sample Y1 , . . . , Y2n of size
2n and split it into two samples Y1 = {Y1 , . . . , Yn } and Y2 = {Yn+1 , . . . , Y2n } of
equal size.
5.2.1 Fourier support estimation We use the first subsample Y1 to construct a
set S̃ that coincides with psupp(θ̂) with high probability. For any j = −bd/2c, . . . , bd/2c,
define,
n
1X d 2
Mj =
|(Yi )j | − σ 2 .
n
i=1

Define the set S̃ by
r
n
log n o
2
S̃ = j ∈ {−bd/2c, . . . , bd/2c} : Mj ≥ Cσ
n
The following proposition shows that S̃ = psupp(θ̂) with high probability.

OPTIMAL RATES FOR MRA

Proposition 5.

15

Fix θ ∈ IRd . Assume that n/(log n) ≥ Cσ 4 . Then
IP[S̃ 6= psupp(θ̂)] ≤ C

σ4
.
n

Proof. This follows from standard concentration arguments. A full proof appears in the Supplementary Materials.
5.2.2 Constrained MLE We use the second sample to construct a constrained
MLE. To that end, for any S ⊂ {1, . . . , bd/2c}, define the projection PS by

 φ̂j if j ∈ S ∪ −S
\
P
(φ)
=
φ̂
if j = 0
S
j
 0
0
otherwise.
The image of PS is a (2|S| + 1)-dimensional real vector space. For convenience,
write φS = PS φ for any vector φ ∈ IRd .
Recall that Y2 = {Yn+1 , . . . , Y2n } denotes the second subsample and define
Y1S , . . . , YnS by
YiS = PS Yn+i + σ(I − PS )ξi

i = 1, . . . , n .

where the ξi ∼ N (0, Id ) are i.i.d, independent from Y2 .
Define the modified MLE θ̌n by
n

1X
θ̌n = argmax
log fφ (YiS̃ ) .
φ∈Im(PS̃ ) n

(5.11)

i=1

Proposition 6. Fix 2 ≤ s ≤ bd/2c, θ ∈ Ts and let S = psupp(θ̂). Then,
there exist positive σ0 , ε0 , c0 that depend on d such that the following hold: for all
σ ≥ σ0 ,
(5.12)

D(θ k φ) ≥ Cσ −4s+2 ρ2 (θ, φ)
−4s+2

(5.13)

D(θ k φ) ≥ Cσ

(5.14)

D(θ k φ) ≥ Cσ −2 ρ2 (θ, φ)

∀ φ ∈ Im(PS ) ρ(θ, φ) ∈ [0, ε0 ) ,
∀ φ ∈ Im(PS ) ρ(θ, φ) ∈ [ε0 , c0 σ) ,
∀ φ ∈ Im(PS ) ρ(θ, φ) ∈ [c0 , ∞) .

A proof appears in Appendix A.
5.3 Proof of upper bound in Theorem 1
Define R = {S̃ = psupp(θ̂)} and observe that
IE[ρ(θ̌n , θ)] ≤ IE[ρ(θ̌n , θ)1IR ] + IE[ρ(θ̌n , θ)1IRc ] .
The first term is controlled by combining Proposition 6 and Theorem 4 to get
σ 2s−1
log n
IE[ρ(θ̌n , θ)1IR ] ≤ C √ + Cσ
,
n
n
where Cσ ≤ Cσ 28s−14 .

16

BANDEIRA, RIGOLLET AND WEED

To bound the second term, we use the Cauchy-Schwarz inequality and Proposition 5 to get
q
σ2
IE[ρ(θ̌n , θ)1IRc ] ≤ C √
IE[ρ(θ̌n , θ)2 ]
n
We now show that IE[ρ(θ̌n , θ)2 ] is bounded uniformly over all choices of S̃ by a
constant multiple of σ 2 using a similar slicing argument as the one employed in
the proof of Lemma B.5 in the Supplementary Materials.
By the triangle inequality,
ρ(θ̌n , θ) ≤ ρ(θ̌n , θS̃ ) + ρ(θS̃ , θ) ≤ ρ(θ̌n , θS̃ ) + 1 .
In view of (5.14), we have
ρ(θ̌n , θS̃ )2 ≤ (c0 σ)2 + Cσ 2 D(θS̃ k θ̌n ) ≤ (C ◦ σ)2 (1 + Gn (θS̃ k θ̌n )) ,
for some constant C ◦ , where
Gn (θS̃ k θ̌n ) = D(θS̃ k θ̌n ) − Dn (θS̃ k θ̌n ) .
For j ≥ 0, define Sj = {φ ∈ IRd : 2j σ ≤ ρ(φ, θS̃ ) ≤ 2j+1 σ} and let J be such
that C ◦ ≤ 2J ≤ 2C ◦ . Observe that
X
IE[ρ(θ̌n , θS̃ )2 ] ≤ 4(C ◦ σ)2 +
IE[ρ(θ̌n , θS̃ )2 1I(θ̌n ∈ Sj )]
j≥J

≤ 4(C ◦ σ)2 + σ 2

X

22j+2 IP[ sup Gn (θS̃ k φ) > C

j>J

≤ 4(C ◦ σ)2 + Cσ 2

X

φ∈Sj

22j
]
n

22j exp(−C22j ) ≤ Cσ 2 ,

j≥0

where we used (B.7) from the Supplementary Materials in the third inequality.
We obtain
IE[ρ(θ̌n , θ)2 ] ≤ 2IE[ρ(θ̌n , θS̃ )2 ] + 2 ≤ Cσ 2 .
We have established that
IE[ρ(θ̌n , θ)] ≤ C

 σ 2s−1
log n σ 3 
√ + σ 26s−13
+
,
n
n
n

which completes the proof of Theorem 1.
6. MINIMAX LOWER BOUNDS
Our minimax lower bounds rely ultimately on Le Cam’s classical two-point
testing method [LeC73]. In particular, the version that we use requires an upper
bound on the KL divergence, which can be obtained using Theorem 3 and a
moment matching argument.

17

OPTIMAL RATES FOR MRA

6.1 Moment matching
The lower bound of Theorem 1 follows from Proposition 7.
Proposition 7.

Fix 2 ≤ s ≤ bd/2c and let θ, φ ∈ IRd satisfy
θ̂m = φ̂m = 0

for m ∈
/ {±(s − 1), ±s}

and
|θ̂m | = |φ̂m |

for m ∈ {±(s − 1), ±s}.

If R is drawn uniformly from F, then for any m = 1, . . . , 2s − 2, it holds
IER [(Rθ)⊗m ] = IER [(Rφ)⊗m ]
Proof. Fix m ≤ 2s − 2. Since IER [(Rθ)⊗m ] and IER [(Rφ)⊗m ] are symmetric
tensors, to show that they are equal it suffices to show that
hIER [(Rθ)⊗m ], u⊗m i = hIER [(Rφ)⊗m ], u⊗m i

∀ u ∈ IRd

or equivalently, that
(6.1)

IER [(u> Rθ)m ] = IER [(u> Rφ)m ]

∀ u ∈ IRd .

Consider the set P = {ζ : |ζ̂j | = |θ̂j | , ∀ j} and note that θ, φ ∈ P. We show
that the function ζ 7→ IER [(u> Rζ)m ] is constant on P, which readily yields (6.1).
For a fixed shift Rz ∈ F, we obtain
bd/2c

X

d
u> Rz ζ = hû, R
z ζi =

û−j ζ̂j z j ,

j=−bd/2c

so

bd/2c

X

(u> Rz ζ)m =

z j1 +···+jm

m
Y

û−jn ζ̂jn .

n=1

j1 ,...,jm =−bd/2c

Taking expectations with respect to a uniform choice of z yields
(6.2)

IER [(u> Rζ)m ] =

X

m
Y

û−jn ζ̂jn ,

j1 +···+jm =0 n=1

where the sums are over all choices of coordinates j1 , . . . , jm ∈ {−bd/2c, . . . , bd/2c}
whose sum is 0.
The Fourier transform of ζ is supported only on coordinates ±(s − 1) and
±s, so we may restrict our attention to sums involving only those coordinates.
Suppose j1 + · · · + jm = 0. Define
α = |{i : ji = s − 1}|

β = |{i : ji = −(s − 1)}|

γ = |{i : ji = s}|

δ = |{i : ji = −s}|

By assumption j1 + · · · + jm = 0, so the tuple (α, β, γ, δ) is a solution to
α(s − 1) + β(−(s − 1)) + γ(s) + δ(−s) = 0

18

BANDEIRA, RIGOLLET AND WEED

or, equivalently,
(α − β)(s − 1) + (γ − δ)s = 0
Since s − 1 and s are coprime, (α − β) and (γ − δ) must be multiples of s and
s − 1, respectively. Since |α − β| + |γ − δ| ≤ m < 2s − 1, in fact α − β = γ − δ = 0.
Therefore the only m-tuples (j1 , . . . , jm ) that appear in the sum on the righthand side of (6.2) are those in which +(s − 1) and −(s − 1) occur an equal
number of times and
Q +s and −s occur an equal number of times. For such mtuples, the product m
n=1 û−jn ζ̂jn can be reduced to a product of terms of the form
û−(s−1) ûs−1 ζ̂s−1 ζ̂−(s−1) and û−s ûs ζ̂s ζ̂−s . Since u and ζ are real vectors, ûj û−j =
|uj |2 and ζ̂j ζ̂−j = |ζ̂j |2 for all j = −bd/2c, . . . , bd/2c, so
m
Y

û−jn ζ̂jn = (|ûs−1 |2 |ζ̂s−1 |2 )a (|ûs |2 |ζ̂s |2 )b ,

n=1

where a and b are the number of occurrences of the pairs ±(s − 1) and ±s,
respectively. This quantity depends only on the moduli |ζ̂s | and |ζ̂s−1 |, hence it is
the same for all ζ ∈ P. This completes the proof of (6.1) and therefore the proof
of the proposition.
6.2 Proof of lower bound in Theorem 1
√
Fix z = eiδ for δ = c1 σ 2s−1 / n for some constant c1 > 0. Let τ be given by

1/2 if |j| = s − 1 or s,
τ̂j =
0
otherwise.
Let φn be given by

1/2 if |j| = s − 1,



z/2
if j = s,
cn =
φ
j
∗
z /2 if j = −s,



0
otherwise.
cn and τ̂ lies in [−s, s] and that kφn − τ k and ρ(φn , τ )
Note that the support of φ
√
are both bounded by c2 σ 2s−1 / n, where c2 can be made arbitrarily small by
taking c1 small enough.
Theorem 3 and Proposition 7 imply that
D(Pnτ k Pnφn ) ≤ Cnσ −4s+2 ρ2 (φn , τ ) ≤ c3 ,
for a positive constant c3 that can be made arbitrarily small by taking c1 small
enough. Using standard minimax lower bounds techniques [Tsy09], we get the
desired result.
APPENDIX A: OMITTED PROOFS
Proof of Lemma 2. We first prove the following simple expression:
1

>

1
IE[e σ2 (θ+σξ) Rθ | ξ]
2
2
,
D(Pθ k Pφ ) = D(θ k φ) = 2 (kφk − kθk ) + IE log
1
>
2σ
IE[e σ2 (θ+σξ) Rφ | ξ]

19

OPTIMAL RATES FOR MRA

where ξ ∼ N (0, Id ).
This claim follows directly from the definition of divergence. Denoting by g(y)
the density of a standard Gaussian random variable with respect to the Lebesgue
measure on IRd , we can write
dPθ
IE[g((y − Rθ)/σ)]
(y) =
dPφ
IE[g((y − Rφ)/σ)]


IE exp − 2σ1 2 (kyk2 − 2y > Rθ + kRθk2 )

= 
IE exp − 2σ1 2 (kyk2 − 2y > R0 φ + kR0 φk2 )


 
IE exp σ12 y > Rθ
1
2
2

 ,
= exp
(kφk − kθk )
2σ 2
IE exp σ12 y > R0 φ
since R is orthogonal. Hence
D(θ k φ) = IEY ∼Pθ log

dPθ
(Y )
dPφ
1

>

IE[e σ2 Y Rθ | Y ]
1
.
= 2 (kφk2 − kθk2 ) + IEY ∼Pθ log
1
>
2σ
IE[e σ2 Y Rφ | Y ]
When Y ∼ Pθ , we can write Y = R0 θ + σξ for a standard Gaussian vector
ξ and rotation R0 independent of R. Since R0 is orthogonal, this has the same
distribution as R0 (θ +σξ). If R and R0 are independent and uniform, then (R0 )> R
has the same distribution as R, so
Y > Rθ =d (θ + σξ)> Rθ .
Therefore
1

(A.1)

>

1
IE[e σ2 (θ+σξ) Rθ | ξ]
2
2
D(θ k φ) = 2 (kφk − kθk ) + IEξ log
.
1
>
2σ
IE[e σ2 (θ+σξ) Rφ | ξ]

We now prove both parts of the Lemma
(i) For convenience write θ̄ = IERθ and φ̄ = IERφ. These vectors satisfy Rθ̄ = θ̄
and Rφ̄ = φ̄ for all R. Hence
(θ + σξ)> Rθ = (θ̄ + ϑ + σξ)> R(θ̄ + ϑ)
= (ϑ + σξ)> Rϑ + (θ̄ + σξ)> θ̄ ,
and similarly
(θ + σξ)> Rφ = (ϑ + σξ)> Rϕ + (θ̄ + σξ)> φ̄ .
Plugging these quantities into the expression for D(θ k φ) yields
1

>

>

1
IE[e σ2 (ϑ+σξ) Rϑ+(θ̄+σξ) θ̄ | ξ]
D(θ k φ) = 2 (kφk2 − kθk2 ) + IEξ log
1
>
>
2σ
IE[e σ2 (ϑ+σξ) Rϕ+(θ̄+σξ) φ̄ | ξ]
1

>

1
1
IE[e σ2 (ϑ+σξ) Rϑ | ξ]
= 2 (kφk2 − kθk2 ) + 2 IEξ (θ̄ + σξ)> (θ̄ − φ̄) + IEξ log
1
>
2σ
σ
IE[e σ2 (ϑ+σξ) Rϕ | ξ]
1
1
= 2 (kφ̄k2 − kθ̄k2 ) + 2 (kθ̄k2 − θ̄> φ̄) + D(ϑ k ϕ)
2σ
σ
1
= 2 kIE[Rθ − Rφ]k2 + D(ϑ k ϕ) .
2σ

20

BANDEIRA, RIGOLLET AND WEED

(ii) The previous claim implies that it suffices to show
1
kIE[(Rθ)⊗2 − (Rφ)⊗2 ]k2 + σ −6 O(ε2 )
4σ 4
for θ and φ satisfying IERθ = IERφ = 0.
We accomplish this by expanding the expression for D(θ k φ) given above as a
power series in σ −1 . Recall
D(θ k φ) =

1

>

1
IE[e σ2 (θ+σξ) Rθ | ξ]
D(θ k φ) = 2 (kφk2 − kθk2 ) + IEξ log
.
1
>
2σ
IE[e σ2 (θ+σξ) Rφ | ξ]
Given a random variable X, write
KX (t) = log IE exp(t> X)
for the cumulant generating function of X. Then
1
>




IE[e σ2 (θ+σξ) Rθ | ξ]
1
1
log
= KRθ
(θ − σξ) − KRφ
(θ − σξ) .
1
>
σ2
σ2
IE[e σ2 (θ+σξ) Rφ | ξ]
Denote by κm and λm the mth cumulant tensors of Rθ and Rφ respectively. Then
X
1
1
D(θ k φ) = 2 (kφk2 − kθk2 ) + IE
hκm − λm , (θ − σξ)⊗m i .
2m
2σ
σ m!
m≥1

Since κ1 = IERθ and λ1 = IERφ,
κ1 − λ1 = IERθ − IERφ = 0 ,
so the first term in the sum vanishes. Moreover, since ξ is a standard Gaussian,
IE[ξ ⊗m ] = 0 for all odd m.
Rearranging to collect powers of σ −1 yields
D(θ k φ) =

1
1
T2 + 4 T4 + o(σ −4 ) ,
σ2
σ

where
1
T2 = (kφk2 − kθk2 + IEhκ2 − λ2 , ξ ⊗2 i)
2
1
1
1
T4 = hκ2 − λ2 , θ⊗2 i + IEhκ3 − λ3 , θ ⊗ ξ ⊗2 i + IEhκ4 − λ4 , ξ ⊗4 i .
2
2
4!
Straightforward calculation yields
IEhκ2 − λ2 , ξ ⊗2 i = kθk2 − kφk2
hκ2 − λ2 , θ⊗2 i = IE[(θ> Rθ)2 − (φ> Rθ)2 ]
IEhκ3 − λ3 , θ ⊗ ξ ⊗2 i = 0
IEhκ4 − λ4 , ξ ⊗4 i = 6IE[(φ> Rφ)2 − 6(θ> Rθ)2 ] .
Combining the above two displays, we conclude
T2 = 0
1
1
T4 = IE[(θ> Rθ)2 − 2(φ> Rθ)2 + (φ> Rφ)] = kIE[(Rθ)⊗2 − (Rφ)⊗2 ]k2 ,
4
4
and the claim follows.

OPTIMAL RATES FOR MRA

21

Proof of Proposition 6. We first suppose ρ(θ, φ) ≤ ε0 .
For any φ ∈ IRd , we have IE[Rφ] = d−1/2 φ̂0 1, where 1 is the all-ones vector
of IRd and φ̂0 is known as the DC component of φ. Therefore ϕ = φ − IE[Rφ]
has Fourier transform ϕ
b = φ̂ − φ̂0 d−1/2 W 1 = φ̂ − φ̂0 e0 , where e0 is the vector of
d
IR indexed over {−bd/2c, . . . , bd/2c} with 1 on the 0 coordinate and 0 on other
coordinates. In other words, the Fourier transform of ϕ is the same as that of φ
except that its DC component ϕ
b0 is set to 0. Similarly, the Fourier transform of
ϑS is the same as that of θ except that its DC component ϑb0 is set to 0.
By Part (i) of Lemma 2,
(A.2)

D(φ) =

1
kIE[Rθ − Rφ]k2 + D(ϑ k ϕ) ,
2σ 2

Write ρ(θ, φ) = ε and suppose first that |θ̂0 − φ̂0 | ≥ 21 ε. Then (A.2) implies
D(φ) ≥

1
ε2
1
kIE[Rθ − Rφ]k2 = 2 (θ̂0 − φ̂0 )2 ≥
≥ Cσ −4s+2 ε2 .
2
2σ
2σ
8dσ 2

Next, if |θ̂0 − φ̂0 | < 21 ε, then
ρ(ϑS , ϕ)2 = ρ(θS , φ)2 − |θ̂0 − φ̂0 |2 ≥ 3ε2 /4 .
Thus, it suffices to show that
D(ϑ k ϕ) ≥ Cσ −4s+2 ρ(ϑ, ϕ)2 .
There are two cases: either ϑS and ϕ have essentially the same power spectrum
(i.e., |ϑ̂k | ≈ |ϕ̂k | for all k) or their power spectra are very different. We will treat
these two cases separately.
Recall that for each j ∈ S, by assumption c−1 ≤ |ϑ̂j | ≤ c. Consider the polar
form ϕ̂j /ϑ̂j = mj eiδj , where mj ≥ 0. Since D(ϕ) = D(Rϕ) for all cyclic shifts,
we may assume that kϑ − ϕk = ε, so that |1 − mj | ≤ Cε for all j.
Suppose first that |1 − mj | ≥ Cε for some j ∈ S. Lemma 2 (ii) yields
1
ε2
⊗2
⊗2 2
kIE[(Rϑ)
−
(Rϕ)
]k
+
C
4σ 4
σ6
2
1
ε
≥ 4 (|ϑ̂j |2 − |ϕj |2 )2 + C 6
4σ
σ
2
ε
c−2
≥ 4 |ϑ̂j |2 (1 − mj )2 + C 6
4σ
σ
2
c−4
ε
≥ 4 (1 − mj )2 − C 6 ≥ Cσ −4 ε2 .
4σ
σ

D(ϕ) =

Hence D(ϕ) ≥ Cσ −4s+2 ε2 .
Next, suppose on the contrary that |1−mj | = o(ε) for all j ∈ S. Since kϑ−ϕk =
ε, we can take the relative phase δj in the polar form to be such that δj ≤ Cε for
all j. By Theorem 3, it is enough to show that there exists an m ≤ 2s − 1 such
that
kIE[(Rϑ)⊗m − (Rϕ)⊗m ]k2 = Ω(ε2 ) .

22

BANDEIRA, RIGOLLET AND WEED

Denote by p the smallest integer in S and observe that
bd/2c

X

ε2 = ρ(ϑ, ϕ)2 = min
z:|z|=1

|1 − mj z j eiδj |2 |ϑ̂j |2 ≤ C

X

|1 − mj e−jδp /p eiδj |2 .

j∈S

j=−bd/2c

Therefore, there exists a coordinate ` ∈ S \ {p} such that
|1 − ei(pδ` −`δp )/p |2 = |1 − m` ei(pδ` −`δp )/p |2 + o(ε2 ) ≥ Cε2 .

(A.3)

Choose m = ` + p. Since `, p ∈ S ⊆ [s] and ` 6= p, the bound m ≤ 2s − 1 holds.
As in the proof of Proposition 7, we have that
kIE[(RϑS )

⊗m

⊗m

− (Rϕ)

m
Y

X

2

]k =

ϑ̂jn −

j1 +···+jm =0 n=1

X

=

2

ϕ̂jn

n=1
m
Y

1−

j1 +···+jm =0

m
Y

mjn e

i

Pm

n=1 δjn

2 m
Y

|ϑ̂jn |2 .

n=1

n=1

Each term in the above sum is positive. One valid solution to the equation j1 +
· · · + jm = 0 is j1 = · · · = j` = −p and j`+1 = · · · = jm = `. We obtain
⊗m

kIE[(RϑS )

⊗m

− (Rϕ)

2

]k ≥ C 1 −
= |1 − e

m
Y

2
i(pδ` −`δp )

mjn e
n=1
i(pδ` −`δp ) 2

| + o(ε2 ) .

Observe that for δ` and δp small enough, it holds
|1 − ei(pδ` −`δp ) |2 ≥ |1 − ei(pδ` −`δp )/p |2 ≥ Cε2
where the last inequality follows from (A.3). Combining the above two displays
proves (5.12).
Now suppose ρ(θ, φ) ≥ ε0 . To show that D(θ k φ) ≥ Cσ −4s+2 , it suffices to
show by the chain rule for divergence that
D(PθN k PφN ) ≥ c for some N ≤ Cσ 4s−2 ,
where c and C are positive constants. We will show the existence of a test which
correctly distinguishes PθN from PφN with probability at least 2/3; the claim will
then follow from the Neyman-Pearson lemma and Pinsker’s inequality.
We first show the existence of such a test when kφk2 ≥ 2. Let N = γσ 4 ≤
γσ 4s−2 for s ≥ 2, with γ > 0 to be chosen later. Let {Yi }N
i=1 be samples from
N
N
d×N
either Pθ or Pφ , and define a test ψ : R
→ {θ, φ} by

P
2
2
θ if N1 N
i=1 kYi k ≤ σ + 1.5
ψ(Y1 , . . . , YN ) =
φ otherwise.
An easy computation shows that
N
1 X
IEθ [
kYi k2 ] = σ 2 + kθk2 ≤ σ 2 + 1
N

varθ [

1
N

i=1
N
X
i=1

kYi k2 ] ≤

4σ 2

+
N

2σ 4

≤

C
γ

N
1 X
IEφ [
kYi k2 ] = σ 2 + kφk2
N

varφ [

1
N

i=1
N
X
i=1

kYi k2 ] ≤

C −2
(σ kφk2 + 1) .
γ

OPTIMAL RATES FOR MRA

23

Together with Chebyshev’s inequality, we get that that for kφk2 ≥ 2,
Pθ (ψ = φ) + Pφ (ψ = θ) ≤

C
≤ 1/3 ,
γ

For γ large enough, as desired.
Next, suppose that kφk2 ≤ 2. For positive integers j, k, denote by (j, k) their
greatest common divisor. Given a vector ζ ∈ Rd , denote by P the following set
of polynomials in the entries of ζ̂:
p0 (ζ) = ζ̂0 ,
pj (ζ) = kζ̂j k2

for 1 ≤ j ≤ bd/2c,

j/(j,k) k/(j,k)
ζ̂j
for

pjk (ζ) = ζ̂−k

1 ≤ j, k ≤ bsc.

If ρ(θ, φ) > 0, then by [KI93, Corollary 2], there exists at least one polynomial
p ∈ P such that
p(θ) 6= p(φ) .
It is easy to see that all the polynomials in P are invariant under the group
action; that is,
p(ζ) = p(Rζ)
for any R ∈ F and p ∈ P. This implies that, given a p ∈ P, the value p(ζ)
can be computed from the moment tensors IE[(Rζ + σξ)⊗k ] for 1 ≤ k ≤ 2s − 1.
Indeed, the entries of IE[(Rζ + σξ)⊗k ] for 1 ≤ k ≤ 2s − 1 generate the ring of
invariant polynomials in the entries of ζ̂ of degree at most 2s − 1, which includes
the set P. Given samples Y1 , . . . , Yn from Pζ , the tensor IE[(Rζ + σξ)⊗k ] can be
consistently estimated by computing the empirical moment tensors. We obtain
that there exists a constant M , depending on s but not on σ, such that for any
p ∈ P there exists an unbiased estimator p̃(ζ) for p(ζ) with variance at most
M σ 4s−2 /N for any ζ ∈ {θ, φ}.
For all θ ∈ Rd , define Bθ,ε0 = {φ : ρ(θ, φ) ≤ ε0 , kφk2 ≤ 2}. It is clear that Bθ,ε0
is compact so that
δ=

inf

inf

min

θ:kθk2 ≤1 φ∈Bθ,ε0 p∈P:p(φ)6=p(θ)

|p(φ) − p(θ)| > 0 .

Note that δ does not depend on θ or φ. Set N = γδ −2 M σ 4s−2 , where γ is to be
chosen later. Since ρ(θ, φ) ≥ ε0 by assumption, there exists a p ∈ P such that
|p(φ) − p(θ)| ≥ δ.
Let ψ : Rd×N → {θ, φ} be the test
ψ(Y1 , . . . , Yn ) = argmin |p̃(ζ) − p(ζ)| .
ζ∈{θ,φ}

By Chebyshev’s inequality, we have
Pθ (ψ = φ) + Pφ (ψ = θ) ≤
for γ sufficiently large.

C
≤ 1/3 ,
γ

24

BANDEIRA, RIGOLLET AND WEED

To conclude the proof, observe that by the chain rule, Pinsker’s inequality, and
the Neyman-Pearson lemma respectively, we get
D(θ k φ) =

1
C
D(PθN k PφN ) ≥
≥ Cσ −4s+2 .
N
N

The desired result (5.13) follows.
√
Finally, let ρ(θ, φ) ≥ c0 σ where c0 = 32 2d. By (A.1), since kθk ≤ 1, we have
1

>

1
IE[e σ2 (θ+σξ) Rθ | ξ]
D(φ) ≥ 2 (kφk2 − 1) + IEξ log
,
1
>
2σ
IE[e σ2 (θ+σξ) Rφ | ξ]
where ξ ∼ N (0, Id ). Next, using the Cauchy-Schwarz inequality and Jensen’s
inequality, we get
1

IEξ log

> Rθ

IE[e σ2 (θ+σξ)
IE[e

1
(θ+σξ)> Rφ
σ2

| ξ]
| ξ]

≥−

1 + kφk
1 + kφk p
IEkθ + σξk ≥ −
1 + dσ 2
2
σ
σ2

√
Hence, for kφk ≥ 16σ 2d, we get D(φ) ≥ kφk2 /(8σ 2 ). Moreover, by the triangle
inequality, kφk/2 ≤ ρ(φ, θ) ≤ 2kφk for all such φ. We obtainthat for all φ ∈ IRd
such that ρ(φ, θ) ≥ c0 σ, it holds D(φ) ≥ kφ − θk2 /(32σ 2 ), which implies the
desired result (5.14).

APPENDIX B: SUPPLEMENTARY MATERIALS
Lemma B.1.

For any δ ∈ (0, 1), γ ∈ [−δ, δ], it holds
1 − 2(1 + γ)m + (1 + 2γ + δ 2 )m ≤ m2 2m−1 δ 2

Proof. Fix δ ∈ (0, 1) and consider the function γ 7→ f (γ, δ) = 1 − 2(1 + γ)m +
(1+2γ +δ 2 )m . Computing the derivative of f , it is easy to see that f is decreasing
on (−δ, −δ 2 ) and increasing on (−δ 2 , δ) so that it achieves its maximum either
at −δ or at δ.
Consider the function δ 7→ g(δ) = f (δ, δ) = 1 − 2(1 + δ)m + (1 + δ)2m and
observe that g(−δ) = f (−δ, δ). Since g(0) = g 0 (0) = 0, and
g 00 (δ) = 2m(1 + δ)m−2 [(2m − 1)(1 + δ)m − (m − 1)] ≤ m2 2m

∀ δ ∈ [−1, 1] ,

it follows from a second order Taylor expansion that
g(δ) ≤ m2 2m−1 δ 2

Lemma B.2. Let ξ ∼ N (0, Id ) be a standard d-dimensional Gaussian. Then
m
for any symmetric order-m tensor and A ∈ IRd there exists constant c = c(A),
independent of m and d, such that
p
√
cm m!kAk ≤ IE|hA, ξ ⊗m i| ≤ (d + m)m kAk .

OPTIMAL RATES FOR MRA

25

Proof. Assume without loss of generality that kAk = 1. The upper bound
follows readily from the Cauchy-Schwarz inequality and the fact that kξk2 ∼ χ2d ,
IEξ |hA, ξ ⊗m i| ≤ IEξ kξ ⊗m k = IEξ kξkm ≤ (d + m)m/2 .
To prove the lower bound, we first show that
IEξ [hA, ξ ⊗m i2 ] ≥ m! .

(B.1)

We do so by expressing the left side in terms of multivariate Hermite polynomials.
Recall that the Hermite polynomials {hk (x)}k≥0 can be normalized to satisfy
the following two properties:
√
1. The function hk (x) is a polynomial with leading term xk / k!,
2. The functions {hk }k≥0 form an orthonormal basis of of L2 (µ), where µ
denotes the standard Gaussian measure on IR.
Given a multi-index α ∈ INd , define the multivariate Hermite polynomial hα
by
hα (x1 , . . . , xd ) =

d
Y

hαi (xi ) .

i=1

The multivariate Hermite polynomials form an orthonormal basis for the space
IR[x1 , . . . , xd ] of d-variate polynomial functions with respect to the inner product
over L2 (µ⊗d ).
Multivariate Hermite polynomials satisfy the following useful property: Given
two multi-indices α, β ∈ INd such that |α| = |β|, we now show that
√
(B.2)
hxβ , hα i = α!δαβ ,
where δ denotes the Kronecker symbol. Indeed, on the one hand if α 6= β, then
since |α| = |β| there exists an index i such that αi > βi . By the definition of
the univariate Hermite polynomials, xβi ∈ span(h1 (x), . . . , hβi (x)), hence orthogonality of the polynomials hαi and hj for 1 ≤ j ≤ βi implies
hxβ , hα i =

d
Y
hxβi , hαi i = 0 .
i=1

On the other hand, if α = β, then
hxα , hα i =

d
d p
Y
Y
√
hxαi , hαi i =
αi ! = α! .
i=1

i=1

The order-m tensor A can be identified with a multilinear—thus polynomial—
map from IRd to IR:
X
X
m!
A(x1 , . . . , xd ) =
Ai1 ...im xi1 . . . xim =
Aα x α ,
α!
d
i1 ,...,im

α∈IN :|α|=m

where in the second equality, we used the fact that A is symmetric. Together
with (B.2), it yields that for any |α| = m,
m!
hA, hα i = √ Aα .
α!

26

BANDEIRA, RIGOLLET AND WEED

Moreover, since kAk = 1, we also have
X
α∈INd :|α|=m

m! 2
A = 1,
α! α

so that
X

hA, hα i2 = m! .

α∈INd :|α|=m

Therefore, using Plancherel’s formula, we get that for any m ≥ 0,
X
X
hA, hα i2 = m! ,
hA, hα i2 ≥
IEξ [hA, ξ ⊗m i2 ] = hA, Ai =
α∈INd :|α|=m

α∈INd

as claimed.
To prove the claimed lower bound on IEξ |hA, ξ ⊗m i|, we employ the following
Khinchine-type inequality:
Theorem B.3 ([Bob00, Theorem 2]). Let f = f (x1 , . . . , xm ) be a degree m
polynomial and let X ∼ N (0, Im ) be a standard Gaussian vector in IRm . Then
there exists a universal constant c such that
p
IE[f (X)2 ] ≤ cm IE[|f (X1 , . . . , Xm )|] .
Applying Theorem B.3 to (B.1) yields the claim.
Lemma B.4. Let H(ζ) and Hn (ζ) be the Hessians of D(φ) and Dn (φ), respectively, evaluated at φ = ζ. If Bε := {φ ∈ IRd : ρ(φ, θ) ≤ ε}, then
IE sup kH(φ) − Hn (φ)k2op ≤ C
φ∈Bε

log n
.
nσ 4

Proof. The matrix Hn (φ) can be written as a sum of independent random
matrices:
n
pθ
1X
Hn (φ) =
Hi (φ) ,
Hi (φ) = ∇2φ log (Yi ) .
n
pφ
i=1

Using symmetrization, we get
(B.3)

IE sup kH(φ) − Hn (φ)k2op ≤
φ∈Bε

n
X
4
IE
sup
εi Hi (φ)
n2 φ∈Bε

2
,
op

i=1

where ε1 , . . . , εn are i.i.d Rademacher random variables that are independent of
the observations Y1 , . . . , Y1 and IEHn (φ) = H(φ).
Using a third order derivative calculation, we find that φ 7→ Hi (φ) satisfies
the following Lipschitz property. For any u, φ, η ∈ IRd , such that kuk = 1 and
φ, η ∈ Bε , we have
|u> Hi (φ)u − u> Hi (η)u| ≤ C
where ξi is the noise in (2.1).

1 + |ξi |3
kφ − ηk
σ3

27

OPTIMAL RATES FOR MRA

Fix γ ∈ (0, ε) and le Z be a γ-net of Bε . In other words, we require that
max min kη − φk ≤ γ .
φ∈Bε η∈Z

We can always choose Z to have cardinality |Z| ≤ (Cε/γ)d for some universal
constant C > 0. Then, by Young’s inequality, we get
(B.4)

sup
φ∈Bε

n
X

2
op

εi Hi (φ)

≤C

i=1

n
n
X
γ 2 X 3 2
|ξ
|
+
C
max
εi Hi (φ)
i
φ∈Z
σ6
i=1

2
op

i=1

The expectation of the first term is controlled using the fact that
(B.5)

IE

n
X

|ξi |3

2

i=1

=

n
X

IE[|ξi ξj |3 ] ≤ Cn2 .

i,j=1

For second term, we employ a standard matrix concentration bound [Tro15, Theorem 4.6.1] to get that
n
n
X




t2 X
k
Hi (φ)2 kop .
IE exp tk
εi Hi (φ)k2op Y1 , . . . , Yn ≤ d exp
2
i=1

i=1

Using standard arguments (see, e.g., [BLM13]), this implies
n
X

IE max
εi Hi (φ)
φ∈Z

2 
op



≤ C(log |Z|)IE max
φ∈Z

i=1

n
X

Hi (φ)2


op

i=1



≤ C(log |Z|)nIE max kH1 (φ)2 kop
φ∈Z

As before, computing a second order derivative, we can show that
kH1 (φ)2

op

2
op

≤ kH1 (φ)

≤ C(1 + |ξ|4 )/σ 4 .

The above two displays yield
n
X

IE max
εi Hi (φ)
φ∈Z

2 
op

≤C

i=1

log(ε/γ)
n
σ4

Combining the last display with (B.3), (B.4), and (B.5), we get
(B.6)

IE sup kH(φ) − Hn (φ)k2op ≤ C
φ∈Bε

 γ2
σ6

+

log(ε/γ) 
log n
≤C
,
4
nσ
nσ 4

for γ = n−1/2 .
Lemma B.5.
satisfies

Assume the conditions of Theorem 4 hold. Then the MLE θ̃n
IE[ρ(θ̃, θ)2 ] ≤ C

σ (8k−4)∧(2`+2)
.
n

28

BANDEIRA, RIGOLLET AND WEED

Proof. As in the proof of Theorem 4, since θ is fixed, we simply write D(φ) =
D(θ k φ) and define
n
fθ
1X
log (Yi ) ,
Dn (φ) =
n
fφ
i=1

where Yi are i.i.d from model (2.1) and we recall that fζ is the density of Pζ , ζ ∈
IRd .
We first establish using Lemma B.7 that the process {Gn (φ)}φ∈IRd defined by
Gn (φ) = D(φ) − Dn (φ) is a subgaussian process with respect to the Euclidean
distance with variance proxy cσ 2 /n for some constant c > 0, i.e., that for any
λ ∈ IR, we have
IE[exp(λ(Gn (φ) − Gn (ζ)))] ≤ exp c


λ2 σ 2
kφ − ζk2 .
n

We then apply the following standard tail bound.
Proposition B.6 ([Ver17, Theorem 8.5.4]). If {Xφ }φ is a (standard) subgaussian process on IRd with respect to the Euclidean metric and Bδ (θ) is a ball
of radius δ around θ, then
2 /δ 2

IP[ sup (Xφ − Xθ ) ≥ Cδ + x] ≤ Ce−Cx

.

φ∈Bδ (θ)

√
The rescaled process σ nGn is standard subgaussian process with respect to
the Euclidean metric, so applying Proposition B.6 and noting that Dn (θ) = 0
yields
"
#
δ
nσ 2 x2 
(B.7)
IP sup Gn (φ) ≥ C √ + x ≤ C exp − C
.
δ2
σ n
φ∈Bδ (θ)
√
For convenience, write vn = n(θ̃ − θ), where θ̃ is a MLE satisfying kθ̃ − θk =
ρ(θ̃, θ). We wish to show that IE[kvn k2 ] is bounded by a constant that depends
on σ but not on n.
Define E to be the compact subset of IRd defined by

E = φ ∈ IRd : ε0 ≤ ρ(θ, φ) ≤ c0 σ ,
where c0 and ε0 are defined in Theorem 4. In particular, for any φ ∈
/ E, it holds
−2k
2
D(φ) ≥ Cσ ρ(θ, φ) .
We employ the so-called slicing (a.k.a peeling) method. Define the sequence
{αj }j≥0 where α0 = 0 and αj = Cσ j(2k−1) for j ≥ 1 for some large enough
√
constant C > 0. For any j ≥ 0, define Sj = {φ ∈ IRd : αj ≤ nρ(φ, θ) ≤
αj+1 } \ E and observe that
X
IE[kvn k2 ] = IE[kvn k2 | θ̃ ∈ E]IP[θ̃ ∈ E] +
IE[kvn k2 | θ̃ ∈ Sj ]IP[θ̃ ∈ Sj ]
j≥0

(B.8)

≤

c20 σ 2 IP[θ̃

∈ E] + Cσ

8k−4

+

X

2
αj+1
IP[θ̃ ∈ Sj ] .

j≥2

We now show that if θ̃ ∈ Sj , j ≥ 2, then Gn (θ̃) = D(θ̃) − Dn (θ̃) is large.
To that end, observe that on the one hand, by definition of the MLE, we have

OPTIMAL RATES FOR MRA

29

Dn (θ̃) ≤ Dn (θ) = 0. On the other hand, D(θ̃) ≥ Cσ −2k ρ(θ̃, θ)2 . Hence, if θ̃ ∈ Sj ,
then Gn (θ̃) ≥ Cσ −2k ρ(θ̃, θ)2 ≥ Cσ −2k αj2 /n. It yields
αj2 
αj2 


≤ IP
sup
Gn (φ) ≥ Cσ −2k
.
IP[θ̃ ∈ Sj ] ≤ IP sup Gn (φ) ≥ Cσ −2k
n
n
φ∈Sj
φ∈B αj+1 (θ)
√
n

Recall αj = Cσ j(2k−1) so that σ −2k αj2 ≥ Cαj+1 /σ, j ≥ 2 and apply (B.7) with
√
δ = αj+1 / n and x = Cσ −2k αj2 /n to get

IP

Gn (φ) ≥ Cσ

sup

−2k

φ∈B αj+1 (θ)
√

αj2 
αj4

≤ C exp − C 2
n
αj+1 σ 4k−2

n


≤ C exp − Cσ 2j(2k−1) .
It yields
(B.9)

X

2
αj+1
IP[θ̃ ∈ Sj ] ≤ C

j≥2

X


σ 4j(2k−1) exp − Cσ 2j(2k−1) ≤ C .

j≥2

When θ̃ ∈ E, we use (5.2) to conclude that if n ≥ σ 2` then



IP[θ̃ ∈ E] ≤ IP sup Gn (φ) ≥ Cσ −` ≤ C exp − Cnσ −2` ,
φ∈E

where in the last inequality, we used (B.7). Together with (B.8) and (B.9), it
implies the desired result.
Lemma B.7. The process {Gn (φ)}φ∈IRd defined by Gn (φ) = D(φ) − Dn (φ) is
a subgaussian process with respect to the `2 distance on IRd with variance proxy
c/(nσ 2 ) for some constant c > 0, i.e., that for any λ ∈ IR, we have
IE[exp(λ(Gn (φ) − Gn (ζ)))] ≤ exp λ2


c
kφ − ζk2 .
2
nσ

Proof. By definition of Gn and the densities fζ and fφ , we have
Gn (φ)−Gn (ζ) = D(φ) − D(ζ) − Dn (φ) + Dn (ζ)
n
1X
= IE[log fζ (Y ) − log fφ (Y )] −
[log fζ (Yi ) − log fφ (Yi )]
n
h
= IE log
h
= IE log

i=1
kY −Rζk2 i
−Rζk2
n
IER exp(− 2σ2 )
IER exp(− kYi 2σ
)
1X
2
−
log
2
kY −Rφk2
kY
−Rφk
i
n
IER exp(− 2σ2 )
IER exp(− 2σ2 )
i=1
n
IER exp(Y > Rζ/σ 2 ) i 1 X
IER exp(Yi> Rζ/σ 2 )
log
−
n
IER exp(Y > Rφ/σ 2 )
IER exp(Yi> Rφ/σ 2 )
i=1
n

= IE[∆(Y )] −

1X
∆(Yi ) ,
n
i=1

where
∆(Y ) = log

IER exp(Y > Rζ/σ 2 )
.
IER exp(Y > Rφ/σ 2 )

30

BANDEIRA, RIGOLLET AND WEED

Next, using a standard symmetrization argument, we get
(B.10)

IE[exp(λ(Gn (φ) − Gn (ζ)))] ≤

n
Y

IE[exp(

i=1

2λ
εi ∆(Yi ))] ,
n

where ε1 , . . . , εn are i.i.d Rademacher random variables that are independent of
the sample Y1 , . . . , Yn . Next, observe that using the Cauchy-Schwarz inequality,
we get that the function ζ 7→ log IER exp(Yi> Rζ/σ 2 ) is kYi k/σ 2 -Lipschitz with
respect to the Euclidean distance. Hence |∆(Yi )| ≤ kYi kkφ − ζk/σ 2 . Moreover,
it follows from (2.1) that conditionally on R`i , the random variable kYi k2 /σ 2 ∼
χ2d (1/σ 2 ) follows a noncentral χ2 distribution with d degrees of freedom and
noncentrality parameter kθk2 /σ 2 = 1/σ 2 ≤ 1. In particular kYi k2 /σ 2 has a finite
moment generating function in the neighborhood of the origin, and we obtain for
some constant c small enough that


σ 2 (εi ∆(Yi ))2 
kYi k2 
IE exp c
≤
IE
exp
c
≤ e.
kφ − ζk2
σ2
It implies (see, e.g., [Ver17] proposition 2.5.2) that


c
2λ
IE exp( εi ∆(Yi )) ≤ exp(λ2 2 2 kφ − ζk2 )
n
n σ
for all λ ∈ IR. Together with (B.10), this yields the desired result.
Lemma B.8. Fix θ ∈ IRd . Assume that n/(log n) ≥ Cσ 4 . For any j =
−bd/2c, . . . , bd/2c, define,
n

1X d 2
|(Yi )j | − σ 2 .
Mj =
n
i=1

Define the set S̃ by
r
n
log n o
S̃ = j ∈ {−bd/2c, . . . , bd/2c} : Mj ≥ Cσ 2
n
Then
IP[S̃ 6= psupp(θ̂)] ≤ C

σ4
.
n

Proof. It is straightforward to check that IE[Mj ] = |θ̂k |2 . To calculate the
variance, we note that
2
2
2
2 d 2
d
\ d
|(Y
i )k | − |θ̂k | − σ = 2σ<((Ri θ)k (ξi )k ) + σ (|(ξi )k | − 1) .
2
2
4
4
d
Hence var[|(Y
i )k | − σ ] ≤ Cσ , and var[Mj ] ≤ Cσ /n.
2
2
2
d
Moreover, the random variable |(Y
i )k | − σ is a σ -subexponential random
variable. Indeed,
2
2
2
2 d 2
2 d 2
2
d
d
|(Y
i )k | − |θ̂k | − σ ≤ 2σ|(ξi )k | + σ (|(ξi )k | − 1) ≤ 2σ (|(ξi )k | − 1) + σ + 1 .

31

OPTIMAL RATES FOR MRA

d
d
The variable (ξ
i )0 is a standard Gaussian, and for k 6= 0 the variable (ξi )k is a
2
2
d
standard complex Gaussian, hence |(ξ
i )k | is a rescaled χ random variable with
2
d
at most 2 degrees of freedom. Hence |(ξ
i )k | is sub-exponential with constant
2
2
2
d
variance proxy, so that |(Y
i )k | − σ is σ -subexponential.
If j ∈ psupp(θ̂), then by Chebyshev’s inequality
r

σ4
log n 
σ4
2
q
≤C
.
≤
C
IP[j ∈
/ S̃] = IP Mj ≤ Cσ
2
n
n
n |θ̂ |2 − Cσ 2 log n
j

n

On the other hand, if j ∈
/ psupp(θ̂), then IE[Mj ] = 0 and by Bernstein’s
inequality
r
log n
C
2
IP[j ∈ S̃] = IP[Mj ≥ Cσ
]≤ .
n
n
The proof follows using a union bound.
REFERENCES
[ABBS14]

[ADBS16]

[APS17]

[BBS16]

[BBV16]

[BCS15]

[BCSZ14]

[BLM13]

[Bob00]

E. Abbe, A. S. Bandeira, A. Bracher, and A. Singer. Decoding
binary node labels from censored edge measurements: Phase transition and efficient recovery. Network Science and Engineering,
IEEE Transactions on, 1(1):10–22, Jan 2014.
Cecilia Aguerrebere, Mauricio Delbracio, Alberto Bartesaghi, and
Guillermo Sapiro. Fundamental limits in multi-image alignment.
Available online at arXiv:1602.01541 [], 2016.
Emmanuel Abbe, Joao Pereira, and Amit Singer. Sample complexity of the boolean multireference alignment problem. Available
online at arXiv:1701.07540 [], 2017.
Afonso S. Bandeira, Nicolas Boumal, and Amit Singer. Tightness
of the maximum likelihood semidefinite relaxation for angular synchronization. Mathematical Programming, pages 1–23, 2016.
A. S. Bandeira, N. Boumal, and V. Voroninski. On the low-rank
approach for semidefinite programs arising in synchronization and
community detection. COLT, 2016.
A. S. Bandeira, Y. Chen, and A. Singer. Non-unique games over
compact groups and orientation estimation in cryo-em. Available
online at arXiv:1505.03840 [], 2015.
Afonso S. Bandeira, Moses Charikar, Amit Singer, and Andy
Zhu. Multireference alignment using semidefinite programming.
In ITCS’14—Proceedings of the 2014 Conference on Innovations
in Theoretical Computer Science, pages 459–470. ACM, New York,
2014.
S. Boucheron, G. Lugosi, and P. Massart. Concentration inequalities. Oxford University Press, Oxford, 2013. A nonasymptotic
theory of independence, with a foreword by Michel Ledoux.
S. G. Bobkov. Remarks on the growth of Lp -norms of polynomials.
In Geometric aspects of functional analysis, volume 1745 of Lecture
Notes in Math., pages 27–35. Springer, Berlin, 2000.

32

[Bou16]
[Bro92]
[BSS11]
[CC16]

[CL11]

[Dia92]

[Edi15]
[FZB02]

[GBM16]

[HK15]
[HUL01]

[JMRT16]

[KI93]

[LeC73]
[LNS99]

[MP00]

[MV10]
[Nog15]

BANDEIRA, RIGOLLET AND WEED

Nicolas Boumal. Nonconvex phase synchronization. SIAM Journal
of Optimization, to appear, 2016.
Lisa Gottesfeld Brown. A survey of image registration techniques.
ACM computing surveys (CSUR), 24(4):325–376, 1992.
A. S. Bandeira, A. Singer, and D. Spielman. The so(3) cheeger
inequality. Unpublished Draft, 2011.
Yuxin Chen and Emmanuel Candes. The projected power method:
An efficient algorithm for joint alignment from pairwise differences.
Available online at arXiv:1609.05820 [], 2016.
T. Tony Cai and Mark G. Low. Testing composite hypotheses,
hermite polynomials and optimal estimation of a nonsmooth functional. Ann. Statist., 39(2):1012–1041, 04 2011.
R. Diamond. On the multiple simultaneous superposition of molecular structures by rigid body transformations. Protein Science,
1(10):1279–1287, October 1992.
Editorial. Method of the year 2015. Nature Methods, 13:1, 2015.
H. Foroosh, J. B. Zerubia, and M. Berthod. Extension of phase
correlation to subpixel registration. IEEE Transactions on Image
Processing, 11(3):188–200, 2002.
Tingran Gao, Jacek Brodzki, and Sayan Mukherjee. The geometry
of synchronization problems and learning group actions. Available
online at arXiv:1610.09051 [], 2016.
Philippe Heinrich and Jonas Kahn. Optimal rates for finite mixture
estimation. arXiv:1507.04313, 2015.
Jean-Baptiste Hiriart-Urruty and Claude Lemaréchal. Fundamentals of convex analysis. Grundlehren Text Editions. SpringerVerlag, Berlin, 2001. Abridged version of ıt Convex analysis and
minimization algorithms. I [Springer, Berlin, 1993; MR1261420
(95m:90001)] and ıt II [ibid.; MR1295240 (95m:90002)].
A. Javanmard, A. Montanari, and F. Ricci-Tersenghi. Phase transitions in semidefinite relaxations. Proceedings of the National
Academy of Sciences of the United States of America (PNAS),
2016.
R. Kakarala and G. J. Iverson. Uniqueness of results for multiple
correlations of periodic functions. Oct. Soc. Am. A, 10:1517–1528,
1993.
L. LeCam. Convergence of estimates under dimensionality restrictions. Ann. Statist., 1:38–53, 1973.
O. Lepski, A. Nemirovski, and V. Spokoiny. On estimation of the
l r norm of a regression function. Probability Theory and Related
Fields, 113(2):221–253, 1999.
Geoffrey McLachlan and David Peel. Finite mixture models. Wiley Series in Probability and Statistics: Applied Probability and
Statistics. Wiley-Interscience, New York, 2000.
Ankur Moitra and Gregory Valiant. Settling the polynomial learnability of mixtures of gaussians. Arxiv:1004.4223v1, 2010.
Eva Nogales. The development of cryo-em into a mainstream structural biology technique. Nature Methods, 13:24–27, 2015.

OPTIMAL RATES FOR MRA

[Pea94]

33

Karl Pearson. Contributions to the mathematical theory of evolution. Philosophical Transactions of the Royal Society of London
A: Mathematical, Physical and Engineering Sciences, 185:71–110,
1894.
[PWB+ 17] A. Perry, J. Weed, A. S. Bandeira, P. Rigollet, and A. Singer. The
sample complexity of multi-reference alignment. Manuscript, 2017.
[PWBM16] A. Perry, A. S. Wein, A. S. Bandeira, and A. Moitra. Messagepassing algorithms for synchronization problems over compact
groups. arXiv:1610.04583, 2016.
[Sad89]
B. M. Sadler. Shift and rotation invariant object reconstruction
using the bispectrum. In Workshop on Higher-Order Spectral Analysis, pages 106–111, Jun 1989.
[Sin11]
A. Singer. Angular synchronization by eigenvectors and semidefinite programming. Appl. Comput. Harmon. Anal., 30(1):20–36,
2011.
[SSK13]
B. Sonday, A. Singer, and I. G. Kevrekidis. Noisy dynamic simulations in the presence of symmetry: Data alignment and model reduction. Computers & Mathematics with Applications, 65(10):1535
– 1557, 2013.
+
[SVN 05]
Sjors H.W. Scheres, Mikel Valle, Rafael Nuñez, Carlos O.S.
Sorzano, Roberto Marabini, Gabor T. Herman, and Jose-Maria
Carazo. Maximum-likelihood multi-reference refinement for electron microscopy images. Journal of Molecular Biology, 348(1):139
– 149, 2005.
[Tro15]
Joel A. Tropp. An introduction to matrix concentration inequalities. Foundations and Trends R in Machine Learning, 8(1-2):1–230,
2015.
[TS12]
D. L. Theobald and P. A. Steindel. Optimal simultaneous superpositioning of multiple structures with missing data. Bioinformatics,
28(15):1972–1979, 2012.
[Tsy09]
Alexandre B. Tsybakov. Introduction to nonparametric estimation.
Springer Series in Statistics. Springer, New York, 2009. Revised and
extended from the 2004 French original, Translated by Vladimir
Zaiats.
[Ver17]
Roman Vershynin. High-Dimensional Probability. Cambridge University Press (to appear), 2017.
[Wal49]
Abraham Wald. Note on the consistency of the maximum likelihood
estimate. Ann. Math. Statistics, 20:595–601, 1949.
[WW84]
E. Weinstein and A. J. Weiss. Fundamental limitations in passive timedelay estimationpart ii: Wide-band systems. IEEE Trans.
Acoust., Speech, Signal Process., 32(5):1064–1078, 1984.
[WY16]
Y. Wu and P. Yang. Minimax rates of entropy estimation on large
alphabets via best polynomial approximation. IEEE Transactions
on Information Theory, 62(6):3702–3720, June 2016.
[ZvdHGG03] J. P. Zwart, R. van der Heiden, S. Gelsema, and F. Groen. Fast
translation invariant classification of hrr range profiles in a zero
phase representation. Radar, Sonar and Navigation, IEE Proceedings, 150(6):411–418, 2003.

34

BANDEIRA, RIGOLLET AND WEED

Afonso S. Bandeira
Department of Mathematics
Courant Institute of Mathematical Sciences
Center for Data Science
New York Univeristy,
New York, NY 10012, USA
(bandeira@cims.nyu.edu)

Philippe Rigollet
Department of Mathematics
Massachusetts Institute of Technology
77 Massachusetts Avenue,
Cambridge, MA 02139-4307, USA
(rigollet@math.mit.edu)

Jonathan Weed
Department of Mathematics
Massachusetts Institute of Technology
77 Massachusetts Avenue,
Cambridge, MA 02139-4307, USA
(jweed@mit.edu)

