Support Recovery in Sparse PCA with Incomplete Data
Authors: Hanbyul Lee
Qifan Song
Jean Honorio
Venue: NeurIPS 2022 Accept (2022)

Abstract: We study a practical algorithm for sparse principal component analysis (PCA) of incomplete and noisy data.
Our algorithm is based on the semidefinite program (SDP) relaxation of the non-convex $l_1$-regularized PCA problem.
We provide theoretical and experimental evidence that SDP enables us to exactly recover the true support of the sparse leading eigenvector of the unknown true matrix, despite only observing an incomplete (missing uniformly at random) and noisy version of it.
We derive sufficient conditions for exact recovery, which involve matrix incoherence, the spectral gap between the largest and second-largest eigenvalues, the observation probability and the noise variance.
We validate our theoretical results with incomplete synthetic data, and show encouraging and meaningful results on a gene expression dataset.

Support Recovery in Sparse PCA with Incomplete Data
Hanbyul Lee Department of Statistics Purdue University
Qifan Song Department of Statistics Purdue University West Lafayette, IN 47906 qf song@purdue.
Jean Honorio
Department of Computer Science
Purdue University
West Lafayette, IN 47906 dee3078@purdue. edu
West Lafayette, IN 47906
jhonorioppurdue. edu
Abstract
We study practical algorithm for sparse principal component analysis (PCA) of incomplete and noisy data. Our algorithm is based on the semidefinite program (SDP) relaxation of the non-convex l-rreuularizdd PCA problem. We provide theoretical and experimental evidence that SDP enables us to exactly recover the true support of the sparse leading eigenvector the unknown only observing an incomplete (missing uniformly at random) and noisy version of it. We derive sufficient conditions for exact recovery, which involve matrix incoherence, the spectral gap between the largest and second- eigenvalues, the observation probability and the noise variance. We validate our theoretical results with incomplete synthetic data, and show encouraging and meaningful results on a gene expression dataset.
matrix, despite
Introduction
Principal component analysis (PCA) is one of the most popular methods to reduce data dimension which is widely used in various applications including genetics, image processing, engineering, and many others. However, standard PCA is usually not preferred when principal components depend on only small number of variables, because it provides dense vectors as solution which degrades interpretability of the result. This can worse especially the high-dimensional setting where solution of standard PCA is inconsistent addressed in several works [Puul 2007, Nadler 2008 Johnstone and Lu 2009]. To solve the inconsistency issue and improve interpretability. sparse PCA has been proposed, which enforces sparsity in the PCA solution that dimension reduction and variable selection can be simultaneously performed. Theoretical and algorithmic researches on sparse PCA have been actively conducted over the past few years [Zou et al. [0066] Amini and Wainwright
2008 Journée et al. 2010 Ma 2013
Lei and Vu 2013 Berk and Bertsimas] 2019 Richtárik et al.
2021
In this paper, we consider a special situation where the data to which sparse PCA is applied are not completely observed, but partially missing. Missing data frequently occurs in wide range of machine learning problems, where sparse PCA is no exception. There are various reasons and situations where data becomes incomplete, such as failures of hardware, high expenses of sampling, and preserving privacy. One concrete example is the analysis of single-cell RNA sequence (scRNA-seq) data [Park and Zhao 2019 where the cells are divided into several distinct types which can be characterized with only a small number of genes among tens of thousands of genes. Sparse PCA can be effectively utilized here reduce the dimension (from numerous cells to few cell types) and to select small number of genes that affect the reduced data. However, since scRNA-seq data usually have many missing values due to technical and sampling issues, the existing sparse PCA theory and method designed for fully observed data cannot be directly applied, and new methodology and theory are in demand.
36th Conference on Neural Information Processing Systems (NcurIPS 2022).

---

Despite the need for theoretical research and algorithmic development of sparse PCA for incomplete data, there have not been many studies yet. Lounici [013] and Kundu al. [2015] considered two different optimization objectives for sparse PCA on incomplete data, which impose / regularization and lo constraint on the classic PCA loss function using (bias-corrected) incomplete matrix, respectively. It was shown that the solution of each problem has a non-trivial error bound under certain conditions, but the optimization problems they considered are either nonconvex or NP-hard, and thus theoretical studies of computational feasible algorithms are still lacking, More recently, Park and Zhao| [2019] proposed a computationally tractable two-step algorithm based on matrix factorization and completion, but its first step is iterative algorithm that requires singular value decomposition in every iteration, which incurs a lot of cost in memory and time under a high- dimensional setting.
With this motivation, we suggest a computational friendly convex optimization problem via a semidefinite relaxation of the I regularized PCA, to solve the sparse PCA on incomplete data. We note that very efficient scalable SDP solvers exist in practice [Yurtsever et al.) 202T]. We assume that the unknown true matrix M E Rddd is symmetric and has a sparse leading eigenvector u1 Our goal is to exactly recover the support of this sparse leading eigenvector, i.c., to find the set correctly where J supp(u1) {i w1,; + 0}. Given noisy observation M for the unknown true matrix M it is intuitive to consider imposing a regularization term on the PCA quadratic loss that aims to find the first principal component. When using the l regularizer, the optimization problem can be written as:
i arg max x Mæ p||=||?
Hence, is estimated with supp()). However, this intuitively appealing objective is nonconvex and very difficult to solve, so the following semidefinite relaxation can be considered as an alternative:
X
arg max
(M P||X
X and {=
By letting X = xx the equivalence of the above two objective functions can be easily justified. Since supp(x) suppl diag we estimate the support J by J pupp(diay) X)) in the semidefinite problem. This kind of relaxation has been studied by d'Apermmoont et al. (2004) and Lei and Vul (2015). but their works were limited to complete data. Surprisingly. without any additional modifications on the relaxation problem such as using matrix factorization or matrix completion, we show that it is possible to exactly recover true support J with the above semidefinite program itself when M is an incomplete observation. Our main contribution is to prove this claim theoretically and experimentally.
In Section 3 we provide theoretical justification (i.e., Theorem I that we can exactly recover the true support J with high probability by obtaining a unique solution of the semidefinite problem, under proper conditions. The conditions involve matrix coherence parameters, the spectral gap between the largest and second- and the noise variance, which are discussed in detail in Corollaries | and 2 Specifically, we show that the sample complexity is related to the matrix coherence parameters well as the matrix dimension d and the support size s. We prove that the observation probability p has the bound of in the worst scenario in terms of the matrix coherence, while has smaller lower
cigenvalues of the true matrix, the observation probability
p=
bound p wi (104 in the best scenario. Section we provide experimental results on incomplete synthetic datasets and gene expression dataset. The experiment on the synthetic datasets validate our theoretical results, and the experiment on the gene expression dataset gives us a consistent result with prior studies. We also show that our SDP algorithm outperforms several other sparse PCA approaches in the synthetic dataset.
2 Preliminaries
2.1 Notation
We first introduce the notations used throughout the paper. Matrices are bold capital, vectors are bold lowercase and scalars entries are not bold. For any positive integer n, denote [n] := (1...., n). For any vector a e dg and index set J C [d]. a, denotes the |J] -ddimennio vector consisting of

---

the entries of a in J. For any matrix A Rdixda and index sets J1 S [di] and J2 C [d2]. A and (A..AA..Ja denote the [Ji] [J2] sub-matrix of A consisting of rows J/1 and columns in J/. and the Jil x d2 (d1 (J2)) sub-matrix of A consisting of rows in J1 (columns in J2), respectively. ||a ||a]|2 and ||a]lcc represent the norm, 2 norm and maximum norm of vector a, respectively. (e; :e e indicates the standard basis of IRd.
variety of norms on matrices will be used: we denote by ||A||2 the spectral norm and by |||||
the Frobenius norm of matrix A. We let ||A||1,1 LiE(d,).EEde |A,p || max |||]]] oc,
ja As.il || max;e(ds] ||A., and || ma ||A,, || represent
the 1,11 norm, the entrywise loc norm, the 12,0 norm and the 1,,0 norm of a matrix A, respectively.
The trace of is denoted tr(),, and the matrix inner product of A and is denoted (A, B).
Also, o;; (A) and (A) represent the ith largest singular value and the ith largest eigenvalue of A, respectively.
The notation C, CI....,,C C1.... denote positive constants whose values may change from line to line. The notation f(z) =(z)=o(g(a))) (z)())) or f(x) < g(x) means lim f(z)/g(z) 0: f(x) = ((g(x))
or f(x) >> g(z) means lim f())/g(r) = oc; f(x) = ((g(z)) or f(z) S g(z) means that there exists constant C' such that f(x) S Cg(a) asymptotically, f(z) ((g(z)) or f(x) g()
means that there exists a constant C' such that f(x) Cg(z) asymptotically, f(x) = O(g(z))
or f(x)= = g(x) means that there exists constants and C" such that Cg(x) S f(x) S C²)(a)
asymptotically.
2.2 Model
We now introduce our model assumption. Suppose that an unknown matrix M' E Raxdd is symmetric. The spectral decomposition of M* is given by
M Σ kE[d]
(M (MM)uu
where X1 (M*) 2... (M are its eigenvalues and UI,... Uh E Rd are the corresponding eigenvectors. We assume that \i (M*) > >²(M*) and the leading eigenvector in of M is sparse, that is, for some set [d].
11,1 + 0 E
U1 =
otherwise.
With a notation supp(a) = (i E [4] az # 0} for any vector a E TR.d, we can write = supp(ui). Also, we denote the size of J by s.
Incomplete and noisy observation M over sampling set S [d] [d]. Specifically, we observe symmetric matrix M Rddd such that
Suppose that we have only noisy observations of the entries of
M.. M₁. 8ū.j (Mi,, + Ei.s)
for 1 S < S d, where 5,, = 1 if (i.j) E Ω and ōt,5 = otherwise, and cg,j is the noise at location (i,j). In this paper, we consider the following assumptions on random sampling and random noise: for 1 < S
* Each (i,j) is included in the sampling set Ω independently with probability p (that is ōig 1.d. Ber(p).)
ōig's is and ei. 's are mutually independent.
E[E,,]] 0 and Var[er.g] o².
| I B almost surely.
Finally, we define the coherence parameters of the sub-matrices M
J. je,y and M M je.ye
3

---

Definition 1 (Coherence parameters), We define the coherence parameters pi Po(M J,y) P1 (M p2(M)"-)) and ps (MMM"' j= as follows:
[MM
J,J I|max
po(M),)) = j,,) (2)))
Pi (M,,
[M j
[ je,
P22 (2 jo.)
Ps j.M : max M je 2ll. || J=l_2,l2,00
|M
We use po, P1, p2 and pi as shorthand for (,3) 11 (,,) je, and (s(M
respectively.
Intuitively, when each coherence parameter is small, all the entries of the corresponding matrix have comparable magnitudes. Note that? S Ho S 1, v8 S f1 S 1, Vald-a) U22 S 2 1s Us $3.11 1.
Remark. We note that our setting different from the case in which data matrix Y is observed with missing entries and the covariance matrix M = =Y Y is analyzed for sparse PCA. Instead, we directly observe a general symmetric incomplete matrix M, which is not necessarily positive semidefinite. Our model can be therefore applied to various data types other than covariance matrices (e.g,, undirected random graphs with missing edges.) When a noisy and incomplete data matrix Y is given and its covariance matrix is considered, our current anulysis cannot be directly applied and one may need some bias correction technique in the algorithm, since E = EE Y is not the true matrix Y where Y EY Furthermore, the independent missingness mechanism over Y does not imply independent missingness over M. We leave this future work.
Main Results
As mentioned in the introduction, we consider the following semidefinite programming (SDP) in order to recover the true support J:
X
arg max
(M, X)
P||XX|1 1|,11
(1)
X (=(
where estimate J by j supp(diag(X)) (X)). We recall that is convex relaxation of the following nonconvex problem:
arg max x T Mæ -|
In Theorem 1, we will show that under appropriate conditions, the solution of IT attains j with high probability. Our main technical tool used in the proof is the primal-dual witness argument Wainwright 2009]. We start with deriving the sufficient conditions for the primal-dual solutions of ①1 to be uniquely determined and satisfy supp(diag(
= J. We then establish a proper candidate solution which meets the derived sufficient conditions, where we make use of the Karush-Kull Kuhn- Tucker (KKT) conditions of (2) to set up a reasonable candidate. We finally develop the conditions under which the established candidate solution satisfies the sufficient conditions from the primal-dual witness argument of (T) with high probability. Detailed proof is given in Appendix B
Theorem 1. Under the model defined in Section
2.2 assume that the following conditions hold:
K1 + ps
2V2. P()((M - (M j,,)
< min
VEJ
> 2 psc ((-p))|M; --,111 + (d- --)ss00 + IM "e,,.,lllll je.j/lmmxx
(K2 +p 104++-1|MM,,,p |M, (1 v3)² < {p- (Ai(M).)) - A2 - K1 2pps
{pp (M1(M)))))))) - A2 - K1 - K3 - pd).
where e 0, and K1, K2 and K3 are defined follows:
Ki (c+ 1). Rí log (2s) + 1) R2 TaVlog(2s)
K₂ = (c+1) Rg log d + / 2(c 1). R4 Ravlog
Kg = (c+ 1). R, log - s)) + 2 1)- Ro log(l2/ s))

---

and
R1 = max ((1 - P)]l J.J mmx + B, JIIMj,yllma
R2 vp(1 (P(1-))|M3 J.J 2,,112,0o + psoo
R3 max ((1 P))lM J. B, p Jr Jr mllmax
R4 max ( p(l1 PMM se, 2 + vpld s)02, vp(1 P)
v /psof)
R2 m1 - (1)) j-
B, pllM "e..--llmaxi
101-ppi Rs = V p jo. + v Vvld - 8302
Then the optimal solution X the problem (7y is unique and satisfies supp(diog(X) J with probability at least 1 - g - d
(2(d - -))=
Remark 1. Our proposed method actually guarantees more than support recovery. Under the
9) where & satisfies
sufficient conditions, the optimal solution X can be expressed as x
|un. - = I minjes (u1ي) which is shown in our proof. This implies consistency of our solution eigenvector estimation.
Remark 2. To deal with missingness, we utilized concentration inequalities such as malrix Bernstein inequality and Chebychev's 's inequality, which is one key technical difference from prior works on sparse PCA with complete data [Lei and Vui 2015y. Beyond this, under the primal-duce dual witness framework, while Lei and Vut (2015) employed implicit constraints in the optimization problem, we made all the constraints explicit and directly derived the sufficient conditions to avoid strong conditions. We note that simply applying concentration inequalities to the approach Lei and Vul (2015) results in stricter requirements on M je.e and p.
We now consider the following two particular scenarios:
(s1) == 02 0, that is, the observation M is noiseless (but still incomplete).
(s2) The rank of M is 1.
to better interpret the conditions of M and p listed in Theorem and understand under what circumstance these conditions hold. For both cases, we set p 0.5 for simplicity. Under the first setting, we can re-express the conditions on M for exact sparse recovery of J in more interpretable way (specifically, in terms of coherence parameters and spectral gap) as well the conditions on In the second setting, we aim to investigate that the maximum level of noise that is allowed by Theorem [ Corollaries [ and 2 include the results of the two settings (sI) and (s2), respectively.
Corollary 1. Assume that B d2 = 0, 2 0.5 and min:e |u1.il = Q (2)) Denote x (M),) (2 (.) by (,)). If the following conditions hold:
(3)
V8 logs
X(M
min P2: s' log }
1 V8
|| M je.j
Ilaxx
(4)
M je Ilmax A(Mi.) min jus, log(d s)
(5)
of min M1 log
P
X(M min (52 s
(6)
M Je,ي mlmxx
X(M X(MG, 3,
1
Vlogld - Viogdd-0)
M', J.. Immxx
lllxx
pk (M)))
0
s2
5

---

Shoorem@hoodd hold asymptotically, that is, when and d are sufficiently large, the optimal solution X to the problem 17 is unique and satisfies supp( (diag( (X) J with probability
then the conditions in Theorem
at least 1 - - d (25) - [2(d -))=- s))-1
Conditions on true matrix M From the conditions in Corollary 1 we can find desirable properties on the matrix M follows:
Incoherence of j.j. and coherence of je.. and je. From the coherence parameter in (3) and those in ( (5) and (6). we see that the sub-matrix M and the sub-matrices M jo and M je_r are expected to be incoherent and coherent, respectively. This is different from other problems involving incomplete matrices, such as matrix completion [Candes and Recht 2009) and standard PCA on incomplete data (Cai et al. 2021], where the entire matrix, not sub-matrix, is required to be incoherent.
We can easily check the need of incoherence of M'..., with an example that the sub-matrix has only one entry with large magnitude while the other entries have relatively small values. Even if the true leading eigenvector of the sub-matrix is not sparse, the sparse PCA algorithm may produce a solution J which has a smaller size than that of the true support J. However, for j and J= coherence is preferable: intuitively speaking, when je.. and M je. are the most coherent, that is, only one entry is nonzero in each sub-matrix, and all other entries are zero, missing the entries in M Jo.j and J does not change the leading eigenvector of M On the other hand, when M Je,J j and M are incoherent, that is, all the entries have comparable magnitudes, missing only a few entries changes the leading eigenvector and its sparsitency, so that sparse PCA is likely to fail to recover J. A simple illustration can be found in the Appendix
Large spectral gup X(M)),) (= X (M, A2 (M,, G.,) )): This can be found in 4. and (6) A sufficiently large spectral gap requirement has been also discussed in the work sparse PCA on the complete matrix (Lei and Vul 2015). It ensures the uniqueness and identifiability of the orthogonal projection matrix with respect to the principal subspace. If the spectral gap of cigenvalues is nearly zero, then the top two eigenvectors are indistinguishable given the observational noise, leading to failure to recover the sparsity of the leading eigenvector.
We also note that λ1 (M (2,))) 2 (M >²(M)) since A (M,)) ) A1 and 22 (M))) < (M". Hence, a large X1 (M*) (M²) implies a large A(M).))
Small magnitudes of M' and M jr This can also be found in ④ (5) and (6. This condition is also natural: the magnitudes are relatively small, missing the entries will not make big impact to the result.
Conditions on p (ratio of missing data) For simplicity, suppose that (M,) O(s) and (log d). We consider two extreme cases where the coherence parameters are maximized and minimized. We discuss the bound of the sample complexity in each case.
The best scenario where the bound of the sample complexity is the lowest: Suppose that ol ( and #3 = 1 (note that when po = ol )) is upper bounded by o 1 Jogs ).) Then the condition (6) can be written as:
p 1+G == o min log sd] M jo, 1 }).
Furthermore, if MM jo Imnax =Ol log sd
then p w( (1o
'This be expressed equivalently as where of min f w / Vlg
J Je.a line
XMG,y)e
viog(d-a) }).
min J w
M Je Iimax

---

The worst scenario where the bound of the sample complexity is the highest: Suppose that and
M1 7s 12
Vold=s)
p3 ~=s" this case, the condition (6) can be written as:
p 62 62 of min { logs
M Je, Imax })-
(sd M Je,a [mmax)
[d]
Suppose that ||M Je,,llmax = 01 (tad)) and || Je, = O +2) for instance, then
p
Next, we consider the second setting (s2) where the rank of M is assumed to be that is, M (M"))uu,, (without loss of generality, we assume λι (M*) > 0). Trivially, M' = M J,J"
J jr J- 0 and Theorem I can be greatly simplified. Here, we focus on analyzing how much noise (parameters B and o²) is allowed.
Corollary 2. Assume that 2 0.5 and the rank of M is 1, that is, M (M yu Let
12+ A=z 48/202- 1/0//22-11---00 -4--=1 - 20 45122
AM > 0. Suppose that and d satisfy J/ S
4/24 4 +122
where 01 = (2 (2- -)
log d (2)))
v /max /max(d-s.s)--log d and 02 = (2 - 8v2 (2g()))
log(2(d- a)) +
(2((=)))))
16s² Va-
16s
If the following conditions hold:
maxijeJ |11,0111
minies |u,
16
max,eJ ju1,s!
p
[1.1 minjeJ
16 2// (2s)
V
< (2p
max LJEJ
|u N1.h
1
2/2. V 5) P
pdd (M")
min ju,a
8v2s
then the optimal solution X to the problem probability least
is unique and satisfies supp(diag(X)) = J with (2 s))
d (2s)
Conditions noise parameters B and o² For simplicity, let (M*) -(ss) and V/uia O Vs! Then the above conditions in Corollary 2 imply that
B BSD p and 0² 33
The condition for is relatively moderate while o² needs to be extremely small satisfy the We comment this is only sufficient condition, and the experimental results
condition in Corollary 2 show that can succeed even with x² larger than the aforementioned bound.
Numerical Results
We perform the SDP algorithm of (1 on synthetic and real data validate theoretic results and show how well the true support of the sparse principal component is exactly recovered. We also compare our SDP algorithm with other sparse PCA methods and show that method performs better. Our experiments were executed on MATLAB and standard CVX code was used, although very efficient scalable SDP solvers exist practice YYurtsever et al. 2021].
4.1 Synthetic Data
We perform two lines of experiments:
With the spectral gap A (M*) d2 A2(M) and the noise parameters B and o² fixed, we compare the results for different and d.
7

---

recovery rales
covery
Exact recovery rate
0.8
Exact recovery rate
0.3
0.4
0.4
5
5
2-20
9-20
0.4
0.6
Rescaled
Exact recovery rates when s 10
Exast recovery rates when 10
Exacc recovery rale
Exaco recovery rale
0.0
0.6
0
d=20
d=20
0.2
d+50
0.2
d+50
d=100
0.4
Rescalee
Figure 1: Results of experiment on synthetic data.
2. With s and d fixed, we compare the results for different spectral gaps and noise parameters.
the leading cigenvector μ1 is set to have number of non-zero entries. A2(M (MM ...... ha(M") are randomly selected from normal distribution with mean 0 and standard deviation 1, and X (M*) ) is set to A2(M") plus the spectral gap. The orthogonal eigenvectors are randomly selected, while the non-zero entries of the leading eigenvector are made to have value of at least 2.5*
In each experiment, we generate the true matrix M as follows:
When generating the observation M. first add to M the entry-wise noise which is randomly selected from a truncated normal distribution with support [-B, B, B). The normal distribution to be truncated is set to have mean O and standard deviation normal- After adding the entry-wise noise, we generate an incomplete matrix M by selecting the observed entries uniformly at random with probability p E (0.1, 0.3, 0.5, 0.7, 0.9}.
In each setting, we run the algorithm ( and verify the solution exactly recovers the true support. We repeat each experiment 30 times with different random seeds, and calculate the rate of exact recovery in each setting.
Experiment In this experiment, fix the spectral gap A (M*) A2 as 20 and the noise parameters B and o² as 5 and 0.01. We use the tuning parameter p = 0.1. We try three different matrix dimensions d E {200, 50, 100) and three different support sizes E {5, 10, 20).
To check whether the bound of the sample complexity obtained in Corollary i is tight, we calculate the coherence parameters and the maximum magnitudes of the sub-matrices at each setting, and calculate the following rescaled parameter:
(MM).))
2s' s/s(d s) } (MM)))P
V I 1-pp min (M) logs, s, [ Jz Illmaxx
min
M ye Ilmax Vlog(d
which is derived from (6). If the exact recovery versus this rescaled parameter is the across different settings, then we empirically justify that the bound the sample complexity we derive is "tight") in the sense that the exact recovery rate is solely determined by this rescaled parameter.
Figure shows the experimental results. The two plots above are the experimental results for different values of when d = 100, and the two plots below are for different values of d when 10.
We note that "tight" does not mean agreement between necessary and sufficient conditions here.

---

recovery
recovery
"gons"
recovery "rorss"
0.8
0.0
0.0
recovery
0.8
0.8
0.5
2 2
3 22 2100
0.21
0.21
0.2
0.1
0.3
0.5
0.7
0.0
0.5
0.5
0.7
0.0
0.3
0.5
0.7
0.9
Figure 2: Results of experiment on synthetic data.
The x-axis of the left graphs represents and the x-axis of the right graphs indicates the rescaled parameter.
We can see from the two graphs on the right that the exact recovery rate versus the rescaled parameter is the same in different settings of d and This means that our bound of the sample complexity is tight.
Another observation we can make is that the exact recovery rate is not necessarily increasing or decreasing as s or increases or decreases. This probably because coherences and maximum magnitudes of sub-matrices are involved in the sample complexity as well.
Here, we fix the matrix dimension d as 100 and the support size as 50. We set B We try three different spectral gaps x (M*) >²(M*) {10, 30, 50) and three different standard deviations of the normal distribution, norral E {0.1, 0.3, 0.5]. We try two different tuning parameters / E {0.1, 0.01] and report the best result.
Experiment 2
Figure demonstrates the experimental results. The three plots show the results when normal is 0.1, 0.3 and 0.5, respectively. The red, green and blue lines indicate the cases where the spectral gap (M)) )- A2 (M") is 50, 30 and 10, respectively. From the plots, we can observe that the exact recovery rate inereases as o² is small and (M*) X²(M*) large, which is consistent with the conditions we have checked in Corollaries and 2
4.2 Gene Expression Data
We analyze gene expression dataset (GSE21385) from the Gene Expression Omnibus website (https://www. http://wwwwwcbi... nlm. nih, gtwlgooo/gggoo/o The dataset examines rheumatoid arthritis synovial fibroblasts, which together with synovial macrophages, are the two leading cell types that invade and degrade cartilage and bone.
The original data set contains 56 subjects and 112 genes. We compute its incomplete covariance matrix, where 87% of the matrix entries are observed since some subject/gene pairs are unobserved. With this incomplete covariance matrix, we solve the semidefinite program in (T) for sparse PCA with p= 2.
By solving T we find that the support of the solution contains 3 genes: beta-1 catenin (CTNNB), hypoxanthine-guann guanine phosphoribosylrrasseena
(HPRTI) and semaphorin IIII (SEMA3F). Our result is consistent with prior studies on rheumatoid arthritis since CTNNB has been found to be upregulated [Iwamoto al. 2018), SEMA3F has been found be downregulated [Tang et al. 2018 and HPRTI is known to be housekeeping gene [Mesko el al. 2013 Additional illustration of this application can be found in Appendix
4.3 Comparison with Other Methods
compare our SDP algorithm with three different methods. First, we consider two sparse PCA algorithms where missing cells are treated as zero: the diagonal thresholding sparse PCA (DTSPCA)
by Johnstone and (009) and the iterative thresholding sparse PCA (ITSPCA) by Ma (0033)
Second, we consider the combination of imputation and our SDP method: we first estimate the

---

Exact recovery rale
0.8
0.6
DTSPCA
impraoo
0.1
0.3
0.5
0.7
0.9
Figure 3: Comparison with other methods.
missing elements of M by using the matrix completion method based on the following Frobenius norm miminization with nuclear norm penalty ||Y|.:
M arg min
(Y) M 1 +||
Y Ysyymmetric
and then implemented the SDP method with the completed matrix M. Details of the simulation
settings and algorithm implementations are provided in Appendix F
F demonstrates demostasse that the exact recovery rate of SDP method is higher than those of the three other methods. The comparison shows that the imputation- baaseddoo---eeddddd strategy fails to work under the small-p regime, and we conjecture the following rationale: it known that the matrix completion algorithm can be successful under the low-rank assumption. However, we do not assume that the true matrix low-rank, but just that its spectral gap ê (M*) >2(M (M*) is large enough, which probably leads to the failure of matrix completion. Over the 20 simulation replications, we found that the average of max,, M, My was 1.0260 while the averages of max;,j |M and min,,; [MG] 2.3200 and 0.0002, respectively. This shows that the matrix completion has not been done very precisely. Given such unsuccessful matrix completion, the imputed cells introduce more noise into the inference and the result of sparse PCA can be even worse than that of simply using zero for the missing entries.
5 Concluding Remarks
We have presented the sufficient conditions exactly recover the true support of the sparse leading eigenvector by solving simple semidefinite programming We have shown that the conditions involve matrix coherence, spectral gap, matrix magnitudes, sample complexity and variance of noise, and provided empirical evidence to justify our theoretical results. To the best of our knowledge, we provide the first theoretical guarantee for exact support recovery with sparse PCA on incomplete data. While currently focus uniformly missing at random setup, an interesting open question is whether possible to provide guarantees for a deterministic pattern of missing entries.
incomplete and noisy observation.
Acknowledgments
This material is based upon work supported by the National Science Foundation under Grant No. 2134209-DMS.
References
Arash A Amini and Martin Wainwright. High High-dimensional analysis of semidefinite relaxations for sparse principal components. In 2008 IEEE international symposium information theory, pages 2454-2458 IEEE, 2008.
Lauren Berk and Dimitris Bertsimas. Certifiably optimal sparse principal component analysis. Mathematical Programming Computation, 11(3):81-4220, 2019.
10

---

Changxiao Cai, Gen Li, Yuejie Chi, H Vincent Poor, and Yuxin Chen. Subspace estimation from unbalanced and incomplete data matrices; f2, so statistical guarantees. The Annals of Statisties, 49 (2):944-667, 2021.
Emmanuel J Candès and Benjamin Recht. Exact matrix completion via convex optimization. Foun- dations of Computational mathematics, 96):717-772, 2009.
Alexandre d' A'Appemoot,, Laurent Ghaoui, Michael Jordan, and Gert Lanckriet. A direct formulation for sparse pea using semidefiiite programming. Advances in neural information processing systems, 17, 2004.
Naoki Iwamoto, Shoichi Fukui, Ayuko Takatani, Toshimasa Shimizu, Masataka Umeda, Ayako Nishino, Takashi Igawa, Tomohiro Koga, Shin-ya Kawashiri, Kunihiro Ichinose, et al. Osteogenic differentiation of fibroblast-like synovial cells in rheumatoid arthritis is induced by microrna-218 through robo/slit pathway. Arthritis research & therapy. 20(1):1-10, 2018.
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):682 (29669,,)8666,,,,,, 2009.
Michel Journée, Yurii Nesterov, Peter Richtárik, and Rodolphe Sepulchre. Generalized power method for sparse principal component analysis. Journal of Machine Learning Research, 11).. 2010.
Abhisek Kundu, Petros Drineas, and Malik Magdon-Ismail. Approximating sparse pea from incom- plete data. Advances in Neural Information Processing Systems, 28, 2015.
Jing Lei and Vincent Q Vu. Sparsistency and agnostic inference in sparse pea. The Annals of Statistics, 43(1):299 322, 2015.
Karim Lounici. Sparse principal component analysis with missing observations. In High dimensional probability VI, pages 327-356. Springer, 2013.
Zongming Ma. Sparse principal component analysis and iterative thresholding. The Annals of Statistics, 81 2013.
Bertalan Mesko, Szilard Poliska, Andrea Váncsa, Zoltan Szekanecz, Karoly Palatka, Zsolt Hollo, Attila Horvath, Laszlo Steiner, Gabor Zahuczky, Janos Podani, et al. Peripheral blood derived gene panels predict response to infliximab in rheumatoid arthritis and crohn's disease. Genome medicine, 5(6):1-10, 2013.
Boaz Nadler. Finite sample approximation results for principal component analysis: matrix perturbation approach. The Annals of Statistics, 36(6):2791-281 2008.
Seyoung Park and Hongyu Zhao. Sparse principal component analysis with missing observations. The Annals of Applied Statistics, 1302):001-11022, 2019.
Debashis Paul. Asymptotics of sample eigenstructure for large dimensional spiked covariance model. Statistica Sinica, pages 1617-1642, 2007.
Peter Richtárik, Majid Jahani, Selin Damla Ahipaşaoglu, and Martin Takáč. Alternating maximization: unifying framework for 8 sparse pea formulations and efficient parallel codes. Optimization and Figgineering 2 1493-1519, 2021.
Man Wai Tang, Beatriz Malvar Fernández, Simon P Newsom, Jaap D Buul, Timothy RDJ Radstake, Dominique Baeten, Paul P Tak, Kris Reedquist, and Samuel García. Class semaphorins modulate the invasive capacity of rheumatoid arthritis fibroblast-like synoviocytes. Rheumatology, 57(5):909-920, 2018.
Joel A Tropp. User-friendly tail bounds for sums of random matrices. Foundations of computational mathematics, 12(4):389 2012.
Martin J Wainwright. Sharp thresholds for high-dimensional and noisy sparsity recovery using E connttrained quadratic programming (lasso). IEFE transactions on information theory, 55(5):
2183-2202, 2009.

---

Alp Yurtsever, Joel A Tropp, Olivier Fercoq, Madeleine Udell, and Volkan Cevher. Scalable semidefi- nite programming. SIAM Journal on Mathematics of Data Science, 3(1):111-200, 2021.
Hui Zou, Trevor Hastie, and Robert Tibshirani. Sparse principal component analysis. Journal of computational and graphical statistics, 15(2):265-2866 2006.
Checklist
1. For all authors...
(a) Do the main claims made in the abstract and introduction accurately reflect the paper's
contributions and scope? [Yes]
(b) Did you describe the limitations of your work? [Yes] See Section 5
(c) Did you discuss any potential negative societal impacts of your work? [No]
(d) Have you read the ethics review guidelines and ensured that your paper conforms to them? [Yes]
2. If you are including theoretical results...
(a) Did you state the full set of assumptions of all theoretical results? [Yes] See Section
(b) Did you include complete proofs of all theoretical results? [Yes] See the supplemental material.
3. If you ran experiments...
(a) Did you include the code, data, and instructions needed to reproduce the main experi- mental results (either in the supplemental material or URL)? [Yes] We provide the code the supplemental material.
(b) Did you specify all the training details (e.g., data splits, hyperparameters, how they were chosen)? [Ys]] See Section
(c) Did you report error bars (e.g., with respect the random seed after running experi- ments multiple times)? [Yes] See Section F
(d) Did you include the total amount of compute and the type of resources used (e.g., type of GPUs, internal cluster, or cloud provider)? [No]
4. If you are using existing assets (e.g., code, data, models) or curating/releasing new assets...
(a) If your work uses existing assets, did you cite the creators? [Yes] See Section 4
(b) Did you mention the license of the assets? [Yes] See Section 4
(c) Did you include any new assets either the supplemental material or as URL? [No]
(d) Did you discuss whether and how consent was obtained from people whose data you're
using/curating? [No]
(e) Did you discuss whether the data you are using/curating contains personally identifiable information or offensive content? [No]
5. If you used crowdsourcing or conducted research with human subjects...
(a) Did you include the full text of instructions given to participants and screenshots, applicable? [N/A]
(b) Did you describe any potential participant risks, with links Institutiona Review Board (IRB) approvals, applicable? [N/A]
(c) Did you include the estimated hourly wage paid to participants and the total amount spent on participant compensation? [N/A]
12


Reviews:
[{"review_id": "OM65Qj2Lf1c", "reviewer": "Reviewer_jWST", "summary": "This paper proposes a SDP relaxation of the non-convex L1-regularized sparse PCA problem under small noise and incomplete data. The resulting convex program can exactly recover the true support of the sparse leading eigenvector with high probability, under certain conditions. Experiments verify the efficacy of the proposed method.", "questions": "1. Given a n*d  data matrix D with missing entries (say uniformly at random), can you please provide details on how to get the d*d matrix M? \n\n2. Since M^* is assumed to be symmetric, why not using M_ij to impute M_ji ? I think this should also be compared in the experiments.", "limitations": "Yes", "rating": 5, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Strengths: the work provides theoretical guarantee in terms of exactly recovering the true support of the sparse leading eigenvector.\n\nWeakness: \n\n1. There is no comparison with any baseline methods. If I understand correctly, the method simply uses 0 for the missing entries to get matrix M. Then, the existing sparse PCA methods without the consideration of missing entries can also be applied on this matrix M, for comparison.\n\n2.  The proposed method can only guarantee recovering the support of the sparse leading eigenvector, not the value. \n", "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": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}, {"review_id": "MBXzXkPK8H", "reviewer": "Reviewer_oQV2", "summary": "This paper studies sparse PCA with incomplete and noisy data, and proposes a SDP relaxation algorithm. This relaxation formulation was previously applied to sparse PCA with complete data, and this paper shows that the same exact formulation can be applied to incomplete observations. Sufficient conditions for exact support recovery are proposed, and experiments on both synthetic and real-world data.", "questions": "- This is likely out of scope for this paper, but I am curious about how the proposed SDP method compares with other possible approaches to solving sparse PCA with incomplete data, such as various imputation methods. Such an experiment that evaluates the performance of the proposed method will further enhance the results of this paper.\n- It may be helpful to expand on the gene expression data experiment, by showing some relevant visualizations of the data and/or results, for example in the Appendix.", "limitations": "Yes.", "rating": 7, "confidence": 4, "soundness": 4, "presentation": 4, "contribution": 4, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "- Overall, this is a well-written paper. It discusses the problem of sparse PCA with incomplete data, and presents a solid and novel theoretical approach to it via SDP relaxation. \n- The experiment on gene expression data is another highlight, though it feels somewhat insufficiently illustrated.\n- I don't see any major weaknesses of this paper.", "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": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}, {"review_id": "LXZrqFcQgTQ", "reviewer": "Reviewer_oWbf", "summary": "## Summary\nThe authors study the semi-definite relaxation of the sparse principal component analysis (PCA) under missingness and establish the sufficient conditions for support recovery of the first principal component. \n\n## Novelty\nSupport recovery of the semi-definite relaxation of the sparse PCA with complete observations was established by Lei and Vu (2015), this work extends the analysis to incomplete observations.  \n\n## Overview of Comments\n\nI really appreciate the clarity and the level of detail, but describing all these in the main paper devoids it of context and seems like a listing of constants and their mathematical form. Due to this, the paper is somewhat unreadable and does not (because of space) address the key difference in analysis and specific challenges it had to overcome (and how) over the previous works. This description is key for the readers to gauge the merits. For instance, the details of how the analysis differs from Lei and Vu (2015) under the missingness is unclear. ", "questions": "## Comments\n\n1. Applicability of the model in case of missing observations. The paper considers decomposition of a symmetric matrix with noisy and potentially missing observations. Arguably in real-world for a data matrix $Y$, this symmetric matrix was derived by performing $M = Y'Y$. But decomposition of $M$ may be ill-posed when columns of data matrix $Y$ have a few observed entries [See Park & Zhao 2019]. It seems that the theoretical results must impose some additional conditions on $Y$ to avoid such a situation. \n   - 1(a) The paper discusses these indirectly in section 3, but implications for original data $Y$ is not discussed. Although it seems these must depend on the incoherence and the sampling probability, these issues must be directly addressed.\n   - 1(b) The bound on p would be more interpretable if directly expressed, than an upper bound on $\\sqrt(1-p/p)$\n   - 1(c) the various notions of incoherence employed must be introduced in the model description since these are inherent assumptions about the matrix $M$. \n\n\n2. In Lei and Vu (2015) the exact support recovery relies on Sparse Principal Subspace (SPS) condition, which includes 1) spectral gap and 2) $supp(diag(\\Pi)) = J$, where $J = Union_{j=1}^{k} supp(u_j)$, where $u_j$ are the eigenvectors. In essence, $J$ is not the support of the principal eigenvector as the proposed work but a union of all eigenvectors. \n\n   - 2 (a) Since the proposed model does not impose sparsity of other principal components and they can be dense, it is unclear if support recovery is even possible if other eigenvectors are dense. What are the additional conditions required even for the complete observation setting. \n\n   - 2(b) This also seems like a special case of sparse PCA when only the learning principal vector is sparse. \n\n3. Writing -- majority of the write-up is dedicated to definitions of constants and results. I appreciate this level of details, but these are best left to the appendix. In the current form these are extremely distracting tending to unreadable. \n\n   - 3(a) Since the papers relies on Landau notation, these need to further simplified. This is also important for the experiments where the paper claims a that the results are \"tight\" but these relationships are not plotted for reference.\n\n   - 3(b) What are the differences in analysis over Lei and Vu (2015)? What do the incoherence conditions imply in practice?   \n", "limitations": "N/A", "rating": 5, "confidence": 4, "soundness": 2, "presentation": 2, "contribution": 2, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "## Strengths\n\n* Precise about the analysis\n\n## Weaknesses\n\n* While the paper is precise about the analysis, it focuses on listing these results and quantities in the main paper, distracting the readers from the key points of the paper. This makes the paper difficult to read. This also makes it difficult to parse the main result and its implications.\n* The paper does not highlight the key differences in analysis and challenges in establishing the results for this case over the prior art.\n* Applicability of the model for real-world applications in unclear (see detailed comments/questions).", "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": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}, {"review_id": "4NnGfLjkF4t", "reviewer": "Reviewer_owug", "summary": "The paper considers the sparse principal component analysis problem for incomplete and noisy data. In particular, it considers its semidefinite programming reformulated problem and develops a new algorithm with theoretical discussions on the uniqueness of the optimal solution with large probability. Numerical justification including the experiments on synthetic and realistic gene data sets have demonstrated the effectiveness. The work looks interesting with insightful theories. ", "questions": "1. In Section 2.2, do the eigenvectors $\\mathbf{u}_i$'s need to be normalized?\n2. In the SDP reformulation in p.2, does it guarantee that the solution is rank-one and even $X=\\mathbf{x}\\mathbf{x}^T$?\n3. In Theorem 1, what is $p$? How is it defined? Also, the existence of the optimal solution is not discussed. ", "limitations": "The sparse principal component problem assumes that the optimal solution has sparse leading eigenvector which may not be the case in real applications and thus limit its applicability. ", "rating": 5, "confidence": 2, "soundness": 3, "presentation": 2, "contribution": 2, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Strengths: In the main result, several conditions have been proposed to guarantee the success of the SDP algorithm. There are some theoretical contributions. The application of gene expression is interesting and useful. \n\nWeaknesses: The proposed algorithm is based on the SDP reformulation of the original sparse PCA problem, and the algorithmic development may be limited. In Section 4.1, there are no other related methods being compared. ", "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": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}]

Author rebuttal:

Decision: Accept