PCA Initialization for Approximate Message Passing in Rotationally Invariant Models
Authors: Marco Mondelli
Ramji Venkataramanan
Venue: NeurIPS 2021 Poster (2021)
Keywords: Statistical Learning Theory
High-dimensional Statistics
Spectral Methods

Abstract: We study the problem of estimating a rank-1 signal in the presence of rotationally invariant noise--a class of perturbations more general than Gaussian noise.  Principal Component Analysis (PCA) provides a natural estimator, and sharp results on its performance have been obtained in the high-dimensional regime. Recently, an Approximate Message Passing (AMP) algorithm has been proposed as an alternative estimator with the potential to improve the accuracy of PCA. However, the existing analysis of AMP requires an initialization that is both correlated with the signal and independent of the noise, which is often unrealistic in practice. In this work, we combine the two methods, and propose to initialize AMP with PCA. Our main result is a rigorous asymptotic characterization of the performance of this estimator. Both the AMP algorithm and its analysis differ from those previously derived in the Gaussian setting: at every iteration, our AMP algorithm requires a specific term to account for PCA initialization, while in the Gaussian case, PCA initialization affects only the first iteration of AMP. The proof is based on a two-phase artificial AMP that first approximates the PCA estimator and then mimics the true AMP. Our numerical simulations show an excellent agreement between AMP results and theoretical predictions, and suggest an interesting open direction on achieving Bayes-optimal performance.

PCA Initialization for Approximate Message Passing in Rotationally Invariant Models
Marco Mondelli IST Austria marco www..mmmmneeeiitt...oooom ac.
Ramji Venkataramanan
University of Cambridge rv2850cam. ac.
Abstract
We study the problem of estimating rank-1 signal in the presence of rotationally invariant noise-a a class of perturbations more general than Gaussian noise. Prin- cipal Component Analysis (PCA) provides a natural estimator, and sharp results its performance have been obtained in the high-dimensional regime. Recently, an Approximate Message Passing (AMP) algorithm has been proposed as an alter- native estimator with the potential to improve the accuracy of PCA. However, the existing analysis AMP requires an initialization that is both correlated with the signal and independent of the noise, which is often unrealistic in practice. In this work, we combine the two methods, and propose to initialize AMP with PCA. Our main result is rigorous asymptotic characterization of the performance of this estimator. Both the AMP algorithm and its analysis differ from those previously derived in the Gaussian setting: at every iteration, AMP algorithm requires specific term to account for PCA initialization. while in the Gaussian case, PCA initialization affects only the first iteration of AMP. The proof is based phase artificial AMP that first approximates the PCA estimator and then mimics the true AMP. Our numerical simulations show an excellent agreement between AMP results and theoretical predictions, and suggest an interesting open direction on achieving Bayes-optimal performance.
two-
Introduction
We consider the problem of estimating rank-1 signal from a noisy data matrix. In the square symmetric case, the data matrix is modeled as
x= uu"n n
W E Ruxn
(1.1)
where u² E 12 is the unknown rank-1 signal, W E Rax* is a symmetric noise matrix, and > 0 captures the signal-to-noise ratio (SNR). In the rectangular case, we observe the data matrix
X = + W m
(1.2)
where 2 E Rhmh and v² E IR" are the unknown signals, and W E mu is rectangular noise matrix.
natural estimator of the signal in the symmetric case is the principal eigenvector of X (singular
vectors, in the rectangular case). The performance of this principal component analysis (PCA) estimator and, more generally, the behavior of eigenvalues and eigenvectors of models like (.22)
has been widely studied in statistics [27 49 and random matrix theory 2 3 10
22 30
If "z",, x² unstructured (e.g., they uniformly distributed
sphere), then it is not generally possible to improve on the PCA estimator. However, a broad range of applications, the unknown signals have some underlying structure, e.g., they may be sparse, their entries may belong certain
35th Conference on Neural Information Processing Systems (NcurIPS 2021).

---

set, or they may be modelled using a prior distribution. Examples of structured matrix estimation problems include sparse PCA
[7
28
6T] non-negative PCA 32 community detection under
and group synchronization 500 Since PCA is ill-equipped to capture the structure of the signal, we aim to improve on it using a family of iterative algorithms known as approximate message passing (AMP). AMP algorithms have two particularly attractive features: (i) they can be tailored to take advantage of prior information on the structure of the signal; and (ii) under suitable model assumptions, their performance in the high-dimensional limit is precisely characterized by a succinct deterministic recursion called state evolution [7] [3 26].
the stochastic block model [ 16 45
AMP algorithms have been applied to wide range of inference problems: estimation in linear models 88] 7 19 31 39), generalized linear models [5]
53 55]. and low-rank matrix The survey 21 provides unified description of AMP for these applications. Using the state evolution analysis, it has been proved that AMP
37 38 42
51
estimation with Gaussian noise 66 [7 23 29 34
44
achieves Bayes-optimal performance in some Gaussian models 17 18 44 and bold conjecture from statistical physics posits that AMP is optimal among polynomial-time algorithms.
We study rank-l matrix estimation in the setting where the noise matrix W is rotationally invariant. This is a much milder assumption than W being Gaussian: only imposes that the orthogonal matrices in the spectral decomposition of W are uniformly random, and allows for arbitrary eigen- values/singular values. Hence, W can capture a more complex correlation structure, which is typical in applications. For the models (I.I (I.)) with rotationally invariant noise, AMP algorithms were
derived in [14 48) and generalized in [20]. In particular, the AMP algorithm of [20] for the problem
produces estimates u¹ E R" as follows:
u' uff'')" f'=
Db..at uu
M 2.
(1.3)
The iteration is initialized with a pilot estimate u'. We can interpret method. Recall that the power method approximates the principal eigenvector of X using the iterative updates i' X / /| 7 Il. For each 1, the function u, can be chosen to exploit any structural information known about the signal (e.g., sparsity). The "memory" coefficients .b..,,...,b,,,
1.3 as a generalized power
have a specific form to ensure that the iterates (f', u u have desirable statistical properties captured by state evolution. A rigorous state evolution result for the iteration (L33) is established in [20], but the algorithm and its analysis require an initialization u² that is correlated with the unknown signal and independent of the noise W. In practice, one typically does not have access to such initialization.
Main contribution. In this paper, we propose an AMP algorithm initialized via the PCA estimator, namely, the principal eigenvector of X for the square case 11.1) and the left singular vector of X for the rectangular case (1.2). Our main technical contribution is a state evolution result for this AMP algorithm, which gives rigorous characterization of its performance in the high-dimensional limit. The challenge is that, as the PCA initialization depends on the noise matrix W, one cannot apply the state evolution machinery of (20). To circumvent this issue, our key idea is to construct and analyze nwo-phase artificial AMP algorithm. In the first phase, the artificial AMP performs a power method approaching the PCA estimator; and in the second phase, it mimics the behavior of the true AMP, We remark that the artificial AMP only serves as a proof technique. Thus, we can initialize it with vector correlated with the signal u* and independent of the noise matrix W, which allows to analyze it using the existing state evolution result.
Our analysis is tight in the sense that our AMP algorithm can be initialized with PCA whenever the PCA estimate has strictly positive correlation with the signal. This requires showing that, when PCA is effective, the state evolution of the first phase of the artificial AMP has a unique fixed point.
To obtain such a result, we exploit free probability tools developed in [T0 I.I. The agreement between the practical performance of AMP and the theoretical predictions of state evolution is demonstrated via numerical results for different spectral distributions of W. Our simulations also show that the performance of AMP-as well as its ability to improve upon the PCA initialization crucially depends on the choice of the denoising functions u, in the algorithm. Thus, the design of a opyes-optimal AMP remains exciting avenue for future research.
Related work. The asymptotic Bayes-optimal error for low-rank matrix estimation has been precisely characterized for Gaussian noise 33 but remains an open problem for rotationally invariant noise. An AMP algorithm with PCA initialization was proposed in 44 for the Gaussian setting, and it was shown to be Bayes-optimal for some signal priors. A recent paper by Zhong et

---

al. [60] shows how AMP with PCA initialization can be used for estimating the top-k principal components applications such as high- nigh-dimensional genomics datasets. The authors use an empirical Bayes method to determine a joint prior distribution for the k principal components, and assuming a Gaussian noise model, employ an AMP algorithm tailored to the prior to improve the PCA estimates of the principal components.
Both our AMP algorithm and proof technique differ significantly from those for Gaussian noise. When W is Gaussian, the PCA initialization affects only the first iteration of AMP. In contrast, for more general noise distributions, the AMP algorithm and its associated state evolution require a correction term at every iteration to account for the PCA initialization. This is due to the fact that, while AMP has single memory term in the Gaussian case, more general noise distributions lead to more involved memory structure, as in (.3)). As regards the proof technique, the argument of consists of decoupling the PCA estimate from the bulk of the spectrum of X. In contrast, our approach is based on two-phase artificial AMP algorithm. This technique has proved successful in the context of generalized linear models 41 42 albeit for Gaussian measurements. Other extensions
44
of AMP beyond the Gaussian setting include Orthogonal AMP [36 [6]. Vector AMP [24 25 52 54),
convolutional AMP [57] and Memory AMP 355. These algorithms have been derived specifically for linear generalized linear models, and extending them (with practical initialization method) to low-rank matrix estimation is an interesting research direction.
Finally, we mention the recent independent work of Zhong et al. [59], which appeared after the original submission of our paper. This work generalizes AMP with PCA initialization to the problem of estimating rank-k matrices in rotationally invariant noise, for & > 1. We remark that, in order to prove a state evolution result for AMP initialized with PCA, in (59 is assumed that the signal strength is sufficiently large. In contrast, our result holds for any signal strength such that the PCA method is effective, but we require the free cumulants of the noise matrix to be non-negative. We also note that, when the signal strength is large, the assumption on the free cumulants can be automatically satisfied (see the footnote on page 6).
2 Preliminaries
Notation and definitions. Given a E IR, we define (a)+ = max(a, 0). Given two integers S j. we define [ (i.....j.. If i > j, then [i, denotes the empty set; products over the empty set are taken to be 1. Given a vector æ F". we denote by ||>||| its Euclidean norm and by (a) its empirical mean, i.e., (ar) Σ = Tiz- The empirical distribution of x (f1,-.., is given by
W Ei=1 82,, where 82, denotes Dirac delta mass on Iy, The notation x - X denotes convergence of the empirical distribution of æ to the random variable X in Wasserstein distance at all orders.
Given a symmetric square matrix A E R=xx, we denote by A1 (A) > )² 2 ... 2 An (A) its eigenvalues sorted in decreasing order. Given rectangular matrix A E Rm with m n, we denote by σ (A) 2 (A) 2 ..... om (m() (A) its singular values sorted in decreasing order.
Rank-1 estimation Symmetric square matrices. Consider the problem of estimating the signal u E R" from the data matrix X in (I.I) We assume that W is rotationally invariant in law, i.c. W = AOA where A- = diag(A) (A) is a diagonal matrix containing the eigenvalues of W and O is a Haar orthogonal matrix independent of A. As + oc, we assume that the empirical distributions of and u² satisfy
A W and u W U.,
(2.1)
where A and U. represent the limiting spectral distribution of the noise and the prior on the signal, respectively. We take ||u||| vn so that E{U?)- = lim ||10"1|2= = 1. We assume that the moment EE(U?++> < Do for some E > We also assume that has compact support, and denote by 6 the supremum of this support. We denote by {<* Jb21 the free cumulants corresponding to the moments {ma J of the empirical eigenvalue distribution of X excluding its largest eigenvalue, i.., my >(x)) )x) (for details, see (A.1) in Appendix A). The assumption (2.
implies that, as n + 00, my + m E{A*} and Ke - 0o where (m J22 and (Ne }A21 are respectively moments and free cumulants of A.
PCA TTA-Symmin -yam-Smmmt square matrices. Let uPca be the principal eigenvector of X, and define 0s. = 1/G(bt), where G(z) = Ef(2 -11 is the Cauchy transform of A, and G(b+) lim,-+6+ G(2).
3

---

Then, for > ag. (X) a.s, (1/a) and (X) a.S, b. where G-1 is the inverse of G;; see Theorem 2.1 in 10 Furthermore, Theorem 2.2 in 10 gives that, for a > Cls.
(UPCA, 2²)² a.s. Pa a G 0.
(2.2)
71
In words, above the spectral threshold dg., the principal eigenvalue of X escapes the bulk of the spectrum and its associated eigenvector becomes strictly correlated with the signal u².
Rank-1 estimation Rectangular matrices. Consider now the problem of estimating the signals "u" E 1n and v² E I given the rectangular data matrix X in (1.2). Without loss of generality, we assume that m n (if m >> one can just exchange the role of u² and ex and consider X in place of X). We assume that W is bi-rotationally invariant in law, i.e., W = =0=AQ AQ,, where = diag(A) A) a m n diagonal matrix containing the singular values of W, and O, Q are Haar orthogonal matrices independent of one another and also of A. As + OC, assume that À " A.
2 U., V- W V, and m/n + 7, for some constant y (0, 1]. We take ||u||| vm and || = v so that [E{U]} = [V(V?} = 1. As before, b < oo is the supremum of the compact support of A, and U.. V,, are assumed to have finite (2 + s)-th moment for some >> 0. To analyze PCA using the framework in [II], we also assume that the entries of u and v" are i.i.d., and their law has zero mean and satisfies log- Sog--oooeev inequality. We denote by (K22) the rectangular free cumulants associated to the even moments {m22)421, with m2k
Li= (X) (for details, 14.II1-A.11 in Appendix A). Furthermore, as n, m co, mgg + mag [{{ } and K2k + 2k where {m cumulants of A.
J J2 and { 222 are respectively even moments and rectangular free
PCA Rectangular matrices. Denote by UPCA and UPCA the left and right principal singular vectors of X, and define ã₂ = 1/ where D(=) d(=)((z), <0=) B==/(=² )²) ((z) >()))= + 2)/2, and D(b))= lim_-b- D(z). Note that the singular value of the rank-one signal u u is ā ê a/v7. Then, for ā > āg. 01 T1(x (x) a.s,
(1/ā²) and o2(X) a.s, b: see Theorem 2.8 in [T]. Furthermore, Theorem 2.9 in [II] gives that, for á dg..
(uPCA. 2² m APCA a²D((D-((1/a3)
-20(D-
0,
(2.3)
(UPCA, a.s. PPCA (2D((D-((1/a))
-20(D=
0.
(2.4)
n
In words, above the spectral threshold â,, the principal singular value escapes from the bulk of the spectrum and the left/right principal singular vectors become correlated with the signal u //"".
3 PCA Initialization for Approximate Message Passing
3.1 Symmetric Square Matrices
We consider family of Approximate Message Passing (AMP) algorithms estimate u² from X 2u'u" T W. We initialize using the PCA estimate UPCA:
u VnuPCA, f f'= b
(3.1)
with b1,1 E 1 Then, for t 2, the algorithm computes
u' = u,(''"=))) f Xuu Σ bu
(3.2)
where the memory coefficients (b., (bEL))) given by b K1, and
b II (
Σ k
be, = 1i+11 II (u((F(-))), for (1- j) E [2, 1).
++1
(3.3)

---

Here, the function u, IR + R is continuously differentiable and Lipschitz, it is applied component- wise to vectors, i.c., u, (f( algorithm in (3.1) is similar the one in [20] Sec. 3.1] (and the ones in [14] 48), with the main differences being the initialization u and the formula for the memory term b We highlight
(r((ft-))), and u denotes its derivative. The AMP ur "o.ff"")......
that the algorithm does not require the knowledge of a or of the noise distribution. In fact, a can be consistently estimated from the principal eigenvalue of X via á = (G(A1 (X))) ) Furthermore, one can compute the moments {m& }221 of the empirical eigenvalue distribution of X (excluding its largest one) and, from these, deduce the free cumulants {fi J:21-
The asymptotic empirical distribution of the iterates f for 2 1, can be succinetly characterized via a deterministic recursion, called state evolution, and expressed via a sequence of mean vectors Mk = (fe (Je))_] and covariance matrices K = (e) For K = 1, set
aPa and 2.2 Then define pk+1 Ex+1 from Pk Σx as follows. Let
01 = a²2 (1 -)2), with Pa given in
,F = UxU. + (Z Zk), where (Z1, Zk) ~ ~O( Σx),
(3.4)
Fi
U. = (F(-n) for 2 S K + and U. = a
for 00 <t t $1.
(3.5)
Then, the entries of μk+1 are given by 14 - aB(U,U,) } for E of Ex+1 can be expressed via the following formula, for s, E [1, K +1]
[1, K + 1]. Furthermore, the entries
(k- (k-+41))++
1
( II B(U((--)))) (II B(u((P-)))) B(Us-JUL--h == i=max(2++1-X)
Oxt
+++
j 0 k=
(3.6)
Our main result, Theorem 1 shows that for 2 1, the empirical joint distribution of the en- tries of (u", f f') converges in Wasserstein distance W2 to the law of the random vector (U., F1.... Fi..... Ft). We provide a proof sketch in Section 5 and the complete proof is deferred to Appendix B This result is stated in terms of Lisuch-Lipschitz test functions. A function ei R + R pseudo- Lipschitz of order 2, i.e., & E PL((), if there a constant C' > such that
|| (r) (v)l C(1 + ||>:||
|lll)) - yll-
(3.7)
The equivalence between convergence in terms of PL(2) functions and convergence in W2 distance follows from [58 Definition 6.7 and Theorem 6.8).
Theorem 1. In the square symmetric model (I.I assumee that a > Ocs and that the free cumulants of order and higher are non-negative, i.e., K& 2 0 for k 2 2. Consider the AMP algorithm with PCA initialization in (3.1) (.2), with continuously differentiable and Lipschitz,fanctions U₁ : - R. (Without loss of generality, assume that (u', .u UCA) 2 0.)
Then, for 2 1 and any PL(2) function 222+22 + R. we almost surely have:
lim
Σ (ut ut,..............
1+1 .......
19) = E (V.UU, UI,... U2t1, Fi...., F.)). (3.8)
70+00 72
1=1
where UI,... U441 and Fi,.... ......, Fe are defined in 3.4
Assumptions of the theorem. The basic assumption that the noise matrix is rotationally invariant is rather mild as it allows for arbitrary cigenvalue distributions. The assumption a > ag ensures that the PCA initialization is correlated with the signal. This condition is necessary and sufficient for PCA to be effective: under the additional requirement that G'(b) = -o, we have that, a < a'g. then the normalized correlation between UPCA and u* vanishes almost surely; see Theorem 2.3 of (TO]. Conversely for a > ag, the asymptotic correlation is strictly non-zero and given by (2.2).
Non-negativity of free cumulants: The assumption that Kg 2 for k: 2 2 appears to be an artifact of the proof technique. As detailed in the proof sketch in Section 5 this assumption needed to show that the state evolution of the artificial AMP in the first phase has unique fixed point. We expect our approach to generalize to any limiting noise distribution A with compact support, and defer such generalization to future work. In support of this view, the simulations of Section 4 verify the
5

---

claim of Theorem setting where the free cumulants of A have alternating signs (corresponding
to an eigenvalue distribution A Uniform[11/2.. 1/2, see Figs. [b] ld and 2 2d). Finally, we
remark that, if W follows a Marcenko- Pastur distribution (W~ = AAA' where A has i.i.d Gaussian entries), then the free cumulants of A are all equal and strictly positive; see (40] Chap. 2. Exercise
111. Thus, the assumption of Theorem T hholds for noise distributions that are sufficiently close to the Marcenko-Pastur one, or for sufficiently large values of the signal-to tortall-t-noisee ratio a I
Continuous differentiability and other technical assumptions: The assumption that u, is continuously differentiable can be weakened to: (i) u, being differentiable almost everywhere, and (ii) satisfying a mild non-degeneracy condition (Assumption 4.2(e) in 2 In this way, we can cover most practically relevant choices of u, such as soft thresholding and ReL.U. Theorem 1 alsoo requires the technical assumptions in (2.T) and the text below it: convergence of the empirical distributions of the signal and of the eigenvalues of the noise matrix; boundedness of the (2 + e)-moment of the signal; and compact support of the spectrum of the noise matrix. We regard these technical assumptions as minor, and remark that they are quite standard in the literature. For the rectangular case, we also need the additional assumption that the law of the signal is zero mean and satisfies a log- Sobolev inequality, which is necessary to apply the framework in [10].
How PCA initialization influences AMP. The form of the memory coefficient b (3.3) reflects the PCA initialization of the AMP iteration. PCA initialization can be interpreted as the result of a first AMP phase with linear denoisers (see the proof sketch in Sec. [ The coefficient b, multiplying the initialization i represents the cumulative effect of this first AMP phase leading to the PCA estimate. The main differences from the AMP algorithm in 20 (where the initialization is independent of W) are the expressions for the coefficient b and the state evolution parameters Txtt (compare (3.3) and (3.6) in this paper with (1.15) and (1.17) in (20). One can interpret the new form of b and Øs.t as memory of the PCA initialization. For the special case of Gaussian noise, the spectral initialization only affects the first iteration of AMP [A.. This is due to the fact that, while in a rotationally invariant model the AMP iterate at step t depends on all previous iterates, in the Gaussian case it depends only on the iterate at step t 1.
Choice of (.). Theorem I hhlld for any choice of denoisers {ur} that are Lipschitz and contin- uously differentiable. Indeed, our analysis shows that by picking (f) = -/a, AMP just returns the PCA estimate; see the proof sketch in Section 5 If some structural information about the signal is available (e-g.. sparsity), denoisers that take advantage of this structure can give substantial im- provements over PCA. Thus, a key question is how to optimally select the u, 's. Theorem I tells us that the empirical distribution of f converges to the law of p,U, + Vo,, Z. for 2 ~ (O 1) and independent of U.. Hence, the quality of the estimate at each iteration ł is governed by the SNR Pe = / Consider running the algorithm for i iterations, and let u T+1 estimate. Then, for each / E [2, J. the Bayes-optimal choice for uy is the one that maximizes pt. i.e., the SNR for the next iteration. In the case of Gaussian noise (44) the maximum is achieved by the posterior mean Un (f) = EQU= | ptU_- + Vor,t2 = f}. For rotationally invariant noise, this choice minimizes the mean- square error 1 l/u' uz 112 (for fixed u1,.... ................................ u but it does not necessarily maximize the SNR pt. We provide an example of this behavior in the simulations reported in Section 4 Therefore, the optimal strategy would be to choose functions Ug...., uy to maximize the SNRs P Pi. and then in the final iteration, to pick u to minimize the desired loss. Note that u, depends on the previously chosen functions #1 ... u in complicated way, due the definition of Fze in (3.6). Thus, finding U, that maximizes the SNR Pt remains an outstanding challenge. Finally, we remark that though only consider one-step denoisers in this paper, Theorem I can be readily extended to cover denoisers with memory, i.e., those of the form u, (f'. ..... 1
( be the final
3.2 Rectangular Matrices
We now present an AMP algorithm to estimate u* and v² from the m x n data matrix X a u u- + W. We initialize the algorithm using the PCA estimate UPCA:
One can add an independent artificial noise matrix with Marcenko- Pastur distribution to the data in order to make the required free cumulants non negative, and the result would hold for a greater than the new spectral threshold.

---

u² = ='=vmuur UPCA. g'=
(1+7)
1
Σ
x Xuu,
v² vi='! - g (3.9)
Kii
a2
im
Then, for & 2 1, we iteratively compute:
f'= Xx' Σ an, = =u4="" g Xu
Σ be w 114 (2
(3.10)
Here, U Ve+1 I + RR are continuously differentiable Lipschitz functions that act component-
wise on vectors. We define a1,1 = a Ei=1 fi2s (222) and for 2 2:
ar,1 = (Mgg II [u "M--
Σ (2(++))
7 (22)
)
(3.11)
-,, - (Wgg)) (M('))
II (M
for j) E [2,1].
(3.12)
imt-j+1
Furthermore, for 1.112 1,
be+1,1 = 7 (f')) MMM (Mar
aa )
7 (2(0+0)
K2r +
(3.13)
II (O))))1("")> K2y - (=t2--)
for (1+1- E [2,4] (3.14)
Similarly to the square case, a can be consistently estimated from the largest singular value of X via a = (o (x)))=I, and the rectangular free cumulants (K2))x>) can obtained from the moments of the empirical distribution of the singular values of X (excluding its largest one).
The asymptotic empirical distributions of the iterates (f* g') can be characterized via a state evo- lution recursion, which specifies sequence of mean vectors M* (Pe)ee Vx = (Vt)rej1,x) These are iteratively defined, starting with the initialization 140 av APCA and 00,0 a²(1 - APCA), where APCA is given by Having defined Pk Σk, / K, K2 let
and covariance matrices Σx = (a..)s, (o.K 52k = (w,,))
(Fo...., = ux U. + (Yo, .., Yk), where (Yo, Yx) NO, Ex),
(3.15)
U. = 1.(FF-1) for 2 < K 11, and U-= Fo for 00 1
(3.16)
a
.......,G)) Pg V. + (Z1, (Z,, Zk), where (Zi,..., [k) ~ ((0,5kk)
(3.17)
Vi v.G))) for < 1 < K +1, and = _ for <151.
(3.18)
Given Hx and Ex. the entries of v are given by V (E(U,,,) (for E [1, K + 1),, and the entries of 22++1 (for s + 1, E [1, + 1)) are given by
D0
15 +) +
s+1
W241,e41
II
x; x
II
X ---1-1)
-2
x2 E U++++1 14441-1-1411 + 22+
T
X
(3.19)
Here, we define x; (Fu( }) if 2, and x; 1/a otherwise: y, [() if and = 7//x otherwise. We note that wi1 is computed by solving the linear equation obtained by setting s = in (3.19) (see (C.6)) Next, given v 1 and Ω for some 2 the entries of #x+1
are /4:- 2E
V.V.) } (for / E 10, K + 1)), and the entries of Σx+1 (for 8, E [0, K + 1)) are
00
U 1)+ 1++ (k-++1)+
Ts11
II
X xv))
X y:) (=
II
/=0 k=
m=ax((,,+111)
N2 B(V (Va-,Ve-k) + xc 2 B B(Us--Ue-R yt-
(3.20)
7

---

Our main result for the rectangular case, Theorem 2 shows that for t 2 1, the empirical joint distribution of the entries of (u', (u" .....f' converges in Wasserstein distance W2 to the law of the random vector (U,,Fi,.... F.). Similarly, the empirical joint distribution of the entries of (v" """""" ..' converges to the law of (V., Gr...., G.). The proof is given in Appendix C As in the square case, we state this result in terms of pseudo-Lisschitz test functions.
Theorem 2. In the rectangular model (1.2) assume that ā > Õs and that High 2 0 for k. 2 1. Consider the AMP algorithm with PCA initialization (39 1.10-.11.. with continuously differentiable and Lipschitz,uuntions functions Ut, Vt : R-- + R. (Assume without loss of generality that (w*, u uCCA) 2 0.)
Then, for & 2 1 and any PL(2) functions vb 122++2 - R and y : R22+1 + R. we almost surely have:
lim 11-+00 m
1
Σ
(u uz
[+1 1. u/
f1) = E (V.U,, UI, U41 Fi...... Fi))
(3.21)
1
lim
p(o ui....,,y, of 9t) = B Vi.... VG, G...,...,GG.). G.)). (3.22)
1=1
where (U1....,II+:) (Fi,..., Ft). (Vi,..., Vi) and (Gi, Gt) are defined as in 3315 3.18
The condition ā > 0g is necessary and sufficient for PCA to be effective: under the additional requirement that o'(bt)= = -oo, if ā < ag. then the normalized correlation between UPCA and u* vanishes almost surely, see [1I] Theorem 2.10]. Comments similar to those at the end of Section 3.] can be made about (i) the requirement that the rectangular free cumulants are non-negative, (ii) the
effect of the PCA initialization on AMP, and (iii) the choice of the denoisers Ut, Vt.
4 Numerical Simulations
We consider the following settings: (i) square model (I.I) with Marcenko-Pasturr noise, i.e., W " LAA" E Rnx", where the entries of A "R" are i.i.d. standard Gaussian, see (a) in the figures;
(ii) square model (1.T) with uniform noise, i.e., O AO E RRax,, where O is a Haar orthogonal matrix and the entries of A are i.i.d. and uniformly distributed in the interval [- /21 1 1/2), see (b) in the figures; (iil rectangular model (1.2) with uniform noise, i.e., W O AQ Rm**, where
O, Q are Haar orthogonal matrices and the entries of A are i.i.d. and uniformly distributed in the interval [0, 1]. see (c)-(d) in the figures.
In the simulations, a is estimated from the largest cigenvaluc/singular value of X. Furthermore, the free cumulants Kg (K22 in the rectangular case) are replaced by their limits ky." (Mo2. resp.), which are obtained as follows. For (a), all the free cumulants of A are equal to e _ p/n, i.e., x² = e for k 2 I, see [40] Chap. 2, Exercise 11]. For (b), the odd free cumulants of A are 0 and the even ones are given by K. = (2n)), where B2n, denotes the 2n-th Bernoulli number. For details, see the derivation of 10 in Appendix A For (c)-(d), the even moments of A are given by m 1/(k+1] 1) and, from these, we numerically compute the rectangular free cumulants via 12 in Appendix A Furthermore, the spectral threshold for the setting in (a) is ag
+ vc: for (b), ag = 0: and for (c)- a_- = O. In (a). we set n = 8000 and c = 2: in (b), we set n = 4000; and in (c)-(d), we set n 8000 and ? 1/2. The signal u² has Rademacher prior, i... its entries are i.i.d. and uniform in (-1,1), 1). In the rectangular case, the signal v² has a Gaussian prior; i.e., it is uniformly distributed on the sphere of radius vn. Given these priors, u, is chosen to be the single-iterate posterior mean denoiser given by u. (x)) tanh(P: 2z/ x1/o1,.). where μ4 and Tt. are the state evolution parameters; these are replaced by consistent estimates in the simulations. For the rectangular case, we choose vx(x) x. Each experiment is repeated for Perials 100 independent runs. We report the average and error bars at 1 standard deviation.
Figure I compares the performance between the proposed AMP algorithm with PCA initialization (PCA+AMP) and the theoretical predictions of state evolution (SE), for two different values of a. On the r-axis, we have the number of iterations of AMP, and on the y-axis the normalized squared correlation between the iterate and the signal. As a reference, we also plot the performance of PCA as a horizontal line. We observe an excellent agreement of AMP with state evolution, even in the settings (b)-(c)-(d) where the free cumulants (resp. rectangular free cumulants) are alternating in sign. This supports our conjecture that Theorems 1 hold for more general noise distributions.
8

---

PONAARI, an
PINTASR
PCh, a-
POL
EK
5K
75A-ME
0.
0.2
Reraaioos
Reeration
Berakion
(a) Square, MP.
(h) Square, uniform.
(c) Rectangular, left.
(d) Rectangular, right.
Figure 1: Comparison between AMP with PCA initialization and the related state evolution (SE). The plots show the normalized squared correlation between iterate and signal, as function the number of iterations.
0.5
0.5
0.
04
34
0.5
6:
Ai
Guaria ADE
PGL Gradta AUD Foll sooert
(leag)
09900
(a) Square, MP.
(b) Square, uniform.
(c) Rectangular, left.
(d) Rectangular, right.
Figure 2: Comparison between AMP with PCA initialization and the PCA method alone. The plots show the normalized squared correlation between signal the estimate (PCA, or AMP+PCA), a function of a.
In Figure 2 we run PCA+AMP until the algorithm converges, and we compare the results with (i) the AMP with PCA initialization developed [44] which assumes that the noise matrix Gaussian (with the correct variance), and (ii) the PCA method alone, as a function of the SNR a. For Marchenko- Pastur noise (setting (a)), PCA+AMP always improves upon the PCA initialization. However, this is not the case when the eigenvalues/singular values of the noise are uniformly distributed (settings (b), (c) and (d)). In fact, we observe a phase transition phenomenon: below a certain critical a, AMP converges to trivial fixed point at O, while PCA shows positive correlation with the signal; above the critical a, PCA+AMP is no worse than PCA. This is due to the sub-optimal choice of uy; recall the discussion on 6 We observe vector (setting (d)), as the prior v² is Gaussian, which case we expect the PCA estimate to be optimal. The interesting behavior demonstrated in Figure motivates the study of the optimal choice for u, V₁ in future work. We also note that, in settings (d). d₂ O, which means that the PCA estimator has non-zero correlation with the signal for all O. However, for 0.1, this correlation remains rather small. Finally, we highlight that our proposed rotationally invariant PCA+AMP always improves upon the Gaussian PCA+AMP. In general, this performance gap will be significant unless the sequence of free cumulants k (N in the rectangular case) decays quickly. For Marchenko- noise, the free cumulants are all equal, and thus the performance gap is significant. If the cigenvalues/singular values of the noise are uniform, then the sequence of free cumulants decays rapidly and the performance gap small.
improvement for the estimation of the right singular
5 Proof Sketch: Symmetric Square Matrices
We consider the following artificial AMP algorithm, whose iterates are denoted by ū " for L 2 1. We initialize with à' Pa /1 and f X ns u' Here, n is standard Gaussian and Pa. is the normalized (limit) correlation of the PCA estimate given in (2.2). We note that this initialization is impractical, as requires the knowledge of the unknown signal u². However, this is not an issue since the artificial AMP serves only as a proof technique. (The true AMP (.2) used for estimation uses the PCA initialization in (3.T)) The subsequent iterates of the artificial AMP are defined in two phases. In the first phase, which lasts up to iteration (T' 1), the functions defining

---

the artificial AMP are chosen so that u T+ is closely aligned with the eigenvector UPCA as T + co. In the second phase, the functions are chosen match those in the true AMP.
The artificial AMP initialization ū' is chosen such that it has non-zero asymptotic correlation with the signal u² Indeed, when the signal prior has zero mean, a random initialization (independent of u*) would be asymptotically uncorrelated with the signal; consequently, the first phase of the artificial AMP would get stuck at a trivial fixed point and the iterates would not be guaranteed to converge to the principal eigenvector. We ensure that this does not happen by defining the initialization üt to be linear combination of the signal and Gaussian noise.
First phase.
For 2 2ts (T + 1), the artificial AMP iterates are
/a: ū' /a; x'u
b
(5.1)
where br,r-j = Kj-1a-1, for (t- j) E [1, We claim that, for sufficiently large T, ū uT+1 approaches the PCA estimate UPCA, that is, lim limn- vz - I
/TupCal || = 0. This
result proved in Lemma B.7 in Appendix [B33 We give a heuristic sanity check here. Assume that T+1 - 0 0.
the iterate ū T++1 converges limit it in the sense that rom limp- /
Then, from 5.I. the limit u satisfies
ū 00 - -x --_-- Σ Ki
a
Kr
Xu
(5.2)
which means that t0 is eigenvector of X. Furthermore, by using known identities in free probability (see (A.4) and (A.6). the eigenvalue a + K; (+) can be re-written as G- (1/aa) Recall that, for a > ag X exhibits a spectral gap and its largest eigenvalue converges to G G- (1/a). Thus, u must be aligned with the principal eigenvector of X, as desired.
key step in our analysis is to show that, T oc, the state evolution of the artificial AMP in the first phase has the unique fixed point (ji = apa, ã = a²2 P2)). This is established in Lemma B.2 proved Appendix B.2 The proof follows the approach developed Section 7 of 20) However, the analysis of (20] requires that is sufficiently large, while our result holds for all > Og. Our idea to exploit the expression of the limit correlation between the PCA estimate and the signal. In particular, we prove that, when the PCA estimate is correlated with the signal, state evolution is close to limit map which is a contraction. The price to pay for this approach the requirement that the free cumulants are non-negative.
Second phase. The second phase of the artificial AMP is designed so that its iterates (üTth, T+ j(t are close to (u", f²), for l 2 2. For t 2 (T 2), the artificial AMP computes:
ū'=u-T(f"
j'= ='=X'
Σ b. u
(5.3)
Here. the functions {ux} 222 are the ones used in the true AMP I z (T + 2) are given by:
The coefficients (b..)} for
(T+1-
bu = K1, br,e-; = Kjt1
II
( )) E 1].
m= (-1+1,7+2)
(5.4)
Since the artificial AMP is initialized with à that is correlated with u² and independent of the noise matrix W, a state evolution result for it can be obtained directly from (20, Theorem 1.1]. We then
show in Lemma [.& in Appendix that the second phase iterates in (5.3) are close to the true AMP iterates in (3.2). and that their state evolution parameters are also close. This result yields Theorem [
as shown in Appendix [B.5] The complete proof of Theorem 2 (rectangular case) given in Appendix
We describe the artificial AMP for this case along with proof sketch in Appendix C.I
10

---

Acknowledgements
M. Mondelli would like to thank László Erdös for helpful discussions. M. Mondelli was partially supported by the 2019 Lopez-Loreta Prize. R. Venkataramanan was partially supported by the Alan Turing Institute under the EPSRC grant EP/N510129/1
References
[1] Emmanuel Abbe. Community detection and stochastic block models: recent developments, The Journal of Machine Learning Research, 653111 2017.
[2] Jinho Baik, Gérard Ben Arous, and Sandrine Péché. Phase transition the largest eigenvalue for nonnull complex sample covariance matrices. Annals of Probability, pages 1643-1697, 2005.
[3] Jinho Baik and Jack Silverstein. Eigenvalues of large sample covariance matrices of spiked population models. Journal of Multivariate Analysis, 9766):1382-1400 2006.
[4] Jean Barbier, Mohamad Dia, Nicolas Macris, Florent Krzakala, Thibault Lesieur, and Lenka Zdeborová. Mutual information for symmetric rank-one matrix estimation: proof of the replica formula. In Neural Information Processing Systems (NeurIPS), pages 424 432, 2016.
[5] Jean Barbier, Florent Krzakala, Nicolas Macris, Léo Miolane, and Lenka Zdeborová. Optimal errors and phase transitions high- high-dimensional generalized linear models. Proceedings of the National Academy of Sciences, 1 (16(1))55551-4600 2019.
[6] Jean Barbier, Nicolas Macris, and Cynthia Rush. All-or- -otthiggg statistical and computational phase transitions sparse spiked matrix estimation. In Neural Information Processing Systems (NeurIPS), 2020.
[7] Mohsen Bayati and Andrea Montanari. The dynamics of message passing on dense graphs, with applications compressed sensing. IEEE Transactions on Information Theory. 57:764-785. 2011.
[8] Mohsen Bayati and Andrea Montanari. The LASSO risk for Gaussian matrices. IEEE Transac- tions on Information Theory, 18:1997-2017, 2012.
[9] Florent Benaych-GGeorges. Rectangular random matrices, related convolution. Probability Theory and Related Fields, 144(3-4):411-555, 2009.
[10] Florent Benaych-GGooges and Raj Rao Nadakuditi. The eigenvalues and eigenvectors of finite, low rank perturbations of large random matrices. Advances in Mathematics, 227(1):494-521, 2011.
[11] Florent Bengych-Georges and Raj Rao Nadakuditi. The singular values and vectors of low rank perturbations of large rectangular random matrices. Journal of Multivariate Analysis, 111:120-135, 2012.
[12] Patrick Billingsley. Probability and measure. John Wiley & Sons, 2008
[13] Erwin Bolthausen. An iterative construction of solutions of the TAP equations for the SiiiingtoonKKkpatrriick model. Communications in Mathematical Physics, 325(1):333-366, 2014.
[14] Burak Çakmak and Manfred Opper. Memory-free dynamics for the Thouless-Anderwwn Palmer equations of Ising models with arbitrary rotation-invarian ensembles of random coupling matrices. Physical Review E, 99(6):062140, 2019.
[15] Mireille Capitaine, Catherine Donati-Martin, and Delphine Féral. The largest eigenvalues of finite rank deformation of large Wigner matrices: convergence and nonuniversality of the fluctuations. The Annals of Probability, 37(1)1147,,2 2009.
[16] Yash Deshpande, Emmanuel Abbe, and Andrea Montanari. Asymptotic mutual information for the balanced binary stochastic block model. Information and Inference, 6. 2016.
11

---

[17] Yash Deshpande and Andrea Montanari. cnoommaiiontthoretiiin optimal sparse PCA. In IEEE International Symposium on Information Theory (ISIT), pages 2197-2201, 2014.
[18] David L. Donoho, Adel Javanmard, and Andrea Montamari. Information-theortiialll optimal compressed sensing via spatial coupling and approximate message passing. IEEE Transactions on Information Theory, (0(11))7434 7464, Nov. 2013.
[19] David L. Donoho, Arian Maleki, and Andrea Montanari. Message Passing Algorithms for Compressed Sensing. Proceedings of the National Academy of Sciences, 106:18914-18919, 2009.
[20] Zhou Fan. Approximate message passing algorithms for rotationally invariant matrices. arXiv:2008. 1892, 2020.
[21] Oliver Y. Feng, Ramji Venkataramanan, Cynthia Rush, and Richard J. Samworth. A unifying tutorial on Approximate Message Passing. arXiv.2105.02180, 2021.
[22] Delphine Féral and Sandrine Péché. The largest eigenvalue of rank one deformation of large Wigner matrices. Communications in mathematical physics, 272(1 ): 111--182 2007.
[23] Alyson K. Fletcher and Sundeep Rangan, Iterative reconstruction of rank-one matrices in noise. Information and Inference: A Journal of the IMA, 7(3):531-562, 2018.
[24] Cédric Gerbelot, Alia Abbara, and Florent Krzakala. Asymptotic errors for high-dimensional convex penalized linear regression beyond Gaussian matrices. In Conference on Learning Theory (COLT), pages 1682-1713, 2020.
[25] Cédric Gerbelot, Alia Abbara, and Florent Krzakala. Asymptotic errors for teacher- student convex generalized linear models (or: How to prove Kabashima's replica formula). arXiv:2006.06581 2020.
[26] Adel Javanmard and Andrea Montanari. State evolution for general approximate message passing algorithms, with applications to spatial coupling. Information and Inference, pages 115-144, 2013.
[27] lain M. Johnstone. On the distribution of the largest eigenvalue in principal components analysis. Annals of Statisties, pages 295-327, 2001.
[28] lain M. Johnstone and Arthur Yu Lu. On consistency and sparsity for principal components analysis in high dimensions. Journal of the American Statistical Association, 104(486), 2009.
[29] Yoshiyuki Kabashima, Florent Krzakala, Mare Mézard, Ayaka Sakata, and Lenka Zdeborová. Phase transitions and sample complexity in Bayes-optimal matrix factorization. IEEE Transac- tions on Information Theory, 62(7):4228-422) 42658,8-446.220 2016.
[30] Antti Knowles and Jun Yin. The isotropic semicircle law and deformation of Wigner matrices. Communications on Pure and Applied Muthematics, 2013.
[31] Florent Krzakala, Mare Mézard, Francois Sausset, Yifan Sun, and Lenka Zdeborová. Probabilis- tic reconstruction in compressed sensing: algorithms, phase diagrams, and threshold achieving matrices. Journal of Statistical Mechanics: Theory and Experiment, 201 (012((0):PP800) 2012.
[32] Daniel D. Lee and H. Sebastian Seung. Learning the parts of objects by non-negative matrix factorization. Nature, 7886555)):788, 1999.
[33] Mare Lelarge and Léo Miolane. Fundamental limits of symmetric low-rank matrix estimation. Probability Theory and Related Fields, 173(3):859-997 929, 2019.
[34] Thibault Lesieur, Florent Krzakala, and Lenka Zdeborová. Constrained low-rank matrix estimation: Phase transitions, approximate message passing and applications. Journal of Statistical Mechanics: Theory and Experiment, 201 2017(7)0073433, 2017.
[35] Lei Liu, Shunqi Huang, and Brian M. Kurkoski. Memory approximate message passing. arX1v.2012.10861, 2020.
12

---

[36] Junjie Ma and Li Ping. Orthogonal AMP. IEEE Access, 5:2020 20330,2033,111 2017.
[37] Junjie Ma, Ji Xu, and Arian Maleki. Optimization-based amp for phase retrieval: The impact of initialization and f2 regularization. IEEE Transactions on Information Theory, 65(6):3600-3629, 2019.
[38] Antoine Maillard, Bruno Loureiro, Florent Krzakala, and Lenka Zdeborová. Phase retrieval in high dimensions: Statistical and computational phase transitions. In Neural Information Processing Systems (NeurIPS), 2020.
[39] Arian Maleki, Laura Anitori, Zai Yang, and Richard G Baraniuk. Asymptotic analysis of complex lasso via complex approximate message passing (CAMP). IEEE Transactions on Information Theory,
2013.
[40] James A. Mingo and Roland Speicher. Free probability and random matrices, volume 35. Springer, 2017.
[41] Marco Mondelli, Christos Thrampoulidis, and Ramji Venkataramanan. Optimal combination of linear and spectral estimators for generalized linear models. arXiv:2008.003 033326, 2020.
[42] Marco Mondelli and Ramji Venkataramanan. Approximate message passing with spectral initialization for generalized linear models. In International Conference on Artificial Intelligence and Statistics (AISTATS), pages 397-405. PMLR, 2021.
[43] Andrea Montanari and Emile Richard. Non-negative principal component analysis: Mes- sage passing algorithms and sharp asymptotics. JEEE Transactions Information Theory, 62(3):1458-1484 2016.
[44] Andrea Montanari and Ramji Venkataramanan. Estimation of low-rank matrices via approximate message passing. Annals of Statistics, 45(1):321-3455 2021.
[45] Cristopher Moore. The computer science and physics of community detection: landscapes, phase transitions, and hardness. arXiv:1702.00461: 2017.
[46] Alexandru Nica and Roland Speicher. Lectures on the combinatorics of free probability, volume 13. Cambridge University Press, 2006.
[47] Jonathan Novak. Three lectures on free probability. Random matrix theory, interacting particle systems, and integrable systems, 65(309-383):13, 2014.
[48] Manfred Opper, Burak Cakmak, and Ole Winther. theory of solving tap equations for Ising models with general invariant random matrices. Journal of Physics A: Mathematical and Theoretical, 49(1 1):1 14002, 2016.
[49] Debashis Paul. Asymptotics of sample cigenstructure for a large dimensional spiked covariance model. Statistica Sinica, 17(4):1617, 2007.
[50] Amelia Perry, Alexander S Wein, Afonso S Bandeira, and Ankur Moitra. Message-passing algorithms for synchronization problems over compact groups. Communications on Pure and Applied Mathematics, 71 71(11):275-2322 2018.
[51] S. Rangan. Generalized Approximate Message Passing for Estimation with Random Linear Mixing. In IEEE International Symposium on Information Theory (ISIT), 2011.
[52] Sundeep Rangan, Philip Schniter, and Alyson K. Fletcher. Vector approximate message passing. IEEE Transactions on Information Theory, 65(10))6664-6684 2019.
[53] Philip Schniter and Sundeep Rangan. Compressive phase retrieval via generalized approximate message passing. IEEE Transactions on Signal Processing, 6344):104331055 2014.
[54] Philip Schniter, Sundeep Rangan, and Alyson K. Fletcher. Vector approximate message passing for the generalized linear model. In 50th Asilomar Conference on Signals, Systems and Computers, pages 1525-1529. IEEE, 2016.
13

---

[55] Pragya Sur and Emmanuel J. Candès. A modern maximum-likelin theory for high- dimensional logistic regression. Proceedings of the National Academy of Sciences, 116129):14516-1452) 2019.
[56] Keigo Takeuchi. Rigorous dynamics of unitarily invariant measurements. IEEE Transactions on Information Theory, 66(1):368-36)) -86, 2020.
-baed signal recovery from
[57] Keigo Takeuchi. Bayes-optimal convolutional AMP. IEEE Transactions Information Theory, 67(7):440544488 2021.
[58] Cédric Villani. Optimal transport: Old and new, volume 338. Springer Science & Business Media, 2008.
[59] Xinyi Zhong, Chang Su, and Zhou Fan. Approximate Message Passing for orthogonally invariant ensembles: Multivariate non-linearities and spectral initialization. LtXivv2110.02318 8 2021.
[60] Xinyi Zhong, Chang Su, and Zhou Fan. 11XV 2021.
Empirical Bayes PCA in high dimensions.
[61] Hui Zou, Trevor Hastie, and Robert Tibshirani. Sparse principal component analysis. Journal of computational and graphical statistics, 15(2):265-28,, 2006.
14


Reviews:
[{"review_id": "zD33zZ8dsG", "reviewer": "Reviewer_aknp", "summary": "This paper considers the problem of estimating a signal from its noisy rank-1 measurement. PCA is a traditional method that can be used for this task. AMP is an alternative that can additionally leverage statistical priors on the unknown signal. AMP algorithms are traditionally analyzed using a scalar iteration known as state evolution (SE).\n\nThe focus of this paper is to theoretically establish SE for AMP under a more general noise matrix starting from an initialization using PCA. Having SE enables the analysis of the error behaviour of AMP in different settings. Empirical results focus on comparing AMP+PCA against PCA on simulated data.\n\nTheoretical contributions: The paper develops a new AMP algorithm initialized using PCA. Theorem 1 and 2 are the main theoretical results in the paper that establish SE for the proposed algorithm for square-symmetric and rectangular matrices.\n\nEmpirical contribution: Validation of the superiority of AMP on simulated data over PCA. ", "questions": "", "limitations": "", "rating": 6, "confidence": 2, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "General comment: While I am familiar with AMP and SE, I am not an AMP researcher and thus have only loose familiarity with corresponding literature. My review thus should be weighted against others who have deeper expertise in AMP/SE. On the other hand, my review provides valid feedback from someone with expertise in optimization and ML.\n\nStrengths: Generalizing the AMP and its analysis to rotationally invariant matrices can lead to new insights into AMP and its usage for statistical estimation from rank-1 measurements. The theorems show that one can establish SE for new AMP algorithms, which is also promising. Simulations highlight the benefit of AMP over PCA when leveraging structure on the signal. \n\nIssues: (1) A persistent issue in the AMP literature is the lack of applications on real data analysis. This paper follows this trend and limits itself to validating the theory on synthetic experiments. In fact, this is a major issue that is often overlooked in theoretical papers, but which limits the significance of the results to data analysis. (2) Empirical validation limits comparisons to PCA. Why not compare to prior AMP methods for estimating from noisy rank-1 matrices (even when there is a mismatch in the noise distribution)? (3) How tolerant is the theoretical analysis to deviations from the assumptions? Which assumptions are more important than others? The assumptions in the manuscript are often presented in a scattered form within other assumptions or paragraphs of text, which reduces clarity.\n\nSummary: For someone not deeply familiar with the theoretical work on AMP, this paper seems to be proposing a novel algorithmic extension and proving that SE is applicable to that extension. As with most prior work on AMP, the method is never used for actual data analysis, but validated on a contrived set of simulations. In short, I think this paper could be of value to AMP researchers, but its potential to broader ML and data science has not been demonstrated. ", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "No issues with this.", "review_text": "", "extra_scores": {"ethical_concerns": "No issues with this.", "needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "2 hours"}, "extra_text": {}}, {"review_id": "TrJv4RmEwLp", "reviewer": "Reviewer_iyWx", "summary": "The submitted paper addresses estimation of rank-1 matrices from observations (1.1) or (1.2) subject to rotationally invariant noise matrix W. The paper proposes approximate message passing (AMP) with long memory Onsager correction (3.2) in the square case or (3.10) in the rectangular case, inspired by an arXiv preprint [20]. AMP uses the conventional principal component analysis (PCA) estimate as the initial message (3.1) or (3.9). To eliminate dependencies between the initial message and the noise matrix W, non-trivial Onsager correction b_{t,1} in (3.3) or a_{t,1}, b_{t+1,1} given in (3.11) and (3.13) is proposed via state evolution while the other correction b_{t,t-j} in (3.3) or (3.12) and (3.14) is the same as in [20]. Theorems 1 and 2 claim state evolution recursion for the square and rectangular cases, respectively. The theoretical predictions are also verified via numerical simulations in artificial setting. ", "questions": "", "limitations": "", "rating": 8, "confidence": 5, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "Originality\n\n The main originality is in the proof strategy to handle dependencies between the initial PCA estimate and the noise matrix W. The strategy is different from in [43], where followed Bolthausen’s conditioning technique to analyze the conditional distribution of the data matrix X for given initial messages. \n\nThe proposed strategy is a two-phase approach: In the first phase, AMP uses a linear denoiser in (5.1) to obtain the conventional PCA estimate after a sufficiently large number of iterations. In the second phase, AMP switches the denoiser to a time-shifted nonlinear denoiser in (5.3). The proof is fully based on state evolution recursion derived in [20] under review (while I believe the correctness of the proof in [20]). \n\nThe first part in the proof is the convergence of AMP in the first phase toward the PCA estimate. The proof utilizes a standard strategy in convergence analysis, such as the existence of a fixed point (Lemma B.4) and a construction of a contractive mapping (Lemma B.5). I guess that the linear denoiser enables such a proof. \n\nThe second part is a proof of the consistency (Lemma B.8) between the actual AMP and the two-phase AMP in the second phase. This is an interesting step to justify the definition of the Onsager correction in the actual AMP. \n\n--Why do Theorems 1 and 2 need the continuous differentiability of the denoiser? The assumption excludes practical functions, such as soft thresholding and ReLU. Conventional state evolution only postulates the Lipschitz continuity, which implies continuous differentiability “almost everywhere.” \n\n--Why do not the authors need the non-negativity of the exponent T+1-(t-j) for 1/\\alpha in (5.4)? In my understanding, this factor originates from the derivatives of the linear denoisers in the first phase for i\\in\\{t-j+1, …, T+1} in the case t-j+1\\leq T+1. Thus, T+1-(t-j) must be non-negative. In the current definition, let j=1 to obtain T+1-(t-j)=T+2-t<0 for t>T+2. \n\n--AMP in the first phase uses the true signal plus noise as the initial message. Why cannot AMP start with only noise, i.e. random initialization?\n\nQuality\n\n As pointed out above, the proof strategy is novel, general, and may be applicable to the other message passing. The mathematical quality in the proof is high while I did not fully read the whole proof in the supplementary material. For a further improvement in terms of references, I would like to introduce a reference [R1]. The reference is a first paper (posted on arXiv 27 March 2020 and published on IEEE TIT 04 May 2021) that presented rigorous state evolution for message passing algorithms with long memory Onsager correction, as cited in both [20] and [34].\n\n[R1] K. Takeuchi, \"Bayes-Optimal Convolutional AMP,\" in IEEE Transactions on Information Theory, vol. 67, no. 7, pp. 4405-4428, July 2021, doi: 10.1109/TIT.2021.3077471.\n\nOptional comment \n\n--It would be better if the authors presented numerical simulations for not artificial but real problems.\n\nClarity\n\nThe paper is well written. Theorems 1 and 2 are presented concisely after defining AMP and state evolution recursion. The main idea in the proof strategy is clear. Nonetheless, a further improvement is possible. \n\n--Clarify the significance of b_{t,1}, (3.11), and (3.13) for readers to understand the main difference between AMP in [20] and the proposed AMP. I reached my understanding mentioned in Originality after I considered several possibilities, such as differences between measurement models or the definitions of free cumulants in [20] and the submitted paper. \n\nSignificance \n\n The submitted paper is a first attempt (i.e. rank-1 case) for estimation of low-rank matrices from rotationally orthogonal measurements. In my understanding, the paper has solved the main challenge in that problem. A generalization to the low rank case is possible. \n", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": " The limitation is that the paper only addressed the rank-1 case. From a practical point of view, the low rank case may be desired. Since real data was not treated, the paper provides no negative social impact.  ", "review_text": "", "extra_scores": {"needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "10"}, "extra_text": {}}, {"review_id": "PeljzHZDNQC", "reviewer": "Reviewer_bZzE", "summary": "The paper studies the asymptotics of the AMP algorithms with PCA initialization.", "questions": "", "limitations": "", "rating": 6, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "I am very sorry that I do not have enough time to spend on this paper especially to go through the long appendix. But I do have concerns about the significance of the contribution.\n1) It is unclear to me why the paper only studies the rank-1 PCA rather than rank-k PCA. It is great to see that the symmetric condition can be relaxed and the results hold for the asymmetric matrix but the rank-1 condition is a very strong condition.\n2) The $u_t$ function in (3.2) is assumed to be continuously differentiable and Lipschitz. This is a standard assumption in AMP literature. However, for some applications of this work like phase retrieval, the $u_t$ function may not be differentiable everywhere depending on the loss function. It would be great to extend the assumption on $u_t$ to tolerate some non-differentiable points.\n3) It would be great to have a section to discuss which application requires a PCA initialization with an AMP algorithm. The only application I know of is phase retrieval. In general, only the non-convex problems require special initialization and hence PCA initialization might be required. Also, we need the AMP algorithm to perform well on the non-convex problem with PCA initialization. These conditions quite limit the potential applications.\n\nIn summary, it seems to be a borderline paper or slightly below average to me.", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "This is a theoretical paper. I think the main limitation is what are the possible applications to the theorems.", "review_text": "", "extra_scores": {"needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "4"}, "extra_text": {}}, {"review_id": "JxdrnrdzBd1", "reviewer": "Reviewer_NNbz", "summary": "The authors consider the problem of estimating a structured rank-1 signal in the presence of rotationally invariant noise, which is a class of perturbations that is more general than Gaussian noise. For this problem, they propose an AMP algorithm that is initialized by PCA, and they derivate a state-evolution characterization that holds in the infinite-dimensional limit. The rationale is that PCA is a near-optimal estimator in the absence of structure, whereas AMP facilitates the exploitation of structure yet allows a rigorous asymptotic analysis. Their work builds upon previous work, such as PCA-initialized AMP with Gaussian noise, but is different in that, to handle rotationally invariant noise, their algorithm uses a PCA correction at every iteration, not just the first one. Their state-evolution analysis uses a two-phase artificial AMP formulation that has previously been used in the context of generalized linear models.  Although their analysis assumes a noise distribution with non-negative free cumulants, they provide empirical simulation results that suggest that derived state evolution is correct with any noise distribution with compact support.", "questions": "", "limitations": "", "rating": 6, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "The contribution is a significant extension of existing results in that it extends PCA-initialized rank-1 AMP from the case of Gaussian noise to more generic noise distributions.  The results are intuitive yet rigorous and mathematically sophisticated. Furthermore, the analysis agrees with simulations conducted outside the stated assumptions (on non-negative free cumulants). The writing is very clear and well organized. \n\nAs pointed out by other reviewers, the lack of experimental results on real data constitutes a major weakness of this paper. Thus, I believe that it would be difficult for this paper to attract a wide audience.  \n\nOverall, this seems like a good contribution to NeurIPS, but not a great one.", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "The reviewer does not see any possibility of negative societal impact.", "review_text": "", "extra_scores": {"ethical_concerns": "The reviewer does not perceive any ethical concerns.", "needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "3 hours"}, "extra_text": {}}, {"review_id": "-diK-HaE_vh", "reviewer": "Reviewer_MuPE", "summary": "This paper considers the rank-1 signal recovery for rotationally invariant noise matrices. It uses PCA to generate an initialization of AMP that is correlated with the signal but independent of the noise. Different from the IID Gaussian noise matrix case, rotationally invariant noise matrix makes the construction and the state evolution of AMP be more difficult. For example, the Onsager term in AMP not only involves the memory in last iteration, but also contains all the memories in preceding iterations. Therefore, the PCA initialization will change the expression of AMP (and its SE) in all the iterations. In a word, this work focuses on a hard rank-1 signal recovery problem with general rotationally invariant noise matrices. The authors propose a modified AMP to solve this problem. Furthermore, the state evolution of the proposed AMP is also rigorously proved. ", "questions": "", "limitations": "", "rating": 10, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "The following are the concerns of the reviewer.\n\n1.\tThe symmetric square matrix case can be treated as an instance of the rectangular matrix (let m=n). Right? If so, it is unnecessary to discuss the symmetric case. The authors can directly focus on the rectangular case and then briefly introduce symmetric case as a special instance. \n\n2.\tThe authors claim that the existing AMP requires an initialization that is correlated with the signal but independent of the noise, which is not practical. Thus, they consider the PCA initialization to solve this issue. However, the reader may be curious that how the existing works solve this issue without PCA. Do you mean that the existing AMP does not work? If yes, how they do the simulation?\n\n3.\tThe reviewer guesses that the AMP structure in (3.1) and (3.2) is not firstly proposed by this paper. The related references should be cited when they introduce them.\n\n4.\tThe main concerns are the two assumptions of Theorem 1. When the assumption “alpha > alpha_S” holds? Is it hard to satisfy this assumption? Are there any examples to show that this assumption holds? \n\n5.\tSome numerical results are provided to support the second (non-negative) assumption of Theorem 2. This assumption rigorously holds for IID Gaussian matrix A. But it is still unclear if this assumption holds for the general rotationally invariant matrix, which is the scenario that this paper focuses. This may limit the contribution of this work.\n\n6.\tIn the simulation results, the proposed AMP always outperforms the existing PCA for IID Gaussian matrix A (see Fig. 1(a) and Fig. 2(a)). However, for other general matrices focused on this paper, the proposed AMP may be worse than the existing PCA in some cases (see Fig. 1(b), Fig. 1(c), Fig. 2(b) and Fig. 2(c)). Sometimes, the proposed AMP is always worse than the existing PCA (see Fig. 1(d) and Fig. 2(d)). This also limits the contribution of this work.\n\n7.\tThe reviewer is looking forward to seeing that the authors can construct a new VAMP/OAMP or memory AMP that always outperforms the existing PCA. \n\n8.  Symbol u_t is confusing. It is used for estimate function and also for function output.", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "yes", "review_text": "", "extra_scores": {"ethical_concerns": "no", "needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "20"}, "extra_text": {}}]

Author rebuttal:

Meta-review: This paper received 5 reviews. The scores/confidences were 6/3, 6/2, 8/5, 6/3, and 10/4, implying that all the reviewers evaluated this paper positively. We can notice that there is a relatively large spread among these review scores, but I think that the review contents are more or less coherent across all the reviews: The major strength of this paper is certainly its high technical quality, as commented by most reviewers. On the other hand, the major weakness is that it does not succeed in demonstrating significance of the theoretical contribution in real applications, so that it would be difficult for this paper to attract a wide audience beyond those who have specific interest in AMP and related iterative inference methods. The five reviews weighed these points differently on the basis of their own expertise, resulting in the spread of the scores. I am thus happy to recommend acceptance of this paper for presentation at the NeurIPS conference. At the same time, I would like to encourage the authors to consider demonstrating the significance of their contribution in real-world settings.

Very minor points:
- Lines 107-109: One would have to assume that $O$ is independent of $\Lambda$.
- Lines 127-129: Similarly, one would have to assume that $O$ and $Q$ are independent and independent of $\Lambda$.
Decision: Accept (Poster)