Conformalized matrix completion
Authors: Yu Gui
Rina Barber
Cong Ma
Venue: NeurIPS 2023 poster (2023)
Keywords: matrix completion
conformal inference
uncertainty quantification
TL;DR: In this work, we propose a distribution-free method for predictive inference in the matrix completion problem which adapts the framework of conformal prediction to offer provable predictive coverage regardless of the accuracy of the low-rank model.

Abstract: Matrix completion aims to estimate missing entries in a data matrix, using the assumption of a low-complexity structure (e.g., low-rankness) so that imputation is possible. While many effective estimation algorithms exist in the literature, uncertainty quantification for this problem has proved to be challenging, and existing methods are extremely sensitive to model misspecification. In this work, we propose a distribution-free method for predictive inference in the matrix completion problem. Our method adapts the framework of conformal prediction, which provides prediction intervals with guaranteed distribution-free validity in the setting of regression, to the problem of matrix completion. Our resulting method, conformalized matrix completion (cmc), offers provable predictive coverage regardless of the accuracy of the low-rank model. Empirical results on simulated and real data demonstrate that cmc is robust to model misspecification while matching the performance of existing model-based methods when the model is correct.

Conformalized matrix completion
Yu Gui
Rina Foygel Barber
Cong Ma
Department of Statistics, University of Chicago {yugui rina, congm/ouchiccgo..ee.
Abstract
Matrix completion aims to estimate missing entries in a data matrix, using the assumption low-complexity structure (e.g., low rank) so that imputation is pos- sible. While many effective estimation algorithms exist in the literature, uncertainty quantification for this problem has proved to be challenging, and existing methods are extremely sensitive to model misspecification. In this work, we propose distribution-fire method for predictive inference in the matrix completion problem. Our method adapts the framework of conformal prediction, which provides confi- dence intervals with guaranteed distribution free validity in the setting of regression, to the problem of matrix completion. Our resulting method, conformalized matrix completion (cmc), offers provable predictive coverage regardless of the accuracy of the low-rank model. Empirical results simulated and real data demonstrate that cmc is robust to model misspecification while matching the performance existing model-based methods when the model is correct.
Introduction
Matrix completion, the task of filling in missing entries in large data matrix, has wide range of applications such as gene-disease association analysis in bioinformatics (Natarajan and Dhillon] 2014), collaborative filtering in recommender systems (Rennie and Srebro 2005). and panel data prediction and inference in econometrics (Amjad et al. [0188] Bai and Ng 2021][Athe et al. (2021). If the underlying signal is assumed to be low-rank, a range of estimation algorithms have been proposed in the literature, including approaches based on convex relaxations of rank (Candes and
Plan 2010 Candès and Tao 2010; Koltchinskii al. 2011 Foygel and Srebro 2011 Negahban and Wainwright [012 Yang and Ma (2022) and nonconvex optimization over the space of low-rank
matrices (Keshavan al. 2010 Jainn Jain al. 2013 Sun and Luo 2016 Ma et al., (2020).
While the problem of estimation has been explored through many different approaches in the literature, the question of uncertainty quantity for matrix completion remains challenging and under-explored.
Cai et al. 2016) and Carppentee et al. (2018) construct confidence intervals for the missing entries using order-wise concentration inequalities, which could be extremely loose. [hen al. (2019) proposes de-biased estimator for matrix completion and characterizes its asymptotic distribution, which then yields entrywise confidence intervals are then derived for the underlying incoherent matrix with i.idd Gaussian noise. See also Xia and Yuan (2021): Farias et al. (022) for related approaches.
The effectiveness of the aforementioned uncertainty quantification methods rely heavily the exact low-rank structure as well as assumptions on the tail of the noise. In practice, the low-rank assumption can only be met approximately, and heavy-tailed noise is quite common in areas such as macroeconomics and genomics (Nair et al.) 2022). Without exact model specifications, the solution the matrix completion algorithm hard to interpret and the asymptotic distribution of the obtained estimator is no longer valid. The dependence of model-based approaches on exact model assumptions motivates to answer the question: can we construct valid and short confidence intervals for the
37th Conference on Neural Information Processing Systems (NcurIPS 2023).

---

missing entries that are free of both model assumptions and of constraints on which estimation algorithms we may use?
1.1 Conformal prediction for distribution-free uncertainty quantification
Our proposed method is built on the idea of conformal prediction. The conformal prediction framework, developed by Vovk al. (2005) and (Shafer and Vovk (2008), has drawn significant attention in recent years in that it enables the construction of distribution-fre confidence intervals that are valid with exchangeable data from any underlying distribution and with any "black-box" algorithm. To reduce computation in full conformal prediction, variants exist based on data splitting (Vovk et al., 20055 Lei al. 2018), leave-one-out (Barber et al.. (02I). cross-validation (Vovk) 2015 2021), etc.
Barber et al.
When distribution shift is present, with covariate shift special case, the exchangeability is violated, and Tibshirani et al. (2019) proposes a more general procedure called weighted conformal prediction that guarantees the validity under the weighted exchangeability of data. [Barber al. (2022) further relaxes the exchangeability assumption on the data and characterizes the robustness of a generalized weighted approach via an upper bound for the coverage gap. Other related works (Lei and Candès 2020 Candès et al. [202T] Jin et al., 2023) adapt the conformal framework to applications in causal inference, survival analysis, etc., and study the robustness of weighted conformal prediction with estimated weights.
1.2 Main contributions
In this paper, we adapt the framework of conformal prediction for the problem of matrix completion, and make the following contributions:
1. We construct distribution- frree confidence intervals for the missing entries in matrix completion via conformal prediction. The validity is free of any assumption on the underlying matrix and holds regardless of the choice of estimation algorithms practitioners use. To achieve this, we prove the (weighted) exchangeability of unobserved and observed units when they are sampled (possibly nonuniformly) without replacement from a finite population.
2. When the sampling mechanism is unknown, we develop provable lower bound for the coverage rate which degrades gracefully as the estimation error of the sampling probability increases.
3. In addition, to improve computational efficiency when faced with a large number of missing entries, we propose a one-shot conformalized matrix completion approach that only computes the weighted quantile once for all missing entries.
2 Problem setup
In this section, we formulate the matrix completion problem and contrast our distribution-free setting with the model-based settings that more common in the existing literature. For a partially-observed dj x d2 matrix, the subset S C [di] x [d2] denotes the set of indices where data is observed-that is, our observations consist of (Maj)it.))ss )s.))sss In much of the matrix completion literature, is common to assume a signal-plus-nois model for the observations,
Mi MG + Esj,
where M' (Mi) (siield) is the "true" underlying signal while (Eij) (E,))s denotes the (typically zero-mean) noise. Frequently, it is assumed that M is low-rank, or follows some other latent low-dimensional structure, so that recovering M from the observed data an identifiable problem. In works of this type, the goal is to construct an estimator M that accurately recovers the underlying M In contrast, in more model-free setting, we may no longer wish hypothesize signal-plus-not model, and can instead assume that we are observing a random subset of entries of deterministic matrix M; in this setting, estimating M itself becomes the goal.
For both frameworks, many existing results focus on the problem of estimation (either of M or of M), with results establishing bounds on estimation errors under various assumptions. For instance, see Candes and Plan (2010); [Candès and Tao (2010); Negahban and Wainwright (2012); [Chen et al.) (2020) for results estimating M (under stronger conditions, e.g., low-rank and incoherent signal,

---

and zero-mean sub-Suuseee Gaussian noise), or Srebro and Shraibman 2005); Foygel and Srebro (2011) for results on estimating M (under milder conditions); note that in the latter case, it common instead refer to predicting the entries of M, e.g., temperature or stock market returns, since we often think of these as noisy observations of some underlying signal. On the other hand, relatively little is known about the problem of uncertainty quantification for these estimates: given some estimator M, can we produce confidence interval around each Mij that is (asymptotically) guaranteed to contain the target with some minimum probability?
The results Chen et al. (019)) address this question under strong assumptions, namely, low-rank and incoherent signal M plus i.i.d. Gaussian noise (t.)))dd Zhao and Udell (2020)
proposes another inferential procedure under the Gaussian copula assumption for the data. In this work, we provide a complementary answer that instead addresses the model-free free setting: without relying on assumptions, we aim to produce a provably valid confidence interval for the entries My of M (not of M', because without assuming a low-rank model, there no meaningful notion of an "underlying signal").
To do so, we will treat the matrix M as deterministic, while the randomness arises purely from the random subset of entries that are observed. More specifically, we assume
Zi j) is observed] Bern(pis); independently for all (i, [di] [d2],
(I)
and let S = ji. [d1] x [d2] Zi, = 1} be the subset of observed locations. We consider the setting where the sample probabilities Piş" s are nonzero. After observing Ms (s)))s
our goal is to provide confidence intervals ((,,))) for all unobserved entries, with coverage over the unobserved portion of the matrix -tht is, we would like our confidence intervals to satisfy
(ngCov(C;M,s)) 2 1 - ax,
where
AvgCov(C:M,S) [8c Euness
MM C(G,))}
measures the "average coverage rate" of the confidence intervals ! an algorithm for constructing confidence intervals (C(i,)))o.ees based on the observed entries Mg. satisfying (2 with no assumptions on M. In particular, given the strong success of various matrix completion algorithms for the problem of estimation, we would like to be able to provide uncertainty quantification around any base astimator-i.ce given any algorithm for producing an estimate M, our method constructs confidence intervals (Õ(t, ((k.))) (r.))ss around this initial estimate.
The goal of this work is to provide
Notation. For matrix A (Asj) E T Xm22, we define || || = (ax [mi] [xmm]] ||||| = (m)))((m) x [ma] |Ayl- We use A., and A to refer to the ith row and jth column of A, respectively. For any t E R, & denotes the distribution given by a point mass at t.
Aijl. and
3 Conformalized matrix completion
We next present our method, conformalized matrix completion (cmc). Our procedure adapts the split conformal prediction method (Vovk et al. 2005) to the problem at hand. As is standard in the conformal prediction framework, the goal of cmc is to provide uncertainty quantification (i.e., confidence intervals) around the output of any existing estimation procedure that the analyst chooses to use use-thi this is a core strength of the conformal methodology, as it allows the analyst to use of the-art estimation procedures without compromising the validity of inference. At high level, the cmc procedure will output a confidence interval of the form
((i.j) My + -
for the matrix value Myy at each unobserved entry j.. Here, after splitting the observed entries into a training set and calibration set, the training set is used to produce Maj (a point estimate for the
'To handle the degenerate case where S [di] [da] (i.e., happen to have observed the entire matrix) and thus S 2, simply detine AvgCov(C; M, = 1.
3

---

target value May) and Sty (a scaling parameter that estimates our uncertainty for this point estimate); then, the calibration set used tune the scalar g, to adjust the width of these confidence intervals and ensure coverage at the desired level 1
Algorithm 1 Conformalized matrix completion (cmc)
1: Input: target coverage level o; data splitting proportion q (0, 1); observed entries Ms.
2: Split the data: draw Wij lid. Bern(y), and define training and calibration sets
Ser = ((t S: Way = 11, and Scal = {(..) S Wey =0}.
3: Using the training data Ms.. indexed by Su S [d1] x [d2], compute:
An initial estimate M using any matrix completion algorithm (with Mi estimating the target Mty);
Optionally, a local uncertainty estimate (with Big estimating our relative uncertainty in the estimate Mty). ). or otherwise set Pij = 1;
An estimate P of the observation probabilities (with Pij estimating Pij. the prohability of entry (i. being observed).
4: Compute normalized residuals on the calibration set, Raj MM, 81g
(i,j) E Scal-
5: Compute estimated odds ratios for the calibration set and test set, Tas and then compute weights for the calibration set and test point,
1-Pu,
(i, Scal US
Piş
max Trg (i-))S-
Tij
Wij heg max hez (L (F
(1,5) E Seal.
Weest
(4)
Σ hey + max hey
(0 (')e)ss-
6: Compute threshold q Quantile (E Wij 8R., + Brst )
7: Output: confidence intervals C(i,j) Mij + ② Big for each unobserved entry (i.j) S
3.1 Exchangeability and weighted exchangeability
We split the set of observed indices S into training and calibration set S Sur U as is shown in Algorithm|| [ We should notice that the sampling without replacement may introduce implicit dependence as well as distribution shift from the i.i.d. sampling. Thus the two sets are dependent, which is different from the split conformal methods in regression problems. Before we present the coverage guarantees for our method, we first examine the role of (weighted) exchangeability this method, to build intuition for how the method is constructed.
3.1.1 Intuition: the uniform sampling case
where the entries are sampled with constant probability, Piy = p. If this the case, then the set of calibration samples and the test point are exchangeable- that is, for j+. denoting the test point that is drawn uniformly from S if we told that the
First, for intuition, consider the simple
combined set Scal U {(...--}} is equal to {(11,11)..... (nna +1)} in no particular order
(where neal (Saai)), () then the test point location (f.,j.) is equally likely to be at any one of these Peat + locations. Consequently, by exchangeability, the calibration residuals {Rij (i,j) E Seat). defined in I
exchangeable with the test residual R4,j. [M M,.j. 1/ as a result, we can construct our confidence interval (iv.. with q determined by the (1 a:)- quantile of the calibration residuals (Rij (i, E Scat).
2 More generally, eme can be implemented with any choice a nonconformity score that determines the shape of the resulting contidence intervals and indeed, this generalization allows for categorical rather than valuud matrix entries Myj. We will address this extension in the Supplementary Material, as well as version of cmc based on full conformal rather than split conformal, in order avoid data splitting.

---

Indeed, this is exactly what cme does in this case: if we are aware that sampling is uniform, it would naturally follow that we estimate Pij = Po by some (potentially data-dependent) scalar Poi consequently, all the weights are given by Miy = I and Ties Pea+t and thus q is the (unweighted) (1 qu-anttllle of the calibration residuals (with small correction term, ie., the term + Wust 0+0c appearing in the definition of to ensure the correct probability of coverage at finite
sample sizes).
3.1.2 Weighted exchangeability for the matrix completion problem
Next, we move to the general case, where now the sampling may be highly nonuniform. First suppose the sampling probabilities Pij are known, and suppose again that we are told that the combined set Sent L [+)) } is equal to (i1:)1),- ..... (fnar +noat +1}} in no particular order. In this case, the test point [+) is not equally likely to at any one of these Real + locations. Instead, we have the following result:
Lemma 3.1. Following the notation defined above, if (le:5a) IS Unif(5"), holds that (((...i)) = (ik, Seal U ((bnn1)) = {111:11..-- (Ineatl, Jnca-1) 1)}, Str = Wixiso
where we define the weights Wiefe
Paggl
for odds ratios given by haj pis
1-Psi_
Evell
heweow
Consequently, while the test point's residual Rigg. = [M,, - Mi 1 // Si. is not exchangeable with the calibration set residuals {Ri (i, E Scalh. the framework of weighted exchangeability (Tibshirani et al. 2019) tells us that places weight Wi on each residual Ricse and in particular,
view R.
draw from the weighted distribution that
MiM.. E Mi.j + F... 81.5.) = S.....
where J. Quantile,-al (E Ebaat+
Wigji 5R ). Noting that this threshold q. may depend on the test location (a): ). to accelerate the computation of prediction sets for all missing entries, we instead propose one-shot weighted conformal approach by upper bounding the threshold
uniformly over all test points: Ph.s- S q² = Quantile (Eiy)es wig 5ney + Wiest 8+00
for (i,j) E Seal.
max hij (i,5)ES=
hij
wi
Wi
(5)
hey + max heg (SS
lest
Σ
hey + max hey ('"
Now the threshold 9"' no longer depends on the location (1+, ) of the test poin If the true probabilities Pij were known, then, we could define "oracle" conformal confidence intervals
(i,j))) My + Rey, E..
As we will see in the proofs below, weighted exchangeability ensures that, if we do have oracle knowledge of the Ptj ss then these confidence intervals satisfy the goal (2) of 1 a average coverage. Since the Piy"s are not known in general, our algorithm simply replaces the oracle weights w" in with their estimates w defined in (4). using Pij as plug-in estimate of the unknown Pij.
3.2 Theoretical guarantee
In the setting where the Piy's are not known but instead are estimated, our cmc procedure simply repeats the procedure described above but with weights calculated using P.y in place of p.j Our theoretical results therefore need to account for the errors in our estimates P.j. to quantify the coverage gap of conformal prediction with the presence of these estimation errors. To do this, define an estimation gap term
A
Σ
Teg
hig
(6)
2
Ee Ses
hey
E E(v hag
((.)))))
'Comparison between the one-shot weighted approach and the exact weighted approach (Algorithm conducted in Section B.]
is
5

---

where (++ denotes the test point. Effectively, A is quantifying the difference between the "oracle" weights, w defined in D and the estimated weights, Wij. defined in ( except that A is defined relative to a specific test point, while the weights in (and wr", for the oracle version) provide "one- shot" procedure that is universal across all possible test points, and is thus slightly more conservative.
Theorem 3.2. Let M, S. and P be estimates constructed using any algorithms, which depend on the data only through the training samples Ms., at locations Str. Then, under the notation and definitions above, conformalized matrix completion ( (cmc) satisfies
E [[AvgCov(C: S)) 2
E(A).
In the homogeneous case, where we are aware that sampling is uniform (i.e., Péj p for some potentially unknown E (0, 1)). our error in estimating the weights is given by A = 0 (since hij constant across all (8,5), and same for and therefore, the true and estimated weights are all equal to FT In this setting, therefore, we would achieve the exact coverage goal 2 with
E AvgCCv (C; S) guaranteed to be at least
3.3 Examples for missingness
To characterize the coverage gap explicitly, we present the following concrete examples for and show how the coverage gap A can be controlled. Technical details are shown in the Supplementary Material.
3.3.1 Logistic missingness
Suppose that the missingness follows a logistic model, with log(Pi// (1 = log (h:z)
U4 + Uy for some E Rd and v E Rd where we assume u 0 for identifiability. This model closely related to logistic regression with a diverging number of covariates (Portnoy) 1988 He and
Shao 2000 Wang 2011). Following Chen et al. (2023), we estimate u. via maximum likelihood.
ú arg max L(u, v) u,v
subject to 1 1T =
Here the log- log-likelihood defined by
C(u, v) =
Σ log (1 + exp( - v3)) + +1 (1))ss; log (1 exp( u,- - 02))}.
(kd))(d)))())
Consequently, we define the estimate of Piy as P.j (1 + exp(-u, 6,)) for which one can
show that [[ S V
(etimax(dy.dd) The proof of this upper bound shown in Section A.3.1
min(dd (d))
3.3.2 Missingness with a general link function
More broadly, we can consider a general link function, where we assume log(p:;/(1 - Ptj)) = - loo (hez) o( Aiz). where & is link function and A Ihuh xda is the model parameter. The logistic model we introduced above corresponds to the case when
the identity function and A u1 + 1 is rank-2 matrix. In this general setup, if rank(A (A) k and || S т, we can apply one-bit matrix completion (Davenport et al. 2014) to estimate {piz]. More precisely, we define the log-likelihood
Ls., (B)
Σ
log(H(B.;)) +
Σ (K.))s)
log(1 (Ba)))
(i.)ES.
where (t()) = q1 e (((()))= (here rescaling by the constant g is necessary to account for subsam- pling the training set Sür from the observed data indices S), and solve the following optimization problem
À = arg max Ls,, (B)
subject to
B < kV/ktddd.
|B ST,
B
where |B|| is the nuclear norm of B. Consequently, we let R.s exp(-oo 0 A. )) which leads to
[|A] < min(d1, d2) -1/4, The proof of this upper bound is shown in Section A.3.2

---

4 Simulation studies
In this section, we conduct numerical experiments to verify the coverage guarantee of conformalized matrix completion (cmc) using both synthetic and real datasets. As the validity of cmc is independent of the choice of base algorithm, we choose alternating least squares (als) (Jain et al. 2013) as the base algorithm (results using convex relaxation (Fazel et al. 2004; Candes and Recht 20 2 Chen et al. 2020) are shown in the Supplementary Material). We use cmc-als to refer to this combination. We use AvgCov(C) (the average coverage, defined in 33), and average confidence interval length, AvgLength(C) = T Eig)es length (en.th((iiji)) to evaluate the performance of different methods. All the results in this section can be replicated with the code available at https //github. com.uujjoommmojiio.a erry/ conf mm
4.1 Synthetic dataset
Throughout the experiments, we set the true rank r* the desired coverage rate 1 0.90, and report the average results over 100 random trials.
Data generation. In the synthetic setting, we generate the data matrix by M M + E, where M' is a rank r" matrix, and E is a noise matrix (with distribution specified below). The low-rank component M is generated by M xU*V*T, where U* E Rdixx" is the orthonormal basis of random d1 x r" matrix consisting of i.i.d. entries from certain distribution Pu.e (specified below) and V Rd is independently generated in the same manner. In addition, the constant is chosen such that the average magnitude of the entries Mayl is 2. Following the observation model (T. we only observe the entries in the random set S.
Implementation. First, when implementing als, we adopt r as the hypothesized rank of the underlying matrix. For model-based methods (Chen et al. 2019), by the asymptotic distributional
characterization My Mi 4 NO0 02)) and Eig 1 (Mij Mt,). one has Maj - Mi S (0, 02 + 2²). Consequently, the confidence interval for the model-based methods can be constructed as
My + 91-0/2
+ 62
where qs the 3th quantile of (V,0 The estimates are obtained via /2 | Ms., Mss.111/1S.. /1// and 07 = (b2 U,, |I )!) with the SVD M = TΣVT, For amc-als, we obtain M, S from Ms. by running als on this training set (and taking 8.j = (0) + 1/)1 The probabilities Piy's 's are estimated via Pij- = (did2g)-21SS) for this setting where the Pıj" 's are homogeneous.
4.1.1 Homogeneous missingness
In the homogeneous case, we consider the following four settings with d, = d2 = 500:
Setting 1: large sample size + Gaussian noise. p = 0.8, Pu,v = =N 1). and Esg ~ 1)..1). Setting 2: small sample size + Gaussian noise. 0.2, P (0, 1), and Eij ~ (0,
Setting 3: large sample size + heavy-tailed noise. 0.8, Pa,. = NO, 1), and Eig 0.2 11.2.
where l, denotes the distribution with degrees of freedom.
Setting 4: violation of incoherence. p 0.8, Pu,e l1.2, and Fig NO,
Figure [ displays the results (i.c., the coverage rate and the interval length) for both cmc-al als and als when we vary the hypothesized rank from 2 to 40. The oracle length is the difference between the (1 a/2 a/2)-th and (a/2)-th quantiles of the underlying distribution for Eij. As can be seen, alc achieves nearly exact coverage regardless of the choice of across all four settings. When the hypothesized rank r is chosen well, cmc-al als often able to achieve an average length similar to that of the oracle.
For the model-based method als, when we underestimate the true rank (r <T"), we may still see adequate coverage since ē tends to over-estimate by including the remaining r" factors noise. However, if we overestimate the rank (r > r*), then als tends to overfit the noise with additional factors and under-estimate regardless of the choice of leading to significant undercoverage
7

---

1.0
15
1.0
Gmc als
cnc ab
AvgCov
0.9
Augl 10
AugCov
Avgl Lenggghh
880
0.8
025
0.7
40
30
10 15 hypothesized rank
typathesized mikk
hypothesized rank
hypothesized rank
(a) Setting large sample size + Gaussian noise
(b) Setting small sample size
Gaussian noise
Augdov
10
AugCov
Augl angth
40
olads
0.8
Augl
20
30
40
nypotheslzed rank
hopothestzed rank
hypothesized rank
hypothesized rank
(c) Sctting 3: large sample size heavy-tailed noise
(d) Setting 4: violation of incoherence
Figure 1: Comparison between cmc-als and als.
Ban
pinoot socres
than
11
11
10
00
-25
(a) Large sample size
Gaussian (b) Large sample size + Gaussian
(c) Violation of incoherence: r
noise: r"
noise: r- 8 < 50 =
6 < Y*
Figure 2: Histogram standardized scores for als and prediction lower and upper bounds for 50 distinct unobserved entries.
in Settings 1, 2 and 4. In Setting with heavy-tailed noise, als tends to be conservative due to overestimating σ. Overall, Figure confirms confres our expectation that the validity of the model-based estimates relies crucially on a well specified model, and it fails to hold when sample size is small (cf. Figure [p), when noise is heavy tailed (cf. Figure Te), and when the underlying matrix is not incoherent (cf. Figure Td).
In Figure we present the histogram of standardized scores (M.y - Ma)// 02 +22 and the plot of the upper and lower bounds for three settings. In Figure 2a when the model assumptions are met and r r", the scores match well with the standard Gaussian and the prediction bounds produced by als and cmc-als als are similar. With the same data generating process, when the rank is overparametrized, the distribution of scores cannot be captured by the standard Gaussian, thus the quantiles are misspecified. As we can see from the confidence intervals, als tends to have smaller intervals which lead to the undercoverage. incoherent. When the rank is underestimated, the r* - r factors will be captured by the noise term and the high heterogeneity the intervals by als are much larger while the conformalized intervals are more adaptive to the magnitude of entries.
the last setting, the underlying matrix is no longer
the entries will further lead to overestimated noise level. As a result,
4.1.2 Heterogeneous missingness
Now we move on to the case with heterogeneous missingness, where the observed entries are no longer sampled uniformly. To simulate this setting, we draw {ax i i di, S k²} i.i.d. from (n, 1) and (biz I < k. S d2) i.i.d. from Unif(-0.5 0.5), and define sampling probabilities by
log(pt/(11 - pt;)) = -i aabi--
8

---

1.00
1.00
and--ae creah
AugCov
0.80
Avg Lenghh
AegCov
0.98
Auglicoo LegLonges
10
0.90
10 15
20
30
35
40
5
10
15 20
25
30 35
40
10
15
20
25
30
35
30 35 Pypothesized rank
10 15 20
25
hypothesized rank
typothesized rank
hypothesized rank
(a) Gaussian noise: k-
1.
(b) Gaussian noise: A
5.
Augl Langth
NogLength
8
cncc-ata crreab
AugCov
0.20
AvgCor
0.05
0.90
0.90
0.85
0.55
20
25
hypothesized rank
hypothesized rank
hyppdhesized
Hypotheeized rank
(c) Adversarial heterogeneous noise: k- 1.
(d) Adversarial heterogeneous noise: k" 5.
1.00
AugCov
AvgLongul Lnnth
AwgCav
AugLo tggLength
20
0.90
20
hypothesized rank
hypothesized mank
hypotesized rank
hypothesized rank
(e) Random heterogeneous noise:
(f) Random heterogeneous noise: k 5.
Figure 3: Comparison between with true and estimated weights under heterogeneous missingness.
We consider two settings with 1-* = and k- 5, respectively. In both settings, Pij is constructed by estimating ãy and bi, via constrained maximum likelihood estimation as is shown Example 3.3.2
To analyze the effect of estimation error predicted by Theorem 3.2 we consider three matrix generation processes for each choice of k", where M* is generated in the same way as in Section 4.1.1 with d1 = d2 = 500, Pa,y = (O,
Gaussian homogeneous noise: Eig NN(0, 1) independently.
Adversarial heterogeneous noise: Egj NN0 ofyy) independently, where we take dig 1 1/2p.-- For this setting, note that the high-noise entries (i.e., those entries that are hardest to predict) occur in locations (i, that are least likely to be observed during training.
Random heterogeneous noise: Eij o(0, independently, where (.)))11))) is drawn from the same distribution as (1))
We write cmc -als to denote the conformalized method (with als as the base algorithm) where we use oracle knowledge of the true Piş' P for weighting, while cmc-als uses estimates Pij. The results are shown in Figures 3 for the settings k- and k* both homogeneous noise and random heterogeneous noise, both cmc -als and cmc-als achieve the correct coverage level, with nearly identical interval length. For adversarial heterogeneous noise, on the other hand, cmc" achieves the correct coverage level as guaranteed by the theory, but cmc-al: ale shows some undercoverage due to the errors in the estimates Pij (since now the variance of the noise is adversarially aligned with these errors); nonetheless, the undercoverage is mild. In Section B.2 in the appendix, the local coverage of cmc-als conditioning Pij = Po for varying Po is presented.
(i.e., the rank of P). Under
4.2 Real data application
We also compare conformalized approaches with model- based approaches using the Rossmann sales datase(") (Farias et al., 2022). This real data provides the underlying fixed matrix M; the missingness pattern is generated artificially, as detailed below. We focus on daily sales (the unit is 1K dollar) of 1115 drug stores workdays from Jan 2013 to July 31, 2015 and hence the underlying matrix has dimension 1115 x 780. The hypothesized rank r is varied from 5 to 60 with step size 5 and we set
attps://wwet http://wwwwglll..oo//////ssaaaaa

---

1.00
1.00
0.66
0112-83
0.65
Ave Lngggn
0.00
0.05
10 15 23 23 25 25 40 45 40 35
10 15 20 25 50 35 40
10 20 50 45 55.60 hypothesized rank
10 19 20 15 40 hypothesized rank
hypothesized rans
hypothesized rank
(a) als vs cmc-als: homogeneous missingness.
(b) als vs cmc-als: logistic missingness.
Figure 4: Comparison between conformalized and model-based matrix completion with sales dataset.
the target level to be 1- a = 0.9. We apply random masking for 100 independent trials and report the average for AvgCov and AvgLength in Figure 4
Homogencous Pij. For each entry, we first apply random masking Zij Bern inde- pendently, In Figure 4al conformalized approach has exact coverage 1 but als tends to be conservative when is small and loses coverage when increases to 50.
Heterogeneous Pij- After drawing piy's from log(pi//(1 ptj)) log (huj) ((Ary) in Section 4.1.2 with = P is obtained via the constrained maximum likelihood estimator in Example 3.3.2 In Figure 4b cmc-als has coverage around and in comparison, when is greater than 40, als fails to guarantee the coverage and AvgCov decreases to 0.75.
Discussion
In this paper. we introduce the conformalized matrix completion method cmc, which offers a finite- sample coverage guarantee for all missing entries average and requires assumptions on the underlying matrix or any specification of the algorithm adopted, relying only on an estimate of the entrywise sampling probabilities. Given an estimate of the entrywise sampling probability, provide an upper bound for the coverage gap and give the explicit form with examples of the logistic well as the general one-bit low-rank missingness. In the implementation, an efficient one-shot weighted conformal approach is proposed with the provable guarantee and achieves nearly exact coverage.
We can compare our findings to few related results in the literature. The work of Chernozhukoy al. (202 applies conformal prediction counterfactual and synthetic control methods, where they include matrix completion as an example with regression-type formulation. However, their approach relies on test statistic that is function of estimated residuals. Consequently, their method requires the assumption of stationarity and weak dependence of errors. Furthermore, the validity of their approach is contingent upon the estimator being either consistent stable. Additionally, Wieczorek/ (2023) studies conformal prediction with samples from deterministic and finite population, but the validity under the sampling without replacement remains an open question in their work.
Our results suggest number of questions for further inquiry. the setting of matrix completion, while the results are assumption-fiee with regard to the matrix M itself, estimating the sampling probabilities that determine which entries are observed remains a key step in the procedure; while our empirical results suggest the method is robust to estimation error in this step, studying the robustness properties more formally is an important open question. More broadly, our algorithm provides example of how the conformal prediction framework can be applied in setting that very different from the usual i.i.d. sampling regime, and may lend insights for how to develop conformal methodologies in other applications with non-i.i.d. sampling structure.
Acknowledgement
R.F.B. was partially supported by the Office of Naval Research via grant N00014-20-1-2337, and by the National Science Foundation via grant DMS-2023109, C.M. was partially supported by the National Science Foundation via grant DMS-23
10

---

References
Amjad, M., Shah, D., and Shen, D. (2018). Robust synthetic control. The Journal of Machine Learning Research, 19(1):02--822
Angelopoulos, A., Bates, S., Malik, J., and Jordan, M. (2020). Uncertainty sets for image classifiers using conformal prediction. arXiv preprint arXiv:2009. 14193.
Athey. S., Bayati. M., Doudchenko, N., Imbens, G., and Khosravi, K. (2021). Matrix comple- tion methods for causal panel data models. Journal of the American Statistical Association, 116(536)1176--733)
Bai, J. and Ng. S. (2021). Matrix completion, counterfactuals, and factor analysis of missing data. Journal of the American Statistical Association, 1 1166536)
Barber, R. F., Candès, E. J., Ramdas, A., and Tibshirani, J. (2021). Predictive inference with the jackknife+. The Annals of Statistics, 49(())486-335 507.
Barber, R. F. Candes, E. J., Ramdas, A., and Tihshirani, R. (2022). Conformal prediction beyond exchangeability. arXiv preprint arXiv:2202. 13415.
Cai, T. T., Liang, T., and Rakhlin, A. (2016). Geometric inference for general high-dimensional linear inverse problems. The Annals of Statistics, 44(4):1536 1563.
Candes, E. and Recht, B. (2012). Exact matrix completion via convex optimization. Communications of the ACM, 55(6):111-119)
Candès, E. J., Lei, L., and Ren, Z. (2021). Conformalized survival analysis. arXiv preprint arXiv:2103.09763.
Candes, J. and Plan, Y. (2010). Matrix completion with noise. Proceedings of the IEEE, 98(6):925- 936.
Candès, E. J. and Tao, T. (2010). The power of convex relaxation: Near-optimal matrix completion. IEEE Transactions on Information Theory, (155):053-2080
Cao, Y. and Xie, Y. (2015). Categorical matrix completion. In 2015 IEEE 6th International Workshop on Computational Advances in Mulii-Sensor Adaptive Processing (CAMSAP), pages 369-372. IEEE.
Carpentier, A., Klopp. O., Löffler, M., and Nickl, R. (2018). Adaptive confidence sets for matrix completion. Bernoulli, 24(4A):2429 2460..
Chen, W., Chun, K.-J,, and Barber, R. F. (2018). Discretized conformal prediction for efficient distribution-free inference. Stat, 7(1):e173.
Chen, W. Wang, Z.. Ha, W., and Barber, R. F. (2016). Trimmed conformal prediction for high- dimensional models. arXiv preprint 1rXiv:1611.0993.
Chen, Y., Chi, Y., Fan, J., Ma, C., and Yan, Y. (2020). Noisy matrix completion: Understanding statis- tical guarantees for convex relaxation via nonconvex optimization. SIAM journal on optimization, 30(4):3098-313
Chen, Y., Fan, J., Ma, C., and Yan, Y. (2019). Inference and uncertainty quantification for noisy matrix completion. Proceedings of the National Academy of Sciences, 116(46)
Chen, Y., Li, C., Ouyang, J., and Xu, G. (2023). Statistical inference for noisy incomplete binary matrix. Journal of Machine Learning Research, 24(95):
Chernozhukov, V., Wüthrich, K., and Zhu, Y. (2021). An exact and robust conformal inference method for counterfactual and synthetic controls. Journal of the American Statistical Association, 116:1849 1864.
Davenport, M. A., Plan, Y., Van Den Berg, E., and Wootters, M. (2014). 1-bit matrix completion. Information and Inference: A Journal of the IMA, 183):189-222

---

Farias, V., L.. A. A.A.A. A., and Peng. T. (2022). Uncertainty quantification for low-rank matrix comple- tion with heterogeneous and sub-exponential noise. In International Conference on Artificial Intelligence and Statistics, pages 1179-1189 PMLR.
Fazel, M., Hindi, H., and Boyd, S. (2004). Rank minimization and applications in system theory. In Proceedings of the 2004 American control conference, volume 4, pages 3273-3278, IEFE.
Foygel, R. and Srebro, N. (2011). Concentration-bass based guarantees for low-rank matrix reconstruction. In Kakade, S. M. and von Luxburg. U., editors, Proceedings of the 24th Annual Conference on Learning Theory, volume 19 of Proceedings of Machine Learning Research, pages 315-340, Budapest, Hungary. PMLR.
He, X. and Shao, Q.-M. (2000). On parameters of increasing dimensions. Journal of multivariate analysis, 73(1):120-135.
Jain, P., Netrapalli, P., and Sanghavi, S. (2013). Low-rank matrix completion using alternating minimization. In Proceedings of the forty-fifth annual ACM symposium Theory of computing. pages 665-674.
Jin, Y., Ren, Z., and Candès, E. J. (2023). Sensitivity analysis of individual treatment effects: A robust conformal inference approach. Proceedings of the National Academy of Sciences, (016)):2211489111
Keshavan, R. H., Montanari, A., and Oh, S. (2010). Matrix complction from a few entries. IEEE transactions on information theory, 56(6): -9988.
Koltchinskii, V., Lounici, K., and Tsybakov, A. B. (2011). Nuclear-norm penalization and optimal rates for noisy low-rank matrix completion. The Annals of Statistics, 39(5::2302 2329.
Lei, J., G'Sell, M., Rinaldo, A., Tibshirani, R. J., and Wasserman, L. (2018). Distribution-fres predictive inference for regression. Journal the American Statistical Association, 113(523):0094- 1111.
Lei, L. and Candès, E. (2020). Conformal inference counterfactuals and individual treatment effects. arXiv preprint 33X1v:2006.06138.
Ma, C., Wang, K., Chi, Y., and Chen, Y. (2020). Implicit regularization nonconvex statistical estimation: Gradient descent converges linearly for phase retrieval, matrix completion, and blind deconvolution. Foundations of Computational Mathematics, 20:451-632.
Nair, J., Wierman, A., and Zwart, B. (2022). The fundamentals of heavy tails: Properties, emergence, and estimation, volume 53. Cambridge University Press.
Natarajan, N. and Dhillon, I. S. (2014). Inductive matrix completion for predicting gene disease associations. Bioinformatics, 30(1 (0:)6)66
Negahban, S. and Wainwright, M. J. (2012). Restricted strong convexity and weighted matrix completion: Optimal bounds with noise. The Journal of Machine Learning Research, 13(1):1665. 1697.
Portnoy, S. (1988). Asymptotic behavior of likelihood methods for exponential families when the number of parameters tends infinity. The Annals of Statistics, pages 356-366.
Rennie, J. D. and Srebro, N. (2005). Fast maximum margin matrix factorization for collaborative prediction. In Proceedings of the 22nd international conference on Machine learning. pages 713-719.
Romano, Y., Sesia, M., and Candes, E. (2020). Classification with valid and adaptive coverage. Advances in Neural Information Processing Systems, 33:3581-3591.
Shafer, G. and Vovk, V. (2008). A tutorial on conformal prediction. Journal of Machine Learning Research, 9(3).
12

---

Srebro, N. and Shraibman, A. (2005). Rank, trace-norm and max-norm In Learning Theory: 18th Annual Conference on Learning Theory, COLT 2005, Bertinoro, Italy, June 27-30, 2005. Proceedinggs 18, pages 545-560. Springer.
Sun, R. and Luo, Z.-Q. (2016). Guaranteed matrix completion via non-convex factorization. IEEE Transactions on Information Theory, 62 (6(11)66535-33 6579.
Tibshirani, R. J., Barber, R. F., Candès, E. J., and Ramdas, A. (2019). Conformal prediction under covariate shift. In NeurIPS.
Vovk, V. (2015). Cross-conformal predictors. Annals of Mathematics and Artificial Intelligence, 74:9-28.
Vovk, V., Gammerman, A., and Shafer, G. (2005). Algorthhiic learning in random world, volume 29. Springer.
Wang, L. (2011). Gee analysis of clustered binary data with diverging number of covariates. The Annals Statistics.
Wieczorek, J. (2023).
based conformal prediction. arXiv preprint arXiv:2303.01422.
Xia, D. and Yuan, M. (2021). Statistical inferences of linear forms for noisy matrix completion. Journal of the Royal Statistical Society Series B: Statistical Methodology, 83(1):58-77.
Yang. Y. and Ma, C. (2022). Optimal tuning-free free convex relaxation for noisy matrix completion. arXiv preprint arXiv:2207.05800
Zhao, and Udell, M. (2020). Matrix completion with quantified uncertainty through low rank gaussian copula. Advances in Neural Information Processing Systems, 33:20977-20988.
13

---

Contents
A Proofs
14
Proof of Lemma
14
Proof of Theorem 3.2
16
A.3
Proofs for examples in Section 3.3
17
A.3.1
Proof for Example
17
A.3.2
Proof of Example 3.3.2
17
A.3.3
Proof of Lemma A.
18
Additional numerical experiments
18
B.
Comparison between one-shot and exact weighted approaches
19
B.2
Local coverage of cmc
19
B.3
Additional results for homogeneous misingness
19
B 4
Estimation error for one-bit bit matrix estimation
19
B,5
Additional results for heterogeneous missingness
20
B.6
Misspecified sampling probability
20
Additional results for the sales dataset
22
B.7.1
Comparison under misspecified missingness
22
B.7.2
Results with convex optimization
22
Additional details of algorithms and extensions
23
Extension to likelihood- scores and categorical matrices
23
C.2
Full conformalized matrix completion
24
Exact split conformalized matrix completion
25
A
Proofs
The section collects the proofs for the technical results in the main text. Section A.1 and A.2
devoted to the proof of Lemma 3.1 and Theorem 3.2 respectively, In the end, Section A.3 provides proofs for the coverage gap bounds in the two examples in Section 3.3
A.1 Proof of Lemma
Fix a set of locations = {(1,, f1), (6nww+1) Imcat+1)) and an index S m S Real + 1. By the definition of conditional probability, we have
P ((i,, j-) = (im. Seal U f(1,, i-)} = B, Snr Su)
P (o)., = (Im Jm). Scal = B (mmr )} Sτr = So, [S] = n)
IP ((Scal U (1i,, [e)] = B Sr = So |8] = n)
P ((i.j..) (Im) Scal B (m((m)))m)) |Str = Sn, |S| = n)
El New ((1,. (11.1)), Scal B ( ji}} Str = So, [S] n
then boils down to computing TP (Seal St, (i,, = (ir, |tr So. |S| n) for any fixed set Sj and fixed location (ir.i) & So U S. To this end, we have the following claim (we defer its
14

---

proof to the end of this section)
P (SSa = S1, (1..ja) = (11 |ttr = So, = n) -u,
II I-Psi
Pii
(8)
did2 n ExEOS [(( 1-Ppy3
Per
where Us.... {AC S [1] x [d2] | = 1Sul: AnS 2). As result, we obtain
P ((+r,+)) = (imimm) Scal U (()r B, Sur = So)
Ia..)BB Psi
I-psj
E1s Meart ILaSEE (((r))) I-Pu5
Pss
(9)
Reatt1 hicf
where big (1 Pis)/ Pis)/Pig- This finishes the proof.
Proof of Equation 8 Denote W the matrix consisting of entries Way defined in Algorithm | The
matrix consists of Zig which are the indicators of non-missingnes, i... Zij = j{ is observed).
Fix any two disjoint subsets So S₁ [d1] [d2] with [So] + [S1] We have
PP (Seal = S1, Sür = So |S| n) =P (supp ()) So S. S supp( (W), supp w) [S| n).
Since W and Z independent, one further has
P (supp(Z) So SI S supp(W), S S supp [|| = n)
-P (Su S supp(W).ct S supp(u))) (supp(Z) So US S n)
| (1 - g) Sil (supp(Z) = Sp U |S| = m)
(1 (1 4) |Su]
II
proj)
II
Ptj
(10)
P[[[| n))
I- Pij
(V' (V dd]]
(iJ)
where the last two identities are based on direct computations.
Based on 10
further compute
(S So |S| = n))
Σ P (SScal = A, Str So
|S| = n.
AERS
q |ool
Σ
II
(1- -Ppg))
Pij
(11)
P([S n)
1 - Pij
(L.J)
where the last identity uses Equation (10).
Now we are ready to prove (8). Recall that the new data point (f.,j.) S is drawn uniformly from Uni Therefore has
P (Seal S1, (i.. (n+1))]++1) Su = So.15 = n)
P (Scall = S1 = Sy, (1r, ) (in+1:In+1) |S] n)
P (Str Su | 12)
((t.:j.) = (n+1) S = U So) (Scal = S1r = So |S| = n)
IP (St So
S =
Io,es I-Psi
Pi.
(12)
Puy di d2 ExeR I( I-PMy
where the last line follows from 10 and II
15

---

A.2
Proof of Theorem 3.2
First, fix any E 00. 1]. Lemma 3.1 together with the weighted conformal prediction framework of Tibshirani et al. 2219), implies that
P M Mr.s. + 2 (a) Be.s. ScrrSall U {(h.) 2 a.
where
J. (a) Quantile
Σ
Wij òR:s
and
haj
Wkg
(((iss))) hey
Indeed, here a can be any function the random variables we are conditioning on that is, a may depend on n., Sue, and on Seal U j.)..
Next define Su. then (and thus also a) can therefore be expressed as a function of n, Str. and Seal U [(),, Therefore, applying the work above, we have
+ A where A is defined in (6) We observe that, since each huj is function of
P (Mr.s. E M. + 2 (a + A) -). Str, and thus, after marginalizing,
U (((m)n
a
P (M... E Misja 1 Fi.j, (a A) Bi.d.) a - E|A].
Next, we verify that holds almost surely-if this is indeed the case, then have shown that
2 -i.j. (a A)
(Mi.s. E Mi... + - Bi.j.) 2 1 - a a[
which establishes the desired result. Thus we only need to show that q 2 F..j. (a+ or equivalently,
Quantile, -a
Beg JR. + Wrest 64000
2 1
Σ
Wij 5R.y
(t.))
((i))S)
Define
hzg
w/s
(13)
Eie Se hez
for all (i.j) E SUU SSU(((i, ((i,))))) Then by definition of t. we see that w'y 2 Tij for (i,j) E S, and therefore,
Quantile
Σ (k
wig
Sn. + Wjest 64
-Qaattile-a
Σ
w/y ōn. + u(B..)) 8400
Quantile -a
Σ
why RR.s
((.))
holds almost surely. Therefore it suffices to show that
Quantile
Σ
wy dr.g Quantile
Σ
wij 8n.
((
(HJ))SSa
holds almost surely. Indeed, we have
dTv Σ
w/y 55
Σ
Wey &
Σ
| A.
(L)
(H))ESen
where dTv denotes the total variation distance. This completes the proof.
16

---

A.3 Proofs for examples in Section 3.3
Recall the definition of A:
haj
haj
A 2 Σ
(U) Exe Pnes ((. hey
Eg (e hep
Note that A & by definition. We start with stating a useful lemma to bound the coverage gap A using the estimation error of haj- Lemma A.1. By the definition of the coverage gap A, we have
( Phis - hisl
(14)
Eieryes heg
A.3.1 Proof for Example 3.3.1
Under the logistic model, one has for any (f,
9 exexp( + uz)
P((i,)) E S..) = Pij 1 + exp(u, u, + vz))
In this case, if i and 8 are the constrained maximum likelihood estimators in Example 3.3.1 Theorem 6 in Chen al. (2023) implies that
u/100 = Or (V log di
|| vllo = =Or -Or log d2
(Vioosd
d2
d1
with the proviso that || + 1101100 T < x., d2 > vil log di and di
>>
(log d2)
Then for
hig exp(-u; vj)) and hij = exp(-iix - ūz), have
fog d2
log d2 d
max Thag heyl max e 1i (B: 1.
u5) - =Op
+
1.)
Further, as ming, bag = min,, exp( - u4 U1)) 2 - then with probability approaching one, for Rugl 2 22 ho- By the upper bound (14)
every (i. [di]] x [d2]. one has hi,i hisj we have
[hij
log di
log d2 d1
log max (d1, d2) min(d1. d2)
d2
Further, as A < 1, we have R(A] S
max (dy
min (d) ,4g
A.3.2 Proof of Example 3.3.2
Define the link function (t) q(1 + +e)(()), where & is monotonic. Applying Theorem in the paper (Davenport et al. 2014), we obtain that with probability at least Ci/(di + dz)
k(d1 + d2)
d | - AIE S V2C,
(15)
did2
with the proviso that dzd2 2 (d1 + d2) log(d1dz). Here Č, = 239/4 9/4//49////// (1 /6) (1
with
(t))
(1) +(t))
Ly = sup
and 3, = sup
(16)
-T (1) (())
(1)(1)22
-TSIST
Denote this high probability event be Eg. Since ||A||] /did2 ( this event Ep, further have
kk(1 + d2)
1/4
A
Al
C
d, d2
did2
17

---

where C, = 5v and ( is a universal constant.
Recall that huj exp(-<(( (Ai)) and || т. On the same event &g, further have
1 H
(k(d1 + d2)
1/4
H S C%
(17)
d,d2
d1d2
where C₁ = 2e V
By the feasibility of the minimizer A and the fact that h.s exp( (A(y)) have his ho e-((r) for all (1,j)- This together with the upper bound (14) implies that
A
E
Ihas R.g honeal
H 1-
hyflea
Define a second high probability event & = {meat 2 (1 c) bound, we have (E1) Co(d1d2) -5 Therefore, on the event & E1, we have
- g) PI
Using the Chernoff
d1d2
[k(d1 + d2)
1/4
A < cho(1 4) PPl
C²
(
did2
Under the assumptions that Pij 1/(1+ ee )(()))), and that for some constant that only depends on T. As result, have to the conclusion that
||A] Al S T, we know that P.y 2 C2
[P]
C2d1d2, which further leads
1
k(11 + d2)
1/4
AS cC2hol1 9) C2
d1d2
on the event Eg n E In addition, on the small probability event (Eg 72)) one trivially has A S 1. Therefore simple combinations of the cases yields the desired bound on E[A].
A.3.3 Proof of Lemma A.1
Reusing the definition of w/ 13 one has
Σ
[ei - w'yl
2
his Ee huy -hu (( doS) hegl
Σ
= 2
E(re)essmm- Tay)) (2 (Crr)))ss harr)
hig Eee Tey Lee heg
1
Σ
<
((b))))
Eirrreesmmu Ray))
(E bazr)
hig higl Eil E( huy hey)
E Eustessuut key)
Σ e2 ((st Key - EesnEsS herl
1 E (hgs. hirl Twy
=
+ 2
2
(VstEssa hwy)) Ea. (hij hig!
(18)
E(r Twy
This completes the proof.
B Additional numerical experiments
In this section, we provide additional simulation results. In Section B.T we compare the proposed one-shot weighted approach with the exact weighted conformal prediction. In Section B.2 the
18

---

0.82
4,0
Conditional coverage Po -o
0.92
3.8
Doverage
length
3,8
3.6
0.90
0.90
0.89
3.6
0.88
3.2
0.88
3.5
3.0
emc-og-ak- ome-1-als approach
cmo-og-ls cmo-1-al approach
0.0
0.2
0.4
0.8
0.8
1.0
0.0
0.2
0.4
0.8
0.8
1,0
Py=
Py
(a) Exact/original (og)
1-shot
(b) Local performance of coverage given Ptj po-
local coverage of cmc-als conditioning on Pij is evaluated. We present the simulation results for the convex relaxation (cvx) and the conformalized convex method (cmc-cvx) in both settings with homogeneous and heterogeneous missingness in Section B.3 and B.5 We further evaluate the
performance of cmc-als when the sampling probability is misspecified in Section B.6
B.1 Comparison between one-shot and exact weighted approaches
We conduct comparison in the setting with heterogeneous missingness and Gaussian noise as shown Section 1.2 with k* = 12. From the result in Figure 5a the performance of two aproaches are essentially the same while the one-shot cmc-als is more computationally efficient.
B.2 Local coverage of cmc
In Figure we evaluate the local performance of emc-als when conditioning on {pij po) i.. subpopulation determined by the value of sampling probability. We use the setting with heterogeneous missingness and Gaussian noise as shown in Section (.1.2 with k 12 and Po ranging from 0.1 to 0.9 with a step size of 0.1. The conditional coverage is approximated via kernel smoothing: with the indicators Aij for coverage, we use the weight Ki coverage by (E (.)) Kej). Here 6p.e is the density function of (N,,, ((), and h 0.05. When increases from 0.2 to 0.8, the conditional coverage increases and stays around 0.9. At Po 0.1, 0.9, the coverage varies much, which is due the small effective sample size for the edge values of Pij. Moreover, for uniform sampling, since the rows and columns are generated i.i.d., there are no meaningful subgroups of the data condition on.
(Pij) and calculate the conditional
B.3 Additional results for homogeneous misingness
In this section, we present the results for synthetic simulation with the convex relaxation cvx and the conformalized convex matrix completion method cmc- cvx. Setting 1, 2. 3, 4 are the same as introduced in Sectior4.1 The true rank is y* = 8 and the hypothesized rank varies from to 24 with the stepsize 4.
The conformalized methods, regardless the based algorithm adopted, have nearly exact coverage around
observe different behaviors between als and since the convex relaxation is free of the choice of r until the projection of Mevx onto the rank-r subspace (Chen et al. 2020). As result, when > r", tends to overestimate the strength of the noise. In Figure 6a 6b 6dd when > r", CVX has coverage rate higher than the target level and the confidence interval more conservative than conformalized methods. Since the accuracy of cvx is based the large sample size, in Figure 6b when the effective sample size insufficient with small Pıjs the residuals from have large deviation from the standard distribution and the intervals the oracle ones. Besides, in Figure bc) when the noise has heavy tails than Gaussian variables, cvx overestimates the noise strength similar to als and is conservative in coverage. When the incoherence condition is violated in Figure 6d] if r < r", both cvx and als the missed factor by overestimating the noise strength and produce extremely large intervals.
But
much larger than
B.4 Estimation error for one-bit matrix estimation
The estimation error in Pij can visualized from the following heatmaps comparing and P. Here the entries sorted by the order of Pty for each row and each column.
19

---

1.0
1.0
16
crre-as
dmc 2l
AvgCov
0.9
AvgCov
Augl eagghh
20
omco cix
omc
0.8
Augl
10
0.7
0.4
16
typothesized mink
hypothesized rank
hypathesized mank
hypothesized rank
(a) Setting I: Large sample size
Gaussian noise
(b) Setting 2: Small sample size+GGussiin noise
60
Augdov
80
AugCov
Aurg angth %
0/X
OX unccle
onedd
3
15
20 hypothesized rank
20 hypothesized rank
hypothestzed rank
hopothesized rank
(c) Setting 3: Large sample size hcavy-tailed noise
(d) Setting 4: Violation of incoherence
Figure Comparison between conformalized and model-based matrix completion approaches.
Tiue
00
Estinated
0.00
True
0.55
Extimoted
066
125
125
0.50
055
055
254
260
0.50
1000
0.4
045
376
0:4
05
250
372
300
125
253
375
500
020
1000
1500
2000
300
1300
150
2003
A="
1: dj
d2
500.
(b) k-
1: dj
d2 2000.
Trus
Estimatedd
True
Extimuled
105
300
500
260
250
1000
1000
04
34
34
372
1801
102
-2
125
200
275
500
125
250
275
200
323
1306
1968
2000
200
1306
1506
2000
(c) k
5: d1
d2 500.
(d) K 5: dr
d2 2000.
Figure 7: Heatmaps for P and P.
B.5 Additional results for heterogeneous missingness
In Figure 8 we present the results for synthetic simulation with the convex relaxation cvx as the base algorithm, where we denote emc-cvx and cmc* -cvx as the conformalized matrix completion method with estimated weights and true weights, respectively. Three settings with heterogeneous missingness are the same as Figure 3
B.6 Misspecified sampling probability
We consider the following four settings:
(a) The underlying missingness follows the rank-one model where Ptj azby, both a; and by are generated i...d. from Unif(0.2, 1). The noise is adversarial. But we estimate Pij via the one-bit matrix completion based the logistic model (working model with hypothesized rank = 5).
(b) The underlying missingness follows the logistic model in Example 3.3.2 with A 5 and we adopt the adversarial noise. But Pıy estimated under the assumption of uniform sampling (working model), i.c. Pij P.
20

---

1.00
1.00
10
crec-aa
AugCov
0.80
WygLengh
omi
BegCov
0.95
Leggeggh
ond cox
emoew
crro-ow
0.90
0.90
10
15
20
10 typothesized rank
20
10
15
20
20 Pypothesized rank
t3
hyppthesized rank
hypothesized rank
(a) Gaussian noise: k*
1.
(b) Gaussian noise: A
5.
AggCov
0.30
Augl Langth
crrc-aai ab
cmel-exs
0.05
LogLength
2
cnd -ls crel-exx
0.90
0.90
0.85
0.55
hypothesized rank
hypothesized rank
hypothesized rank
Hypotheeized rank
(c) Adversarial heterogeneous noise: k-
(d) Adversarial heterogeneous noise: k" 5.
1.00
AggCov
AvgLongli Lnnth
AugCov
TogLongth
2
0.90
0.85
hypothesized rank
hypothesized rank
hypothesized
hypothesized
(e) Random heterogeneous noise:
(f) Random heterogeneous noise: k- 5.
Figure 8: Comparison under heterogenous missingness.
1.00
1.00
cro--lb-Tra arokk
aro--la-mTa
10.0
mmosls
AvgCov
AvgLengli
Avv
0.90
0.90
5.0
0.80
0.80
0.0
hypotheeized rank
typotheeized rank
hepoeheelzed
hypotheeized
(a) Estimating rank-1 model with logistic model (ad- versarial noise).
(b) Estimating logistic model with uniform sampling (adversarial noisc).
1:00
25
croosberTex anoels
1.00
AvgLength
0.95
10.0
amooetss
AvgCCv
0.90
0.90
5.0
0.85
0.85
0.80
0.80
0.8
hypothesized rank
hypothesized rank
hypothesized
hypotheeized
(c) Estimating rank-1 model with logistic model (ran- dom noise).
(d) Estimating logistic model with uniform sampling (random noise).
(c) Same as (a) except that we use the random heterogeneous noise.
(d) Same as (b) except that we use the random heterogeneous noise.
From the results, we can see that when the noise is generated following the random heterogeneous model, where the values of entries are independent of the sampling probabilities, the misspecification of the sampling model only slightly affects the coverage. Moreover, when the noise generated in the adversarial way, where the values of entries depend a misspecified sampling model lower than the target level, but above 0.85 in practice, which depends the divergence between the true and the working sampling models.
Pij" is, we can see that the coverage with
21

---

B.7 Additional results for the sales dataset
Denote M the underlying matrix in the sales dataset. In Figure 10a we plot singular values of M and top-5 singular values contain large proportion of the information. In Figure [Ob] we plot the histogram of entries My's 's of the underlying matrix, and the sales dataset has the range from 0 to over 20 thousand with heavy tail.
10
0.15
10
10
10'
0.05
0.00
200
400
600
800
10
20
30
40
(a) Singular values of the underlying matrix M.
(b) Histogram of entries in M.
Figure 10: Descriptive plots for the underlying matrix M.
B.7.1 Comparison under misspecified missingness
We consider the sales dataset where the missingness occurs the way that: for weekdays, each entry is observed with 0.8. On weekends, as the stores are likely to be operated by less experienced interns are less frequent to report the sales data, each entry observed with a lower probability, e.g. 0.8/3. Moreover, as there could be subgroup of stores that are less frequent in sales reporting, 200 stores are randomly sampled, for which p 0.8/3. We use the logistic model with k 5 to estimate Pij. Comparison between als and cmc-ale is shown in Figure 11
B.7.2 Results with convex optimization
In Figure CmC-Cvx has nearly exact coverage at ax cvx tends to have higher coverage than the target level. Besides, the convex approach has much larger intervals when r is large, which can caused by the overfitting of the observed entries. conformalized approach leaves out proportion of observed entries as the training set, intervals produced by cmc-cvx are less accurate than cvx due to the poorly behaved base algorithm.
cmb-als
db
AveCov
0.90
AveLo
0.85
0.80
15
25
35
45
55
15
25
35
45
55
hypothesized rank
hypothesized rank
Figure 11: Sales dataset with heterogeneous and misspecified missingness.
22

---

1.90
1.00
20
0.95
355
0.35
Anglangthh
15
OFrO-CI CIX.
0.90
Avgl
OYX
0.80
0.75
0.75
25
25
25
hypothessie eeizeeeeeeeddddee rank
Iypothesized rank
hypothesized rank
hypothesized rank
(a) Homogeneous missingness.
(b) Heterogencous missingness with k"
Figure 12: Comparison between conformalized and model-based matrix completion with sales dataset.
C Additional details of algorithms and extensions
C.1 Extension to likelihood -based scores and categorical matrices
We will show in this section that cmc can also be applied to categorical matrix completion (Cao and Xic 2015) or more general setting in the validity of which is also guaranteed by the presented theorem.
Setup To formulate the problem, consider an underlying parameter matrix M' E the observations {My E [di].] [d2]} are drawn from the distribution
[d1] x [dg] and
Mij Mij PM;;'
(19)
where {Po can be a family of parametric distributions with probability density pe- The cat- egorical matrix completion is specific example where the support of My is finite or countable.
For example, a Poisson matrix is generated by Myy ~ Pois(Mi;), where Mīj is the Poisson mean.
Similar to the previous setup, we treat M as deterministic, and entries the subset S C [d1] x [d2] available. Here S is sampled in the same manner before with the matrix P [d1] [d2].
Split conformal approach Consider the split approach with the partition S S,, U Seat in Algorithm [ that
With the training set MEr, we obtain an estimated likelihood function ((m; i, such
A(m; Part (m),
which an estimate for the true likelihood of Mij at m. The estimation can be feasible given certain low-complexity structures. For example, if hypothesized distribution family (Qo)eee with probability density qo is given and the underlying mean matrix M is assumed to be low-rank. Then, M can be viewed perturbation of M and with entries Mt.. Denote M as the estimate for M then we have the estimated likelihood
estimate M via matrix completion algorithms
im;; i, (. 9M (m).
The odds ratios are also estimated from the training set, i.e. hijs from which we compute the weights
heg
mij
(i,5) Seal, Wrest
Σ key max Tug (E
Σ
Tuy + max hez
(R
For each (i,5) E Seal, calculate the likelihood-based nonconformity score
Rij - - {MMyi (.).
Then, for any test point (in,ja) S", construct the confidence interval
Ō(i..i)) {m [K] 7(m; j.) S 7
where à is the weighted quantile
g Quantile
Σ
Wej 8R. + Brea 8400
(KJ)
More examples of conformal methods for classification are shown in Romano et (2020), An- gelopoulos al. 2020), etc.
23

---

C.2 Full conformalized matrix completion
the procedure of the full conformalized matrix completion (full-cmc) is presented. This full conformal version of cmc offers the coverage guarantee as given in Theorem 3.2 for the split version of cmc (except with the entire observed set S in place of Seal, when defining the formal proof of this bound for full conformal very similar to the proof of Theorem 3.2 using an analogous weighted exchangeability argument as in Lemma 3.1 and so we omit it here.
In Algorithm2 2
To define this algorithm, we need some notation: given the observed data Ms, plus test point location (i.,j.) and hypothesized value m for the test point Mi.j. define matrix M(m) with entries
Mij-
(i,J) E
M
(m)
m,
(i.j)=(i+.33 J++.
(20)
otherwise,
where, abusing notation, "Mij = 0" denotes that no information is observed this entry.
Algorithm 2 full-cme: full conformalized matrix completion
1: Input: target level oz partially observed matrix Ms.
2: Using the training data Ms, compute an estimate P of the observation probabilities (with Pig estimating Pij. the probability of entry (i,j) being observed).
3: for (i., J+) in S do
4: for m E M do
5: Augment Ms with one additional hypothesized entry, {Mi.... = m.. to obtain M(m)
defined as in (20).
6:
Using the imputed matrix _ compute:
An initial estimate ((m) using any matrix completion algorithm (with M. ( estimating
the target M,y);
Optionally, a local uncertainty estimate g(m)) (with 8/) estimating our relative uncertainty
in the estimate M; or otherwise set <(m) = 1;
An estimate P of the observation probabilities (with Pig estimating Pij. the probability of entry (i.j) being observed).
7:
Compute normalized residuals for (i,j) E S U f(i.,.))
Mij - (
(m) Rij
((m)
8ij
8:
Compute weights
Tag
1- Dij
iij
hij =
(1.)) e SU(....))
Eve (ES hez
Dij
9: Compute the weighted quantile
(i = Quantile
Σ
Wind Tind + (ir... ò
(21)
(1,))S
10: end for
11: end for
Output: (i..)))) - { m M: R(T)) (inn)))) (i.,))) E
We note that, when M IR (or an infinite subset of R), the computation of the prediction set impossible in most of the cases. In that case, our algorithm can be modified via a trimmed or discretized approximation; these extensions presented for the regression setting in the work of

---

Chen et al. 2016 2018 and can be extended to the matrix completion setting in a straightforward
way.
C.3 Exact split conformalized matrix completion
In Algorithm we present the exact split approach, which is less conservative than our one-shot approach given in Algorithm [ but may be less computationally efficient. In this version of the algorithm, the quantile (()..=>) needs to be computed for each missing entry since the weight vector depends on the value of P.....
Algorithm 3 split-cme: split conformalized matrix completion
Input: target coverage level a; data splitting proportion g E (0. 1); observed entries Ms.
2: Split the data: draw Wij lid. Bern(g), and define training and calibration sets,
Sur = ((i,)) E Waj 1}, Seal {(i.)) E S: Wij 0}. 3: Using the training data Ms.. indexed by S. S [di] x [d2], compute:
An initial estimate M using any matrix completion algorithm (with Mi estimating the target Mij);
Optionally, a local uncertainty estimate (with Fiy estimating our relative uncertainty the estimate Mij). or otherwise set 8ij = 1;
An estimate P of the observation probabilities (with P.y estimating Pij. the probability of entry (i.j) being observed).
4: Compute normalized residuals on the calibration set,
Mi Maj (i,g) E Seal-
Raj-
8ig
5: Compute estimated odds ratios for the calibration set and set,
Teg
Dij
(i, Seal U
P.j"
6:
for (ix, E S do
7:
Compute weights for the calibration set and test point,
h Tey Pri.jo
Way
(1,5) Seal; Brest
Σ
Teg Teis
(V (Vid))ESe
8: Compute threshold
((e.j.)) Quantile
Wij 8R. + Trest dero0)
(
where ò, denotes the point mass at t.
9: end for
10:
Output: confidence intervals
Ö(i.)j.) Mr.5. ((insaa) Biag-
for each unobserved entry (ir,j.) S".
25


Reviews:
[{"review_id": "jeQeOIJECd", "reviewer": "Reviewer_bWVG", "summary": "This paper proposes to use conformal prediction for uncertainty quantification of matrix completion. The proposed conformalized matrix completion offers provable predictive coverage regardless of the accuracy of the low-rank model. Empirical results on simulated and real data demonstrate that cmc is robust to model misspecification while matching the performance of existing model-based methods when the model is correct. ", "questions": "Please refer to the \"Weaknesses\".\n", "limitations": "The limitations have not been adequately addressed in the main text. ", "rating": 5, "confidence": 2, "soundness": 3, "presentation": 2, "contribution": 2, "strengths": "Strength 1. This paper studies conformal prediction for matric completion with interesting theoretical and algorithmatical findings. \n1.1） Through conformal prediction, this paper proposes distribution-free confidence intervals for the completed entries in matrix completion. Benefiting from the use of conformal prediction, the validity is free of any assumption on the underlying matrixs and holds regardless of the choice of estimation algorithms. By proving the (weighted) exchangeabiity of unobserved and observed units when they are (non-uniformly) sampled without replacement from a finite population.  \n1.2） A provable lower bound for the coverage rate is provided when the sampling mechanism is unknown. \n1.3） A one-shot conformalized matrix completion approach is proposed for higher computational efficiency.\n\nStrength 2. Experiments on both synthetic and real data suggest the effectivenss of the proposal. ", "weaknesses": "Weakness 1: One weakness of the proposed uncertainty quantification method for model calibration is that the probability bounds it provides may not be very precise for specific observation models and optimization algorithms. While the method is advantageous in that it is not limited to any particular model or algorithm, there is a possibility that the uncertainty estimates can be too broad for certain scenarios. This potential issue should be acknowledged and addressed in the research.\n\nWeakness 2: Another weakness is that the proposed approach does not assume any specific structure on the underlying matrix. As a result, it may not effectively utilize any existing structure in the matrix. To overcome this limitation, the proposed method can be extended to vector completion problems, which could exploit the matrix structure more efficiently. The authors should consider discussing this point and its potential implications.\n\nWeakness 3: The paper could benefit from improved clarity in its explanations. It would be more reader-friendly if the authors provided more intuitive explanations for the newly introduced quantities, such as the quantity labeled as (6) in the paper. Enhancing clarity in the presentation of the research would greatly improve its accessibility to readers.\n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you for your helpful comments and feedback! Our replies below address the individual points raised in your review.\n\n- In response to Weakness 1 and Weakness 2, \nThank you for these points. These questions are inherent to the conformal prediction framework, and are not specific to our specific method or to the specific setting of matrix completion. Our replies are therefore general, and apply to any implementation of the conformal framework to any prediction problem:\n\n    - To reply to Weakness 1, the conformal prediction method is designed to provide uncertainty quantification to ANY base estimation algorithm – it is indeed true that the resulting intervals may be very wide if the base algorithm is a poor estimate, but this is not a drawback of conformal. A useful analogy is using cross validation to estimate the error of a model. If we fit a model that actually has very high out-of-sample error, and then cross validation estimates a very high risk for this model, that is a success of CV (i.e., we have correctly identified that this model does not generalize well), not a failure. Analogously, if conformal prediction is applied to a poor estimation model and consequently gives a very wide prediction interval, that is a success of conformal (i.e., we have correctly identified that, when building a valid prediction interval around this particular \\hat{M}, we need to make it very wide if we want to achieve coverage at a certain level), not a failure. This view of conformal can be captured by describing conformal as a “wrapper method” – its job is to provide guard rails around any base algorithm (i.e., any procedure that we use to produce \\hat{M}), and it’s up to the analyst to choose a good base algorithm.\n\n    - To turn to Weakness 2, now we can see that actually the nature of conformal as a wrapper method, addresses this concern from the referee. While conformal itself does not use any model or property for the matrix, the base algorithm that produces \\hat{M} is free to use ANY prior knowledge about structure in the matrix (and we do not need to worry if our assumptions are exactly correct!). For example, in our implementation of cmc-als, the base algorithm to produce \\hat{M} is the als algorithm, which (implicitly) assumes that M follows a signal-plus-noise structure, M* + (iid noise), with incoherent and low rank M*. If this assumption is correct or approximately correct, our estimator \\hat{M} will be excellent and conformal will provide a precise prediction interval – this is why cmc-als is competitive with als in terms of interval width, in the setting where the assumptions of the model hold. The benefit of cmc (and conformal in general) is that it maintains validity even if the assumptions of the base algorithm are wrong (as we see in our simulations for misspecified settings). We did not discuss these ideas extensively in the paper because this interpretation of the conformal framework is broadly applicable to all the conformal prediction literature and is not specific to our setting. However, we will add a bit more about these ideas into our revised paper to give more context to readers who are less familiar with the conformal prediction framework.\n\n\n- In response to Weakness 3, \nThanks for the suggestions on presentation. Lines 179 - 182 explain intuitions for the newly defined $\\Delta$. In words, $\\Delta$ quantifies the estimation error in $\\hat{w}$ w.r.t.~the oracle weight $w^*$. The smaller $\\Delta$, the better the estimation of the weights, and hence the smaller the coverage gap from the nominal level $1-\\alpha$. \nWe will go through the paper and add more explanations for other quantities to improve clarity. \n"}, {"review_id": "5hSHRqNp3J", "reviewer": "Reviewer_1TvL", "summary": "This paper presents a distribution-free method for constructing prediction intervals in the matrix completion problem, where randomness only arises from the sampling of observed entries. The approach utilizes weighted conformal prediction and establishes a lower bound on the probability of each unobserved entry being included in the prediction interval. The lower bound relies on the estimation error of sampling probabilities, which can be negligible if accurately estimated.", "questions": "* Are there any experimental results when the sampling probabilities are completely misspecified, e.g., when the heterogeneous missingness is misspecified as uniform sampling?\n* I'm curious about the impact of the estimated values of $\\hat{s_{ij}}$'s. Does misspecification of these values have a substantial negative impact on the results? Are there any empirical findings regarding this?\n* Is there a rationale for the estimation of theta at Line 237? If there is a reference concerning this, could you provide it?\n* Since the constructed prediction interval is for M, rather than M*, wouldn't it be appropriate to employ an approach for approximately low-rank matrices in experiments? With an approximately low-rank matrix completion approach, how does the estimation of theta change?\n* What does Line 277 mean?", "limitations": "Yes.", "rating": 5, "confidence": 2, "soundness": 2, "presentation": 3, "contribution": 3, "strengths": "* Overall, the paper is clearly written.\n* The utilization of weighted conformal prediction in the matrix completion problem is quite novel.\n* The study shows that the resulting prediction interval performs well with well-estimated sampling probabilities.", "weaknesses": "* Incorrect estimation of sampling probabilities may significantly degrade the quality of the prediction interval.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you for your helpful comments and feedback! Our replies below address the individual points raised in your review.\n\n- In response to the Weakness and to the 1st bullet point under Questions (“Are there any experimental results when the sampling probabilities are completely misspecified…”),\nIn the attached pdf, we conducted simulations in the following 4 settings:\n\n    - Setting [a] [Figure 1(a)] The underlying missingness follows the rank-one model where $p_{ij} = a_i b_j$, both $a_i$ and $b_j$ are generated i.i.d. from Unif(0.2,1). The noise is adversarial (line 273). But we estimate p_{ij} via the one-bit matrix completion based on the logistic model (working model with hypothesized rank k=5).\n\n    - Setting [b] [Figure 1(b)] The underlying missingness follows the logistic model as in Example 3.3.2 with $k^*=5$ and we adopt the adversarial noise on line 273. But $p_{ij}$ is estimated under the assumption of uniform sampling (working model), i.e. $\\hat{p}_{ij}$ = $\\hat{p}$ (line 239).\n\n    - Setting [c] [Figure 1(c)] Same as Setting [a] except that we use the random noise on line 276.\n\n    - Setting [d] [Figure 1(d)] Same as Setting [b] except that we use the random noise on line 276.\n\n\tFrom the results, we can see that when the noise is generated following line 276 (i.e., random noise model), where the values of entries are independent of the sampling probabilities, the misspecification of the sampling model only slightly affects the coverage. Moreover, when the noise is generated in the adversarial way, where the values of entries depend on $p_{ij}$’s, we can see that the coverage with a misspecified sampling model is lower than the target level, but is above 0.85 in practice, which depends on the divergence between the true and the working sampling models. We will add this experiment to the appendix. \n\n- In response to the 2nd bullet point under Questions\nIn fact, $\\hat{s}_{ij}$ doesn’t affect the validity of the prediction intervals. Similar to the estimate $\\hat{M}$, they are both initial estimates. Regardless of their accuracy, our conformal methods help correct them and offer provable predictive coverage. \n\n- In response to the 3rd bullet point under Questions (“Is there a rationale for the estimation of theta at Line 237?”):\nEq [19] in the paper “Inference and uncertainty quantification for noisy matrix completion” published in PNAS provides the rationale for this $\\hat{\\theta}$. In words, it can be viewed as the asymptotic variance of the model-based estimate when the model is correctly specified. \n\n- In response to the 4th bullet point under Questions (“Since the constructed prediction interval is for $M$, rather than $M^*$...”):\nWhen we use model-based estimate ALS, we used the prediction intervals for M (i.e., low-rank matrix plus noise), instead of confidence intervals for $M^*$. So our comparisons are fair in this sense. \n\n- In response to the 5th bullet point under Questions (“What does Line 277 mean?”)\nWe apologize for any confusion. We will reword this in the revision. To clarify:\nThe matrix ($p_{ij}$) is drawn by setting $log[ p_{ij}/(1-p_{ij}) ] = \\langle a_i,b_j \\rangle$ where vectors $a_i$, $b_j$ are drawn as in lines 264–266. To draw the sigma’s, we can think of first drawing an independent copy of p – say, ($p^\\prime_{ij}$) – and then setting $\\sigma_{ij} = 1/2p^\\prime_{ij}$ for each (i, j). Equivalently, we are drawing new $a_i$, $b_j$ vectors, and then defining $\\sigma_{ij}$ by taking $log[(1/2\\sigma_{ij}) / (1-1/2\\sigma_{ij}) ] = \\langle a_i,b_j\\rangle$."}, {"review_id": "lp2xsinXL6", "reviewer": "Reviewer_1PHH", "summary": "This paper utilizes conformal inference techniques to address uncertainty quantification in the matrix completion problem. The authors present a novel method for constructing prediction sets for the missing entries estimated by any given matrix completion algorithm, employing a simple data hold-out strategy. The underlying assumption is that the data matrix $M$ is fixed, while randomness arises from the missingness of the matrix entries.\n\nSpecifically, each index is assumed to have a distinct probability of being observed, independent of the others. The authors demonstrate that when the missingness probabilities are uniform across all entries, the data become exchangeable, allowing for the direct application of standard conformal inference methods.\n\nIn cases with heterogeneous missing probabilities, the problem becomes more intricate and aligns with the conformal inference under the covariate-shift framework introduced by Tibshirani et al. (2019). To tackle this scenario, the authors propose applying the method developed by Tibshirani et al. (2019) and employ data-driven estimates of the missingness probabilities derived from parametric models.\n\nThe validity of the proposed method is established through both theoretical arguments leveraging existing results and empirical evaluation using synthetic and real data experiments. These experiments effectively demonstrate the efficacy and practical applicability of the proposed approach.", "questions": "1)    Usefulness of the one-shot shortcut. The paper introduces a one-shot weighted conformal approach to reduce computational costs. However, it is important to consider potential issues with this relaxation. Firstly, the presence of extremely small probabilities $p_{ij}$ for some test points may lead to excessively large odd ratios $h_{ij}$, resulting in overly conservative predictions that may not be useful for other test points. Secondly, it is not entirely clear why the proposed algorithm without the one-shot relaxation would be prohibitively computationally expensive. Since the weight $w_{ij}$ has a simple form, it should not significantly increase evaluation time in practice. \n\n2)    Improvement of technical notation. Some technical details require clarification. For example, in line 157 on page 5, the odd ratio is defined as $(1-p_{ij})/p_{ij}$, which poses a problem when $p_{ij}=0$ for some index $(i,j)$. Additionally, in the Supplement Material, lines 464-467, the variables $\\mathbf{Z}$ and $\\mathbf{W}$ are not explicitly defined, making it more difficult to verify the proof.\n\n3)    Missing citation. In line 233-237, the benchmark prediction sets are mentioned, which are believed to follow from the asymptotic results in Chen et al. (2019). However, the citation for this reference is missing. Including the appropriate citation will provide proper attribution and give readers an opportunity to explore the referenced work.\n\n4)    Consideration of alternative methods. There are other matrix completion algorithms available that provide uncertainty quantification, such as matrix completion with Gaussian Copula (Zhao et al, 2020). It would be beneficial for the authors to implement additional benchmarks using these alternative methods to strengthen the validity and robustness of the data experiments.", "limitations": "The main limitations of this paper concern its technical novelty compared to prior work on conformal inference and the depth of its empirical evaluations.", "rating": 6, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 2, "strengths": "1)    This paper is both interesting and original, as it contributes to connection between conformal inference and matrix completion. The authors successfully present a reasonable missingness model and introduce a principled framework to apply weighted split conformal methods.\n\n2)    The practical usefulness and potential impact of this paper are noteworthy. Until recently, uncertainty estimation in matrix completion was relatively unexplored, and prior methods were limited due to their heavy reliance on assumptions that often do not hold in practice. Therefore, the assumption-lean approach presented in this paper holds substantial promise and has the potential to significantly impact related fields.\n\n3)    The paper is well-written, and the proposed method is clearly and effectively explained.", "weaknesses": "1)    Technical novelty and originality of theoretical contributions. While this paper seems to rely heavily on the results of Tibshirani et al. (2019) and Barber et al. (2022), the relationship between these works could be explained more clearly. In particular, it would be beneficial to clarify which aspects of the proofs are novel and which can be considered as special instances of prior work. \n\n2)    Thoroughness of the numerical experiments. While the proposed method appears promising in theory, the numerical experiments lack convincing evidence. For example, in the four settings depicted in Figure 1, it would be desirable to observe under-coverage in $\\texttt{als}$ when the signal-to-noise ratio is low or when the incoherence condition is violated, even with an oracle rank. However, both $\\texttt{cmc}$ and $\\texttt{als}$ exhibit similar performance when the rank is correctly chosen, with the advantages of $\\texttt{cmc}$ primarily stemming from tuning the hypothesized rank. Further investigation and more comprehensive experiments would be useful to establish the illustrate of the proposed method.\n\n3)    Realism of real-data experiments. The real data application may not provide significantly more information than a synthetic experiment since the authors begin with the full matrix $M$ and manually sample missing entries using a logistic missingness model. However, it remains unclear whether this logistic model accurately represents practical scenarios. Designing experiments based on more realistic missingness patterns would enhance the informativeness of the results. Additionally, it is expected that the benchmark $\\texttt{als}$ would exhibit under-coverage in some heterogeneous settings, but the authors did not compare it to $\\texttt{als}$. By neglecting this comparison, the practical relevance of the work is not made as clear as it could be.\n\n4)    Limited empirical evaluation. The authors only consider the ``average coverage rate\" as defined in Equation (3), which provides a relatively weak coverage guarantee. It would be helpful if the numerical experiments also included appropriate conditional coverage metrics to provide a more comprehensive evaluation.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you for your helpful comments and feedback! Our replies address the individual points raised in your review.\n- [Weakness 1] Thank you for the question. The main novelty of our work is in reformulating the matrix completion problem as an instance of weighted exchangeability, for which weighted conformal prediction (WCP) can be applied (Lemma 3.1). Secondly, many existing theoretical results using WCP either assume that the true weights w are known, or require ||\\hat{w}-w|| to be small with high probability; this can be quite loose, and we have been able to work with an error term that measures only the expected value of the gap. We will explain the novelty and contributions more clearly in our revision.\n- [Weakness 2] We set the hypothesized rank r=8 in these two settings. [Figure 1(j)] We modify Setting 2 (line 243) with d=400 and p=0.1. The coverage rate of als drops to 0.74 while cmc-als has the desired guarantee. [Figure 1(k)] We modify Setting 4 (line 246) by replacing t_{1.2} with Cauchy(0,1) for the entries. als does not have coverage guarantee while the coverage of cmc-als is exact.\n- [Weakness 3] Real data: In the pdf [Figure 1(e)], we consider the sales dataset where the missingness occurs in the way that: For weekdays, each entry is observed with pr=0.8. On weekends, as the stores are likely to be operated by less experienced interns or are less frequent to report the sales data, each entry is observed with a lower probability, e.g. 0.8/3. Moreover, as there could be a subgroup of stores that are less frequent in sales reporting, 200 stores are randomly sampled, for which pr=0.8/3. We use the logistic model with k=5 to estimate p_{ij}.\nIn the existing literature with theoretical validity guarantee, the implementation of als is not designed for heterogeneous sampling. To compare with als, we estimate p_{ij} via the uniform sampling model.\nIn the pdf [Figure 1(h)], we use the setting same as Figure 2(f) with hypothesized r=8 and observe a coverage gap of 1% for als. When the noise is adversarial [Figure 1(i)], the cmc-als has a coverage gap of 1% due to the estimation error of w while the coverage of als drops to 0.86.\n- [Weakness 4] In the pdf [Figure 1(f)], we evaluate the local performance of cmc-als when conditioning on (p_{ij} = p_0), i.e. a subpopulation determined by the value of sampling probability. We use the setting same as Figure 2(b) with the hypothesized r=12 and p_0 ranging from 0.1 to 0.9 with a step size of 0.1. The conditional coverage is approximated via kernel smoothing: with the indicators A_{ij} for coverage, we use the weight K_{ij} = \\phi_{p_0,h}(p_{ij}) and calculate the conditional coverage by (\\sum A_{ij} K_{ij})/(\\sum K_{ij}). Here \\phi_{\\mu,\\sigma} is the density function of \\calN(\\mu,\\sigma^2) and h=0.05.\nWhen p_0 increases from 0.2 to 0.8, the conditional coverage increases and stays around 0.9. At p_0=0.1 and 0.9, the coverage varies much, which is due to the small effective sample size for the edge values of p_{ij}. Moreover, for uniform sampling, since the rows & columns are generated i.i.d., there are no meaningful subgroups of the data to condition on.\n- [Question 1] Thanks for this question – we realize now that our explanation in the paper was not clear. While there is a small computational benefit, indeed as the reviewer points out, it is not really significant. The real benefit is that the one-shot method produces a much more interpretable answer: the final output of the method is of the form on line 168 for a single value \\hat{q} – that is, the same scaling or inflation factor is applied across the entire matrix.\nWe would like to pause to clarify why we would NOT want to have different rescaling factors \\hat{q}_{ij} at different entries. In fact, it may appear initially that this output would be more meaningful – we might have different levels of uncertainty at different locations (i,j), so would it not be beneficial to choose \\hat{q}_{ij} adapting to the inflation of an entry (i,j)?\nThe answer to this is that, without the one-shot simplification, the different values of \\hat{q}_{ij} are actually unrelated to the issue of higher or lower uncertainty about our estimates of \\hat{M} and \\hat{s} at this entry. They only have to do with differences in the weights vectors - an entry (i,j) would have a higher \\hat{q}_{ij} based solely on its (estimated) weight \\hat{h}_{ij} that reflects its probability of being sampled. Thus, the \\hat{q}_{ij}'s are not actually reflecting a meaningful notion of local noise/uncertainty and in our opinion the method is more interpretable with the simplification.\nWith that said, of course we do not want to propose a method that is needlessly conservative. To that end, we have carried out a simulation to verify that the two versions give essentially the same performance in the pdf [Figure 1(g)]. The setting is the same as Figure 2(b) with hypothesized r=12.\n- [Question 2] Thanks for pointing this out. We assume that the p_{ij}'s are nonzero in order for WCP to be well defined. We apologize for omitting to state this explicitly and will add this assumption to our revision.\nNotations Z & W – we apologize for the confusing notations in the appendix. The matrix W consists of entries W_{ij} defined on line 2 in Algorithm 1. The matrix Z consists of Z_{ij} which are the indicators of non-missingness, i.e. Z_{ij} = 1{(i,j) is observed}.\n- [Question 3] The reviewer is correct. For model-based inferential results, we use results from Chen et al. (2019). We will add this reference in place.  \n- [Question 4] Thanks for this reference. This paper gives theoretical guarantees under the Gaussian copula assumption for the data, and thus does not provide an alternative mechanism for distribution-free theory. We will add the paper to our discussion. Unfortunately, due to space constraints, we cannot add an additional experiment to compare, but we expect to see a loss of coverage when the model assumption is strongly violated."}, {"review_id": "bQtsOHJ43A", "reviewer": "Reviewer_u41g", "summary": "This paper addresses the problem of uncertainty quantification in matrix completion by developing a distribution-free method for predictive inference.  The authors propose a novel approach based on the conformal prediction framework, aiming to overcome the limitations imposed by stringent model assumptions such as low-rankness of the underlying matrix and the light-tailed noise.  Their method, referred to as Conformalized Matrix Completion (CMC; Algorithm 1), can be combined with any arbitrary matrix completion algorithm to provide confidence intervals for the imputed matrix entries.  Here, the “confidence” is measured with respect to the randomness inherent in the measurements of the matrix entries, assuming a randomized probabilistic measurement model.  \n\nTo support their proposed method, the authors also present an intuitive explanation of the method and a theoretical guarantee.  They start by giving an intuitive exposition of how the (weighted) exchangeability plays a crucial role in constructing confidence intervals within the conformal prediction framework (Section 3.1).  Thereafter, they present Theorem 3.2, which establishes a theoretical guarantee on the average coverage of the confidence intervals, which serves as a key theoretical result of this paper.  Additionally, the authors illustrate their analysis by describing two examples of the missingness patterns (i.e., the random observation models) in Section 3.3.  Through this illustration, they explicitly evaluate the rate of the expected “weight gap” as expressed in Eq. (6), ultimately demonstrating satisfactory coverage of the constructed confidence intervals. \n\nIn Section 4, the authors report the results of their numerical simulations. These simulations serve two purposes: (1) comparing the performance of the proposed conformalized method against a model-based baseline (Sections 4.1.1 & 4.2), and (2) investigating the stability of the proposed method when the oracle knowledge of the measurement probability for each entries is not available (Section 4.1.2).", "questions": "1.  I suggest the authors consider providing a clarification of the assumption on the measurement model, as mentioned in the first point of the \"Weakness\" section. This clarification would help readers better understand the specific assumptions and limitations of the proposed method.\n\n2.  It would be valuable if the authors could compare their method with the approaches mentioned in the second paragraph (lines 24-30) by conducting simulations or other appropriate means. Such comparisons would enable a clearer assessment of how the proposed method outperforms or differs from existing methods.\n\n3.  I am curious about the degree to which the proposed CMC method relies on the uncertainty estimate \\hat{s}. Many matrix estimation algorithms only provide point estimates for the entries without accompanying uncertainty estimates. In Algorithm 1, when uncertainty estimates are unavailable, the authors set \\hat{s}_{ij} = 1 by default (line 3). However, this default choice can be problematic. For instance, consider two scenarios: (1) estimating M and (2) estimating 10M. In scenario 2, the uncertainty should be ten times larger than in scenario 1, but the default choice cannot account for this. It would be beneficial if the authors could address this issue and discuss potential alternatives or adjustments to handle situations where uncertainty estimates are not provided.\n\n4.  In Section 3.3, it might be helpful for readers to have either the upper bound for the expected weight gap or the resulting lower bound for average coverage presented as explicit corollaries.\n\n5.  The caption of Figure 2 does not provide sufficient information to understand what is being compared. It would be helpful to revise the caption to clearly indicate the elements being compared or provide a brief description of the comparison being made.", "limitations": "While the authors outline potential future research directions, they do not explicitly address the limitations of their proposed method and approach. It would be beneficial if the authors could include a brief discussion of the technical limitations, such as the assumptions made in their approach. Additionally, providing insights into potential negative impacts when applying the method in real-world applications would further enhance the paper's practical relevance. However, it should be noted that as this paper primarily focuses on theoretical aspects, addressing these limitations in detail may not be deemed critical.", "rating": 6, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 2, "strengths": "One of the primary strengths of this paper lies in its adaptation of the conformal prediction framework to address the challenge of uncertainty quantification in matrix completion. Unlike previous approaches in the literature, this method does not heavily rely on specific model assumptions. This combination of conformal prediction and matrix completion is a significant and original contribution, making the paper stand out in its field.\n\nAnother notable strength is the careful organization of the paper, which effectively presents the essential components of the study to support the authors' main claims. Section 2 provides a concise overview of the problem setup and evaluation metric. The proposed method is then described in detail, accompanied by an intuitive explanation based on the concept of exchangeability in Section 3.1. The main theoretical result is presented neatly in Section 3.2, and to further enhance the understanding, concrete examples of the missingness patterns are provided in Section 3.3. The paper is further strengthened by the inclusion of numerical study results in Section 4 and a comprehensive discussion in Section 5. As a result, the paper forms a clear and self-contained report of the study, ensuring that the readers can grasp the main ideas and findings with ease.\n\nOverall, the adaptation of conformal prediction, the clear organization of the paper, and its ability to convey the study's main contributions effectively are key strengths of this research work.", "weaknesses": "While the paper exhibits notable strengths, there are a few areas that could benefit from further improvement and clarification. I have identified three main concerns that warrant attention.\n\nFirstly, the authors frequently emphasize that their proposed method is \"distribution-free\" and does not rely on “any” assumptions about the underlying matrix (e.g., in line 60). However, it is important to explicitly acknowledge that the proposed method and analysis do depend on the probabilistic random observation model, which is necessary for the conformal prediction framework. While this limitation does not appear to be critical, it would be beneficial to provide a careful clarification of the assumptions and limitations to avoid potential confusion.\n\nSecondly, although the authors highlight the limitations of existing uncertainty quantification approaches (lines 24 - 30), there is a lack of comparisons between the proposed method and these previous approaches. Consequently, it remains unclear whether and how the proposed method surpasses these existing methods. Including such comparisons would enhance the clarity and strengthen the argument for the superiority of the proposed approach. \n\nAdditionally, the issue of heavy-tailed noise is briefly mentioned by the authors (lines 32 - 33), who suggest that their method is less sensitive to the noise tail (lines 32-33, lines 257-261). However, this point is not adequately addressed beyond the comparison of the vanilla ALS method to the CMC-ALS method in Figure 1-(c), where the latter seems to be slightly better than the former, but still remains overly conservative; the difference in performance between the ALS and the CMC-ALS does not appear to be substantial (approximately 9x worse vs 6x worse than the oracle). To improve the cohesiveness of the paper, it would be beneficial to either elaborate further on this point or consider removing it if it does not significantly contribute to the main findings.\n\nIn summary, while the paper possesses strengths, it would benefit from addressing these weaknesses. Clearer clarification of assumptions and limitations, comparative analyses with existing methods, and a more comprehensive discussion of the performance in the presence of heavy-tailed noise would enhance the overall quality and cohesiveness of the paper.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you for your helpful comments and feedback! Our replies below address the individual points raised in your review.\n\n- In response to the 1st Weakness (“Firstly, the authors frequently emphasize that their proposed method is \"distribution-free\" “),\nThanks for the suggestions. We will add clarifications about the assumptions and limitations in the revision. \n\n- In response to the 2nd Weakness (“Secondly, although the authors highlight the limitations of existing uncertainty quantification approaches…“),\nIn fact, in our experiment, we did compare our conformal methods with the existing uncertainty quantification approaches (lines 24 - 30). But we forgot to add references. In the updated manuscript, we will add corresponding references to existing approaches. \n\n- In response to the 3rd Weakness (“Additionally, the issue of heavy-tailed noise is briefly mentioned…“),\nThanks for the suggestion on this heavy-tailed setting. In the appended pdf [Figure 1(l)], we further show the ratio between the obtained length and the oracle length for both als and cmc-als. Although the heavy-tailed problem is indeed a hard problem in terms of statistical inference, we can see that cmc-als actually corrects the conservativeness a lot over als. Moreover, in this heavy-tailed setting, we can see from the error bars that als is much more variable than cmc-als.\n\n- Question 1 is addressed in the response to the 1st Weakness, above.\n- Question 2 is addressed in the response to the 2nd Weakness, above.\n- In response to Question 3,\nThank you for this question - this is an important point and we will be sure to explain more clearly in our revision. The values $\\hat{s}\\_{ij}$ capture the estimated relative, rather than absolute, noise levels across the different entries. Rescaling the entire matrix $\\hat{s}$ by a constant (i.e., replacing $\\hat{s}\\_{ij}$ with $\\hat{s}^\\prime_{ij} = c\\cdot \\hat{s}\\_{ij}$, for the same value c for all (i,j) does not change the outcome of the procedure; running with $\\hat{s}^\\prime\\_{ij}$'s instead of $\\hat{s}\\_{ij}$'s would result in a quantile value $\\hat{q}^\\prime = \\hat{q}/c$, and the resulting prediction intervals would therefore be identical. This means that, if we set $s_{ij}=1$ for all (i,j), this ensures that all the entries will be given prediction intervals of equal width, and can therefore be interpreted as assuming that the entries all have equal variance (but we are not assuming that the entries have variance equal to 1).\n\n\n- In response to Question 4,\nThank you for the suggestion! The detailed result for Example 3.3.1 is shown in Section A.3.1 in the appendix and the detailed result for Example 3.3.2 is shown in Section A.3.2. To improve clarity, we will write them as explicit theorems and add pointers to them in the main text.\n- In response to Question 5,\nThanks for the suggestions. We will add more details into the figure captions for Fig 1, 2,3 "}]

Author rebuttal:
We are grateful to all the reviewers for their helpful feedback, comments, and suggestions on our manuscript. In the comments below, we have replied to each reviewer’s points individually.
The attached pdf contains additional simulation results and details of each setting are stated in the rebuttal.


Meta-review: The reviews seem high quality and agree that this is a nice contribution that advances the literature on matrix completion using ideas from conformal inference to establish confidence bounds on imputations.
Decision: Accept (poster)