Matrix completion

Matrix completion is the task of filling in the missing entries of a partially observed matrix, which is equivalent to performing data imputation in statistics. A wide range of datasets are naturally organized in matrix form. One example is the movie-ratings matrix, as appears in the Netflix problem: Given a ratings matrix in which each entry (i,j) represents the rating of movie j by customer i , if customer i has watched movie j and is otherwise missing, we would like to predict the remaining entries in order to make good recommendations to customers on what to watch next. Another example is the document-term matrix: The frequencies of words used in a collection of documents can be represented as a matrix, where each entry corresponds to the number of times the associated term appears in the indicated document. Without any restrictions on the number of degrees of freedom in the completed matrix, this problem is underdetermined since the hidden entries could be assigned arbitrary values. Thus, we require some assumption on the matrix to create a well-posed problem, such as assuming it has maximal determinant, is positive definite, or is low-rank. For example, one may assume the matrix has low-rank structure, and then seek to find the lowest rank matrix or, if the rank of the completed matrix is known, a matrix of rank r that matches the known entries.…

Matrix completion is the task of filling in the missing entries of a partially observed matrix, which is equivalent to performing data imputation in statistics. A wide range of datasets are naturally organized in matrix form. One example is the movie-ratings matrix, as appears in the Netflix problem: Given a ratings matrix in which each entry (i,j) represents the rating of movie j by customer i , if customer i has watched movie j and is otherwise missing, we would like to predict the remaining entries in order to make good recommendations to customers on what to watch next. Another example is the document-term matrix: The frequencies of words used in a collection of documents can be represented as a matrix, where each entry corresponds to the number of times the associated term appears in the indicated document. Without any restrictions on the number of degrees of freedom in the completed matrix, this problem is underdetermined since the hidden entries could be assigned arbitrary values. Thus, we require some assumption on the matrix to create a well-posed problem, such as assuming it has maximal determinant, is positive definite, or is low-rank. For example, one may assume the matrix has low-rank structure, and then seek to find the lowest rank matrix or, if the rank of the completed matrix is known, a matrix of rank r that matches the known entries. The illustration shows that a partially revealed rank-1 matrix (on the left) can be completed with zero-error (on the right) since all the rows with missing entries should be the same as the third row. In the case of the Netflix problem the ratings matrix is expected to be low-rank since user preferences can often be described by a few factors, such as the movie genre and time of release. Other applications include computer vision, where missing pixels in images need to be reconstructed, detecting the global positioning of sensors in a network from partial distance information, and multiclass learning. The matrix completion problem is in general NP-hard, but under additional assumptions there are efficient algorithms that achieve exact reconstruction with high probability. In statistical learning point of view, the matrix completion problem is an application of matrix regularization which is a generalization of vector regularization. For example, in the low-rank matrix completion problem one may apply the regularization penalty taking the form of a nuclear norm R(X) = \lambda\|X\|_* & \text{subject to} & X_{ij} = M_{ij} & \;\; \forall i,j \in E\\ \end{align} Candès and Recht proved that with assumptions on the sampling of the observed entries and sufficiently many sampled entries this problem has a unique solution with high probability. An equivalent formulation, given that the matrix M to be recovered is known to be of rank r , is to solve for X where X_{ij} = M_{ij} \;\; \forall i,j \in E Candès and Recht define the coherence of a matrix U with column space an r- dimensional subspace of \mathbb{R}^n as \mu (U) = \frac{n}{r} \max_{i , where P_U is the orthogonal projection onto U . Incoherence then asserts that given the singular value decomposition U\Sigma V^\dagger of the m by n matrix M , \mu (U), \; \mu (V) \leq \mu_0 The entries of \sum_k u_k v_k ^\dagger have magnitudes upper bounded by \mu_1 \sqrt{\frac{r}{mn}} for some \mu_0, \; \mu_1 . Y_{ij} = M_{ij} + Z_{ij}, (i,j) \in \Omega, where {Z_{ij}:(i,j) \in \Omega} is a noise term. Note that the noise can be either stochastic or deterministic. Alternatively the model can be expressed as P_\Omega(Y) = P_\Omega(M) + P_\Omega(Z), where Z is an n \times n matrix with entries Z_{ij} for (i,j) \in \Omega assuming that \|P_\Omega(Z)\|_F\leq\delta for some \delta > 0 .To recover the incomplete matrix, we try to solve the following optimization problem: \begin{align} & \underset{X}{\text{min}} & \|X\|_* \\ & \text{subject to} & \|P_\Omega(X-Y)\|_F \leq \delta\\ \end{align} Among all matrices consistent with the data, find the one with minimum nuclear norm. Candès and Plan have shown that this reconstruction is accurate. They have proved that when perfect noiseless recovery occurs, then matrix completion is stable vis a vis perturbations. The error is proportional to the noise level \delta . Therefore, when the noise level is small, the error is small. Here the matrix completion problem does not obey the restricted isometry property (RIP). For matrices, the RIP would assume that the sampling operator obeys (1-\delta)\|X\|^2_F \leq \frac{1}{p}\|P_\Omega(X)\|^2_F \leq (1+\delta)\|X\|^2_F for all matrices X with sufficiently small rank and \delta sufficiently small. The methods are also applicable to sparse signal recovery problems in which the RIP does not hold. & \text{subject to} & & X_{ij} = M_{ij} \;\; \forall i,j \in E\\ & & & \begin{bmatrix} W_1 & X \\ X^T & W_2 \end{bmatrix} \succeq 0 \end{align} The complexity of using SDP to solve the convex relaxation is O(\text{max}(m,n)^4) . State of the art solvers like SDPT3 can only handle matrices of size up to 100 by 100. An alternative first order method that approximately solves the convex relaxation is the Singular Value Thresholding Algorithm introduced by Cai, Candès and Shen. Candès and Recht show, using the study of random variables on Banach spaces, that if the number of observed entries is on the order of \max{\{\mu_1^2, \sqrt{\mu_0}\mu_1, \mu_0 n^{0.25}\}}nr \log n (assume without loss of generality m ), the rank minimization problem has a unique solution which also happens to be the solution of its convex relaxation with probability 1-\frac{c}{n^3} for some constant c . If the rank of M is small ( r \leq \frac{n^{0.2}}{\mu_0} ), the size of the set of observations reduces to the order of \mu_0 n^{1.2} r \log n . These results are near optimal, since the minimum number of entries that must be observed for the matrix completion problem to not be underdetermined is on the order of nr \log n . This result has been improved by Candès and Tao. They achieve bounds that differ from the optimal bounds only by polylogarithmic factors by strengthening the assumptions. Instead of the incoherence property, they assume the strong incoherence property with parameter \mu_3 . This property states that: | \langle e_a, P_U e_{a'} \rangle - \frac{r}{m} 1_{a = a'} | \leq \mu_3 \frac{\sqrt{r}}{m} for a, a' \leq m and | \langle e_b, P_U e_{b'} \rangle - \frac{r}{n} 1_{b = b'} | \leq \mu_3 \frac{\sqrt{r}}{n} for b, b' \leq n The entries of \sum_i u_i v_i^\dagger are bounded in magnitude by \mu_3 \sqrt{\frac{r}{mn}} Intuitively, strong incoherence of a matrix U asserts that the orthogonal projections of standard basis vectors to U has magnitudes that have high likelihood if the singular vectors were distributed randomly. Candès and Tao find that when r is O(1) and the number of observed entries is on the order of \mu_3^4 n(\log n)^2 , the rank minimization problem has a unique solution which also happens to be the solution of its convex relaxation with probability 1-\frac{c}{n^3} for some constant c . For arbitrary r , the number of observed entries sufficient for this assertion hold true is on the order of \mu_3^2 nr (\log n)^6 Another convex relaxation approach is to minimize the Frobenius squared norm under a rank constraint. This is equivalent to solving \begin{align} & \min\limits_{X} & & \Vert X \Vert_F^2 \\ & \text{subject to} & & X_{ij} = M_{ij} \;\; \forall i,j \in E\\ & & & \operatorname{Rank}(X) \leq k. \end{align} By introducing an orthogonal projection matrix Y (meaning Y^2=Y, Y=Y' ) to model the rank of X via X=YX, \text{trace}(Y)\leq k and taking this problem's convex relaxation, we obtain the following semidefinite program \begin{align} & \min\limits_{X, Y, \theta} & & \text{trace}(\theta) \\ & \text{subject to} & & X_{ij} = M_{ij} \;\; \forall i,j \in E\\ & & & \operatorname{trace}(Y) \leq k, 0 \preceq Y \preceq I\\ & & & \begin{pmatrix} Y & X \\ X^\top & \theta \end{pmatrix}\succeq 0. \end{align} If Y is a projection matrix (i.e., has binary eigenvalues) in this relaxation, then the relaxation is tight. Otherwise, it gives a valid lower bound on the overall objective. Moreover, it can be converted into a feasible solution with a (slightly) larger objective by rounding the eigenvalues of Y greedily. Remarkably, this convex relaxation can be solved by alternating minimization on X and Y without solving any SDPs, and thus it scales beyond the typical numerical limits of state-of-the-art SDP solvers like SDPT3 or Mosek. This approach is a special case of a more general reformulation technique, which can be applied to obtain a valid lower bound on any low-rank problem with a trace-matrix-convex objective. Solve \min_{X,Y} \min_{S \in \mathbb{R}^{r \times r}} \frac{1}{2} \sum_{i,j \in E} (M_{ij} - (XSY^\dagger)_{ij})^2 + \rho G(X,Y) where G(X,Y) is some regularization function by gradient descent with line search. Initialize X,\;Y at X_0,\;Y_0 where \text{Tr}(M_E) = X_0 S_0 Y_0^\dagger . Set G(X,Y) as some function forcing X, \; Y to remain incoherent throughout gradient descent if X_0 and Y_0 are incoherent. Return the matrix XSY^\dagger . Steps 1 and 2 of the algorithm yield a matrix \text{Tr}(M^E) very close to the true matrix M (as measured by the root mean square error (RMSE)) with high probability. In particular, with probability 1-\frac{1}{n^3} , \frac{1}{mnM_{\text{max}}^2} \| M - \text{Tr}(M^E) \|_F^2 \leq C \frac{r}{m|E|} \sqrt{\frac{m}{n}} for some constant C . \| \cdot \|_F denotes the Frobenius norm. Note that the full suite of assumptions is not needed for this result to hold. The incoherence condition, for example, only comes into play in exact reconstruction. Finally, although trimming may seem counter intuitive as it involves throwing out information, it ensures projecting M^E onto its first r principal components gives more information about the underlying matrix M than about the observed entries. In Step 3, the space of candidate matrices X,\;Y can be reduced by noticing that the inner minimization problem has the same solution for (X,Y) as for (XQ,YR) where Q and R are orthonormal r by r matrices. Then gradient descent can be performed over the cross product of two Grassmann manifolds. If r \ll m,\;n and the observed entry set is in the order of nr\log n , the matrix returned by Step 3 is exactly M . Then the algorithm is order optimal, since we know that for the matrix completion problem to not be underdetermined the number of entries must be in the order of nr\log n . X= UV^T ; the algorithm then alternates between finding the best U and the best V . While the overall problem is non-convex, each sub-problem is typically convex and can be solved efficiently. Jain, Netrapalli and Sanghavi have given one of the first guarantees for performance of alternating minimization for both matrix completion and matrix sensing. The alternating minimization algorithm can be viewed as an approximate way to solve the following non-convex problem: \begin{align} & \underset{U, V \in \mathbb{R}^{n\times k}}{\text{min}} & \|P_\Omega(UV^T)-P_\Omega(M)\|^2_F \\ \end{align} The AltMinComplete Algorithm proposed by Jain, Netrapalli and Sanghavi is listed here: Input: observed set \Omega , values P_\Omega(M) Partition \Omega into 2T+1 subsets \Omega_0,\cdots,\Omega_{2T} with each element of \Omega belonging to one of the \Omega_t with equal probability (sampling with replacement) \hat{U}^0 = SVD(\frac{1}{p}P_{\Omega_0}(M), k) i.e., top- k left singular vectors of \frac{1}{p}P_{\Omega_0}(M) Clipping: Set all elements of \hat{U}^0 that have magnitude greater than \frac{2\mu\sqrt{k}}{\sqrt{n}} to zero and orthonormalize the columns of \hat{U}^0 for t = 0, \cdots, T-1 do \quad \hat{V}^{t+1}\leftarrow \text{argmin}_{V\in \mathbb{R}^{n\times k}}\|P_{\Omega_{t+1}}(\hat{U}V^T-M)\|^2_F \quad \hat{U}^{t+1}\leftarrow \text{argmin}_{U\in \mathbb{R}^{m\times k}}\|P_{\Omega_{T+t+1}}(U(\hat{V}^{t+1})^T-M)\|^2_F end for Return X= \hat{U}^T(\hat{V}^T)^T They showed that by observing |\Omega| = O((\frac{\sigma_1^*}{\sigma_k^*})^6k^7\log n \log (k \|M\|_F/\epsilon)) random entries of an incoherent matrix M , AltMinComplete algorithm can recover M in O(\log(1/\epsilon)) steps. In terms of sample complexity ( |\Omega| ), theoretically, Alternating Minimization may require a bigger \Omega than Convex Relaxation. However empirically it seems not the case which implies that the sample complexity bounds can be further tightened. In terms of time complexity, they showed that AltMinComplete needs time O(|\Omega|k^2\log(1/\epsilon)) . It is worth noting that, although convex relaxation based methods have rigorous analysis, alternating minimization based algorithms are more successful in practice. (\boldsymbol{X}) = \frac{1}{2} \left\| P_{\Omega}(\boldsymbol{X} - \boldsymbol{O}) \right\|_F^2 ensures fidelity to the observed entries, with P_{\Omega} as the projection onto the observed set \Omega and \boldsymbol{O} as the observed matrix. g(\boldsymbol{X}) = \|\boldsymbol{X}\|_* is the nuclear norm to enforce low-rank structure. r(\boldsymbol{X} \mid 0) = \sum_{k=1}^ \left\| \operatorname{vec}_{\overline{\Omega}}(\boldsymbol{X}) - a_k \mathbf{1} \right\|_0 is the discrete-space regularizer, with \mathcal{A} being the discrete alphabet (e.g., {1, 2, 3, 4, 5}) and \overline{\Omega} the set of unobserved entries. To solve this non-convex problem, the \ell_0 -norm is approximated by a continuous function. This approximation is convexized using fractional programming, transforming the problem into a series of convex subproblems. The algorithm iteratively updates the matrix estimate by applying proximal operations to the discrete-space regularizer and singular value thresholding to enforce the low-rank constraint. Initializing the process with the solution from the \ell_1 -norm-based method can accelerate convergence. Simulation results, tested on datasets like MovieLens-100k, demonstrate that this method outperforms both its \ell_1 -norm-based predecessor and other state-of-the-art techniques, particularly when the ratio of observed entries is low (e.g., 20% to 60%). ===Collaborative filtering=== Collaborative filtering is the task of making automatic predictions about the interests of a user by collecting taste information from many users. Companies like Apple, Amazon, Barnes and Noble, and Netflix are trying to predict their user preferences from partial knowledge. In these kind of matrix completion problem, the unknown full matrix is often considered low rank because only a few factors typically contribute to an individual's tastes or preference. ===System identification=== In control, one would like to fit a discrete-time linear time-invariant state-space model \begin{align} x(t+1)&=Ax(t)+Bu(t)\\ y(t)&=Cx(t)+Du(t) \end{align} to a sequence of inputs u(t) \in \mathbb{R}^m and outputs y(t) \in \mathbb{R}^p, t = 0, \ldots, N . The vector x(t) \in \mathbb{R}^n is the state of the system at time t and n is the order of the system model. From the input/output pair, one would like to recover the matrices A,B,C,D and the initial state x(0) . This problem can also be viewed as a low-rank matrix completion problem. ===Internet of things (IoT) localization=== The localization (or global positioning) problem emerges naturally in IoT sensor networks. The problem is to recover the sensor map in Euclidean space from a local or partial set of pairwise distances. Thus it is a matrix completion problem with rank two if the sensors are located in a 2D plane and three if they are in a 3D space. ===Social networks recovery=== Most of the real-world social networks have low-rank distance matrices. When we are not able to measure the complete network, which can be due to reasons such as private nodes, limited storage or compute resources, we only have a fraction of distance entries known. Criminal networks are a good example of such networks. Low-rank Matrix Completion can be used to recover these unobserved distances. ==See also== Matrix regularization Netflix Prize Collaborative filtering System identification Convex optimization Imputation

Text derived from https://en.wikipedia.org/w/index.php?oldid=1377560297 under CC BY-SA 4.0; modifications have been made.