#  Exact Error in Matrix Completion: Approximately Low-Rank Structures and
Missing Blocks

Agostino Capponi Email<EMAIL_ADDRESS>Department of Industrial
Engineering and Operations Research, Columbia University, 500 West 120th
Street, New York, NY 10027. Mijailo Stojnic  Email<EMAIL_ADDRESS>Department of Industrial Engineering and Operations Research, Columbia
University, 500 West 120th Street, New York, NY 10027.

###### Abstract

We study the completion of approximately low rank matrices with entries
_missing not at random_ (MNAR). In the context of typical large-dimensional
statistical settings, we establish a framework for the performance analysis of
the nuclear norm minimization ($\ell_{1}^{*}$) algorithm. Our framework
produces _exact_ estimates of the _worst-case_ residual root mean squared
error (RMSE) and the associated _phase transitions_ (PT), with both exhibiting
remarkably simple characterizations. Our results enable to precisely quantify
the impact of key system parameters, including data heterogeneity, size of the
missing block, and deviation from ideal low rankness, on the accuracy of
$\ell_{1}^{*}$-based matrix completion. To validate our theoretical worst-case
RMSE estimates, we conduct numerical simulations, demonstrating close
agreement with their numerical counterparts.

62B10,

94A16,

62D10,

Matrix Completion,

Approximately low rank,

Phase transitions,

Nuclear norm,

###### keywords:

[class=MSC]

###### keywords:

##  1 Introduction

Since the inception of matrix completion (MC), a large body of research (see,
e.g., [14, 36, 13, 43]) has focused on the completion when entries are missing
at randomly (MAR) – the most popular example being the completion of the
movie-ratings matrix from the Netflix problem. In this paper, we study a
variant of this problem which has received significant attention in recent
years (see, e.g., [3, 5, 50, 8, 11, 22, 25, 4, 32]), namely the low rank
completion (recovery) with entries missing not at random (MNAR). Specifically,
we consider a missingness pattern which follows a _block structure_ , common
in many social and health science applications. For instance, suppose we want
to estimate whether liver transplant surgery increases the life expectancy of
individuals with cirrhosis. Then the liver transplant is an irreversible
treatment, and estimating (_causally inferring_) the counterfactual entries of
the treated group (i.e., individuals subject to liver transplant surgery)
would amount to estimating the block of entries where rows correspond to
treated patients and columns to times since the treatment started. Another
important example of block recovery arises in financial time series analysis,
where many accounting based variables are missing after having been previously
observed, leading to missing time blocks (see [10]).

The non-randomness of the missing entries invalidates the use of standard MC
analytical methods. Exploring the block structure on the other hand, enables
us to develop a novel analytical method and use it to characterize the
performance of the nuclear norm minimization heuristic (which we denote by
$\ell_{1}^{*}$) when employed for matrix completion and data imputation.
However, some major challenges are faced along the way. A bit of
contextualization helps to better understand them.

_Main challenges:_ The performance analyses of the MC algorithms are usually
done on a _qualitative/scaling_ level. Such analyses provide useful
descriptive characterizations that unfortunately lack a full precision. The
alternative, _precise_ ones, are much harder to obtain and typically fall
under the umbrella of the so-called _phase-transition_ (PT) considerations.
Over the past couple of decades, considerable effort has been dedicated to
obtaining the PTs of MC and their more general counterparts in low-rank
recovery (LRR). This effort has yielded substantial success, as evidenced by
studies such as [21, 18, 35, 15]. However, (with the exception of [15] which,
on the other hand, is a bit more limited in the range of covered parametric
space) all these methods typically relate to, analytically slightly easier,
LRR problems and usually rely on an underlying randomness in acquiring (or
missing to acquire) the data. In other words, they relate to the LRR with the
above mentioned MAR type of data acquisition paradigm. Despite the significant
progress made in these directions, four major challenges have persisted as
formidable hurdles:

(_i_) Providing precise PT performance characterizations, visibly lacking in
existing research, for _truly_ MNAR in the ideal low rank MC; (_ii_) Bridging
the gap between MNAR low-rank recovery (LRR) results and their counterparts in
MC; (_iii_) Establishing a connection between ideal and approximately low rank
MNAR MC; and (_iv_) Fully characterizing the MNAR MC scenarios with matrices
whose internal structure slightly deviates from the ideal low rankness.

_Main contributions:_ We develop a novel mathematical machinery that helps
address all of the four above mentioned challenges. To precisely state our
contributions, a few technical preliminaries are needed. In particular, we
assume the existence of a partially observed (with unobserved entries forming
a block) _approximately low rank_ ground-truth matrix $X_{sol}$. This implies
that $X_{sol}$ is characterized by a limited number of dominant singular
values, which primarily define/describe its heterogeneity, and a larger number
of smaller singular values. Under mild technical assumptions stated in Section
3, we consider the use of the $\ell_{1}^{*}$ heuristic for completion/recovery
of $X_{sol}$. Assuming that $\hat{X}$ is the $X_{sol}$’s estimate produced by
$\ell_{1}^{*}$, our main focus is studying the behavior of
$\textbf{RMSE}_{(ns)}\triangleq\|\hat{X}-X_{sol}\|_{F}$. The main
contributions of our study are as follows: _(i)_ We consider a statistical
approach with minimal distributional assumptions and leverage results from
random matrix theory (in particular, free probability theory (FPT) [46, 47,
48]), to solve the problem of characterizing the performance of $\ell_{1}^{*}$
in the large dimensional approximately low rank MC scenario. _(ii)_ We show
that $\ell_{1}^{*}$ exhibits a sharp _phase-transitioning_ (PT) behavior,
characterized by a boundary that separates regions of system parameters where
the finiteness of the RMSE is guaranteed or not. Consequently, we present a
precise phase transition type of analysis (where the qualitative/scaling types
of estimates are not allowed) and determine the _exact_ location of the
regions separating PT curve. _(iii)_ We show that the allowable _data
heterogeneity_ , which ensures the finiteness of the RMSE, grows linearly with
matrix size. Furthermore, our analysis precisely determines the involved
proportionality constants. _(iv)_ We obtain _explicit exact worst case_
estimates for the (appropriately scaled) residual RMSE and precisely explain
the connections/relations between the obtained RMSEs and PTs.

We conduct a set of numerical experiments which demonstrate a strong agreement
between the theoretical and simulated values, even for rather small matrix
dimensions (of order of a few tens) .

###  1.1 Related prior work

We position our work with respect to three streams of MC literature: (_i_)
general nuclear norm MC; (_ii_) MC with observations _truly_ missing not at
random (MNAR); and (_iii_) performance analysis of nuclear norm optimization
in _true_ MNAR scenarios.

(_i_) The pioneering studies [37, 14, 43] introduced the nuclear norm, in the
ideal low rank context, as a computationally efficient alternative to the
highly nonconvex lowest rank minimization, and provided its generic
performance analysis. Follow-up studies on algorithmic and performance
characterization include, e.g., [16, 29, 26, 38, 31, 33]. Many generic MC
studies, however, rely on various forms of evenly balanced missingness most
efficiently characterized through the so-called missing completely at random
(MCAR) observations. On the other hand, many practically relevant
observational scenarios are not of MCAR type. One, then, typically
distinguishes two settings: 1) observations are still missing at random (MAR),
but not completely at random and potentially dependently on other
observations; and 2) the missing patterns are either fully deterministic (we
say, slightly differently from [39], of _true_ MNAR type) or strongly resemble
a fully deterministic structure. The MAR type of observations have been
extensively studied recently with strong algorithmic and theoretical results
regarding consistency and residual errors. See [24, 28, 12, 44], where the
effects of the nonuniformly distributed missingness are studied on standard
algorithms, and [51, 49, 9, 42, 40, 30], where various algorithmic adaptations
along the lines of IPW, primePCA, modified USVT and Nuclear norm, EM/FISTA,
inverse propensity scoring (IPS), are respectively introduced and their
performances analyzed to account for particular missingness types.

(_ii_) The _true_ MNAR scenarios are analytically very hard and the available
results are scarce. To the best of our knowledge (with the exception of our
own [18]111Some of the results appearing here, particularly those related to
the ideal low rank scenario, were initially presented in a preliminary form in
the conference proceedings [18].), we are the first to analyze a true MNAR
instance on a _phase transition_ level of precision (where scaling types of
estimates, $O(\cdot)$, are not allowed). Nonetheless, it is worth mentioning
several works that study related problems and obtain different types of
results. For example, [3] considers various MNAR scenarios, introduces
specific algorithmic designs, and theoretically analyzes the consistency of
the individual matrix entries estimators (similar results were obtained in
[25, 4, 32]). Recently, methods akin to MC have been employed for factor
modeling in finance. In particular, [50] discusses the role of PCA and the
resulting consistency estimates in several MNAR settings, including both the
block and the staggered ones (see also [22] for further extensions in these
directions). Related lines of work [8] and [11] introduce different types of
algorithms, called Tall-wide and Tall-project, and prove their consistency
estimates.

(_iii_) The most closely related studies to our work are [5] and [18]. In [5],
a noisy emulation of an approximately low rank scenario is used to obtain
qualitative/scaling estimates of the nuclear norm’s RMSE upper bounds. It is
not an MNAR setting, because it uses randomness in treatment selection but
allows for a staggered scenario analysis (for the relevance of such a
scenario, see, e.g., [7, 41, 6]). However, the study also discusses the use of
the nuclear norm for block missing entries. In [18], we consider an idealized
setup where ground truth matrices have low rank (but _not_ approximately low
rank) and provide phase transition results.

_Key Highlights:_ We highlight two major departures of our work relative to
existing studies in the area: _(i)_ (_“noisy vs approximately low rank”_) The
majority of surveyed studies uses low rank matrix plus noise to emulate
approximately low rank settings. While in a qualitative/scaling type of
analysis such an emulation can be appropriate, in phase transition types of
analysis such an emulation may lead to incorrect conclusions. On the other
hand, our model is designed to reflect real-word datasets, such as the
California tobacco data widely used in many studies (e.g.[2]). In such data
sets, deviations from ideal low rankness are directly attributable to the
data’s inherent structure, rather than to inaccuracies stemming from erroneous
recording. _(ii)_ (_“homogeneous vs heterogeneous”_) Most, if not all, of the
above mentioned works either assume fixed rank ($k$) of the dominating part of
the matrix or allow for a weak dependence with the matrix size ($n$), i.e.
allow for $k$ such that $\beta=\lim_{n\rightarrow\infty}\frac{k}{n}\rightarrow
0$. This means that their setups can only capture extremely homogeneous data
settings. In many application contexts, including finance and health sciences,
the large amount of available data is highly heterogeneous. Our paper shows
that the methods we analyze can handle such data settings as well.
Specifically, we not only establish the recoverability of ranks linearly
proportional to the data size but also accurately determine the precise
constants governing this proportionality. Our findings indicate that these
methods can handle high data heterogeneity on a much larger scale than what
assumed in previous studies.

##  2 Preliminaries of mathematical models

We start with the basic description of matrix completion (MC). As it is well
known (see, e.g., [14, 43]), MC assumes that there exists a _“ground truth”_
(presumably approximately low rank) data matrix $X\in{\mathbb{R}}^{n\times n}$
and the so-called masking matrix $M\in{\mathbb{R}}^{n\times n}$

$M_{i,j}=\begin{cases}1,&\mbox{component $(i,j)$ observed}\\\
0,&\mbox{otherwise}.\end{cases}.$ (1)

The following is at the heart of matrix completion

Find $\displaystyle X$ subject to $\displaystyle Y=M\circ X,$ (2)

where $\circ$ stands for the component-wise multiplication and $X$ is “not far
away” from $X_{sol}$ in a metric of choice. Intuitively, the above means the
following: one observes via mask $M$ a collection of elements of $X_{sol}$,
stores the observations in $Y$ and them, inversely, attempts to recover $X$
that fits the observations and is not far away from $X_{sol}$. To provide any
notion and quantification of “not far away” concept the structure of $X_{sol}$
must be known and properly taken into account.

###  2.1 Choice of $X_{sol}$

Various forms of $X_{sol}$’s structure can be considered. We, here, assume a
scenario typically characterized as

$\emph{Approximately low rank:}\qquad Y=M\circ
X_{sol},X_{sol}=U\mbox{diag}(\sigma+\epsilon_{\sigma})V^{T},\|\sigma\|_{0}=k,$

where $U^{T}U=I_{n\times n},V^{T}V=I_{n\times n}$, and the vector of singular
values imperfections, $\epsilon_{\sigma}$, can be small with respect to the
the vector of dominating part of singular values, $\sigma$, in a variety of
metrics. Above, $I_{n\times n}$ denotes the $n\times n$ identity matrix, and
$\mbox{diag}(\cdot)$ the operator that creates a column vector of the diagonal
entries of its matrix argument (or reversely, the diagonal matrix of the
column vector entries).

The standard foundational matrix completion literature (see, e.g., [14, 43,
36, 13]) recognizes two additional forms as of potentially equally important
interest. One of them considers observations imperfections as consequences of
a noisy process $\epsilon_{Y}$

$\begin{array}[]{lcl}\emph{Noisy low rank:}&&Y=M\circ
X_{sol}+\epsilon_{Y},X_{sol}=U\mbox{diag}(\sigma)V^{T},\|\sigma\|_{0}=k.\end{array}$

The other approximates the deviations from low rankness by actually removing
them and assuming the

$\emph{Idealized low rank sceanrio:}\qquad Y=M\circ
X_{sol},X_{sol}=U\mbox{diag}(\sigma)V^{T},\|\sigma\|_{0}=k<n.$

The main advantage of such an idealized simplification is that it usually
allows for mathematical tractability. The potential drawback is a loss of
modeling and estimation accuracy as the real data rarely perfectly fit such an
idealized description.

Which of the three forms is to be used heavily depends on the concrete problem
at hand. For example, although the first two may look as more complete and
advantageous, it is not clear that they would be the right choice in scenarios
where the accuracy of the idealized one is within acceptable tolerance. Also,
choosing between the first and the second needs to be done very carefully and
with respect to the true source of the imperfections (deviations from the
ideal low rankness) which, again, almost exclusively depends on the problem at
hand. As a concrete real data example where the importance of such a choice is
especially important, it is instructive to look at the California tobacco
study from [2], which is a golden standard in econometric literature on
synthetic control. In such an example, one postulates that it is likely that
over years there is a parallel trend in tobacco consumption across different
states. That implies that a matrix containing annual consumption data of each
state is indeed approximately low rank. This intuition is confirmed by
choosing the tobacco consumption matrix that contains the annual per capita
cigarette consumption in the United States and plotting its spectrum (see the
Appendix). Observe that the approximate low rankness comes as a consequence of
the internal structure of the data rather than as a consequence of incorrectly
recorded data. This clearly indicates that the first approach should be
preferable over the second. Moreover, it provides a generic conclusion, that
the first option is preferable whenever the deviation from the ideal low
rankness is caused by the internal structure rather than brought on externally
(as in, say, incorrectly recorded observations scenarios).

It is worth adding that the first two forms can be (and have often been
throughout the literature) interchangeably used as each other’s emulations
(particularly so as facilitating tools in heavy mathematical analyses). The
analysis that we conduct is of the phase transition precision level type,
which implies that the two models _cannot_ be used as replacements of each
other. Moreover, given the econometrics motivating applications, the
approximately low rank one has to be considered. For the concreteness, we
choose $X_{sol}=U\bar{\sigma}V^{T}$, where
$\bar{\sigma}=\mbox{diag}(\sigma,\epsilon_{\sigma})$ and
$\sigma\in{\mathbb{R}}^{k},\epsilon_{\sigma}\in{\mathbb{R}}^{n-k},\sigma_{\epsilon}\in{\mathbb{R}}$,
$\sigma\gg\sigma_{\epsilon}{\bf 1},\epsilon_{\sigma}\leq\sigma_{\epsilon}{\bf
1}$, with $\sigma$, $\epsilon_{\sigma}$, and $\sigma_{\epsilon}$ not changing
as $n$ grows (${\bf 1}$’s are (here and everywhere in the paper) vectors of
all ones of appropriate dimensions; notation $a\gg b$ means that there is a
sufficiently large _constant_ $C$ such that $a\geq Cb$).

###  2.2 Choice of $M$

The problem in (2) (with a random $M$) fits within a standard MC setup. To
capture the non-random block structure of missing entries (visualized in
Figure 1), we consider the masking matrix

$\displaystyle M^{(l_{1},l_{2})}$ $\displaystyle\triangleq$ $\displaystyle{\bf
1}_{n\times 1}{\bf 1}_{n\times 1}^{T}-I^{(l_{1})}(I^{(l_{1})})^{T}{\bf
1}_{n\times 1}{\bf 1}_{n\times 1}^{T}I^{(l_{2})}(I^{(l_{2})})^{T}$
$\displaystyle I^{(l)}$ $\displaystyle\triangleq$
$\displaystyle\begin{bmatrix}{\bf 0}_{l\times(n-l)}\\\
I_{(n-l)\times(n-l)}\end{bmatrix},$ (3)

where ${\bf 1}/{\bf 0}$ are the column vectors and matrices of all ones/zeros
with the dimensions specified in their subscripts.

Figure 1: Matrix $M\triangleq M^{(l_{1},l_{2})}$ – block causal inference
setup

We study the most challenging large $n$ linear regime and, for simplicity of
the presentation, assume the so called square block case with $l\triangleq
l_{1}=l_{2}$ and $\beta\triangleq\lim_{n\rightarrow\infty}\frac{k}{n}$ and
$\eta\triangleq\lim_{n\rightarrow\infty}\frac{l}{n}$. Such an assumption makes
the analysis somewhat easier while imposing no substantial limitations on the
generality of the final results. Namely, the obtained results can always be
used as bounds on the non-square ones (simply trimming rectangles so that they
become squares results in using less available data and obtaining larger
residual RMSE’s).

##  3 Mathematical analysis

We start by defining $\ell_{p}^{*}(X)$ as the so-called $\ell_{p}$ (quasi)
norm of $\sigma(X)$ (the vector of singular values of $X$),

$\displaystyle\ell_{p}^{*}(X)\triangleq\ell_{p}(\sigma(X)),p\in{\mathbb{R}}_{+},\quad\mbox{with}\quad\ell_{0}^{*}(X)\triangleq\lim_{p\rightarrow
0}\ell_{p}^{*}(X).$ (4)

Instead of solving (2), the MC literature often opts for its tightest convex
relaxation

$\displaystyle\hat{X}\triangleq\min_{X}$ $\displaystyle\ell_{1}^{*}(X)$
subject to $\displaystyle Y=M\circ X,Y\triangleq M\circ
X_{sol}.\vspace{-.0in}$ (5)

We refer to the introductory papers [37, 14, 43, 23] for additional background
on the relation between $\ell_{0}^{*}$ and $\ell_{1}^{*}$ and connections to
low rank recovery (LRR) and MC (related considerations, can be found in a long
line of excellent follow-up work, e.g., [15, 16, 17, 26, 27, 28, 29, 33, 34,
38, 31]). It is well known in the structured recovery literature that, in the
idealized scenario, $\sigma_{\epsilon}=0$, the above heuristics often produce
the true solution, i.e., one often has $\hat{X}=X_{sol}$. What is somewhat
hidden in such a folklore knowledge, is that the masking matrix $M$ is
typically viewed as well balanced or comprised of randomly located zeros and
ones. The problem we study, though, presents three key differences: 1)
scenario of interest is _non-idealized_ , i.e., $\sigma_{\epsilon}>0$; 2) $M$
is not random; and 3) the locations of zeros and ones form a very particular
_block_ structure. These differences make it impossible to follow
standard/known CS/MC analytical paths. Instead, we take the differences to our
advantage, avoid the known analytical paths, and create a new one to
characterize the $\ell_{1}^{*}$’s performance. Since the considered scenario
is non-idealized, there will be a nonzero residual estimation error
$\hat{W}\triangleq\hat{X}-X_{sol}$. Any $\ell_{1}^{*}$ performance analysis
then must account for a characterization of $\hat{W}$. Our approach is a no
exception in that sense. However, as particularly relevant to us, we highlight
below the _precision_ of the conducted analysis.

Ultimate performance analysis goal:
$\displaystyle\mbox{\text@underline{{\emph{Precise}}} characterization
of}\hskip
2.168pt\mathbf{RMSE}_{(ns)}\triangleq\|\hat{W}\|_{F}\triangleq\|\hat{X}-X_{sol}\|_{F}.$

Below, we provide such a characterization relying on typical statistical
viewing of the problem and utilizing advanced sophisticated probabilistic
concepts that we introduce along the way. We begin with the following
algebraic result.

######  Theorem 3.1.

(Algebraic characterization of $W$) Consider a
$\bar{U}\in{\mathbb{R}}^{n\times k}$ such that $\bar{U}^{T}\bar{U}=I_{k\times
k}$ and a $\bar{V}\in{\mathbb{R}}^{n\times k}$ such that
$\bar{V}^{T}\bar{V}=I_{k\times k}$ and an approximately rank $k$ matrix
$X_{sol}=X\in{\mathbb{R}}^{n\times n}$, such that
$X_{sol}=U\bar{\sigma}V^{T},\bar{\sigma}=\mbox{diag}(\sigma,\epsilon_{\sigma}),\sigma\in{\mathbb{R}}^{k},\epsilon_{\sigma}\in{\mathbb{R}}^{n-k},\sigma\gg\sigma_{\epsilon},\epsilon_{\sigma}\leq\sigma_{\epsilon}$.
Also, let the orthogonal spans $\bar{U}^{\perp}\in{\mathbb{R}}^{n\times(n-k)}$
and $\bar{V}^{\perp}\in{\mathbb{R}}^{n\times(n-k)}$ be such that
$U\triangleq\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix}$ and
$V\triangleq\begin{bmatrix}\bar{V}&\bar{V}^{\perp}\end{bmatrix}$ and

$U^{T}U\triangleq\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix}^{T}\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix}=I_{n\times
n}\quad\mbox{and}\quad
V^{T}V\triangleq\begin{bmatrix}\bar{V}&\bar{V}^{\perp}\end{bmatrix}^{T}\begin{bmatrix}\bar{V}&\bar{V}^{\perp}\end{bmatrix}=I_{n\times
n}.$ (6)

With $M\in{\mathbb{R}}^{n\times n}$ as in (2.2), assume that $Y=M\circ
X_{sol}$ and let $\hat{X}$ be the solution of (3). Set

$\hat{W}=\arg\min_{W,M\circ W={\bf 0}_{n\times n}}{\rm
tr}\,(\bar{U}^{T}W\bar{V})+\ell_{1}^{*}(\mbox{diag}(\epsilon_{\sigma})+(\bar{U}^{\perp})^{T}W\bar{V}^{\perp}).$
(7)

Then for any $X_{sol}$ that satisfies the above setup one has

$\textbf{\emph{RMSE}}_{(ns)}=\|\hat{X}-X_{sol}\|_{F}\leq\|\hat{W}\|_{F}.$ (8)

Moreover, there exists an $X_{sol}$ such that

$\textbf{\emph{RMSE}}_{(ns)}=\|\hat{X}-X_{sol}\|_{F}=\|\hat{W}\|_{F}.$ (9)

###### Proof.

Presented in the Appendix. ∎

To analyze the optimization problem presented in (115), we adopt an approach
that relies on using a generic Lagrangian mechanism. The initial step in this
process involves explicitly formulating the optimization problem from (115)

$\displaystyle f_{pr}(M;U,V)\triangleq\min_{W}$ $\displaystyle{\rm
tr}\,(\bar{U}^{T}W\bar{V})+\ell_{1}^{*}(\mbox{diag}(\epsilon_{\sigma})+(\bar{U}^{\perp})^{T}W\bar{V}^{\perp})$
subject to $\displaystyle M\circ W={\bf 0}_{n\times n}.$ (10)

One can then write the corresponding Lagrangian and the Lagrange dual
function,

$\displaystyle{\cal L}(W,\Lambda)$ $\displaystyle\triangleq$
$\displaystyle{\rm
tr}\,(\bar{U}^{T}W\bar{V})+\ell_{1}^{*}(\mbox{diag}(\epsilon_{\sigma})+(\bar{U}^{\perp})^{T}W\bar{V}^{\perp})+\Theta(M\circ
W)$ (11) $\displaystyle=$ $\displaystyle\max_{\Lambda,\Lambda^{T}\Lambda\leq
I}\left({\rm tr}\,(\bar{U}^{T}W\bar{V})+{\rm
tr}\,(\Lambda(\mbox{diag}(\epsilon_{\sigma})+(\bar{U}^{\perp})^{T}W\bar{V}^{\perp}))+\Theta(M\circ
W)\right)$ $\displaystyle=$ $\displaystyle\max_{\Lambda,\Lambda^{T}\Lambda\leq
I}\left({\rm
tr}\,\left((\bar{V}\bar{U}^{T}+\bar{V}^{\perp}\Lambda(\bar{U}^{\perp})^{T}+\Theta\circ
M)W\right)+{\rm tr}\,(\Lambda\mbox{diag}(\epsilon_{\sigma}))\right),$

and

$\displaystyle g(\Theta)$ $\displaystyle\triangleq$
$\displaystyle\min_{W}{\cal L}(W,\Lambda,\Theta,\gamma)$ (12) $\displaystyle=$
$\displaystyle\min_{W}\max_{\Theta,\Lambda,\Lambda^{T}\Lambda\leq I}\left({\rm
tr}\,\left((\bar{V}\bar{U}^{T}+\bar{V}^{\perp}\Lambda(\bar{U}^{\perp})^{T}+\Theta\circ
M)W\right)+{\rm tr}\,(\Lambda\mbox{diag}(\epsilon_{\sigma})\right).$

Utilizing the Lagrangian duality one then further has

$\displaystyle g(\Theta)$ $\displaystyle=$
$\displaystyle\min_{W}\max_{\Lambda,\Lambda^{T}\Lambda\leq I}\left({\rm
tr}\,\left((\bar{V}\bar{U}^{T}+\bar{V}^{\perp}\Lambda(\bar{U}^{\perp})^{T}+\Theta\circ
M)W\right)+{\rm tr}\,(\Lambda\mbox{diag}(\epsilon_{\sigma})\right)$ (13)
$\displaystyle\geq$ $\displaystyle\max_{\Lambda,\Lambda^{T}\Lambda\leq
I}\min_{W}\left({\rm
tr}\,\left((\bar{V}\bar{U}^{T}+\bar{V}^{\perp}\Lambda(\bar{U}^{\perp})^{T}+\Theta\circ
M)W\right)+{\rm tr}\,(\Lambda\mbox{diag}(\epsilon_{\sigma})\right).$

Moreover, using the definition of $M$ from (2.2) we have

$\displaystyle f_{pr}(M;U,V)$ $\displaystyle\geq$
$\displaystyle\max_{\Theta}g(\Theta)$ (14) $\displaystyle=$
$\displaystyle\max_{\Lambda,\Lambda^{T}\Lambda\leq
I,(I^{(l)})^{T}(\bar{V}\bar{U}^{T}+\bar{V}^{\perp}\Lambda(\bar{U}^{\perp})^{T}I^{(l)}=0}{\rm
tr}\,(\Lambda\mbox{diag}(\epsilon_{\sigma})).$

It is now critically important to observe that the above holds as long as the
following optimizing condition is feasible

$\exists\Lambda|\Lambda^{T}\Lambda\leq
I\quad\mbox{and}\quad(I^{(l)})^{T}\left(\bar{V}\bar{U}^{T}+\bar{V}^{\perp}\Lambda(\bar{U}^{\perp})^{T}\right)I^{(l)}=0.$
(15)

On the other hand, as soon as the system parameters are such that it becomes
infeasible, one has that both optimization in (13) and the residual error
norm, $\|W\|_{F}$, are unbounded. Clearly, to be able to conduct the analysis
of the residual error, one must first determine the range or parameters where
such error is bounded, i.e., where the optimizing condition in (14) (or (117))
is feasible.

###  3.1 Feasibility of the optimizing condition in (14)

For the time being we assume $k\leq l$ (later on this assumption will be
rigorously justified). From (117), one then easily has

$\Lambda=((I^{(l)})^{T}\bar{V}^{\perp})^{-1}(I^{(l)})^{T}\bar{V}\bar{U}^{T}I^{(l)}((\bar{U}^{\perp})^{T}I^{(l)})^{-1}\Longrightarrow(I^{(l)})^{T}\left(\bar{V}\bar{U}^{T}+\bar{V}^{\perp}\Lambda(\bar{U}^{\perp})^{T}\right)I^{(l)}=0,$
(16)

where $(\cdot)^{-1}$ stands for the pseudo-inverse. To make writing a bit
easier we can set

$\displaystyle\Lambda_{opt}$ $\displaystyle\triangleq$
$\displaystyle\Lambda_{V}\Lambda_{U}^{T}$ $\displaystyle\Lambda_{V}$
$\displaystyle\triangleq$
$\displaystyle((I^{(l)})^{T}\bar{V}^{\perp})^{-1}(I^{(l)})^{T}\bar{V}$
$\displaystyle\Lambda_{U}$ $\displaystyle\triangleq$
$\displaystyle((I^{(l)})^{T}\bar{U}^{\perp})^{-1}(I^{(l)})^{T}\bar{U}.$ (17)

Let $\lambda_{max}(\cdot)$ be the maximum eigenvalue of its symmetric matrix
argument. After combining (117) and Theorem 9, we conclude that if

$\displaystyle\lambda_{max}(\Lambda_{opt}^{T}\Lambda_{opt})\leq 1,$ (18)

then (117) will be satisfied. After basic algebraic transformations, (18) can
also be rewritten as

$\lambda_{max}(\Lambda_{opt}^{T}\Lambda_{opt})=\lambda_{max}(\Lambda_{opt}\Lambda_{opt}^{T})=\lambda_{max}(\Lambda_{V}\Lambda_{U}^{T}\Lambda_{U}\Lambda_{V}^{T})=\lambda_{max}(\Lambda_{V}^{T}\Lambda_{V}\Lambda_{U}^{T}\Lambda_{U})\leq
1.$ (19)

From (19), it is clear that the spectrum of
$\Lambda_{V}^{T}\Lambda_{V}\Lambda_{U}^{T}\Lambda_{U}$ as well as the spectra
of $\Lambda_{V}^{T}\Lambda_{V}$ and $\Lambda_{U}^{T}\Lambda_{U}$ play an
important role in the $\ell_{0}^{*}-\ell_{1}^{*}$-equivalence. We first
observe a _worst case_ bound. Namely, since

$\displaystyle\lambda_{max}(\Lambda_{opt}^{T}\Lambda_{opt})=\lambda_{max}(\Lambda_{V}^{T}\Lambda_{V}\Lambda_{U}^{T}\Lambda_{U})\leq\lambda_{max}(\Lambda_{V}^{T}\Lambda_{V})\lambda_{max}(\Lambda_{U}^{T}\Lambda_{U}),$
(20)

one has that if the individual spectra of $\Lambda_{V}^{T}\Lambda_{V}$ and
$\Lambda_{U}^{T}\Lambda_{U}$ do not exceed one then the
$\ell_{0}^{*}-\ell_{1}^{*}$-equivalence holds. Given the obvious importance of
these spectra, we below look at them in more detail. Clearly, due to symmetry,
it is sufficient to focus on only one of them. To that end we start by
observing

$\displaystyle((I^{(l)})^{T}\bar{V}^{\perp})^{-1}$ $\displaystyle=$
$\displaystyle(\bar{V}^{\perp})^{T}I^{(l)}\left((I^{(l)})^{T}\bar{V}^{\perp}(\bar{V}^{\perp})^{T}I^{(l)}\right)^{-1}.$
(21)

From (21), one quickly finds

$\displaystyle\left(((I^{(l)})^{T}\bar{V}^{\perp})^{-1}\right)^{T}((I^{(l)})^{T}\bar{V}^{\perp})^{-1}$
$\displaystyle=$
$\displaystyle\left((I^{(l)})^{T}\bar{V}^{\perp}(\bar{V}^{\perp})^{T}I^{(l)}\right)^{-1}.$
(22)

We can then write

$\displaystyle(I^{(l)})^{T}\bar{V}\bar{V}^{T}I^{(l)}$ $\displaystyle=$
$\displaystyle(I^{(l)})^{T}\left(I-\bar{V}^{\perp}(\bar{V}^{\perp})^{T}\right)I^{(l)}.$
(23)

From (17), we also have

$\displaystyle Q_{1}\triangleq\Lambda_{V}^{T}\Lambda_{V}$ $\displaystyle=$
$\displaystyle\left(\left(((I^{(l)})^{T}\bar{V}^{\perp})^{-1}(I^{(l)})^{T}\bar{V}\right)^{T}\left(((I^{(l)})^{T}\bar{V}^{\perp})^{-1}(I^{(l)})^{T}\bar{V}\right)\right)$
(24) $\displaystyle=$
$\displaystyle\bar{V}^{T}I^{(l)}\left(((I^{(l)})^{T}\bar{V}^{\perp})^{-1}\right)^{T}((I^{(l)})^{T}\bar{V}^{\perp})^{-1}\left((I^{(l)})^{T}\bar{V}\right).$

Now, we find it more convenient to work with a slightly changed version of
matrix $Q_{1}$. Namely, after combining (17), (22), and (23), we obtain

$\displaystyle Q$ $\displaystyle\triangleq$
$\displaystyle\left(\left(((I^{(l)})^{T}\bar{V}^{\perp})^{-1}\right)^{T}((I^{(l)})^{T}\bar{V}^{\perp})^{-1}\left((I^{(l)})^{T}\bar{V}\bar{V}^{T}I^{(l)}\right)\right)$
(25) $\displaystyle=$
$\displaystyle\left(\left((I^{(l)})^{T}\bar{V}^{\perp}(\bar{V}^{\perp})^{T}I^{(l)}\right)^{-1}\left((I^{(l)})^{T}\left(I-\bar{V}^{\perp}(\bar{V}^{\perp})^{T}\right)I^{(l)}\right)\right)$
$\displaystyle=$
$\displaystyle\left((I^{(l)})^{T}\bar{V}^{\perp}(\bar{V}^{\perp})^{T}I^{(l)}\right)^{-1}-I.$

Clearly, all the nonzero eigenvalues of $Q_{1}$ and $Q$ are identical. When
$k\leq n-l$ then $Q$ has all the eigenvalues of $Q_{1}$ plus $n-l-k$ extra
zeros. On the other hand, when $k\geq n-l$ then $Q_{1}$ has all the
eigenvalues of $Q$ plus $k-(n-l)$ extra zeros. Since adding or removing zeros
from the spectra does not change any of their features of our interests here,
instead of working directly with $Q_{1}$, we can work with $Q$. In particular,
we have

$\lambda_{max}(Q_{1})=\lambda_{max}(\Lambda_{V}^{T}\Lambda_{V})=\lambda_{max}\left(\left((I^{(l)})^{T}\bar{V}^{\perp}(\bar{V}^{\perp})^{T}I^{(l)}\right)^{-1}\right)-1=\lambda_{max}(Q).$
(26)

We are now in position to establish a spectral alternative to Corollary 2 from
[19].

######  Proposition 1.

(Idealized low rank scenario — $\ell_{0}^{*}-\ell_{1}^{*}$-equivalence
condition via mask-modified bases spectra) Assume the setup of Proposition 1
and Corollaries 1 and 2 from [19] with $k\leq l$. Let $\lambda_{V}$ and
$\lambda_{U}$ be defined as in (17). Then

$\ell_{0}^{*}\Longleftrightarrow\ell_{1}^{*}\quad\mbox{and}\quad\textbf{\emph{RMSE}}=\|\hat{X}-X_{sol}\|_{F}=0$
$\mbox{If and only
if}\quad\lambda_{max}(\Lambda_{V}^{T}\Lambda_{V}\Lambda_{U}^{T}\Lambda_{U})\leq
1.$ (27)

Moreover, if

$\left(\lambda_{max}\left(\left((I^{(l)})^{T}\bar{V}^{\perp}(\bar{V}^{\perp})^{T}I^{(l)}\right)^{-1}\right)-1\right)\left(\lambda_{max}\left(\left((I^{(l)})^{T}\bar{U}^{\perp}(\bar{U}^{\perp})^{T}I^{(l)}\right)^{-1}\right)-1\right)\leq
1,$ (28)

then again $\ell_{0}^{*}\Longleftrightarrow\ell_{1}^{*}$ and
$\textbf{\emph{RMSE}}=\|\hat{X}-X_{sol}\|_{F}=0$.

###### Proof.

The first part follows from Corollaries 1 and 2 from [19], (16), (17), (19),
the above discussion and some additional considerations while the second part
follows by additionally taking into account (20) and (26). We below present
all the details split into three parts: the first two relate to the
equivalence condition (equation (27)) while the third one relates to (28).

_1) $\Longrightarrow$ – The “if part” of condition (27):_ Choosing
$\Lambda=\Lambda_{opt}$

$\displaystyle\Lambda_{opt}$ $\displaystyle\triangleq$
$\displaystyle\Lambda_{V}\Lambda_{U}^{T}=-((I^{(l)})^{T}\bar{V}^{\perp})^{-1}(I^{(l)})^{T}\bar{V}\bar{U}^{T}I^{(l)}((\bar{U}^{\perp})^{T}I^{(l)})^{-1},$
(29)

(where, we recall that, $(\cdot)^{-1}$ stands for the pseudo-inverse) one
ensures

$(I^{(l)})^{T}\left(\bar{V}\bar{U}^{T}+\bar{V}^{\perp}\Lambda(\bar{U}^{\perp})^{T}\right)I^{(l)}=0.$
(30)

Let $\lambda_{max}(\cdot)$ be the maximum eigenvalue of its symmetric matrix
argument. A combination of (117)-(30) ensures that if

$\displaystyle\lambda_{max}(\Lambda_{opt}^{T}\Lambda_{opt})\leq 1,$ (31)

then $\Lambda_{opt}$ satisfies (117). (31) is implied by (27) since

$\lambda_{max}(\Lambda_{opt}^{T}\Lambda_{opt})=\lambda_{max}(\Lambda_{U}\Lambda_{V}^{T}\Lambda_{V}\Lambda_{U}^{T})=\lambda_{max}(\Lambda_{V}^{T}\Lambda_{V}\Lambda_{U}^{T}\Lambda_{U})\leq
1,$

which suffices to complete the proof of the “if part”.

_2) $\Longleftarrow$ – The “only if part” of condition (27):_ Consider SVDs

$B\triangleq(I^{(l)})^{T}V^{\perp}=U_{B}\Sigma_{B}V_{B}^{T},\quad
C\triangleq(I^{(l)})^{T}U^{\perp}=U_{C}\Sigma_{C}V_{C}^{T}$ (32)

with unitary $U_{B},V_{B},U_{C},V_{C}$ and diagonal (with no zeros on the main
diagonal) $\Sigma_{B},\Sigma_{C}$. Any $\Lambda$ can be parameterized as

$\Lambda=V_{B}H^{T}+V_{B}^{\perp}D^{T},\quad H\triangleq
V_{C}E+V_{C}^{\perp}F,$ (33)

for some $E,F,D$ and unitary $V_{B}^{\perp}$ and $V_{C}^{\perp}$ such that
$V_{B}^{T}V_{B}^{\perp}=V_{C}^{T}V_{C}^{\perp}=0$. Also, one can set
$\Lambda_{*}$ and write the SVD of $E$

$\Lambda_{*}\triangleq V_{B}E^{T}V_{C}^{T},\quad E=U_{E}\Sigma_{E}V_{E}^{T},$
(34)

where $U_{E},V_{E}$ are unitary and $\Sigma_{E}$ is diagonal with entries on
the main diagonal being nonzero and in the ascending order. Let ${\bf u}_{e}$
be the last column of $U_{E}$ (i.e. the eigenvector of $EE^{T}$ that
corresponds to its largest eigenvalue). Since $\|V_{C}{\bf u}_{e}\|_{2}=1$,

$\displaystyle\lambda_{max}(\Lambda^{T}\Lambda)$ $\displaystyle\geq$
$\displaystyle{\bf u}_{e}^{T}V_{C}^{T}\Lambda^{T}\Lambda V_{C}{\bf u}_{e}$
(35) $\displaystyle=$ $\displaystyle{\bf u}_{e}^{T}V_{C}^{T}HH^{T}V_{C}{\bf
u}_{e}+{\bf u}_{e}^{T}V_{C}^{T}DD^{T}V_{C}{\bf u}_{e}$ $\displaystyle\geq$
$\displaystyle{\bf
u}_{e}^{T}V_{C}^{T}(V_{C}E+V_{C}^{\perp}F)(V_{C}E+V_{C}^{\perp}F)^{T}V_{C}{\bf
u}_{e}$ $\displaystyle=$ $\displaystyle{\bf u}_{e}^{T}EE^{T}{\bf
u}_{e}=\lambda_{max}(EE^{T})$ $\displaystyle=$
$\displaystyle\lambda_{max}(V_{C}EE^{T}V_{C}^{T})=\lambda_{max}(\Lambda_{*}^{T}\Lambda_{*}).$

If $\Lambda$ satisfies the condition of (117) then a combination of (117) and
(32)-(34) gives

$(I^{(l)})^{T}\bar{V}\bar{U}I^{(l)}+B\Lambda_{*}C^{T}=0,$ (36)

and a combination of (17), (32), and (36) gives

$\Lambda_{V}\Lambda_{U}^{T}=-B^{-1}(I^{(l)})^{T}\bar{V}\bar{U}I^{(l)}(C^{T})^{-1}=\Lambda_{*}.$
(37)

Finally, for $\Lambda$ that fits (117), from (35) and (37) one has

$\displaystyle 1$ $\displaystyle\geq$
$\displaystyle\lambda_{max}(\Lambda^{T}\Lambda)>\lambda_{max}(\Lambda_{*}^{T}\Lambda_{*})=\lambda_{max}(\Lambda_{*}\Lambda_{*}^{T})$
(38) $\displaystyle=$
$\displaystyle\lambda_{max}(\Lambda_{V}\Lambda_{U}^{T}\Lambda_{U}\Lambda_{V}^{T})=\lambda_{max}(\Lambda_{V}^{T}\Lambda_{V}\Lambda_{U}^{T}\Lambda_{U}),$

which completes the proof of the “only if part”.

_3) Suffciency of the condition ( 28):_ Since

$\lambda_{max}(\Lambda_{V}^{T}\Lambda_{V}\Lambda_{U}^{T}\Lambda_{U})\leq\lambda_{max}(\Lambda_{V}^{T}\Lambda_{V})\lambda_{max}(\Lambda_{U}^{T}\Lambda_{U}),$
(39)

one has that if the individual spectra of $\Lambda_{V}^{T}\Lambda_{V}$ and
$\Lambda_{U}^{T}\Lambda_{U}$ do not exceed one then the
$\ell_{0}^{*}-\ell_{1}^{*}$-equivalence holds. Repeating (21)-(25), with $V$
replaced by $U$, $Q_{1}$ by $Q_{1}^{\perp}$, and $Q$ by $Q^{\perp}$ one
arrives at the following analogue of (26)

$\lambda_{max}(\Lambda_{U}^{T}\Lambda_{U})=\lambda_{max}(Q_{1}^{\perp})=\lambda_{max}(Q^{\perp}).$
(40)

A combination of (39) and (26) - (40) completes the proof of the condition
(28)’s sufficiency for the $\ell_{0}^{*}-\ell_{1}^{*}$-equivalence. ∎

The feasibility of the optimizing condition is precisely what determines the
so-called phase transition in the idealized scenario where
$\epsilon_{\sigma}=0$. The above machinery is actually powerful enough to
precisely determine, in a statistical context, the location of the PT curve.
We assume a generic _typical_ statistical scenario, where the orthogonal spans
of interests are Haar distributed (this is a bit less restrictive than the
usual i.i.d. Gaussian assumptions from, say, factor models literature, see,
e.g., [50], and it does not impose _per se_ the so-called _strong_
factor/loading assumptions see, e.g. [8, 11] (note that Gaussians are only a
special case of the Haar distributed objects)).222The assumed rotational
invariance is effectively analogous to the assumption of randomly positioning
the sparse portion within the unknown vectors in basic compressed sensing
(with the exception of results of [20], this assumption is almost always in
place in the theoretical compressed sensing _phase transition_
characterizations). In the absence of any _a priori_ statistically available
knowledge about the data, the above is the minimum that realistically _must_
be assumed. In other words, if, beyond the degree of heterogeneity, nothing
more is known about the data, no modeling concept can be proven as better than
_uniform_ randomness.

We first need the following lemma.

######  Lemma 3.2.

Assume large $n$ linear regime with

$\displaystyle\beta\triangleq\lim_{n\rightarrow\infty}\frac{k}{n}\quad\mbox{and}\quad\eta\triangleq\lim_{n\rightarrow\infty}\frac{l}{n}.$
(41)

Let $\bar{V}^{\perp}\in{\mathbb{R}}^{n\times(n-k)}$ be a Haar distributed
basis of an $n-k$-dimensional subspace of $R^{n}$. Analogously, let
$\bar{U}_{D}^{\perp}\in{\mathbb{R}}^{n\times(n-k)}$ be a Haar distributed
basis of an $n-l$-dimensional subspace of $R^{n}$. Moreover, let
$\bar{V}^{\perp}\in{\mathbb{R}}^{n\times(n-k)}$ and
$\bar{U}_{D}^{\perp}\in{\mathbb{R}}^{n\times(n-k)}$ be independent of each
other. Also, let ${\cal V}$, ${\cal U}$, and $\tilde{D}$ be defined as

${\cal V}\triangleq\bar{V}^{\perp}(\bar{V}^{\perp})^{T},\quad{\cal
U}\triangleq\bar{U}_{D}^{\perp}(\bar{U}_{D}^{\perp})^{T},\quad\tilde{D}\triangleq{\cal
V}{\cal U}.$ (42)

Set $x_{l}$ and $x_{u}$ as

$\displaystyle x_{l}$ $\displaystyle\triangleq$
$\displaystyle\beta+\eta-2\beta\eta-\sqrt{(\beta+\eta-2\beta\eta)^{2}-(\beta-\eta)^{2}}$
$\displaystyle x_{u}$ $\displaystyle\triangleq$
$\displaystyle\beta+\eta-2\beta\eta+\sqrt{(\beta+\eta-2\beta\eta)^{2}-(\beta-\eta)^{2}},$
(43)

Then the limiting spectral distribution of $\tilde{D}$, $f_{\tilde{D}}(x)$, is

$\displaystyle f_{\tilde{D}}(x)$ $\displaystyle=$ $\displaystyle
f_{0}\delta(x)+f_{\tilde{D}}^{(b)}(x)+f_{1}f_{0}\delta(x-1)$ (44)
$\displaystyle=$
$\displaystyle\max(\beta,\eta)\delta(x)+f_{\tilde{D}}^{(b)}(x)+\max(1-(\beta+\eta),0)\delta(x-1),$

with

$\displaystyle f_{\tilde{D}}^{(b)}(x)$ $\displaystyle=$
$\displaystyle\begin{cases}\frac{\sqrt{-(x-(\beta+\eta))^{2}-4\beta\eta(x-1)}}{2\pi(x-x^{2})},&\mbox{if
}x_{l}\leq x\leq x_{u}.\\\ 0,&\mbox{otherwise}.\end{cases}$ (45)

###### Proof.

Presented in the Appendix. ∎

Now we write

$\ell_{0}^{*}-\ell_{1}^{*}-\mbox{equivalence}\quad\Longleftrightarrow\quad\frac{1}{x_{l}}-1\leq
1.$ (46)

Recalling on (198)

$\displaystyle x_{l}$ $\displaystyle\triangleq$
$\displaystyle\beta+\eta-2\beta\eta-\sqrt{(\beta+\eta-2\beta\eta)^{2}-(\beta-\eta)^{2}},$
(47)

and combining further with (46), we obtain

$\displaystyle\begin{array}[]{rrcl}&\frac{1}{2}&\leq&x_{l}\\\
\Longleftrightarrow&\frac{1}{2}&=&\beta+\eta-2\beta\eta-\sqrt{(\beta+\eta-2\beta\eta)^{2}-(\beta-\eta)^{2}}\\\
\Longleftrightarrow&(\beta+\eta-2\beta\eta)^{2}-(\beta-\eta)^{2}&\leq&\left(\beta+\eta-2\beta\eta-\frac{1}{2}\right)^{2}\\\
\Longleftrightarrow&(\beta+\eta-2\beta\eta)^{2}-(\beta-\eta)^{2}&\leq&\left(\beta+\eta-2\beta\eta\right)^{2}-\left(\beta+\eta-2\beta\eta\right)+\frac{1}{4}\\\
\Longleftrightarrow&-(\beta-\eta)^{2}&\leq&-\left(\beta+\eta-2\beta\eta\right)+\frac{1}{4}\\\
\Longleftrightarrow&0&\leq&\beta^{2}+\eta^{2}-\left(\beta+\eta\right)+\frac{1}{4}\\\
\Longleftrightarrow&\eta-\eta^{2}&\leq&\left(\frac{1}{2}-\beta\right)^{2}\\\
\Longleftrightarrow&\beta&\leq&\frac{1}{2}-\sqrt{\eta-\eta^{2}}.\end{array}$
(56)

From (46) and (56), we finally have

$\ell_{0}^{*}-\ell_{1}^{*}-\mbox{equivalence}\quad\Longleftrightarrow\quad\beta\leq\frac{1}{2}-\sqrt{\eta-\eta^{2}}.$
(57)

We are now in position to formalize the above results and establish the so-
called phase-transition phenomenon as well as its a precise _worst case_
location in a _typical_ statistical scenario.

######  Theorem 3.3.

($\ell_{1}^{*}$ – phase transition – (typical worst case)) Consider a rank-$k$
matrix $X_{sol}=X\in{\mathbb{R}}^{n\times n}$ with the Haar distributed (_not
necessarily independent_) bases of its orthogonal row and column spans
$\bar{U}^{\perp}\in{\mathbb{R}}^{n\times(n-k)}$ and
$\bar{V}^{\perp}\in{\mathbb{R}}^{n\times(n-k)}$
($X_{sol}^{T}\bar{U}^{\perp}=X_{sol}\bar{V}^{\perp}={\bf 0}_{n\times(n-k)}$).
Let $M\triangleq M^{(l)}\in{\mathbb{R}}^{n\times n}$ be as defined in (2.2).
Assume a large $n$ linear regime with
$\beta\triangleq\lim_{n\rightarrow\infty}\frac{k}{n}$ and
$\eta\triangleq\lim_{n\rightarrow\infty}\frac{l}{n}$ and let $\beta_{wc}$ and
$\eta$ satisfy the following

$\ell_{1}^{*}$ worst case phase transition (PT) characterization
$\xi_{\eta}^{(wc)}(\beta)\triangleq\beta-\frac{1}{2}+\sqrt{\eta-\eta^{2}}=0.$
(58)

If and only if $\beta\leq\beta_{wc}$

$\lim_{n\rightarrow\infty}{\mathbb{P}}(\ell_{0}^{*}\Longleftrightarrow\ell_{1}^{*})=\
\lim_{n\rightarrow\infty}{\mathbb{P}}(\mathbf{RMSE}=0)=1,$ (59)

and the solutions of (2) and (3) coincide with overwhelming probability.

###### Proof.

The “if” part follows from Proposition 1, Lemma 3.2, and the above discussion.
The “only if” part additionally assumes $\bar{U}^{\perp}=\bar{V}^{\perp}$
which, as discussed above, ensures that the results are, in the worst case,
achievable. ∎

We visually illustrate the results from the above theorem in Figure 2. As can
be seen from the figure, the phase transition curve splits the entire
$(\beta,\eta)$ region into two subregions. The first of the subregions is
below (or to the right of) the curve and in that region the
$\ell_{0}^{*}-\ell_{1}^{*}$-equivalence phenomenon occurs. This means that one
can recover an ideally low rank $X_{sol}$ masked by $M$ as in (2) via the
$\ell_{1}^{*}$ heuristic from (3) with the residual mean square error (RMSE)
equal to zero. In other words, for the system parameters $(\beta,\eta)$ that
belong to the subregion below the curve one has a perfect recovery with
$X_{sol}$ and $\hat{X}$ (the respective solutions of (2) and (3)) being equal
to each other and consequently with $\textbf{RMSE}=\|\hat{X}-X_{sol}\|_{F}=0$.
On the other hand, in the subregion above the curve, the $\ell_{1}^{*}$
heuristic fails and one can even find an $X_{sol}$ for which
$\textbf{RMSE}\rightarrow\infty$.

Figure 2: Typical _worst case_ $\ell_{1}^{*}$ phase transition (ideal low rank
context (block causal inference – C-inf))

As mentioned after (14), one should observe that (14) holds as long as the
optimizing condition is feasible. The above theorem precisely identifies the
relation between the system parameters $\beta$ and $\eta$ when that happens.
Visually speaking, that happens whenever $(\eta,\beta)$ are below (or to the
right of) the curve in Figure 2. As soon as the system parameters are above
the curve, the optimizing condition becomes infeasible, and one has that both
the optimization in (13) and the residual error norm, $\|W\|_{F}$, are, in
general, unbounded. Moreover, as mentioned immediately after (14), the
feasibility of the optimizing condition is precisely what determines the phase
transition in the idealized scenario (where $\epsilon_{\sigma}=0$) and what is
fully characterized by Theorem 3.3.

In regimes below the phase transitions, the residual error norm, denoted as
$\|W\|_{F}$, can be effectively characterized. We provide a detailed
characterization of this in the subsequent section.

###  3.2 Operating in the regimes below phase transitions

As mentioned above, one obtains the worst case scenario for $V=U$. In such a
scenario one also has $V_{B}=V_{C}$, $\Lambda$ becomes symmetric and its
parametrization from (32) and (33) becomes

$\Lambda=V_{B}EV_{B}^{T}+V_{B}F(V_{B}^{\perp})^{T}+V_{B}^{\perp}D.$ (60)

Moreover, although it is not strictly needed for our considerations here, we
mention that from (117), (32) and (33) one can explicitly determine

$E=-\Sigma_{B}^{-1}V_{B}^{T}GV_{B}\Sigma_{B}^{-1},\quad\mbox{where}\quad
G=(I^{(l)})^{T}\bar{U}^{\perp}(\bar{U}^{\perp})^{T}I^{(l)}.$ (61)

Recalling (14) and utilizing Proposition 1 together with the above
parametrization of $\Lambda$ from (60), we further have

$\displaystyle f_{pr}(M;U,V)$ $\displaystyle=$
$\displaystyle\max_{\Lambda,\Lambda^{T}\Lambda\leq
I,(I^{(l)})^{T}(\bar{V}\bar{U}^{T}+\bar{V}^{\perp}\Lambda(\bar{U}^{\perp})^{T}I^{(l)}=0}{\rm
tr}\,(\Lambda\mbox{diag}(\epsilon_{\sigma}))$ $\displaystyle=$
$\displaystyle\max_{\Lambda=V_{B}EV_{B}^{T}+V_{B}F(V_{B}^{\perp})^{T}+V_{B}^{\perp}D,\Lambda\Lambda\leq
I}{\rm
tr}\,((V_{B}EV_{B}^{T}+V_{B}F(V_{B}^{\perp})^{T}+V_{B}^{\perp}D)\mbox{diag}(\epsilon_{\sigma})).$

We continue to consider the worst case and take
$\epsilon_{\sigma}=\max(\epsilon_{\sigma}){\bf 1}$ which implies
$\mbox{diag}(\epsilon_{\sigma})=\max(\epsilon_{\sigma})I$. Moreover, since it
will turn out that, in the worst case, the particular scaling plays no role,
to facilitate writing we will assume $\mbox{diag}(\epsilon_{\sigma})=I$. We
then from (LABEL:eq:supp3) have

$\displaystyle f_{pr}(M;U,V)$ $\displaystyle=$
$\displaystyle\max_{\Lambda=V_{B}EV_{B}^{T}+V_{B}^{\perp}D_{S}(V_{B}^{\perp})^{T},\Lambda\Lambda\leq
I}{\rm
tr}\,((V_{B}EV_{B}^{T}+V_{B}^{\perp}D_{S}(V_{B}^{\perp})^{T})\mbox{diag}(\epsilon_{\sigma})).$
(63) $\displaystyle=$ $\displaystyle{\rm tr}\,(V_{B}EV_{B}^{T})+n-k-l,$

where the optimum is achieved for $D_{S}=I$ and ultimately

$\tilde{\Lambda}=V_{B}EV_{B}^{T}+V_{B}^{\perp}D_{S}(V_{B}^{\perp})^{T}=V_{B}EV_{B}^{T}+V_{B}^{\perp}(V_{B}^{\perp})^{T}.$
(64)

We then have the following key result.

######  Theorem 3.4.

(Residual optimization characterization) Assume the setup of Theorem 9. One
then has the following

Worst case RMSE (deterministic) estimate:
$\|\bar{W}\|_{F}=\sqrt{-\min_{\bar{\Lambda},\bar{\Lambda}(I-(((I^{(l)})^{T}\bar{V}^{\perp})^{T}(I^{(l)})^{T}\bar{V}^{\perp}))=0}\left(\frac{1}{4}\|(I^{(l)})^{T}(\bar{V}^{\perp}\bar{\Lambda}(\bar{V}^{\perp})^{T}I^{(l)}\|_{F}^{2}-{\rm
tr}\,(\bar{\Lambda}\mbox{diag}(\epsilon_{\sigma}))\right)}.$ (65)

###### Proof.

We start by considering the following optimization

$\displaystyle\|\hat{W}\|_{F}^{2}=\min_{W}$ $\displaystyle\|W\|_{F}^{2}.$
subject to $\displaystyle{\rm
tr}\,(\bar{U}^{T}W\bar{V})+\ell_{1}^{*}(\mbox{diag}(\epsilon_{\sigma})+(\bar{U}^{\perp})^{T}W\bar{V}^{\perp})\leq
f_{pr}(M;U,V)={\rm tr}\,(\tilde{\Lambda})$ (66) $\displaystyle M\circ W=0.$

Following again the Lagrangian mechanism discussed above, we write

$\displaystyle\mathcal{L}(W,\lambda,\Theta)$ $\displaystyle=$
$\displaystyle\|W\|_{F}^{2}+\lambda{\rm
tr}\,(\bar{U}^{T}W\bar{V})+\lambda\ell_{1}^{*}(\mbox{diag}(\epsilon_{\sigma})+(\bar{U}^{\perp})^{T}W\bar{V}^{\perp})$
(67) $\displaystyle+\Theta(M\circ W)-\lambda{\rm tr}\,(\tilde{\Lambda})$
$\displaystyle=$
$\displaystyle\max_{\Lambda,\Lambda^{T}\Lambda\leq\lambda^{2}I}\Bigg{.}\Bigg{(}\|W\|_{F}^{2}+\lambda{\rm
tr}\,(\bar{U}^{T}W\bar{V})+{\rm
tr}\,(\Lambda(I+(\bar{U}^{\perp})^{T}W\bar{V}^{\perp}))$
$\displaystyle+\Theta(M\circ W)-\lambda{\rm
tr}\,(\tilde{\Lambda})\Bigg{.}\Bigg{)}$ $\displaystyle=$
$\displaystyle\max_{\Lambda,\Lambda^{T}\Lambda\leq\lambda^{2}I}\Bigg{(}\Bigg{.}\|W\|_{F}^{2}+{\rm
tr}\,\left((\lambda\bar{V}\bar{U}^{T}+\bar{V}^{\perp}\Lambda(\bar{U}^{\perp})^{T}+\Theta\circ
M)W\right)$ $\displaystyle+{\rm
tr}\,(\Lambda-\lambda\tilde{\Lambda})\Bigg{.}\Bigg{)}.$

We then also have for the Lagrange dual function

$\displaystyle g(\lambda,\Theta)$ $\displaystyle=$
$\displaystyle\min_{W}\mathcal{L}(W,\lambda,\Theta)$ $\displaystyle=$
$\displaystyle\min_{W}\max_{\Lambda,\Lambda^{T}\Lambda\leq\lambda^{2}I}\Bigg{(}\Bigg{.}\|W\|_{F}^{2}+{\rm
tr}\,\left((\lambda\bar{V}\bar{U}^{T}+\bar{V}^{\perp}\Lambda(\bar{U}^{\perp})^{T}+\Theta\circ
M)W\right)$ $\displaystyle+{\rm
tr}\,(\Lambda-\lambda\tilde{\Lambda})\Bigg{.}\Bigg{)}.$

Connecting (66)-(LABEL:eq:supp8) and utilizing the Lagrange duality gives

$\displaystyle\|\hat{W}\|_{F}^{2}$ $\displaystyle\geq$
$\displaystyle\max_{\Lambda,\Lambda^{T}\Lambda\leq\lambda^{2}I}\min_{W}\Bigg{(}\Bigg{.}\|W\|_{F}^{2}+{\rm
tr}\,\left((\lambda\bar{V}\bar{U}^{T}+\bar{V}^{\perp}\Lambda(\bar{U}^{\perp})^{T}+\Theta\circ
M)W\right)$ $\displaystyle+{\rm
tr}\,(\Lambda-\lambda\tilde{\Lambda})\Bigg{.}\Bigg{)}$ $\displaystyle=$
$\displaystyle\max_{\Lambda,\Lambda^{T}\Lambda\leq\lambda^{2}I}\left(-\frac{1}{4}\|(I^{(l)})^{T}(\lambda\bar{V}\bar{U}^{T}+\bar{V}^{\perp}\Lambda(\bar{U}^{\perp})^{T})I^{(l)}\|_{F}^{2}+{\rm
tr}\,(\Lambda-\lambda\tilde{\Lambda})\right)$ $\displaystyle=$
$\displaystyle-\min_{\Lambda,\Lambda^{T}\Lambda\leq\lambda^{2}I}\left(\frac{1}{4}\|(I^{(l)})^{T}(\lambda\bar{V}\bar{U}^{T}+\bar{V}^{\perp}\Lambda(\bar{U}^{\perp})^{T})I^{(l)}\|_{F}^{2}-{\rm
tr}\,(\Lambda-\lambda\tilde{\Lambda})\right).$ (69)

After setting

$\displaystyle\bar{\Lambda}=\Lambda-\lambda\tilde{\Lambda},$ (70)

one can rewrite (69) as

$\displaystyle\|\hat{W}\|_{F}^{2}$ $\displaystyle=$
$\displaystyle-\min_{\Lambda,\Lambda^{T}\Lambda\leq\lambda^{2}I}\left(\frac{1}{4}\|(I^{(l)})^{T}(\lambda\bar{V}\bar{U}^{T}+\bar{V}^{\perp}\Lambda(\bar{U}^{\perp})^{T})I^{(l)}\|_{F}^{2}-{\rm
tr}\,(\Lambda-\lambda\tilde{\Lambda})\right)$ $\displaystyle=$
$\displaystyle-\min_{\Lambda,(\bar{\Lambda}+\lambda\tilde{\Lambda})^{T}(\bar{\Lambda}+\lambda\tilde{\Lambda})\leq\lambda^{2}I}\left(\frac{1}{4}\|(I^{(l)})^{T}(\lambda\bar{V}\bar{U}^{T}+\bar{V}^{\perp}(\bar{\Lambda}+\lambda\Lambda)(\bar{U}^{\perp})^{T})I^{(l)}\|_{F}^{2}-{\rm
tr}\,(\bar{\Lambda})\right)$ $\displaystyle=$
$\displaystyle-\min_{\Lambda,(\bar{\Lambda}+\lambda\tilde{\Lambda})^{T}(\bar{\Lambda}+\lambda\tilde{\Lambda})\leq\lambda^{2}I}\left(\frac{1}{4}\|(I^{(l)})^{T}\bar{V}^{\perp}\bar{\Lambda}(\bar{U}^{\perp})^{T}I^{(l)}\|_{F}^{2}-{\rm
tr}\,(\bar{\Lambda})\right).$

As mentioned above, one has in the worst case $V=U$ and $\bar{V}=\bar{U}$ and
the following parametrization of $\bar{\Lambda}$

$\displaystyle\bar{\Lambda}=\begin{bmatrix}V_{B}^{\perp}&V_{B}\end{bmatrix}\begin{bmatrix}Z_{1}&Z_{2}\\\
Z_{2}&Z_{3}\end{bmatrix}\begin{bmatrix}V_{B}^{\perp}&V_{B}\end{bmatrix}^{T}.$
(72)

After recalling on parametrization of $\tilde{\Lambda}$ from (64) we find

$\displaystyle\tilde{\Lambda}=\begin{bmatrix}V_{B}^{\perp}&V_{B}\end{bmatrix}\begin{bmatrix}I&0\\\
0&E\end{bmatrix}\begin{bmatrix}V_{B}^{\perp}&V_{B}\end{bmatrix}^{T},$ (73)

and

$\displaystyle\bar{\Lambda}+\lambda\tilde{\Lambda}=\begin{bmatrix}V_{B}^{\perp}&V_{B}\end{bmatrix}\begin{bmatrix}\lambda
I+Z_{1}&Z_{2}\\\ Z_{2}&\lambda
E+Z_{3}\end{bmatrix}\begin{bmatrix}V_{B}^{\perp}&V_{B}\end{bmatrix}^{T}.$ (74)

Moreover, we then also have

$\displaystyle\lambda^{2}I-(\bar{\Lambda}+\tilde{\Lambda})^{2}=\begin{bmatrix}V_{B}^{\perp}&V_{B}\end{bmatrix}\begin{bmatrix}\lambda^{2}I-(\lambda
I+Z_{1})^{2}-Z_{2}^{T}Z_{2}&Z_{4}\\\
Z_{4}&Z_{5}\end{bmatrix}\begin{bmatrix}V_{B}^{\perp}&V_{B}\end{bmatrix}^{T}.$
(75)

For $\lambda^{2}I-(\bar{\Lambda}+\tilde{\Lambda})^{2}$ to be positive semi-
definite one would need $\lambda^{2}I-(\lambda I+Z_{1})^{2}-Z_{2}^{T}Z_{2}\geq
0$ which means that $Z_{1}\leq 0$. Recalling from (32),
$(I^{(l)})^{T}V^{\perp}=U_{B}\Sigma_{B}V_{B}^{T}$, one then from
(LABEL:eq:supp11) has that $\bar{\Lambda}$ does not belong to the span of
$V_{B}^{\perp}$, i.e. it does not belong to the span orthogonal to
$(I^{(l)})^{T}V^{\perp}$. Together with (LABEL:eq:supp11) and the strong
duality (which, due to the underlying convexity, ensures that the inequality
in (69) can be replaced by an equality) this gives the equation in the
statement of Theorem 3.4 and ultimately completes the prof of the theorem. ∎

A combination of Theorem 3.4 and Lemma 3.2 gives the following statistical
RMSE characterization.

######  Theorem 3.5.

(Worst case RMSE) Assume the setup of Theorem 9 with $U=V$ Haar distributed
and a large $n$ linear regime with
$\beta\triangleq\lim_{n\rightarrow\infty}\frac{k}{n}$ and
$\eta\triangleq\lim_{n\rightarrow\infty}\frac{l}{n}$. Assume that $\hat{X}$ is
the solution of (3). One then has the following

Worst case RMSE (statistical) estimate: $\displaystyle x_{l/u}$
$\displaystyle\triangleq$
$\displaystyle\beta+\eta-2\beta\eta\pm\sqrt{(\beta+\eta-2\beta\eta)^{2}-(\beta-\eta)^{2}},$
$\displaystyle\xi$ $\displaystyle\triangleq$
$\displaystyle\sigma_{\epsilon}\frac{\sqrt{\int_{x_{l}}^{x_{u}}\frac{\sqrt{-(x-(\beta+\eta))^{2}-4\beta\eta(x-1)}}{2\pi
x^{2}(x-x^{2})}dx+\max(1-(\beta+\eta),0)}}{\sqrt{1-\beta}(1-\eta)},$
$\displaystyle\forall X_{sol}$
$\displaystyle\lim_{n\rightarrow\infty}{\mathbb{P}}\left(\textbf{RMSE}=\frac{n\|\hat{X}-X_{sol}\|_{F}}{\sqrt{n-k}(n-l)}\leq\xi\right)\longrightarrow
1,$ $\displaystyle\exists X_{sol}$
$\displaystyle\lim_{n\rightarrow\infty}{\mathbb{P}}\left(\textbf{RMSE}=\frac{n\|\hat{X}-X_{sol}\|_{F}}{\sqrt{n-k}(n-l)}\in((1-\epsilon)\xi,(1+\epsilon)\xi)\right)\longrightarrow
1.\hskip 14.45377pt$

###### Proof.

From considerations presented in (LABEL:eq:supp11)-(75) we have

$\|\hat{W}\|_{F}^{2}=-\min_{Z_{3},\bar{\Lambda}=V_{B}Z_{3}V_{B}^{T},\lambda,(\bar{\Lambda}+\lambda\tilde{\Lambda})^{T}(\bar{\Lambda}+\lambda\tilde{\Lambda})\leq\lambda^{2}I}\left(\frac{1}{4}\|(I^{(l)})^{T}\bar{U}^{\perp}\bar{\Lambda}(\bar{U}^{\perp})^{T}I^{(l)}\|_{F}^{2}-{\rm
tr}\,(\bar{\Lambda})\right).$ (77)

Substituting $(I^{(l)})^{T}\bar{U}^{\perp}=U_{B}\Sigma_{B}B_{B}^{T}$ we obtain

$\|\hat{W}\|_{F}^{2}=-\min_{Z_{3},\lambda,(Z_{3}+\lambda E)^{T}(Z_{3}+\lambda
E)\leq\lambda^{2}I}\left(\frac{1}{4}{\rm
tr}\,(\Sigma_{B}Z_{3}\Sigma_{B}^{2}Z_{3}\Sigma_{B})-{\rm tr}\,(Z_{3})\right).$
(78)

Elementary algebraic considerations give that optimal $Z_{3}$ is diagonal and
since $\lambda$ can be arbitrarily large we have from (78)

$\|\hat{W}\|_{F}^{2}=-\min_{\bar{z}}\left(\frac{1}{4}\sum_{i=1}^{n-l}\bar{z}_{i}^{2}\bar{\sigma}_{i}^{4}-\sum_{i=1}^{n-l}\bar{z}_{i}\right)=\sum_{i=1}^{n-l}\frac{1}{\bar{\sigma}_{i}^{4}},\quad\mbox{where}\quad\bar{\sigma}=\mbox{diag}(\Sigma_{B}).$
(79)

Keeping in mind that the distribution of $\bar{\sigma}_{i}^{2}$ is as given in
(199) and (200), one finally has

$\displaystyle{\mathbb{E}}\|\hat{W}\|_{F}^{2}$ $\displaystyle=$
$\displaystyle{\mathbb{E}}\sum_{i=1}^{n-l}\frac{1}{\bar{\sigma}_{i}^{4}}$
$\displaystyle=$
$\displaystyle\left(\int_{x_{l}}^{x_{u}}\frac{\sqrt{-(x-(\beta+\eta))^{2}-4\beta\eta(x-1)}}{2\pi
x^{2}(x-x^{2})}+\max(1-(\beta+\eta),0)\delta(x-1)\right).$

Concentrations, scalings by $\sqrt{n-k}(n-l)/n$, and adjusting the maximum
magnitude of $\epsilon_{\sigma}$, $\max(\epsilon_{\sigma})$, to be any
quantity $\sigma_{\epsilon}$ (instead of 1) completes the proof of the
theorem. ∎

###  3.3 A summarized interpretation of the methodology and obtained
theoretical results

We find it useful to provide a more concrete interpretation of the results
given in Theorem 3.5 within the causal inference context (see, e.g., [5] for a
pioneering work in establishing matrix completion $\longleftrightarrow$ block
causal inference (MC$\longleftrightarrow$ block C-inf) connection). We first
recall that within the MC $\longleftrightarrow$ block C-inf context, the rows
and columns of $X_{sol}$ correspond to the units and time periods (at which
the units can potentially be exposed to treatments), respectively. Assuming
the so-called irreversible simultaneous treatments during $n-l$ consecutive
time periods on a subset of $n-l$ units, the block missing part of $X_{sol}$
corresponds to the so-called _counterfactuals_ (values of the units had the
treatment not been applied). The approximate recovery of the counterfactual
part of $X_{sol}$ is then the ultimate block C-inf’s goal. The optimization
problem in (2) is posed so that its optimal solution hopefully closely matches
$X_{sol}$. However, solving (2) in an idealized scenario is NP hard and in an
non-idealized one even uniquely impossible. As is a common practice in matrix
completion literature, we consider $\ell_{1}^{*}$ heuristic as the tightest
convex relaxation of $\ell_{0}^{*}$ (which, in the idealized scenario, is
known to often be successful) hoping that its solution, $\hat{X}$, also
closely matches $X_{sol}$. Under statistical assumptions on the origin of
$X_{sol}$, Theorem 3.5 effectively provides the exact relation among three key
system and algorithmic parameters: (_i_) data heterogeneity
($\text{rank}(X_{sol})$ or in a scaled form
$\beta=\lim_{n\rightarrow\infty}\frac{\text{rank}(X_{sol})}{n}=\lim_{n\rightarrow\infty}\frac{k}{n}$);
(_ii_) the number of the untreated units and the treatment starting time ($l$
or in a scaled form $\eta=\lim_{n\rightarrow\infty}\frac{l}{n}$); (_iii_) the
absolute spectral deviation from the ideal low rankness ($\sigma_{\epsilon}$);
and (_iv_) the residual RMSE error, $\|\hat{X}-X_{sol}\|_{F}$.

As particularly relevant, we also emphasize the following two observed
features of the obtained theoretical results: (_i_) (_Simplicity_) Both the
phase transition characterization from Theorem 3.3 and the corresponding
residual RMSE one from Theorem 3.5 are remarkably simple, single line
functional equations where the key system parameters
$(\beta,\eta,\sigma_{\epsilon})$ explicitly appear. (_ii_) (_Generality_) The
characterization from Theorem 3.5 exhibits a particularly generic behavior.
Namely, the residual error is characterized in the worst case sense, ensuring
that no matter what the relation between $\sigma$ and $\sigma_{\epsilon}$ is,
the error remains bounded by the quantity given in the theorem. In other
words, no matter what the level of relative deviation from the ideal low
rankness is, the error is always bounded by the quantity given in the theorem.
Moreover, for $\sigma\gg\sigma_{\epsilon}$, the bound is achievable, which
implies that better (more generic) results cannot be proven.

##  4 Numerical results

To complement the above theoretical findings, we conduct a set of numerical
experiments in several different relevant regimes. Our focus is on comparing
our theoretical $\ell_{1}^{*}$ performance analysis results with the
corresponding $\ell_{1}^{*}$ numerical estimates.

_1) Worst case regime:_ As the theory suggests, somewhat counter-intuitively,
the worst case RMSE estimates are achieved for matrices with a large gap
between dominant singular values and the negligible ones. To emulate such a
regime we set $\sigma_{\epsilon}=1$ and magnitude of each of the large
$\sigma$’s components to be 50. We find the ratio 50:1 as sufficiently large
to emulate the theoretically needed infinite one. The obtained results are
shown in Table 1. It is easily observed that there is a close agreement
between theoretical and simulated results across the range of allowed $\beta$.
One additionally notes that for $\beta=0.2$ a slight deviation from what the
theory predicts can be seen. The reason, of course, is the fact that at
$\beta=0.2$, one is exactly on the phase transition and above such a $\beta$
the residual error can even be infinite.

Table 1: Block missingness matrix completion via nuclear norm minimization – _worst case_ RMSE; $\eta=0.9$ Method | $\beta$
---|---
(to estimate RMSE) | $\mathbf{0.05}$ | $\mathbf{0.1}$ | $\mathbf{0.15}$ | $\mathbf{0.2}$
theoretical $\frac{n{\mathbb{E}}\|\hat{X}-X_{sol}\|_{F}}{\sigma_{\epsilon}\sqrt{n-k}(n-l)}$ | $\mathbf{3.446}$ | $\mathbf{3.755}$ | $\mathbf{4.161}$ | $\mathbf{4.617}$
simulated $\frac{n{\mathbb{E}}\|\hat{X}-X_{sol}\|_{F}}{\sigma_{\epsilon}\sqrt{n-k}(n-l)}$ | $\mathbf{3.46}$ | $\mathbf{3.82}$ | $\mathbf{4.24}$ | $\mathbf{5.00}$
asymmetric: simulated $\frac{n{\mathbb{E}}\|\hat{X}-X_{sol}\|_{F}}{\sigma_{\epsilon}\sqrt{n-k}(n-l)}$ | $\mathbf{2.33}$ | $\mathbf{2.66}$ | $\mathbf{3.18}$ | $\mathbf{3.73}$

_2) Asymmetric case:_ Differently from $U=V$ which achieves the worst case but
is highly unlikely to be observed in real data, we consider its an asymmetric
anti-pod, where not only aren’t $U$ and $V$ equal to each other, they are
actually completely independent of each other. This scenario is expected to
more faithfully capture real data. In the absence of any other _a priori_
available knowledge about the underlying parts of the model, such a scenario
might in fact be the most natural one to assume. We show the results for this
scenario with identical system parameters as in the worst case. As the results
in Table 1 show, the worst case results are indeed upper bounds on the ones
obtained in the asymmetric scenario.

_3) Varying $\epsilon_{\sigma}$:_ For the worst case scenario, we set the
magnitude of $\epsilon_{\sigma}$ equal to $\sigma_{\epsilon}=1$. We also
simulated the scenario where the magnitudes are standard normals and uniformly
random while we maintained the overall norm of $\epsilon_{\sigma}=\sqrt{n-k}$
to ensure fairness in comparison with the worst case scenario. The results are
shown in Table 2. As it can be seen from the table, the results are dominated
by the worst case ones as expected by the theoretical predictions. One should
though observe that the results obtained under the Gaussian assumption are
actually fairly close to the flat ones (where all the magnitudes of
$\epsilon_{\sigma}$ are equal to $\sigma_{\epsilon}=1$).

Table 2: Block missingness matrix completion via nuclear norm minimization – _worst case_ RMSE; $\eta=0.9$; different statistics of $\epsilon_{\sigma}$ Method | $\beta$
---|---
(to estimate RMSE) | $\mathbf{0.05}$ | $\mathbf{0.1}$ | $\mathbf{0.15}$ | $\mathbf{0.2}$
theoretical $\frac{n{\mathbb{E}}\|\hat{X}-X_{sol}\|_{F}}{\sigma_{\epsilon}\sqrt{n-k}(n-l)}$ | $\mathbf{3.446}$ | $\mathbf{3.755}$ | $\mathbf{4.161}$ | $\mathbf{4.617}$
Flat $\epsilon_{\sigma}$: simulated $\frac{n{\mathbb{E}}\|\hat{X}-X_{sol}\|_{F}}{\sigma_{\epsilon}\sqrt{n-k}(n-l)}$ | $\mathbf{3.46}$ | $\mathbf{3.82}$ | $\mathbf{4.24}$ | $\mathbf{5.00}$
Gaussian $\epsilon_{\sigma}$: simulated $\frac{n{\mathbb{E}}\|\hat{X}-X_{sol}\|_{F}}{\sigma_{\epsilon}\sqrt{n-k}(n-l)}$ | $\mathbf{3.40}$ | $\mathbf{3.54}$ | $\mathbf{3.92}$ | $\mathbf{5.00}$
Uniform $\epsilon_{\sigma}$: simulated $\frac{n{\mathbb{E}}\|\hat{X}-X_{sol}\|_{F}}{\sigma_{\epsilon}\sqrt{n-k}(n-l)}$ | $\mathbf{3.15}$ | $\mathbf{3.44}$ | $\mathbf{3.83}$ | $\mathbf{5.03}$

_4) Varying magnitudes of $\sigma$:_ As discussed above, to achieve the worst
case behavior, one needs to choose a rather large magnitude of dominating
singular values $\sigma$. A natural question is then, what happens when the
magnitudes of $\sigma$, $\sigma_{mag}$, are not necessarily as large. We, in
Table 3, show the obtained results as magnitudes of $\sigma$ change (to ensure
a fair comparison, we continue to keep, as earlier, $\epsilon_{\sigma}$ equal
to $\sigma_{\epsilon}=1$). As it can be seen from the table, the results are
converging in a rather fast manner. The ratios
$\frac{\sigma_{mag}}{\sigma_{\epsilon}}=\frac{50}{1}$ that we, out of
precaution, chose above is not even needed. The worst case behavior is fairly
closely approached, already for ratios around $6:1$.

Table 3: Block missingness matrix completion via nuclear norm minimization; $\eta=0.9$; RMSE as a function of $\frac{\sigma_{mag}}{\sigma_{\epsilon}}$ Method | $\frac{\sigma_{mag}}{\sigma_{\epsilon}}$ |  $\frac{\sigma_{mag}}{\sigma_{\epsilon}}$ (theory)
---|---|---
(to estimate RMSE) | $\mathbf{1}$ | $\mathbf{2}$ | $\mathbf{3}$ | $\mathbf{4}$ | $\mathbf{6}$ | $\mathbf{\infty}$
Flat $\epsilon_{\sigma}$: simulated $\frac{n{\mathbb{E}}\|\hat{X}-X_{sol}\|_{F}}{\sigma_{\epsilon}\sqrt{n-k}(n-l)}$ | $\mathbf{3.43}$ | $\mathbf{3.71}$ | $\mathbf{3.87}$ | $\mathbf{3.93}$ | $\mathbf{4.01}$ | $\mathbf{4.16}$

We choose $n=40$, as a rather small matrix dimension. Despite the fact that
the proofs are obtained assuming large dimensional settings, we observe a
strong agreement between theoretical and simulated results even for fairly
small sizes. Also, to show an agreement with the theoretical results, we,
consider the so-called _typical_ behavior with the singular values of the
unknown targeted matrices equal to one. Choices of this type are known to
serve as the worst case examples in establishing the reversal
$\ell_{0}^{*}-\ell_{1}^{*}$-equivalence conditions in idealized scenarios (for
completeness, we also ran simulations where the singular values were randomly
chosen and the obtained results were either identically matching or improving
upon the ones shown in Table 1).

##  5 Conclusion

In this paper, we studied the nuclear norm minimization heuristic and
theoretically analyzed its performance for completing matrices with _missing
not at random_ (MNAR) entries. We place special emphasis on _approximately low
rank_ matrices and on the _block missing patterns_ encountered in the so-
called irreversible simultaneous treatment scenarios often seen in various
health and social science applications. In such scenarios, an exact recovery
of the incomplete matrices is not attainable and the residual errors are
present as well. To assess the quality of the underlying nuclear norm based
algorithms, we analytically characterized the residual errors.

Relying on a typical statistical view of the problem and utilizing random
matrix theory related to free probability theory (FPT), we first introduced a
generic performance analysis framework. We then conducted a concrete,
statistical phase transition type of analysis, and obtained _explicit exact
worst case_ estimates for the residual RMSE which hold below the so-called
(worst case) phase-transition line. The resulting characterizations are
simple, single line functional equations where all three key system parameters
— (_i_) the data heterogeneity ($\beta$), (_ii_) the size of the missing
blocks ($\eta$), and (_iii_) the magnitude of the spectral deviations from the
ideal low rankness ($\sigma_{\epsilon}$)— are clearly visible.

To complement our theoretical findings, we conducted several numerical
experiments as well. The following are the main takeaways: _(i)_ A very strong
agreement between the theoretical predictions and numerical simulations is
observed. _(ii)_ Despite the fact that our theoretical predictions are
obtained in a large dimensional asymptotic regime, we observe a strong match
between theory and simulations, even for rather small problems with dimensions
on the order of few tens.

##### Appendix

This appendix contains proofs of technical results stated in the main body of
the paper. We also demonstrate that our analysis can handle data datasets
encountered in real applications, such as the California tobacco data, widely
used in many studies (e.g. [2]). Our analysis, particular examining the
spectrum of this data set, reveals that its approximately low rank structure
is directly attributable to the data’s inherent structure.

##  Appendix A Proof of Theorem 9

The theorem follows from the proofs of much stronger fundamental propositions
and corollaries that relate to a more generic concept of the low rank recovery
(LRR). In such a concept, the fairly particular linear constraints of (2) and
(3) are replaced by a more generic ${\bf y}=\mbox{vec}(Y)=A\mbox{vec}(X)$,
where $A\in{\mathbb{R}}^{m\times n^{2}}$ is a generic measurement matrix.

######  Proposition 2.

($\ell_{0}^{*}-\ell_{1}^{*}$-equivalence condition (LRR) – symmetric $X$)
Consider a $\bar{U}\in{\mathbb{R}}^{n\times k}$ such that
$\bar{U}^{T}\bar{U}=I_{k\times k}$ and a $\mbox{rank}-k$ a priori known to be
symmetric matrix $X_{sol}=X\in{\mathbb{R}}^{n\times n}$ with all of its
columns belonging to the span of $\bar{U}$. For concreteness, and without loss
of generality, assume that $X$ has only positive nonzero eigenvalues. For a
given matrix $A\in{\mathbb{R}}^{m\times n^{2}}$ ($m\leq n^{2}$) assume that
${\bf y}=A\mbox{vec}(X)\in{\mathbb{R}}^{m}$. If

$(\forall W\in{\mathbb{R}}^{n\times n}|A\mbox{vec}(W)={\bf 0}_{m\times
1},W=W^{T}\neq{\bf 0}_{n\times n})\quad-{\rm
tr}\,(\bar{U}^{T}W\bar{U})<\ell_{1}^{*}((\bar{U}^{\perp})^{T}W\bar{U}^{\perp}),$
(81)

then the solutions of (2) and (3) coincide. Moreover, if

$(\exists W\in{\mathbb{R}}^{n\times n}|A\mbox{vec}(W)={\bf 0}_{m\times
1},W=W^{T}\neq{\bf 0}_{n\times n})\quad-{\rm
tr}\,(\bar{U}^{T}W\bar{U})\geq\ell_{1}^{*}((\bar{U}^{\perp})^{T}W\bar{U}^{\perp}),$
(82)

then there is an $X$ from the above set of the symmetric matrices with columns
belonging to the span of $\bar{U}$ such that the solutions of (2)) and (3) are
different.

###### Proof.

The proof can be obtained by a step-by-step adaptation of Theorem 2 in [45] to
matrices or their singular/eigenvalues. For concreteness and without loss of
generality we also assume that the eigen-decomposition of $X$ is

$\displaystyle X=U\Lambda
U^{T}=\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix}\begin{bmatrix}\bar{\Lambda}_{X}&{\bf
0}_{k\times(n-k)}\\\ {\bf 0}_{(n-k)\times
k}&\bar{\Lambda}_{X}^{\perp}\end{bmatrix}\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix}^{T}.$
(83)

(i) $\Longrightarrow$ (the if part): Following step-by-step the proof of
Theorem $2$ in [45], we start by assuming that $\hat{X}$ is the solution of
(3). Then we want to show that if (81) holds then $\hat{X}=X$. As usual, we
instead of that, assume opposite, i.e. we assume that (81) holds but
$\hat{X}\neq X$. Then since ${\bf y}=A\mbox{vec}(\hat{X})$ and ${\bf
y}=A\mbox{vec}(X)$ must hold simultaneously there must exist $W$ such that
$\hat{X}=X+W$ with $W\neq 0$, $A\mbox{vec}(W)=0$. Moreover, since $\hat{X}$ is
the solution of (3) one must also have

$\displaystyle\begin{array}[]{r r r l@{\
}}&\ell_{1}^{*}(X+W)=\ell_{1}^{*}(\hat{X})&\leq&\ell_{1}^{*}(X)\\\
\Longleftrightarrow\hskip 21.68121pt$
$&\ell_{1}^{*}(\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix}^{T}(X+W)\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix})&\leq&\ell_{1}^{*}(X)\\\
\Longrightarrow\hskip 21.68121pt$
$&\ell_{1}^{*}(\bar{U}^{T}(X+W)\bar{U})+\ell_{1}^{*}((\bar{U}^{\perp})^{T}(X+W)\bar{U}^{\perp})&\leq&\ell_{1}^{*}(X).\end{array}$
(87) (88)

The last implication follows after one trivially notes

$\displaystyle\ell_{1}^{*}(\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix}^{T}(X+W)\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix})$
$\displaystyle=$ $\displaystyle\max_{\Lambda_{*}=\Lambda_{*}^{T}\in{\cal
L}_{*}}{\rm
tr}\,(\Lambda_{*}\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix}^{T}(X+W)\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix})$
(89) $\displaystyle\geq$
$\displaystyle\max_{\Lambda_{*}=\Lambda_{*}^{T}\in{\cal L}_{*}^{0}}{\rm
tr}\,(\Lambda_{*}\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix}^{T}(X+W)\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix})$
$\displaystyle=$
$\displaystyle\ell_{1}^{*}(\bar{U}^{T}(X+W)\bar{U})+\ell_{1}^{*}((\bar{U}^{\perp})^{T}(X+W)\bar{U}^{\perp}),$

where

$\displaystyle{\cal L}_{*}^{0}$ $\displaystyle\triangleq$
$\displaystyle\left\\{\Lambda_{*}\in{\mathbb{R}}^{n\times
n}|\Lambda_{*}=\Lambda_{*}^{T},\Lambda_{*}\Lambda_{*}^{T}\leq
I,\Lambda_{*}=\begin{bmatrix}\Lambda_{*,1}&0_{k\times(n-k)}\\\ 0_{(n-k)\times
k}&\Lambda_{*,2}\end{bmatrix}\right\\}$ (90) $\displaystyle\subseteq$
$\displaystyle\left\\{\Lambda_{*}\in{\mathbb{R}}^{n\times
n}|\Lambda_{*}=\Lambda_{*}^{T},\Lambda_{*}\Lambda_{*}^{T}\leq
I\right\\}\triangleq{\cal L}_{*}.$

Next, the key observation is to note that the absolute values can be removed
in the nonzero part and that the $\ell_{1}^{*}(\cdot)$ can be “ _replaced_ ”
by ${\rm tr}\,(\cdot)$. Such a simple observation is the most fundamental
reason for all the success of the RDT when used for the exact performance
characterization of the structured objects’ recovery. From (88) we then have

$\displaystyle\begin{array}[]{r r r l@{\
}}&\ell_{1}^{*}(\bar{U}^{T}(X+W)\bar{U})+\ell_{1}^{*}((\bar{U}^{\perp})^{T}(X+W)\bar{U}^{\perp})&\leq&\ell_{1}^{*}(X)\\\
\Longrightarrow\hskip 21.68121pt$ $&{\rm
tr}\,(\bar{U}^{T}(X+W)\bar{U})+\ell_{1}^{*}((\bar{U}^{\perp})^{T}(W)\bar{U}^{\perp})&\leq&\ell_{1}^{*}(X)\\\
\Longleftrightarrow\hskip 21.68121pt$ $&{\rm
tr}\,(\bar{U}^{T}W\bar{U})+\ell_{1}^{*}((\bar{U}^{\perp})^{T}W\bar{U}^{\perp})&\leq&0.\end{array}$
(94)

We have arrived at a contradiction as the last inequality in (94) is exactly
the opposite of (81). This implies that our initial assumption $\hat{X}\neq X$
cannot hold and we therefore must have $\hat{X}=X$. This is precisely the
claim of the first part of the theorem.

(ii) $\Longleftarrow$ (the only if part): We now assume that (82) holds, i.e.

$(\exists W\in{\mathbb{R}}^{n\times n}|A\mbox{vec}(W)={\bf 0}_{m\times
1},W\neq{\bf 0}_{n\times n})\quad-{\rm
tr}\,((\bar{U})^{T}W\bar{U})\geq\ell_{1}^{*}((\bar{U}^{\perp})^{T}W\bar{U}^{\perp})$
(95)

and would like to show that for such a $W$ there is a symmetric rank-$k$
matrix $X$ with the columns belonging to the span of $\bar{U}$ such that ${\bf
y}=A\mbox{vec}(X)$, and the following holds

$\ell_{1}^{*}(X+W)<\ell_{1}^{*}(X).$ (96)

Existence of such an $X$ would ensure that it both, satisfies all the
constraints in (3) and is not the solution of (3) . One can reverse all the
above steps from (95) to (88) with strict inequalities and arrive at the first
inequality in (88) which is exactly (96). There are two implications that
cause problems in such a reversal process, the one in (95) and the one in
(88). If these implications were equivalences everything would be fine. We
address these two implications separately.

1) the implication in (94) – particular $X$ to “overwhelm” $W$: Assume
$X=\bar{U}\Lambda_{x}\bar{U}^{T}$ with $\Lambda_{x}>0$ being a diagonal matrix
with arbitrarily large elements on the main diagonal (here it is sufficient
even to choose diagonal of $\Lambda_{x}$ so that its smallest element is
larger than the maximum eigenvalue of $\bar{U}^{T}W\bar{U}$). Now one, of
course, sees the main idea behind the “removing the absolute values” concept
from [45]. Namely, for such an $X$ one has that
$\ell_{1}^{*}(\bar{U}^{T}X+W)\bar{U})=tr(\ell_{1}^{*}(\bar{U}^{T}X+W)\bar{U}))$
since for symmetric matrices the $\ell_{1}^{*}(\cdot)$ (as the sum of the
argument’s _absolute_ eigenvalues) and ${\rm tr}\,(\cdot)$ (as the sum of the
argument’s eigenvalues) are equal. That basically means that when going
backwards the second inequality in (94) not only follows from the first one
but also implies it as well. In other words, for
$X=\bar{U}\Lambda_{x}\bar{U}^{T}$ (with $\Lambda_{x}>0$ and arbitrarily large)

$\displaystyle\begin{array}[]{r r r l@{\ }}&{\rm
tr}\,(\bar{U}^{T}W\bar{U})+\ell_{1}^{*}((\bar{U}^{\perp})^{T}W\bar{U}^{\perp})&\leq&0\\\
\Longleftrightarrow\hskip 21.68121pt$ $&{\rm
tr}\,(\bar{U}^{T}(X+W)\bar{U}^{)}+\ell_{1}^{*}((\bar{U}^{\perp})^{T}(W)\bar{U}^{\perp})&\leq&\ell_{1}^{*}(X)\\\
\Longleftrightarrow\hskip 21.68121pt$
$&\ell_{1}^{*}(\bar{U}^{T}(X+W)\bar{U})+\ell_{1}^{*}((\bar{U}^{\perp})^{T}(X+W)\bar{U}^{\perp})&\leq&\ell_{1}^{*}(X),\end{array}$
(100)

which basically mans that there is an $X$ that can “overwhelm” $W$ (in the
span of $\bar{U}$) and ensures that the “ _removing the absolute values_ ” is
not only a _sufficient_ but also a _necessary_ concept for creating the
relaxation equivalence condition.

2) the implication in (88): One would now need to somehow show that the third
inequality in (88) not only follows from the second one but also implies it as
well. This boils down to showing that inequality in (89) can be replaced with
an equality or, alternatively, that ${\cal L}^{0}$ and ${\cal L}$ are
provisionally equivalent. Neither of these statements is generically true.
However, since we have a set of $X$ at our disposal there might be an $X$ for
which they actually hold. We continue to assume
$X=\bar{U}\Lambda_{x}\bar{U}^{T}$ with $\Lambda_{x}>0$ being a diagonal matrix
with arbitrarily large entries on the main diagonal. Then the last equality in
(89) gives

$\displaystyle\begin{array}[]{r r r l@{\ }}$
$&\ell_{1}^{*}(\bar{U}^{T}(X+W)\bar{U})+\ell_{1}^{*}((\bar{U}^{\perp})^{T}(X+W)\bar{U}^{\perp})&\leq&\ell_{1}^{*}(X)\\\
\Longleftrightarrow\hskip 21.68121pt$
$&\max_{\Lambda_{*}=\Lambda_{*}^{T}\in{\cal L}_{*}^{0}}{\rm
tr}\,(\Lambda_{*}\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix}^{T}(X+W)\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix})&\leq&\ell_{1}^{*}(X).\end{array}$
(103)

Also, one has

$\displaystyle\begin{array}[]{r r r l@{\ }}&\hskip
21.68121pt\max_{\Lambda_{*}=\Lambda_{*}^{T}\in{\cal L}_{*}^{0}}{\rm
tr}\,(\Lambda_{*}\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix}^{T}(X+W)\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix})&\leq&\ell_{1}^{*}(X)\\\
\Longleftrightarrow$
$&\max_{\Lambda_{*,i}=\Lambda_{*,i}^{T},\Lambda_{*,i}\Lambda_{*,i}^{T}\leq
I,i\in\\{1,2\\}}{\rm
tr}\,(\Lambda_{*,1}\bar{U}^{T}X\bar{U}+\Lambda_{*,2}(\bar{U}^{\perp})^{T}W\bar{U}^{\perp})&\leq&\ell_{1}^{*}(X)\\\
\Longleftrightarrow$
$&\max_{\Lambda_{*,i}=\Lambda_{*,i}^{T},\Lambda_{*,i}\Lambda_{*,i}^{T}\leq
I,i\in\\{1,2\\}}{\rm
tr}\,(\Lambda_{*,1}\Lambda_{x}+\Lambda_{*,2}(\bar{U}^{\perp})^{T}W\bar{U}^{\perp})&\leq&{\rm
tr}\,(\Lambda_{x}).\end{array}$ (107) (108)

Now, if at least one of the elements on the main diagonal of $\Lambda_{*,1}$,
$\mbox{diag}(\Lambda_{*,1})$, is smaller than 1, then the corresponding
element on the diagonal of $\Lambda_{x}$ can be made arbitrarily large
compared to the other elements of $\Lambda_{x}$ and one would have

$\displaystyle\begin{array}[]{r r r l@{\
}}&\max_{\Lambda_{*,i}=\Lambda_{*,i}^{T},\Lambda_{*,i}\Lambda_{*,i}^{T}\leq
I,i\in\\{1,2\\}}{\rm
tr}\,(\Lambda_{*,1}\Lambda_{x}+\Lambda_{*,2}(\bar{U}^{\perp})^{T}W\bar{U}^{\perp})&<&{\rm
tr}\,(\Lambda_{x})\\\ \Longleftrightarrow\hskip 21.68121pt$
$&\max_{\Lambda_{*}=\Lambda_{*}^{T}\in{\cal L}_{*}^{0}}{\rm
tr}\,(\Lambda_{*}\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix}^{T}(X+W)\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix})&<&\ell_{1}^{*}(X)\\\
\Longleftrightarrow\hskip 21.68121pt$
$&\max_{\Lambda_{*}=\Lambda_{*}^{T}\in{\cal L}_{*}}{\rm
tr}\,(\Lambda_{*}\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix}^{T}(X+W)\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix})&<&\ell_{1}^{*}(X),\end{array}$
(112) (113)

where the last equivalence holds since the difference of the terms on the
left-hand side in the last two inequalities is bounded independently of $X$.
Also, the last inequality in (113) together with the first equality in (89)
and the first inequality in (88) produces (96). Therefore the only scenario
that is left as potentially not producing (96) is when all the elements on the
main diagonal are larger than or equal to 1. However, the two lemmas preceding
the theorem show that in such a scenario ${\cal L}^{0}={\cal L}$ and one
consequently has an equality instead of the inequality in (89) which then,
together with (88), implies (96). This completes the proof of the second (“the
only if”) part of the theorem and therefore of the entire theorem. ∎

The above derivations can be repeated for asymmetric matrices $X_{sol}$ as
well. The final result differs in a cosmetic change. Namely, the right
$\bar{U}$ should be replaced by a corresponding $\bar{V}$. Also, the condition
in the theorem relates matrix $W$ to the null-space of matrix $A$ and as such
belongs to the class of so-called vector matrix terminology (VMT) based
conditions. In the MC and C-inf cases that are of our interest here, it is
more convenient to deal with its a masking matrix terminology (MMT) analogue.
Recalling on the proof of Theorem 9 and the role of matrix $W$ within that
proof, one has that stating that $W$ belongs to the null-space of $A$ is
basically equivalent to stating that $M\circ W={\bf 0}_{n\times n}$. In other
words, one has the equivalence between the following two sets

$(W\in{\mathbb{R}}^{n\times n}|A\mbox{vec}(W)={\bf 0}_{m\times 1},W\neq{\bf
0}_{n\times n})\Longleftrightarrow(W\in{\mathbb{R}}^{n\times n}|M\circ W={\bf
0}_{n\times n},W\neq{\bf 0}_{n\times n}).$ (114)

Continuing further in the spirit of the RDT, the following corollary of the
above theorem can be established as well.

######  Corollary 1.

($\ell_{0}^{*}-\ell_{1}^{*}$-equivalence condition – general $X$) Assume the
setup of Proposition 2 with $X_{sol}$ being the unique solution of (2). Let
the masking matrix $M\in{\mathbb{R}}^{n\times n}$ have $m$ ones and
$(n^{2}-m)$ zeros and let $A$ be generated via $M$, i.e. let $A$ be the matrix
obtained after removing all the zero rows from
$\mbox{diag}^{-1}(\mbox{vec}(M))I_{n^{2}\times n^{2}}$. If and only if

$\min_{W,W^{T}W=1,M\circ W={\bf 0}_{n\times n}}{\rm
tr}\,(\bar{U}^{T}W\bar{V})+\ell_{1}^{*}((\bar{U}^{\perp})^{T}W\bar{V}^{\perp})\geq
0,$ (115)

then

$\ell_{0}^{*}\Longleftrightarrow\ell_{1}^{*}\quad\mbox{and}\quad\textbf{\emph{RMSE}}=\|\hat{X}-X_{sol}\|_{F}=0,$
(116)

and the solutions of (2) and (3) coincide.

###### Proof.

Follows immediately as a combination of 81, 82, and (114). ∎

Remark: Carefully comparing the conditions in (81) and (115) one can observe
that a strict inequality is loosened up a bit at the expense of the uniqueness
assumption. With a little bit of extra effort one may avoid this. However, to
make writings below substantially easier we will work with a non-strict
inequality.

One also has the following alternative to the above corollary.

######  Corollary 2.

($\ell_{0}^{*}-\ell_{1}^{*}$-equivalence condition via masking matrix –
general $X$) Assume the setup of Proposition 2 and Corollary 1. Let $I^{(l)}$
be as in (2.2). Then

$\ell_{0}^{*}\Longleftrightarrow\ell_{1}^{*}\quad\mbox{and}\quad\textbf{\emph{RMSE}}=\|\hat{X}-X_{sol}\|_{F}=0$
$\mbox{ If and only if}\quad\exists\Lambda|\Lambda^{T}\Lambda\leq
I\quad\mbox{and}\quad(I^{(l)})^{T}\left(\bar{V}\bar{U}^{T}+\bar{V}^{\perp}\Lambda(\bar{U}^{\perp})^{T}\right)I^{(l)}=0.$
(117)

###### Proof.

The “if” part follows from Corollary 1 and (14). The “only if” part follows
after noting that all the inequalities in (10)-(14) are written for generic
instructional purposes. Due to the underlying convexity and the strong duality
they all actually can be replaced with equalities as well. ∎

Proposition 2 and Corollaries 1 and 2 consider the so-called idealized low
rank scenario. In other words, they consider scenarios where the underlying
matrix $X_{sol}$ is of rank $k<n$. Corollary 1 particulary states that in such
a scenario one has that the nuclear norm optimization in (3) (when it comes to
its equivalence with the corresponding $\ell_{0}^{*}$ from (2)) effectively
boils down to the residual optimization problem given in (115). Also, in the
idealized context, it is possible to find scenarios when the optimal $W$ in
(115) is zero. In the non-idealized ones, i.e. when $X_{sol}$ has full rank
that is not possible and one has the residual error $W=\hat{X}-X_{sol}\neq 0$.
The following theorem (basically Theorem 9) is the key upgrade of the above
corollaries from idealized low rank to approximately low rank scenario, that
ultimately characterizes the residual error.

######  Theorem A.1.

(Algebraic characterization of $W$) Consider a
$\bar{U}\in{\mathbb{R}}^{n\times k}$ such that $\bar{U}^{T}\bar{U}=I_{k\times
k}$ and a $\bar{V}\in{\mathbb{R}}^{n\times k}$ such that
$\bar{V}^{T}\bar{V}=I_{k\times k}$ and an approximately rank $k$ matrix
$X_{sol}=X\in{\mathbb{R}}^{n\times n}$, such that
$X_{sol}=U\bar{\sigma}V^{T},\bar{\sigma}=\mbox{diag}(\sigma,\epsilon_{\sigma}),\sigma\in{\mathbb{R}}^{k},\epsilon_{\sigma}\in{\mathbb{R}}^{n-k},\sigma\gg\sigma_{\epsilon},\epsilon_{\sigma}\leq\sigma_{\epsilon}$.
Also, let the orthogonal spans $\bar{U}^{\perp}\in{\mathbb{R}}^{n\times(n-k)}$
and $\bar{V}^{\perp}\in{\mathbb{R}}^{n\times(n-k)}$ be such that
$U\triangleq\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix}$ and
$V\triangleq\begin{bmatrix}\bar{V}&\bar{V}^{\perp}\end{bmatrix}$ and

$U^{T}U\triangleq\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix}^{T}\begin{bmatrix}\bar{U}&\bar{U}^{\perp}\end{bmatrix}=I_{n\times
n}\quad\mbox{and}\quad
V^{T}V\triangleq\begin{bmatrix}\bar{V}&\bar{V}^{\perp}\end{bmatrix}^{T}\begin{bmatrix}\bar{V}&\bar{V}^{\perp}\end{bmatrix}=I_{n\times
n}.$ (118)

With $M\in{\mathbb{R}}^{n\times n}$ as in (2.2), assume that $Y=M\circ
X_{sol}$ and let $\hat{X}$ be the solution of (3). Set

$\hat{W}=\arg\min_{W,M\circ W={\bf 0}_{n\times n}}{\rm
tr}\,(\bar{U}^{T}W\bar{V})+\ell_{1}^{*}(\mbox{diag}(\epsilon_{\sigma})+(\bar{U}^{\perp})^{T}W\bar{V}^{\perp}).$
(119)

Then for any $X_{sol}$ that satisfies the above setup one has

$\textbf{\emph{RMSE}}_{(ns)}=\|\hat{X}-X_{sol}\|_{F}\leq\|\hat{W}\|_{F}.$
(120)

Moreover, there exists an $X_{sol}$ such that

$\textbf{\emph{RMSE}}_{(ns)}=\|\hat{X}-X_{sol}\|_{F}=\|\hat{W}\|_{F}.$ (121)

###### Proof.

Follows directly from Corollary 1 by repeating all arguments and accounting
for an additional term $\bar{U}^{\perp}\epsilon_{\sigma}(\bar{V}^{\perp})^{T}$
in the representation of $X_{sol}$ (that term is propagated as an addition to
$W$ and eventually appears only in the second term of (119)). Effectively,
(120) states that the residual error will be no worse (larger in
$\|\cdot\|_{F}$ norm) than the solution of the optimization in (119).
Moreover, the tightness of the analyses leading to Proposition 2 and Corollary
1 ensures that such an error is also achievable (as stated in (121)). ∎

##  Appendix B Proof of Lemma 3.2

###  B.1 Basics of free probability theory (FPT) – random matrix variables

We find it useful to first recall on some FPT basics. To start things off we
consider two symmetric matrices, $A=A^{T}\in{\mathbb{R}}^{n\times n}$ and
$B=B^{T}\in{\mathbb{R}}^{n\times n}$, as random variables. We also assume
large $n$ regime and that the eigenspaces of these matrices are Haar
distributed. Moreover, we will assume that their individual respective
spectral laws are $f_{A}(\cdot)$ and $f_{B}(\cdot)$. Furthermore, we recall on
the so-called Stieltjes transform (or as we will often call it G-transform) of
a pdf $f(\cdot)$

$\displaystyle G(z)$ $\displaystyle\triangleq$
$\displaystyle\int_{I_{f}}\frac{f(x)}{z-x}dx,\quad z\in{\mathbb{C}}\setminus
I_{f},$ (122)

where $I_{f}$ is the domain of $f(\cdot)$. One then also has the inverse
relation (somewhat analogous to the above relation between the inverse Fourier
and the underlying pdf of the sum of scalar random variables)

$\displaystyle f(x)=\lim_{\epsilon\rightarrow
0^{+}}\frac{G(x-i\epsilon)-G(x+i\epsilon)}{2i\pi}\quad\mbox{or}\quad
f(x)=-\lim_{\epsilon\rightarrow
0^{+}}\frac{\mbox{imag}(G(x+i\epsilon))}{\pi}.$ (123)

For the above to hold it makes things easier to implicitly assume that $f(x)$
is continuous. We will, however, utilize it even in discrete (or semi-
discrete) scenarios since the obvious asymptotic translation to continuity
would make it fully rigorous. A bit later though, when we see some concrete
examples where things of this nature may appear, we will say a few more words
and explain more thoroughly what exactly can be discrete and how one can deal
with such a discreteness. In the meantime we proceed with general principles
not necessarily worrying about all the underlying technicalities that may
appear in scenarios deviating from the typically seen ones and potentially
requiring additional separate addressing. To that end we continue by
considering the $R(\cdot)$\- and $S(\cdot)$-transforms that satisfy the
following

$\displaystyle R(G(z))+\frac{1}{G(z)}=z,$ (124)

and

$\displaystyle S(z)=\frac{1}{R(zS(z))}\quad\mbox{and}\quad
R(z)=\frac{1}{S(zR(z))}.$ (125)

Let $f_{A}(\cdot)$ and $f_{B}(\cdot)$ be the spectral distributions of $A$ and
$B$ and let $R_{A}(z)$/$S_{A}(z)$ and $R_{B}(z)$/$S_{B}(z)$ be their
associated $R(\cdot)$-/$S(\cdot)$-transforms. One then has the following

Key Voiculescu’s FPT concepts [46, 47]: $\displaystyle\begin{array}[]{r c l l
r c l}C&=&A+B&\quad\Longrightarrow\quad$ $&R_{C}(z)&=&R_{A}(z)+R_{B}(z)\\\
C&=&AB&\quad\Longrightarrow\quad$ $&S_{C}(z)&=&S_{A}(z)S_{B}(z).\end{array}$
(128)

Now it is relatively easy to see that (122)-(128) are sufficient to determine
the spectral distribution of the sum or the product of two independent
matrices with given spectral densities and the Haar distributed bases of
eigenspaces. The above is, of course, a generic principle. It can be applied
pretty much always as long as one has access to the statistics of the
underlying matrices $A$ and $B$. In the following section, we raise the bar a
bit higher and show that, in the cases of our interest in this paper, one can
use all of the above in such a manner that eventually all the relevant
quantities are explicitly determined. Moreover, although the methodology may,
on occasion, seem a bit involved, the final results will turn out to be
presentable in fairly neat and elegant closed forms.

###  B.2 Spectral considerations

We start with a trivial observation. Let $f_{{\cal V}}(\cdot)$ be the spectral
distribution of ${\cal V}$. Then

$\displaystyle f_{\cal V}(x)=(1-\beta)\delta(1-x)+\beta\delta(x),$ (129)

where $\delta(\cdot)$ stands for the standard delta function with nonzero
value only when its argument takes value zero. Using the definition of the
$G$-transform from (122) we can find

$\displaystyle G_{\cal V}(z)=\int\frac{f_{\cal
V}(x)}{z-x}dx=\int\frac{(1-\beta)\delta(1-x)+\beta\delta(x)}{z-x}dx=\frac{1-\beta}{z-1}+\frac{\beta}{z}=\frac{z-\beta}{z^{2}-z}.$
(130)

Also, we have

$\displaystyle R_{\cal V}(y)=z-\frac{1}{y}\quad\mbox{with}\quad y=G_{\cal
V}(z)\quad\mbox{and}\quad z=G_{\cal V}^{-1}(y).$ (131)

From (130) and (131) we further find

$\displaystyle G_{\cal
V}(z)=y\quad\Longleftrightarrow\quad\frac{z-\beta}{z^{2}-z}=y\quad\Longleftrightarrow\quad
z^{2}y-z(y+1)+\beta=0.$ (132)

Solving for $z$ gives

$\displaystyle z=\frac{y+1\pm\sqrt{(y+1)^{2}-4\beta y}}{2y}.$ (133)

Combining (131) and (133) we obtain for the $R$-transform

$\displaystyle R_{\cal V}(y)=z-\frac{1}{y}=\frac{y-1\pm\sqrt{(y+1)^{2}-4\beta
y}}{2y},$ (134)

where we for the completeness adopt the strategy to keep both $\pm$ signs. To
determine the $S$-transform we start by combining (125) and (134)

$\displaystyle S_{\cal V}(z)=\frac{1}{R_{\cal V}(zS_{\cal
V}(z))}=\frac{1}{\frac{zS_{\cal V}(z)-1\pm\sqrt{(zS_{\cal V}(z)+1)^{2}-4\beta
zS_{\cal V}(z)}}{2zS_{\cal V}(z)}}.$ (135)

After a bit of algebraic transformations we have

$\displaystyle\begin{array}[]{c r c l}&zS_{\cal
V}(z)-1-2z&=&\mp\sqrt{(zS_{\cal V}(z)+1)^{2}-4\beta zS_{\cal V}(z)}\\\
\Longleftrightarrow&(zS_{\cal V}(z)-1-2z)^{2}&=&(zS_{\cal V}(z)+1)^{2}-4\beta
zS_{\cal V}(z)\\\ \Longleftrightarrow&(zS_{\cal V}(z))^{2}-2(2z+1)zS_{\cal
V}(z)+(2z+1)^{2}&=&(zS_{\cal V}(z))^{2}+2zS_{\cal V}+1-4\beta zS_{\cal
V}(z)\\\ \Longleftrightarrow&-2(2z+1)zS_{\cal V}(z)+4z^{2}+4z&=&2zS_{\cal
V}(z)-4\beta zS_{\cal V}(z)\\\
\Longleftrightarrow&4z^{2}+4z&=&(4z^{2}+4z)S_{\cal V}(z)-4\beta zS_{\cal
V}(z)\\\ \Longleftrightarrow&z+1&=&S_{\cal V}(z)(z+1-\beta).\end{array}$ (142)
(143)

From (143) we finally have

$\displaystyle S_{\cal V}(z)=\frac{z+1}{z+1-\beta}.$ (144)

As this is a very generic result it is useful to have it formalized in the
following lemma.

######  Lemma B.1.

Let $\bar{V}^{\perp}\in{\mathbb{R}}^{n\times(n-k)}$ be Haar distributed
unitary basis of an $(n-k)$-dimensional subspace of ${\mathbb{R}}^{n}$. Let
${\cal V}$ be defined as

$\displaystyle{\cal V}\triangleq\bar{V}^{\perp}(\bar{V}^{\perp})^{T}.$ (145)

In the large $n$ linear regime, with
$\beta\triangleq\lim_{n\rightarrow\infty}\frac{k}{n}$, the $S$-transform of
the spectral density of ${\cal V}$, $f_{\cal V}(\cdot)$, is

$\displaystyle S_{\cal V}(z)=\frac{z+1}{z+1-\beta}.$ (146)

###### Proof.

Follows from the above discussion. ∎

Since ${\cal V}$ and ${\cal U}$ are structurally identical (with the only
difference being one of their dimensions) we easily have

$\displaystyle S_{\cal U}(z)=\frac{z+1}{z+1-\eta}.$ (147)

A combination of (128), (144), and (147) gives

$\displaystyle S_{\tilde{D}}(z)=\frac{(z+1)^{2}}{(z+1-\beta)(z+1-\eta)}.$
(148)

From (125) we also have

$R_{\tilde{D}}(z)=\frac{1}{S_{\tilde{D}}(zR_{\tilde{D}}(z))}=\frac{1}{\frac{(zR_{\tilde{D}}(z)+1)^{2}}{(zR_{\tilde{D}}(z)+1-\beta)(zR_{\tilde{D}}(z)+1-\eta)}}=\frac{(zR_{\tilde{D}}(z)+1-\beta)(zR_{\tilde{D}}(z)+1-\eta)}{(zR_{\tilde{D}}(z)+1)^{2}}.$
(149)

Moreover, (123) gives

$\displaystyle R_{\tilde{D}}(G_{\tilde{D}}(z))+\frac{1}{G_{\tilde{D}}(z)}=z,$
(150)

and

$\displaystyle
G_{\tilde{D}}(z)R_{\tilde{D}}(G_{\tilde{D}}(z))=zG_{\tilde{D}}(z)-1,$ (151)

From (149) one further finds

$\displaystyle
R_{\tilde{D}}(G_{\tilde{D}}(z))=\frac{(G_{\tilde{D}}(z)R_{\tilde{D}}(G_{\tilde{D}}(z))+1-\beta)(G_{\tilde{D}}(z)R_{\tilde{D}}(G_{\tilde{D}}(z))+1-\eta)}{(G_{\tilde{D}}(z)R_{\tilde{D}}(G_{\tilde{D}}(z))+1)^{2}}.$
(152)

After plugging (151) in (152), we have

$\displaystyle R_{\tilde{D}}(G_{\tilde{D}}(z))$ $\displaystyle=$
$\displaystyle\frac{(zG_{\tilde{D}}(z)-1+1-\beta)(zG_{\tilde{D}}(z)-1+1-\eta)}{(zG_{\tilde{D}}(z)-1+1)^{2}}$
(153) $\displaystyle=$
$\displaystyle\frac{(zG_{\tilde{D}}(z)-\beta)(zG_{\tilde{D}}(z)-\eta)}{(zG_{\tilde{D}}(z))^{2}}.$

A combination of (150) and (153) further gives

$\displaystyle
z-\frac{1}{G_{\tilde{D}}(z)}=\frac{(zG_{\tilde{D}}(z)-\beta)(zG_{\tilde{D}}(z)-\eta)}{(zG_{\tilde{D}}(z))^{2}}.$
(154)

From (154) we quickly find

$\displaystyle
z^{3}(G_{\tilde{D}}(z))^{2}-z^{2}G_{\tilde{D}}(z)=z^{2}(G_{\tilde{D}}(z))^{2}-(\beta+\eta)zG_{\tilde{D}}(z)+\beta\eta,$
(155)

and

$\displaystyle(G_{\tilde{D}}(z))^{2}(z^{3}-z^{2})-G_{\tilde{D}}(z)(z^{2}-z(\beta+\eta))-\beta\eta=0.$
(156)

Solving for $G_{\tilde{D}}(z)$ finally gives

$\displaystyle
G_{\tilde{D}}^{\pm}(z)=\frac{z^{2}-z(\beta+\eta)\pm\sqrt{(z^{2}-z(\beta+\eta))^{2}+4\beta\eta(z^{3}-z^{2})}}{2(z^{3}-z^{2})},$
(157)

or

$\displaystyle
G_{\tilde{D}}^{\pm}(z)=\frac{z-(\beta+\eta)\pm\sqrt{(z-(\beta+\eta))^{2}+4\beta\eta(z-1)}}{2(z^{2}-z)}.$
(158)

The above is sufficient to establish the following lemma.

######  Lemma B.2.

Let $\bar{V}^{\perp}\in{\mathbb{R}}^{n\times(n-k)}$ and
$\bar{U}_{D}^{\perp}\in{\mathbb{R}}^{n\times(n-k)}$ be Haar distributed
unitary bases of $(n-k)$-dimensional subspaces of ${\mathbb{R}}^{n}$. Let
${\cal V}$, ${\cal U}$, and $\tilde{D}$ be defined as

$\displaystyle{\cal V}$ $\displaystyle\triangleq$
$\displaystyle\bar{V}^{\perp}(\bar{V}^{\perp})^{T}$ $\displaystyle{\cal U}$
$\displaystyle\triangleq$
$\displaystyle\bar{U}_{D}^{\perp}(\bar{U}_{D}^{\perp})^{T}$
$\displaystyle\tilde{D}$ $\displaystyle\triangleq$ $\displaystyle{\cal V}{\cal
U}.$ (159)

In the large $n$ linear regime, with
$\beta\triangleq\lim_{n\rightarrow\infty}\frac{k}{n}$, the $G$-transform of
the spectral density of $\tilde{D}$, $f_{\tilde{D}}(\cdot)$, is

$\displaystyle
G_{\tilde{D}}^{\pm}(z)=\frac{z-(\beta+\eta)\pm\sqrt{(z-(\beta+\eta))^{2}+4\beta\eta(z-1)}}{2(z^{2}-z)}.$
(160)

###### Proof.

Follows from the above discussion. The “$+$/$-$” signs are taken for
negative/positive imaginary part under the root. ∎

One then relies on (123) to determine $f_{\tilde{D}}(x)$ as

$\displaystyle f_{\tilde{D}}(x)=-\lim_{\epsilon\rightarrow
0^{+}}\frac{\mbox{imag}(G_{\tilde{D}}(x+i\epsilon))}{\pi}.$ (161)

The above is a generic procedure and we in Figure 4 show the results that one
can get for two concrete values $\beta=0.2$ and $\eta=0.6$. One should note
that it is not clear _a priori_ which of the two $\pm$ signs should be used.
As Figure 4 indicates one most definitely has to be fairly careful and account
for both signs. From Figure 4 one further observes that there are four
critical points in the spectrum itself: the locations of the two delta
functions, zero and one, and two edges of the spectrum’s bulk, $x_{l}$ and
$x_{u}$. The values of these points are shown in the plots on the right hand
side. In general one can actually determine their closed forms as well.
Moreover, it turns out that one can determine the closed form of the entire
spectral function. The section that follows analyzes the spectrum of
$\tilde{D}$ in more details and eventually provides the closed form
expressions for all the relevant spectral features.

Figure 3: Both $G_{\tilde{D}}^{+}(z)$ and $G_{\tilde{D}}^{-}(z)$ need to be
taken into account

Figure 4: Both $G_{\tilde{D}}^{+}(z)$ and $G_{\tilde{D}}^{-}(z)$ need to be
taken into account

####  B.2.1 The spectrum of $\tilde{D}$ – closed form expressions

As one of our main concerns in this paper is the utilization of the final
results that we will get in this section and not necessarily the presentation
of the tiny details needed to get them, we will sketch all the key arguments
and leave out all the unnecessary minute details. However, we do emphasize
that the sketch will contain all the key pointers so that with a little bit of
effort one, if in a need, can fill in all the missing pieces of the overall
mosaic.

As mentioned above, looking at the denominator of (158) and keeping in mind
the $f\longleftrightarrow G$ connection from (161) one observes that the pdf
of interest, $f_{\tilde{D}}(x)$, potentially has two delta functions, one at
zero and the other one at one. Moreover, the bulk of the spectrum will be in
the range where the real part under the root is negative. It also goes almost
without saying that the entire spectrum will be located between zero and one.
Finally, the breaking point, $x_{c}$, where one needs to switch from
$G_{\tilde{D}}^{+}(z)$ to $G_{\tilde{D}}^{-}(z)$ in (161) is determined as the
value where the imaginary part under the root changes its sign. Equipped with
these observations one can then proceed to actually concretely determine some
of the relevant quantities.

Based on what we have just observed above, we first express $f_{\tilde{D}}(x)$
as the sum of its three key constitutive parts (two delta functions and the
bulk)

$\displaystyle
f_{\tilde{D}}(x)=f_{0}\delta(x-0)+f_{\tilde{D}}^{(b)}(x)+f_{1}\delta(x-1).$
(162)

From (162) one has that $f_{\tilde{D}}(x)$ will be fully specified if one can
determine the delta multipliers $f_{0}$ and $f_{1}$, and the bulk pdf
$f_{\tilde{D}}^{(b)}(\cdot)$.

1) Finding $f_{0}$: To determine $f_{0}$ we start by observing from (161) for
$x=0$

$\displaystyle f_{\tilde{D}}(0)=-\lim_{\epsilon\rightarrow
0^{+}}\frac{\mbox{imag}(G_{\tilde{D}}(i\epsilon))}{\pi}.$ (163)

Utilizing (158) we further have

$\displaystyle f_{\tilde{D}}(0)$ $\displaystyle=$
$\displaystyle-\frac{1}{\pi}\lim_{\epsilon\rightarrow
0^{+}}\mbox{imag}\left(\frac{i\epsilon-(\beta+\eta)\pm\sqrt{(i\epsilon-(\beta+\eta))^{2}+4\beta\eta(i\epsilon-1)}}{2((i\epsilon)^{2}-i\epsilon)}\right)$
(164) $\displaystyle=$ $\displaystyle-\frac{1}{\pi}\lim_{\epsilon\rightarrow
0^{+}}\mbox{imag}\left(\frac{-(\beta+\eta)\pm\sqrt{-\epsilon^{2}+(\beta+\eta)^{2}-4\beta\eta-2i\epsilon(\beta+\eta-2\beta\eta)}}{2(-\epsilon^{2}-i\epsilon)}\right)$
$\displaystyle=$ $\displaystyle-\frac{1}{\pi}\lim_{\epsilon\rightarrow
0^{+}}\mbox{imag}\left(\frac{-(\beta+\eta)\pm\sqrt{(\beta-\eta)^{2}}}{2(-i\epsilon)}\right)$
$\displaystyle=$ $\displaystyle-\frac{1}{\pi}\lim_{\epsilon\rightarrow
0^{+}}\mbox{imag}\left(\frac{-(\beta+\eta)-|\beta-\eta|}{-2i\epsilon}\right)$
$\displaystyle=$
$\displaystyle\left(\beta+\eta+|\beta-\eta|\right)\left(-\frac{1}{\pi}\lim_{\epsilon\rightarrow
0^{+}}\mbox{imag}\left(\frac{1}{i\epsilon}\right)\right)$ $\displaystyle=$
$\displaystyle\max(\beta,\eta)\delta(0),$

where the fourth equality (the choice of the “$-$” sign in $\pm$) follows
since $0\leq x_{c}$ (the spectrum belongs to the interval $[0,1]$ and $x_{c}$
must be in the spectrum) and the last equality follows since by convention

$\displaystyle\delta(0)=\left(-\frac{1}{\pi}\lim_{\epsilon\rightarrow
0^{+}}\mbox{imag}\left(\frac{1}{i\epsilon}\right)\right).$ (165)

To see the rationale behind (165) we briefly digress and start with

$g(x)=\delta(x).$ (166)

Then from (122)

$G(z)=\int_{x}\frac{\delta(x)dx}{z-x}=\frac{1}{z},$ (167)

and from (123)

$\delta(x)=-\lim_{\epsilon\rightarrow
0^{+}}\frac{\mbox{imag}(G(x+i\epsilon))}{\pi}.$ (168)

For $x=0$ then

$\delta(0)=-\lim_{\epsilon\rightarrow
0^{+}}\frac{\mbox{imag}(G(i\epsilon))}{\pi}=-\lim_{\epsilon\rightarrow
0^{+}}\mbox{imag}\left(\frac{1}{\pi
i\epsilon}\right)=-\frac{1}{\pi}\lim_{\epsilon\rightarrow
0^{+}}\mbox{imag}\left(\frac{1}{i\epsilon}\right),$ (169)

which is identical to (165). The above description of the delta function may
not necessarily be the most adequate one. However, for what we need here it is
conceptually sufficient. Namely, we are here interested in determining the
proportionality constants that multiply the delta functions rather than the
functions’ expressions themselves. One way to make everything more adequate
would be to translate everything into the continuous domain by choosing a
continuous function as an asymptotic replacement for $\delta(x)$. For example,
one can use the Gaussian continual approximation

$\delta(x)=\lim_{\sigma\rightarrow
0^{+}}\frac{e^{-\frac{x^{2}}{2\sigma^{2}}}}{\sqrt{2\pi\sigma^{2}}}.$ (170)

Then all the above holds for small $\sigma=\epsilon\sqrt{\pi}/2$ and

$\delta(x)\rightarrow\lim_{\sigma\rightarrow
0^{+}}\frac{e^{-\frac{x^{2}}{2\sigma^{2}}}}{\sqrt{2\pi\sigma^{2}}}\rightarrow\lim_{\sigma\rightarrow
0^{+}}\frac{e^{-\frac{x^{2}}{\pi\epsilon^{2}}}}{\pi\epsilon}\quad\mbox{and}\quad\delta(0)\rightarrow\lim_{\sigma\rightarrow
0^{+}}\frac{1}{\pi\epsilon}.$ (171)

While computing and manipulating the transforms mentioned, it is necessary to
consider the minor $\epsilon$-differences they induce. These differences are
conceptually insignificant, and our results remain valid for small values of
$\sigma$ or $\epsilon$. Including these $\epsilon$-modifications in our
writing would lead to a more cumbersome presentation, demanding the addition
of many minor details to demonstrate their marginal impact. Since, on the
other hand, they contribute exactly nothing to the essence of the arguments
and final results we chose to operate in a semi-discrete domain with the delta
functions. As a consequence one has the expressions given in (165) and (169).
We believe that a little bit of conventional inadequacy is better than to
overwhelm the presentation with details which would make the overall content
less accessible.

2) Finding $f_{1}$: To determine $f_{1}$ we follow the above methodology and
start by observing from (161) for $x=1$

$\displaystyle f_{\tilde{D}}(1)=-\lim_{\epsilon\rightarrow
0^{+}}\frac{\mbox{imag}(G_{\tilde{D}}(1+i\epsilon))}{\pi}.$ (172)

Further utilization of (158) gives

$\displaystyle f_{\tilde{D}}(1)$ $\displaystyle=$
$\displaystyle\\-\frac{1}{\pi}\lim_{\epsilon\rightarrow
0^{+}}\mbox{imag}\left(\frac{1+i\epsilon-(\beta+\eta)\pm\sqrt{(1+i\epsilon-(\beta+\eta))^{2}+4\beta\eta
i\epsilon}}{2((1+i\epsilon)^{2}-1-i\epsilon)}\right)$ (173) $\displaystyle=$
$\displaystyle-\frac{1}{\pi}\lim_{\epsilon\rightarrow
0^{+}}\mbox{imag}\left(\frac{1-(\beta+\eta)\pm\sqrt{-\epsilon^{2}+(1-(\beta+\eta))^{2}-2i\epsilon(-1+\beta+\eta-2\beta\eta)}}{2(-\epsilon^{2}+i\epsilon)}\right)$
$\displaystyle=$ $\displaystyle-\frac{1}{\pi}\lim_{\epsilon\rightarrow
0^{+}}\mbox{imag}\left(\frac{1-(\beta+\eta)\pm\sqrt{(1-(\beta-\eta))^{2}}}{2(i\epsilon)}\right)$
$\displaystyle=$ $\displaystyle-\frac{1}{\pi}\lim_{\epsilon\rightarrow
0^{+}}\mbox{imag}\left(\frac{1-(\beta+\eta)+|1-(\beta+\eta)|}{2i\epsilon}\right)$
$\displaystyle=$
$\displaystyle\left(1-(\beta+\eta)+|1-(\beta+\eta)|\right)\left(-\frac{1}{\pi}\lim_{\epsilon\rightarrow
0^{+}}\mbox{imag}\left(\frac{1}{i\epsilon}\right)\right)$ $\displaystyle=$
$\displaystyle\max(1-(\beta+\eta),0)\delta(0),$

where the fourth equality (the choice of the “$+$” sign in $\pm$) follows
since now $x_{c}\leq 1$ and the last equality follows by the above discussed
$\delta(0)$ convention.

3) Finding $f_{\tilde{D}}^{(b)}(x)$: To determine $f_{\tilde{D}}^{(b)}(x)$ for
$x\notin\\{0,1\\}$ we again start with (161) and, for a general $x$, we write
the following

$\displaystyle f_{\tilde{D}}^{(b)}(x)=-\lim_{\epsilon\rightarrow
0^{+}}\frac{\mbox{imag}(G_{\tilde{D}}(x+i\epsilon))}{\pi}.$ (174)

Relying once again on (158) we, for $x\notin\\{0,1\\}$, have

$\displaystyle f_{\tilde{D}}^{(b)}(x)$ $\displaystyle=$
$\displaystyle-\frac{1}{\pi}\lim_{\epsilon\rightarrow
0^{+}}\mbox{imag}\left(\frac{x+i\epsilon-(\beta+\eta)\pm\sqrt{(x+i\epsilon-(\beta+\eta))^{2}+4\beta\eta(x+i\epsilon-1}}{2((x+i\epsilon)^{2}-x-i\epsilon)}\right)$
(175) $\displaystyle=$ $\displaystyle-\frac{1}{\pi}\lim_{\epsilon\rightarrow
0^{+}}\mbox{imag}\left(\frac{x-(\beta+\eta)\pm\sqrt{-\epsilon^{2}+(x-(\beta+\eta))^{2}+4\beta\eta(x-1)-2i\epsilon(-x+\beta+\eta-2\beta\eta)}}{2(x^{2}-x-\epsilon^{2}+i\epsilon(2x-1))}\right)$
$\displaystyle=$ $\displaystyle-\frac{1}{\pi}\lim_{\epsilon\rightarrow
0^{+}}\mbox{imag}\left(\frac{x-(\beta+\eta)\pm\sqrt{(x-(\beta+\eta))^{2}+4\beta\eta(x-1)-2i\epsilon(-x+\beta+\eta-2\beta\eta)}}{2(x^{2}-x)}\right)$
$\displaystyle=$ $\displaystyle-\frac{1}{\pi}\lim_{\epsilon\rightarrow
0^{+}}\mbox{imag}\left(\frac{\pm\sqrt{(x-(\beta+\eta))^{2}+4\beta\eta(x-1)-2i\epsilon(-x+\beta+\eta-2\beta\eta)}}{2(x^{2}-x)}\right).$

Now, since one is interested in the imaginary part of interest is the region
of $x$ where the real part under the root is negative (outside that region,
i.e in the region of $x$ where the real part under the root is nonnegative
$f_{\tilde{D}}^{(b)}(x)$ is zero). To determine the region of interest we
start by setting

$\displaystyle
T_{\tilde{D}}\triangleq\\{x\in{\mathbb{R}}|(x-(\beta+\eta))^{2}+4\beta\eta(x-1)\leq
0\\}\quad\mbox{and}\quad x_{c}\triangleq\beta+\eta-2\beta\eta.$ (176)

To explicitly characterize $T_{\tilde{D}}$ we look at the following

$\displaystyle\begin{array}[]{r r c l}$
$&(x-(\beta+\eta))^{2}+4\beta\eta(x-1)&=&0\\\
\Longleftrightarrow&x^{2}-2x(\beta+\eta-2\beta\eta)+(\beta+\eta)^{2}-4\beta\eta&=&0\\\
\Longleftrightarrow&x^{2}-2x(\beta+\eta-2\beta\eta)+(\beta-\eta)^{2}&=&0.\end{array}$
(180)

Solving for $x$ one finds

$\displaystyle x$ $\displaystyle=$
$\displaystyle\frac{2(\beta+\eta-2\beta\eta)\pm\sqrt{(2(\beta+\eta-2\beta\eta))^{2}-4(\beta-\eta)^{2}}}{2}$
(181) $\displaystyle=$
$\displaystyle\beta+\eta-2\beta\eta\pm\sqrt{(\beta+\eta-2\beta\eta)^{2}-(\beta-\eta)^{2}}.$

Setting

$\displaystyle x_{l}$ $\displaystyle\triangleq$
$\displaystyle\beta+\eta-2\beta\eta-\sqrt{(\beta+\eta-2\beta\eta)^{2}-(\beta-\eta)^{2}}$
$\displaystyle x_{u}$ $\displaystyle\triangleq$
$\displaystyle\beta+\eta-2\beta\eta+\sqrt{(\beta+\eta-2\beta\eta)^{2}-(\beta-\eta)^{2}},$
(182)

one has

$\displaystyle T_{\tilde{D}}=\\{x\in{\mathbb{R}}|x\in[x_{l},x_{u}]\\}.$ (183)

Moreover, from (182), one also has

$\displaystyle 0\leq x_{l}\leq x_{u}\leq 1,$ (184)

with

$\displaystyle x_{l}=0\quad\mbox{if}\quad\beta=\eta\quad\mbox{and}\quad
x_{u}=1\quad\mbox{if}\quad\beta=\eta=0.5.$ (185)

The first two inequalities in (182) are trivial, whereas the third one follows
after noting

$\displaystyle\beta+\eta-2\beta\eta\leq\max(\beta,1-\beta)\leq 1,$ (186)

and observing the following sequence

$\displaystyle\begin{array}[]{r r c
l}&\beta+\eta-2\beta\eta+\sqrt{(\beta+\eta-2\beta\eta)^{2}-(\beta-\eta)^{2}}&\leq&1\\\
\Longleftrightarrow&(\beta+\eta-2\beta\eta)^{2}-(\beta-\eta)^{2}&\leq&(1-(\beta+\eta-2\beta\eta))^{2}\\\
\Longleftrightarrow&-(\beta-\eta)^{2}&\leq&1-2(\beta+\eta-2\beta\eta)\\\
\Longleftrightarrow&-(\beta-\eta)^{2}+2(\beta+\eta-2\beta\eta)-1&\leq&0\\\
\Longleftrightarrow&-(\beta+\eta)^{2}+2(\beta+\eta)-1&\leq&0\\\
\Longleftrightarrow&-(1-(\beta+\eta))^{2}&\leq&0.\end{array}$ (193)

Returning to (175) we further have for $x\in T_{\tilde{D}}=[x_{l},x_{u}]$

$\displaystyle f_{\tilde{D}}^{(b)}(x)$ $\displaystyle=$
$\displaystyle-\frac{1}{\pi}\lim_{\epsilon\rightarrow
0^{+}}\mbox{imag}\left(\frac{\pm\sqrt{(x-(\beta+\eta))^{2}+4\beta\eta(x-1)-2i\epsilon(-x+\beta+\eta-2\beta\eta)}}{2(x^{2}-x)}\right)$
(194) $\displaystyle=$ $\displaystyle-\frac{1}{\pi}\lim_{\epsilon\rightarrow
0^{+}}\mbox{imag}\left(\frac{\pm\sqrt{(x-(\beta+\eta))^{2}+4\beta\eta(x-1)-2i\epsilon(x_{c}-x)}}{2(x^{2}-x)}\right)$
$\displaystyle=$ $\displaystyle-\frac{1}{\pi}\lim_{\epsilon\rightarrow
0^{+}}\mbox{imag}\left(\frac{i\sqrt{-(x-(\beta+\eta))^{2}-4\beta\eta(x-1)}}{2(x^{2}-x)}\right),$

where the “$+$” plus sign is chosen if $x_{c}\leq x\leq x_{u}$ and the “$-$”
sign is chosen if $x_{l}\leq x\leq x_{c}$. Finally, from (194) one easily
finds

$\displaystyle f_{\tilde{D}}^{(b)}(x)$ $\displaystyle=$
$\displaystyle\frac{\sqrt{-(x-(\beta+\eta))^{2}-4\beta\eta(x-1)}}{2\pi(x-x^{2})}\quad\mbox{
if }\quad x_{l}\leq x\leq x_{u}.$ (195)

The above is then sufficient to completely characterize the spectral
distribution $f_{\tilde{D}}(x)$. We summarize the results in the following
lemma (basically a mirrored version of Lemma 3.2).

######  Lemma B.3.

Assume large $n$ linear regime with

$\displaystyle\beta\triangleq\lim_{n\rightarrow\infty}\frac{k}{n}\quad\mbox{and}\quad\eta\triangleq\lim_{n\rightarrow\infty}\frac{l}{n}.$
(196)

Let $\bar{V}^{\perp}\in{\mathbb{R}}^{n\times(n-k)}$ be a Haar distributed
basis of an $n-k$-dimensional subspace of $R^{n}$. Analogously, let
$\bar{U}_{D}^{\perp}\in{\mathbb{R}}^{n\times(n-k)}$ be a Haar distributed
basis of an $n-l$-dimensional subspace of $R^{n}$. Moreover, let
$\bar{V}^{\perp}\in{\mathbb{R}}^{n\times(n-k)}$ and
$\bar{U}_{D}^{\perp}\in{\mathbb{R}}^{n\times(n-k)}$ be independent of each
other. Also, let ${\cal V}$, ${\cal U}$, and $\tilde{D}$ be defined as

$\displaystyle{\cal V}$ $\displaystyle\triangleq$
$\displaystyle\bar{V}^{\perp}(\bar{V}^{\perp})^{T}$ $\displaystyle{\cal U}$
$\displaystyle\triangleq$
$\displaystyle\bar{U}_{D}^{\perp}(\bar{U}_{D}^{\perp})^{T}$
$\displaystyle\tilde{D}$ $\displaystyle\triangleq$ $\displaystyle{\cal V}{\cal
U}.$ (197)

Set $x_{l}$ and $x_{u}$ as in (182), i.e.

$\displaystyle x_{l}$ $\displaystyle\triangleq$
$\displaystyle\beta+\eta-2\beta\eta-\sqrt{(\beta+\eta-2\beta\eta)^{2}-(\beta-\eta)^{2}}$
$\displaystyle x_{u}$ $\displaystyle\triangleq$
$\displaystyle\beta+\eta-2\beta\eta+\sqrt{(\beta+\eta-2\beta\eta)^{2}-(\beta-\eta)^{2}},$
(198)

Then the limiting spectral distribution of $\tilde{D}$, $f_{\tilde{D}}(x)$, is

$\displaystyle f_{\tilde{D}}(x)$ $\displaystyle=$ $\displaystyle
f_{0}\delta(x)+f_{\tilde{D}}^{(b)}(x)+f_{1}f_{0}\delta(x-1)$ (199)
$\displaystyle=$
$\displaystyle\max(\beta,\eta)\delta(x)+f_{\tilde{D}}^{(b)}(x)+\max(1-(\beta+\eta),0)\delta(x-1),$

with

$\displaystyle f_{\tilde{D}}^{(b)}(x)$ $\displaystyle=$
$\displaystyle\begin{cases}\frac{\sqrt{-(x-(\beta+\eta))^{2}-4\beta\eta(x-1)}}{2\pi(x-x^{2})},&\mbox{if
}x_{l}\leq x\leq x_{u}.\\\ 0,&\mbox{otherwise}.\end{cases}$ (200)

###### Proof.

Follows through a combination of (162), (164), (173), (195), and the above
discussion. ∎

In Figure 5 we show the spectral function obtained based on the above lemma
for $\beta=0.1$ and $\eta=0.8$. We observe a very strong agreement between the
simulated results and the above theoretical predictions. Simulation results
were obtained using moderately large $n=4000$.

Figure 5: $f_{\tilde{D}}(x)$ – spectral function of $\tilde{D}$; $\beta=0.1$
and $\eta=0.8$

In Figure 6 we show the spectral function obtained based on the above lemma
for $\beta=0.2$ and $\eta=0.9$. Due to a remarkable property of the underlying
functions the spectrum is identical as in Figure 5 apart from the fact that
the multiplier of the delta function at zero is increased from $0.8$ to $0.9$
at the expense of removing the delta function at one. We also again observe a
very strong agreement between the simulated results and the theoretical
predictions. As in Figure 5, Simulation results were again obtained for
$n=4000$.

Figure 6: $f_{\tilde{D}}(x)$ – spectral function of $\tilde{D}$; $\beta=0.2$
and $\eta=0.9$

Various other proofs of the lemma are possible. We, however, found the above
one as very instructive.

##  Appendix C Real data properties

As we have mentioned in the introduction, it is very important to properly
choose the model for the analyzed system. Since the algorithmic methodology
that we analyzed applies to a wide range of causal inference settings (and
especially so in the block causal inference), we want to ensure that the
underlying model is properly selected so to represent a good fit to real data.
The approximately low rank model that we consider is typically the preferred
choice when the deviation from the ideal low rankness is caused by internal
data properties rather by the errors in observations. A typical example where
that happens is the so-called California tobacco study from [2, 1] (a golden
standard, most often used in real data studies in recent years within the
synthetic control methods).

Namely, as we discussed in the paper, in this example one postulates that it
is likely that over years there is a similar trend in tobacco consumption
across a set of different states in North America. This in turn implies that a
matrix containing (in each row) annual consumption data of each state should
be approximately low rank. In Figure 7, we plot the spectrum of the resulting
tobacco consumption matrix. As it can be seen from the figure, the spectrum
clearly indicates the approximately low rank properties. Magnitudes of several
largest singular values are substantially larger than magnitudes of the
remaining ones (the matrix also happens to have a strong spiked component that
doesn’t change anything regarding the overall approximate low rankness; in
fact it actually strengthens it).

For completeness, we do mention that the matrix of interest is of size
$38\times 31$ and that it contains the annual average cigarette consumption
for $38$ states over a span of $31$ years (this basically implies that the
size of $X_{sol}$ is $38\times 31$ and the total number of singular values is
$31$ as the figure indicates). We once again point out that one should observe
that the approximate low rankness comes as a consequence of the internal
structure of the data rather than as a consequence of incorrectly recorded
data.

Figure 7: California smoking data set. We plot the spectrum of $X_{sol}$
(without and with the largest component) – (upper part) without rank one
spike; (lower part) with rank one spike. One of the singular values is
substantially larger than any of the remaining ones (the largest singular
value is above $4000$, whereas all others are smaller than $300$).

## References

  * [1] [author] Abadie, A.A. and Cattaneo, M. D.M. D. (2018). Econometric Methods for Program Evaluation. Annual Review of Economics 10 465-503.
  * [2] [author] Abadie, A.A., Diamond, A.A. and Hainmueller, J.J. (2010). Synthetic Control Methods for Comparative Case Studies: Estimating the Effect of California’s Tobacco Control Program. Journal of the American Statistical Association 105 493-505.
  * [3] [author] Agarwal, A.A., Dahleh, M.M., Shah, D.D. and Shen, D.D. (2021). Causal Matrix Completion. available online at http://arxiv.org/abs/2109.15154.
  * [4] [author] Agarwal, A.A., Shah, D.D. and Shen, D.D. (2020). Synthetic interventions. available online at http://arxiv.org/abs/2006.07691.
  * [5] [author] Athey, S.S., Bayati, M.M., Doudchenko, N.N., Imbens, G.G. and Khosravi, K.K. (2021). Matrix Completion Methods for Causal Panel Data Models. Journal of the American Statistical Association 116 1716-1730.
  * [6] [author] Athey, S.S. and Imbens, G. W.G. W. (2018). Design-Based Analysis in Differencein- Differences Settings With Staggered Adoption. Technical Report, National Bureau of Economic Research.
  * [7] [author] Athey, S.S. and Stren, S.S. (2002). The Impact of Information Technology on Emergency Health Care Outcomes. The RAND Journal of Economics 33 399-432.
  * [8] [author] Bai, J.J. and Ng, S.S. (2019). Matrix completion, counterfactuals, and factor analysis of missing data. available online at http://arxiv.org/abs/1910.06677.
  * [9] [author] Bhattacharya, S.S. and Chatterjee, S.S. (2021). Matrix completion with data-dependent missingness probabilities. available online at http://arxiv.org/abs/2106.02290.
  * [10] [author] Bryzgalova, S.S., Lerner, S.S., Lettau, M.M. and Pelger, M.M. (2023). Missing Financial Data. Working Paper. Available online at https://papers.ssrn.com/sol3/papers.cfm?abstract_id=4106794.
  * [11] [author] Cahan, E.E., Bai, J.J. and Ng, S.S. (2021). Factor-Based Imputation of Missing Values and Covariances in Panel Data of Large Dimensions. available online at http://arxiv.org/abs/2103.03045.
  * [12] [author] Cai, T. T.T. T. and Zhou, W. X.W. X. (2016). Matrix completion via max-norm constrained optimization. Electronic Journal of Statistics 10 1493-1525.
  * [13] [author] Candes, E. J.E. J. and Plan, Y.Y. (2010). Matrix Completion With Noise. Proceedings of the IEEE 98 925-936.
  * [14] [author] Candes, E. J.E. J. and Recht, B.B. (2009). Exact Matrix Completion via Convex Optimization. Foundations of Computational Mathematics 9 717-772.
  * [15] [author] Candes, E. J.E. J. and Recht, B.B. (2011). Simple bounds for recovering low-complexity models. Math Program 1-11.
  * [16] [author] Candes, E. J.E. J. and Tao, T.T. (2010). The Power of Convex Relaxation: Near-Optimal Matrix Completion. IEEE Transactions on Information Theory 56 2053-2080.
  * [17] [author] Candes E. J., Ma Y. Wright J.M. Y. W. J. Li X (2009). Robust principal component analysis? JACM 58\.
  * [18] [author] Capponi, A.A. and Stojnic, M.M. (2023). Causal Inference – closed form expressions for worst case typical phase transitions. Proceedings of the 2023 IEEE International Symposium on Information Theory.
  * [19] [author] Capponi, A.A. and Stojnic, M.M. (2023). Supplement to Recovery of Approximately Low Rank Matrices with Block Missing Patterns.
  * [20] [author] Donoho, D.D. (2006). High-dimensional centrally symmetric polytopes with neighborlines proportional to dimension. Disc. Comput. Geometry 35 617-652.
  * [21] [author] Donoho, D. L.D. L., Gavish, M.M. and Montanari, A.A. (2013). The phase transition of matrix recovery from Gaussian measurements matches the minimax MSE of matrix denoising. Proceedings of the National Academy of Sciences (PNAS) 110 8405-8410.
  * [22] [author] Duan, J.J., Pelger, M.M. and Xiong, R.R. (2023). Target PCA: Transfer Learning Large Dimensional Panel Data. available online at http://arxiv.org/abs/2308.15627.
  * [23] [author] Fazel, M.M., Hindi, H.H. and Boyd, S. P.S. P. (2001). A rank minimization heuristic with application to minimum order system approximation. Proceedings of the IEEE American Control Conference 4734-4739.
  * [24] [author] Foygel, R.R. and Srebro, N.N. (2011). Concentration-based guarantees for low-rank matrix reconstruction. In Proceedings of the 24th Annual Conference on Learning Theory 315-340.
  * [25] [author] Kallus, N.N., Mao, X.X. and Udell, M.M. (2018). Causal inference with noisy and missing covariates via matrix factorization. Proceedings of the 32rd International Conference on Neural Information Processing 6921-6932.
  * [26] [author] Keshavan, R. H.R. H., Montanari, A.A. and Oh, S.S. (2010). Matrix Completion From a Few Entries. IEEE Transactions on Information Theory 56 2980-2998.
  * [27] [author] Keshavan, R. H.R. H., Montanari, A.A. and Oh, S.S. (2010). Matrix Completion From Noisy Entries. Journal of Machine Learning Research 11 2057-2078.
  * [28] [author] Klopp, O.O. (2014). Noisy Low-Rank Matrix Completion With General Sampling Distribution. Bernoulli 20 282-303.
  * [29] [author] Koltchinskii, V.V., Lounici, K.K. and Tsybakov, A. B.A. B. (2011). Nuclear-Norm Penalization and Optimal Rates for Noisy Low-Rank Matrix Completion. The Annals of Statistics 39 2302-2329.
  * [30] [author] Ma, W.W. and Chen, G. H.G. H. (2019). Missing Not at Random in Matrix Completion: The Effectiveness of Estimating Missingness Probabilities Under a Low Nuclear Norm Assumption. Proceedings of the 33rd Conference on Neural Information Processing Systems (NeurIPS), Vancouver, Canada. available online at http://arxiv.org/abs/1910.12774.
  * [31] [author] Mazumder, R.R., Hastie, T.T. and Tibshirani, R.R. (2010). Spectral Regularization Algorithms for Learning Large Incomplete Matrices. Journal of Machine Learning Research 11 2287-2322.
  * [32] [author] Muhammad, A.A., Shah, D.D. and Shen, D.D. (2018). Robust synthetic control. The Journal of Machine Learning Research 19 802-852.
  * [33] [author] Negahban, S. N.S. N. and Wainwright, M. J.M. J. (2011). Estimation of (Near) Low-Rank MatricesWith Noise and High-Dimensional Scaling. The Annals of Statistics 39 1069-1097.
  * [34] [author] Negahban, S. N.S. N. and Wainwright, M. J.M. J. (2012). Restricted Strong Convexity and Weighted Matrix Completion: Optimal Bounds With Noise. Journal of Machine Learning Research 13 1685-1697.
  * [35] [author] Oymak, S.S. and Hassibi., B.B. (2010). New null space results and recovery thresholds for matrix rank minimization. available online at http://arxiv.org/abs/1011.6326.
  * [36] [author] Recht, B.B. (2011). A Simpler Approach to Matrix Completion. Journal of Machine Learning Research 12 3413-3430.
  * [37] [author] Recht, B.B., Fazel, M.M. and Parrilo, P. A.P. A. (2010). Guaranteed minimum-rank solutions of linear matrix equations via nuclear norm minimization. SIAM Rev 52 471-501.
  * [38] [author] Rohde, A.A. and Tsybakov, A. B.A. B. (2011). Estimation of High-Dimensional Low-Rank Matrices. The Annals of Statistics 39 887-930.
  * [39] [author] Rubin, D. B.D. B. (1976). Biometrika. 63 58-592.
  * [40] [author] Schnabel, T.T., Swaminathan, A.A., Singh, A.A., Chandak, N.N. and Joachims, T.T. (2016). Target PCA: Transfer Learning Large Dimensional Panel Data. Proceedings of the 33 rd International Conference on Machine Learning, New York, NY, USA, 2016. JMLR: W&CP volume 48. available online at http://arxiv.org/abs/1602.05352.
  * [41] [author] Shaikh, A.A. and Toulis, P.P. (2019). Randomization Tests in Observational Studies With Staggered Adoption of Treatment. University of Chicago, Becker Friedman Institute for Economics Working Paper 144\.
  * [42] [author] Sportisse, A.A., Boyera, C.C. and Josse, J.J. (2018). Imputation and low-rank estimation with Missing Not At Random data. available online at http://arxiv.org/abs/1812.11409.
  * [43] [author] Srebro, N.N., Alon, N.N. and Jaakkola, T. S.T. S. (2005). Generalization Error Bounds for Collaborative Prediction With Low-Rank Matrices. in Advances in Neural Information Processing Systems 17 1321-1328. eds. L. K. Saul, Y.Weiss, and L. Bottou.
  * [44] [author] Srebro, N.N. and Shraibman, A.A. (2005). Rank, trace-norm and maxnorm. International Conference on Computational Learning Theory 543-560.
  * [45] [author] Stojnic, M.M. (2009). A simple performance analysis of $\ell_{1}$-optimization in compressed sensing. ICASSP, International Conference on Acoustics, Signal and Speech Processing 3021-3024. Taipei, Taiwan.
  * [46] [author] Voiculescu, D.D. (1986). Addition of certain non-commuting random variables. J. Funct. Anal. 66 323-346.
  * [47] [author] Voiculescu, D.D. (1987). Multiplication of certain noncommuting random variables. J. Operator Theory 18 2223-2235.
  * [48] [author] Voiculescu, D.D. (1991). Limit laws for random matrices and free products. Invent. Math. 104 201-220.
  * [49] [author] Wang, J.J., Wong, R. K. W.R. K. W., Mao, X.X., Chuen, K.K. and Chan, G.G. (2021). Matrix Completion with Model-free Weighting. Proceedings of the 38th International Conference on Machine Learning, PMLR 139\.
  * [50] [author] Xiong, R.R. and Pelger, M.M. (2018). Large Dimensional Latent Factor Modeling with Missing Observations and Applications to Causal Inference. Electronic copy available at: https://ssrn.com/abstract=3465357.
  * [51] [author] Zhu, Z.Z., Wang, T.T. and Samworth, R.R. (2019). High-dimensional principal component analysis with heterogeneous missingness. available online at http://arxiv.org/abs/1906.12125.