## 1 Introduction  
In many real-world applications, the input data is usually high-dimensional. On one hand, this is a serious challenge for storage and computation. On the other hand, it makes a lot of machine learning algorithms unworkable due to the curse of dimensionality [1]. Means of obtaining a concise and informative data representation for high-dimensional data has become a highly significant focus. Matrix factorization is one kind of popular and effective model of data representation, and finds two or more low-rank matrix factors and their product can well approximate the data matrix. Various matrix factorization methods have been proposed, adopting different constraints on matrix factors. The classical matrix factorization models include Principal Component Analysis (PCA) [2], Singular Value Decomposition (SVD), QR decomposition, vector quantization. Among the various matrix factorization approaches, Non-negative Matrix Factorization (NMF) [3] is a promising one. In NMF, data matrix \(\mathbf{X}\) is decomposed into a non-negative basic matrix \(\mathbf{U}\) which reveals the latent semantic structure, and a non-negative encoding matrix \(\mathbf{V}\), which denotes a new representation with respect to the basis matrix. Because of the non-negative constraints, NMF only allows pure additive combinations, and leads to a parts-based representation. Due to its psychological and physiological interpretation, NMF and its variants have been widely used in computer vision [3], pattern recognition [4], image processing [5], document analysis [6]. Standard NMF performs factorization in Euclidean space. It is unable to discover geometrical structures in data space, which is critical in real-world applications. Therefore, lots of recent work has focused on preserving the intrinsic geometry of the data space by adding different constraints to the objective function of NMF. [PERSON] et al. [7] proposed graph regularized NMF (GNMF) by constructing a nearest neighbor graph while preserving the local geometrical information of the data space. [PERSON] et al. [8] proposed Manifold Regularized Sparse NMF for hyperspectral unmixing, in which manifold regularization was introduced into sparsity-constrained NMF for unmixing. [PERSON] et al. [9] proposed Neighborhood-Preserving Non-negative Matrix Factorization, which imposed an additional constraint on NMF that each item be able to be represented as a linear combination of its neighbors. All the mentioned graph regularized NMF methods construct a graph to encode the geometrical information and use graph Laplacian as a smooth operator. Despite the successful application of graph Laplacian in semi-supervised and unsupervised learning, it still suffers from the problems that the solution is biased towards a constant, as well as its lack of extrapolating power [10]. Sparsity regularization methods that focus on selecting the input variables that best describe the output have been widely investigated. [PERSON] [11] proposed a sparse constraint NMF and added the \(\ell_{1}\) norm constraint on the basis and encoding matrices, which were able to discover sparse representations better than those given by standard NMF. [PERSON] et al. [12] proposed Unified Sparse Subspace Learning (USSL) for learning sparse projections by using a \(\ell_{1}\) norm regularizer. The limitation of the \(\ell_{1}\) norm penalty is that it is unable to guarantee successful models in cases of categorical predictors, for the reason that each dummy variable is selected independently [13]. So \(\ell_{1}\) norm is not feasible for conducting feature selection. To settle this issue, [PERSON] et al. [14] proposed a robust feature selection approach by imposing \(\ell_{2,1}\) norm on loss functions. [PERSON] et al. [15] proposed \(\ell_{2,1}\) norm regularized discriminative feature selection for unsupervised learning. [PERSON] et al.