Robust Learning from Noisy Side-information by Semidefinite Programming

1 Introduction

Semidefinite programming (SDP) studies optimization problems with a convex objective function over semidefinite constraints [Vandenberghe and Boyd, 1996; Boyd and Vandenberghe, 2004]. Many machine learning problems can be reduced as SDPs [Lemon et al., 2016]. Prominent examples include embedding and clustering [Kulis et al., 2007; Royer, 2017], sparse PCA [D'aspremont et al., 2007; Zou and Xue, 2018], maximum variance unfolding (MVU) [Weinberger et al., 2004; Song et al., 2008], non-parametric kernel learning (NPKL) [Li et al., 2008; Zhuang et al., 2011]. Generally, the SDP optimization problem is formulated as

where {Qτ , tτ}mτ=1 comes from the training data (such as side information), A is a symmetric matrix to regularize Z (depending on applications), S+ is the cone of positive semidefinite (PSD) matrices, and γ > 0 is a hyper-parameter.

1 Quanming Yao is the corresponding author.

The PSD constraint, i.e., Z ∈ S+, is the most challenging part in solving (1) [Lemon et al., 2016]. For example, it costs O(n3 ) time at each iteration using the interior point algorithm [Helmberg et al., 1996] if Z is of size n × n. Another example is the projection gradient descent algorithm, in which the projection to PSD cone will also cost O(n3 ) time [Jaggi, 2013]. To drop the PSD constraint, there is an (efficient) matrix factorization method to go about fitting a low rank model [Burer and Monteiro, 2003; Lemon et al., 2016]. Namely, factorizing Z to XX⊤, then (1) can be converted into as below

While the problem is nonconvex, instead of optimizing w.r.t. Z ∈ S+, we only need to solve an unconstrained optimization problemwith variable X ∈ Rn×r . Moreover, it is theoretically shown that the factorized problem (2) is equivalent to (1) when the rank of solution is deficient [Srinadh et al., 2016; Zheng and Lafferty, 2015].

Many algorithms have been proposed to solve (2) and are all much more efficient than interior method and projection gradient descent for (1). When the objective F is linear, L-BFGS is introduced in [Burer and Monteiro, 2003; Nocedal and Wright, 2006] for optimization. However, the convergence properties of L-BFGS are unclear for the nonconvex problem here. When F is convex and smooth, blockcyclic coordinate minimization has been used in [Hu et al., 2011] to solve a special nonconvex program of SDP, but a closed-form solution is preferred in each block coordinate update, which might be overly restrictive. More recently, gradient descent based methods [Srinadh et al., 2016; Zheng and Lafferty, 2015] have been developed as the stateof-the-art for low-rank SDP. These algorithms have a convergence guarantee, and linear/sub-linear convergence rate are also established for some low-rank SDP formulations [Srinadh et al., 2016; Pumir et al., 2018].

In above applications, the squared loss is used in F to encourage the learned Z to be consistent with give sideinformation. However, since the squared loss is sensitive to outliers, all existing SDP algorithms are not robust. Moreover, robustness is of a real demand. The side-information utilized in SDP may not be accurate for real applications, e.g., samples can corrupted in MVU [Dekel et al., 2010] and links collected for kernel learning can come from spammer or attacker [Raykar et al., 2010]. Such corruptions and noise can significantly deteriorate performance of learning models [Raykar et al., 2010].

Motivated by the success of making matrix factorization robust by replacing the squared loss with ℓ1 loss [Lin et al., 2017; Yao and Kwok, 2018], we also proposed to use the ℓ1 loss in (2) for SDP, and illustrate such need with three applications, i.e., robust NPKL, robust colored MVU and sparse PCA. However, the resulting optimization problem is neither convex nor smooth, and none of existing SDP algorithms can be used for optimization. To solve the new objective, we propose a new optimization algorithm based on MajorizationMinimization (MM), of which the crux is constructing a good surrogate. Besides, while MM generally only guarantees producing limit points, we prove that by iteratively optimizing the constructed surrogate, the proposed algorithm ensures a convergence to some critical points. Finally, we demonstrate the efficiency and robustness of the proposed algorithm using above three applications. Results show that the proposed algorithm is not only faster but also better in recovery over state-of-the-arts. As a summary, we hight-light our contributions as follows:

• We are the first to introduce the robust loss, i.e., ℓ1 loss, into SDP. To show its necessity, we further illustrate the usage of the new objective with three applications, i.e., robust NPKL, robust colored MVU and sparse PCA;
• As no existing SDP algorithms can be applied, we propose a novel optimization algorithm to solve the new objective, which is not only efficient but also guaranteed converging to some critical points;
• Finally, the robustness of the new objective and the effectiveness of the proposed algorithm over state-of-the-arts SDP algorithms are empirically verified on above three applications.

Notation

We use an uppercase letter to indicate a matrix, and a lowercase letters to a scalar. The transpose of vector or matrix is denoted by the superscript (·)⊤. The identity matrix is denoted by I; for a matrix A = [aij], tr(A) is the trace of a square matrix, ‖A‖F = (∑ij a2ij )1/2 is its Frobenius norm. |a| is the absolute of scaler a, and |S| is cardinal number of set S.

2 Robust Semidefinite Programming

Most SDP learning algorithms assume perfect side information (i.e., perfect triplet constraints in our case). This however is not always the case in practice because in many real-world applications, the constraints are derived from the side information such as users implicit feedbacks and citations among articles. As a result, these constraints are usually noisy and consist of many mistakes. We refer to the problem of learning SDP matrix from noisy side information as robust SDP learning. Feeding the noisy constraints directly into a SDP learning algorithm will inevitably degrade its performance, and more seriously.

2.1 Proposed Formulation

Inspired by the recent success of using ℓ1-loss instead of the squared loss in making matrix factorization robust [Lin et al., 2017; Yao and Kwok, 2018], we also propose to replace the squared loss in (1) by the more robust ℓ1 loss, which makes low-rank SDP less sensitive to corruptions in the training data. This leads to the objective of robust SDP as

where λ > 0 and the last term is to further prevent overfitting. This new objective is neither convex nor smooth. As a result, none of existing algorithms for low-rank SDP, e.g., L-BFGS [Burer and Monteiro, 2003], gradient descent [Srinadh et al., 2016; Zheng and Lafferty, 2015], and coordinate descent [Hu et al., 2011], can be applied. In next Section, we will design an efficient algorithm for (3) based on MM, which also has a convergence guarantee.

2.2 Application Examples

However, before that, we illustrate the usage and importance of the new formulation using three examples.

Example 1: Robust NPKL

Given n patterns, let M be the must-link set containing pairs that should belong to the same class, and C be the cannot-link set containing pairs that should not belong to the same class. Denote T =M∪C. Non-parametric kernel learning (NPKL) [Hoi et al., 2007] tries to build a kernel matrix utilizing above side information. We adopt the formulation in [Li et al., 2008; Zhuang et al., 2011], which learns a kernel matrix using the following SDP problem

where Z is the target kernel matrix, L is the graph Laplacian matrix of the data, M and C are encoded into {(Qτ , tτ )}|T | τ=1 . Let Qτ = I(:, j) (I(:, i)) ⊤, then tr(QτZ) = Zij. Thus tτ = 1 if (i, j) ∈ M and 0 if (i, j) ∈ C. The first term of objective in (4) measures the difference between Zij and tτ , and the second term tr(ZL) encourages smoothness on the data manifold by aligning Z with L.

The side information in this application is those "can" and "cannot" links. These links are usually provided by humans, e.g., labeled by experts or crowdsourced from the web. As human may not be reliable and there can be spammers and attackers in the crowdsourcing platform, errors and noise can exist in these links [Raykar et al., 2010]. These again inspire a more robust formulation of NPKL as

Examples 2: Robust CMVU

Maximum variance unfolding (MVU) [Weinberger et al., 2004; Weinberger and Saul, 2006] is an effective method for dimensionality reduction. It produces a low-dimensional representation of the data by simultaneously maximizing the variance of their embeddings and preserving the local distances of the original data. MVU can be viewed as a nonlinear generalization of principal component analysis. The colored maximum variance unfolding (CMVU) is a "colored" variants of MVU [Song et al., 2008], subjected to class labels information. In [Song et al., 2008], it is formulated as a low-rank SDP problem as

where Eij = I(:, i) - I(:, j), Qτ = EijEij ⊤, dτ = dij denotes the Euclidean distance between the i-th and j-th objects in primal space, N denotes the set of neighbor pairs, whose distances are to be preserved in the embedding, T is a kernel matrix of the labels, Hij = δij- 1n centers the data and the labels in the feature space, and λ controls the tradeoff between dependence maximization and distance preservation.

In this example, the side information is the local distance dτ from the original data. However, during the data collection, outliers or corrupted samples can be introduced into feature space [Dekel et al., 2010]. This motivates our robust formulation as:

Example 3: Sparse PCA

Finally, for the last example, we consider sparse PCA [Zou and Xue, 2018]. Here, we are not making sparse PCA more robust, but show sparse PCA also fall into our objective (3), and thus can be solved with proposed algorithm. For a given covariance matrix Σ ∈ Rn×n , sparse PCA tries to find a sparse vector x that maximizes x⊤Σx, i.e. a sparse principal component of Σ. Following [D'aspremont et al., 2007], sparse PCA can be relaxed into SDP problem as:

We propose to factorize Z = XX⊤, and solve the following approximated objective instead

3 Optimization Algorithm

Majorization minimization (MM) is a general technique to make difficult optimization problems easier [Hunter and Lange, 2004; Lange et al., 2000]. Recently, it has been applied in robust matrix factorization (RMF) [Lin et al., 2017; Yao and Kwok, 2018]. Inspired by such success and the fact that no existing SDP algorithms can be applied here, we give a MM algorithm to solve (3) in the sequel.

3.1 Majorization-Minimization (MM)

Consider a function g(X), which is hard to optimize. Let the iterate at the kth MM iteration be Xk. The next iterate is generated asXk+1 = Xk + arg min

where hk is a surrogate that is being optimized instead of g. A good surrogate should have the following properties [Lange et al., 2000]:

The first two conditions (a) and (b) ensure {g(Xk)} generated from MM is a non-increasing sequence, and (c) encourages sub-problems hk ( ˜X; Xk) can be easily solved.

However, MM only guarantees that the objectives obtained in successive iterations are non-increasing, but does not guarantee convergence of the sequence {Xk} [Hunter and Lange, 2004; Lange et al., 2000; Lin et al., 2017].

3.2 Constructing the Convex Surrogate

Here, we show how a surrogate can be constructed from R, which can meet the above three conditions of MM. First, we upper bound R (see (3)) in following Lemma 1 based on (7).

However, the upper bound in Lemma 1 is not convex, as the term |tr( ˜X⊤Qτ ˜X)| is convex only when Qτ ∈ S+ [Boyd and Vandenberghe, 2004], which is not guaranteed. Let (γi, vi)'s be the eigen-pairs of a symmetric square matrix M, we use (·)+ and (·)- to denote positive and negative eigen values of M, i.e., M+ = ∑i max(γi, 0)vivi⊤ and M- = ∑i min(γi, 0)vivi⊤, thus M = M+ +M-. To address this issue, we make use of following Lemma 2.

Combining Lemma 1 and 2, a convex surrogate is constructed as follow:

Proposition 1. Let B = Q+ 12 (λI + γA+), C = A+ λγ I, Q = ∑mτ=1 ¯Qτ , bkτ = 12 (tr(Xk ⊤QτXk) - tτ), and ck = γ2 tr(Xk ⊤(A + λγ I)Xk), then R( ˜X + Xk) ≤ Hk ( ˜X,Xk)

and the equality holds iff ˜X = 0.

Obviously, Hk ( ˜X,X) is convex w.r.t. ˜X, so it is a convex surrogate of R. Besides, from Proposition 1, it can also been seen that R( ˜X + Xk) ≤ Hk ( ˜X;Xk) for any ˜X, 0 = argmin ˜X (Hk ( ˜X;Xk) - R( ˜X + Xk)) and R(Xk) = Hk (0;Xk). Thus, all three desired properties on surrogate of MM (Section 3.1) are satisfied.

Remark 1. MM algorithm, recently, has also been considered in RMF-MM [Lin et al., 2017] and RMFNL [Yao and Kwok, 2018] for RMF. While our algorithm is also based on MM and adopt matrix factorization, since our objective comes from low-rank SDP, the way to construct the surrogate is significantly different. Specifically, in [Lin et al., 2017; Yao and Kwok, 2018], Z is factorized as XY ⊤, X and Y are bounded separately; we need to find other ways to bound X here, which is enabled by Lemma 2.

3.3 Solving the Surrogate

According to the framework of MM algorithm, i.e., (7), at each iteration we need to update Xk by Xk+1 = Xk + argmin ˜X Hk ( ˜X,Xk). Previously, ADMM (Alternating Direction Method of Multipliers) [Boyd et al., 2011] is used in [Lin et al., 2017] for RMF, and APG (Accelerated Proximal Gradient) [Beck and Teboulle, 2009] is later proposed in [Yao and Kwok, 2018] to further explore data sparsity in RMF. Here, we show both ADMM and APG can still be applied to solve the optimization problem argmin ˜X Hk ( ˜X,Xk). In Table 1, we can see APG can need less space than ADMM when E is sparse.

Using ADMM

Since terms in the ℓ1 loss is complex, we reformulate Hk as

Then, we can introduce dual parameter pτ for each linear constraint, and use ADMM algorithm [Boyd et al., 2011; Lin et al., 2017] to solve the augmented Lagrangian of (8). It can be easily checked that updates of ˜X, eτ and pτ have closed-form solutions.

Using APG

In the other way, using ‖x‖1 = infz x⊤z : ‖z‖∞ ≤ 1 [Boyd and Vandenberghe, 2004]. We can derive the dual form Dk of Hk from

which is given by

where P(zτ) = (∑mτ=1 zτQτ + γ2 C)Xk. Since the dual problem is a smooth and convex optimization problem with simple box constraints. Thus, as [Yao and Kwok, 2018], we also solve this problem by APG, and recover ˜X = -B-1 (∑mτ=1 zτQτ + γ2 C)Xk. Let Q = ∑mτ=1 Qτ , it only requires nnz(Q) entries to store dual variables zτ 's when using APG.

3.4 Complete Algorithm

Based on the above analysis, we list the complete steps for solving (3) in Algorithm 1.

Algorithm 1 RSDP: Robust semi-definite programming by majorization-minimization.

• 1: Initialization: X1 = 0.
• 2: for k = 1, . . . ,K do
• 3: ˜Xk =argmin ˜X H( ˜X,Xk) via ADMM or APG;
• 4: update Xk+1 = ˜Xk +Xk;
• 5: end for
• 6: return XK+1.

The convergence guarantee is in Theorem 1. Note that, as in Section 3.1, MM generally only guarantees the convergence of {R(Xk)} not {Xk}. Besides, the proofs in RMFMM and RMFNL cannot be directly applied neither, due to the difference in Remark 1.

Theorem 1. If lim‖X‖F→∞R(X) = ∞ and infX R(X) > -∞, then for Algorithm 1, we have

• (a). there exists a constant α > 0 such that R(Xk) - R(Xk+1) ≥ α2 ‖Xk+1 - Xk‖2 ;
• (b). the sequence {Xk} is bounded;
• (c). any limit points of {Xk} are also critical points of R.

An overall comparison of the proposed Algorithm 1 and other algorithms used in Section 4 is summarized in Table 1.

4 Empirical Study

In this section, we perform experiments on three applications of SDP, namely, robust NPKL (Section 4.1), robust MVU (Section 4.2), and sparse PCA (Section 4.3). These are also three applications we discussed in Section 2.2. The following algorithms based on squared loss will be compared:

1. FW [Laue, 2012]: an application of Frank-Wolf algorithm [Jaggi, 2013] in SDP problem (1);
2. L-BFGS [Nocedal andWright, 2006]: solve (2) with the most commonly used quasi-Newton solver for smooth minimization problem;
3. nmAPG [Li and Lin, 2015]: an application of state-ofthe-art accelerated gradient descent algorithm for problem (2);

The following algorithms based on ℓ1-loss are compared:

1. SADMM [Boyd et al., 2011]: replace the squared loss in (1) by ℓ1 loss, and solve the resulting nonsmooth but convex problem with ADMM;
2. SDPLR [Burer and Monteiro, 2003]: solve the same problem as SADMM, but the SDPLR package is used;
3. SDPNAL [Toh et al., 2015]: solve the same problem as SADMM, but SDPNAL package (a newton-CG augmented Lagrangian method ) is used;
4. RSDP(ADMM): the proposed Algorithm 1 for the robust objective (3); and ADMM is used as solver for the convex surrogate;
5. RSDP(APG): same as RSDP(ADMM) but APG is used as the solver instead of ADMM.

All algorithm is stopped when the relative change of objective values in successive iterations is smaller than 10-5 or when the number of iterations reaches 2000. As for the rank r of initial solution X, in Sections 4.1 and 4.2 we follow [Burer and Monteiro, 2003] and set its value to be the largest r satisfying r(r + 1) ≤ m, where m is the total number of observed data (i.e, m is the number of must-link and cannotlink pairs in Section 4.1, the number of given neighbor pairs in Section 4.2 respectively). In Section 4.3 we set r = 10. Finally, all algorithms are implemented in Matlab run on a PC with a 3.07GHz CPU and 24GB RAM. To reduce statistical variability, all results are averaged over five repetitions. Availability of codes and data sets are in Appendix.B.

4.1 Robust NPKL

Experiments are performed on the adult data sets that has been commonly used as benchmark data about NPKL learning [Zhuang et al., 2011]. Let the number of training samples be ¯n, we randomly sample 6¯n pairs and construct set T = {(Qτ , tτ)}, i.e., |T | = 6¯n. We randomly sample 20% pairs from T for training, 20% for validation, and the rest for testing. For performance evaluation, we follow [Lin et al., 2017; Yao and Kwok, 2018] and use the (i) testing root mean square error, RMSE = (∑¯nt τ=1 (tr( ¯X⊤Qτ ¯X) - tτ )2 /¯nt)1/2 , where ¯X is the output of the algorithm, ¯nt is a the number of the testing pairs; and (ii) CPU time (sec).

Robustness Against Gaussian Noise

Gaussian noise is the most natural noise type, here, to test the robustness of RSDP against such noise. Specifically, we randomly sample respectively 5% or 10% pairs from training pairs; and for selected pairs, all their labels tτ 's are added with Gaussian noise N(0, 5). Table 2 2 shows the performance of the all compared algorithms. As can be seen, while squared loss based algorithms, i.e., FW, nmAPG and L-BFGS are very fast, they produce much higher testing RMSEs than those based on the ℓ1-loss, i.e, ADMM, SDPNAL, SDPLR and RSDP. Among algorithms for the ℓ1-loss, RSDP(APG) is the fastest and RSDP(ADMM) is the second fastest due to the usage of factorization. These demonstrate the robustness of the proposed formulation and the efficiency of the proposed robust SDP algorithms.

2 For all tables in the sequel, the best and comparable results according to the pair-wise 95% significance test are high-lighted.

Robustness Against Flipping Labels

In real applications, some attackers want to deteriorate the system's learning performance by directly flipping the labels [Raykar et al., 2010]. This is also the worst case of label noise. To further test the robustness our method, we consider such scenario here. Specifically, we take 5% or 10% pairs from training pairs; for selected pairs, all their labels tτ 's are reverse (i.e., making tτ = 1-tτ ). Table 2 shows performance of the all compared algorithms. As can be seen, all algorithms based on ℓ1-loss again produces lower testing RMSE than the squared loss based algorithms; and RSDP is faster than other ℓ1-loss based algorithms, i.e., SADMM, SDPNAL and SDPLR.

4.2 Robust CMVU

Newsgroups 20, is used here. As in [Song et al., 2008], we construct the set N by considering the 1% nearest neighbor pairs of each point. We first construct 1%n square Euclid distances {dτ}'s as reference 'labels' using noiseless data, then add random or outlier to data in MVU modle to output embdding ¯X by all compared algorithms. The tradeoff parameter γ is set to 0.01 as a default. For performance evaluation, same as in robust NPLK, we also use the testing RMSE and CPU time (sec).

Robustness Again Small Deviations

We add Gaussian noise N (0, 0.01¯x), where ¯x is a vector contains means of each features, to Newsgroups 20 data in each dimension. Tabel 3 shows the performance. The observations are the same as that of Section 4.1. We can see that RSDP produces much lower testing RMSE and is the most efficient among algorithms working with the ℓ1 loss.

Robustness Against Large Outliers

Except small deviations, there also can be outliers in the data. Here, we randomly smaple 5% and 10% data points respectively, for the selected each data points xi, we convert it into outliers by adding large random noise, N(0, 5˜x) where ˜x is a vector made from largest elements in the absolute value among each feature. Tabel 3 shows the performance against outliers. Again, by adopting factorization and the ℓ1-loss, RSDP is not only efficient but also achieves the lowest testing RMSE.

4.3 Sparse PCA

We use the colon cancer data set which contains 2000 microarray readings from 62 subjects. In order to vary the problem dimension n, we randomly sampled readings. We set γ = 10 to obtain sparse solution. Results are reported in Table 4 (including the running CPU time (sec), the objective value at convergence, the sparsity of the solution and the captured variance, and these measurements have been used in [Laue, 2012; D'aspremont et al., 2007]). As can been seen, all compared algorithms produce close solutions. The RSDP(ADMM) is comparable with the state-of-the-arts, and RSDP(APG) is much faster than the others when n is large.

5 Conclusion

In this paper, we propose a robust formulation of semidefinite programing (SDP) by replacing commonly used squared loss with ℓ1 loss in standard SDP applications. As the resulting optimization problem is neither convex nor smooth, where existing SDP algorithms cannot be applied, we propose a new algorithm based on majorization-minimization. The algorithm is not efficient, but also guaranteed converging to some critical points. Finally, we demonstrate the efficiency and robustness over state-of-the-arts SDP solvers using three applications as kernel learning, matrix variance unfolding and sparse PCA. As a future work, we plan to automatically select the choice of loss function by automated machine learning [Yao et al., 2018].

Acknowledgments

This research was supported in part by the National Natural Science Foundation of China (No. 61663049, 61165012 and 2019010108203008).