# Constrained Dual Graph Regularized Orthogonal Nonnegative Matrix Tri-Factorization for Co-Clustering

Authors: Shaodi Ge; Hongjun Li; Liuhong Luo
Journal: Mathematical Problems in Engineering
 (2019)

---

## Abstract

Coclustering approaches for grouping data points and features have recently been receiving extensive attention. In this paper, we propose a constrained dual graph regularized orthogonal nonnegative matrix trifactorization (CDONMTF) algorithm to solve the coclustering problems. The new method improves the clustering performance obviously by employing hard constraints to retain the priori label information of samples, establishing two nearest neighbor graphs to encode the geometric structure of data manifold and feature manifold, and combining with biorthogonal constraints as well. In addition, we have also derived the iterative optimization scheme of CDONMTF and proved its convergence. Clustering experiments on 5 UCI machine-learning data sets and 7 image benchmark data sets show that the achievement of the proposed algorithm is superior to that of some existing clustering algorithms.


---

## Body

## 1. Introduction

Clustering analysis divides data sample points into clusters according to their similarity, and it has been commonly applied in data mining [1–4], machine learning [5], and image processing [6]. And coclustering can cluster rows and columns of the data matrix simultaneously, thereby capturing the binary relationship between data and features and obtaining better clustering results [7, 8].The effectiveness of nonnegative matrix factorization (NMF) [9] has been proved to be effective in solving clustering and coclustering problems [10–13]. A nonnegative data matrix is decomposed into two or three nonnegative factors, and their product approximately linearizes the original matrix. Since the size of each factor matrix is much smaller than the original matrix, NMF provides a parts-based new representation to the original data [9, 14, 15]. Due to the nonnegative constraints of factor matrices, NMF with the additive combination can approximately match the procedure of processing raw data by human cognition and usually produces sparse data representation to effectively mine the potential structure of data [16, 17].According to the needs of practical clustering tasks, several improved methods based on NMF are proposed. Deng et al. [18] put forward an unsupervised learning algorithm called graph-regularized NMF (GNMF) that creates a nearest neighbor graph to retain the neighborhood information of high-dimensional data in the low-dimensional representation space. It is worthy to be mentioned that graph Laplacian regularizer has been included in various NMF frameworks [19–21] after its first application in GNMF. Peng et al. [22] proposed NMF with integrated graph and feature learning (NMF2L) by introducing an error matrix to distinguish features and constructing a data graph to integrate feature learning and manifold learning. Shang et al. [23] presented an algorithm named as local discriminative based sparse subspace learning for feature selection, which introduces a local discriminant model into the feature selection framework of subspace learning. In addition, many researchers have found that semisupervised algorithms can enhance learning accuracy even with only a small amount of label information [24–26]. On this basis, Liu et al. introduced a constrained NMF (CNMF) algorithm for the scenario of limited label data [27], in which existing label information could serve as an additional hard constraint. Meng et al. [28] proposed a semisupervised graph regularized deep NMF with biorthogonal constraints (SGDNMF). SGDNMF can study the hidden layer representation according to the unknown and various attributes of the original high-dimensional data to improve the effect of dimensionality reduction and the clustering performance simultaneously.However, the aforementioned approaches focus mainly on one-side clustering, i.e., grouping similar objects. Coclustering methods involve the relationship between features and data points to group the data of columns (data) and rows (features) at the same time [29–33]. Del Buono and Pio [34] improved the average intercluster similarity by optimizing block matrix. Shang et al. [35] developed a graph dual regularization NMTF (DNMTF) algorithm that constructs two associated graphs under the objective function of NMF and achieves better clustering performance. Meng et al. [36] proposed a dual-graph regularized NMF with sparse and orthogonal constraints (SODNMF) and obtained encouraging clustering results. Wang and Huang [7] suggested the penalized NMTF algorithm (PNMT), in which three penalty terms are introduced to ensure that the clustering indicator matrix is close to orthogonality. Salehani et al. [19] combined graph regular constraints with L2,1 norm as additional conditions for revealing the intrinsic geometry and discriminant form of feature and data spaces and proposed a robust graph-based nonnegative block value decomposition (GSR-NBVD).Considering the advantages of the algorithms mentioned above, we in this paper propose a new algorithm, constrained dual graph regularized orthogonal nonnegative matrix trifactorization (CDONMTF). The main highlights of the proposed method are as follows:(i)
The priori label information of original data samples, intrinsic geometric structures in both data manifold and feature manifold, and biorthogonal constraints are integrated into a unified framework(ii)
An improved optimization problem is constructed by adopting the priori label information, geometric structures of data and feature, and biorthogonal constraints as well(iii)
The multiplicative update rules are derived to solve the optimization problem, and it is theoretically proved that the objective function does not increase under the update rules(iv)
The performance of the proposed algorithm is verified by clustering experiments on some open data setsThe following sections of this paper are presented as follows. Section2 briefly reviews the related works and explicitly describes the proposed CDONMTF algorithm, including the multiplicative update rules and convergence analysis. Section 3 includes the clustering experiments and the result analysis. Finally, the summary is made in Section 4.


## 2. CDONMTF

This section reviews several NMF-based algorithms and describes the proposed CDONMTF algorithms in detail. We then derive the multiplicative iterative rules for the corresponding optimization problem and analyze its convergence.


### 2.1. Related Works

Given a nonnegative data matrixX∈R+m×n, in which xi=X:,i∈R+m×1,i=1,2,…,n, are m-dimensional nonnegative column vectors, typically representing data points. yi=Xi,:T∈R+n×1,i=1,2,…,m, are n-dimensional nonnegative row vectors, typically representing features. The coclustering method can group the data space into p clusters, as well as the feature space into q clusters [7].A NMF algorithm can find two nonnegative matricesU∈R+m×r and V∈R+n×r, whose product UVT can be approximate to X. It can be further understood that individual original data point xi can be explicated as a linear combination of column vectors in U, which can be thought of as a base matrix, and matrix VT is an encoding matrix.The error between the matrixX and the product UVT is defined as a metric function DNMFXUVT, which should be minimized and usually called an objective function. A commonly used objective function is based on the square of the Frobenius norm [9]:(1)DNMFXUVT=X−UVTF2s.t.U≥0,V≥0.Shang et al. [35] considered the geometric structure of both data manifold and feature manifold by establishing two associated graphs embedded in the objective function of NMF and proposed a clustering method called DNMTF. The objective function is as follows:(2)DDNMTFXUSVT=X−USVTF2+αTrVTLVV+βTrUTLUUs.t.U≥0,S≥0,V≥0,where S is an auxiliary matrix, Tr⋅ denotes the trace of the matrix, LV and LU are graph Laplacians of the data graph and the feature graph, respectively, and α and β are regularized parameters.It is usually infeasible to obtain labels for all samples in a data set, but the cost of obtaining a small portion of data labels is relatively low. Under this consideration, Liu et al. [27] extended NMF to a semisupervised learning algorithm called constrained NMF (CNMF).The objective function of CNMF is(3)DCNMFXUVT=X−UZTATF2s.t.U≥0,Z≥0,where A represents a hard constraints matrix and Z is an auxiliary matrix.


### 2.2. Problem Formulation

DNMTF maintains the geometric structure of data manifold and feature manifold through two nearest neighbor graphs, and CNMF makes use of priori label information of original data effectively through hard constraints matrix. The clustering results of these two methods are much better than that of classical NMF methods. We try to integrate the idea of graph embedding and semisupervised learning and utilize biorthogonal constraints as an independent penalty to build the CDONMTF algorithm.Firstly, we construct a data graph and a feature graph to simulate the geometric structures of data space and feature space [37]. Following the graph embedding formula, the inner compactness of data is characterized by the nearest neighbor graph. Sk is used to measure the smoothness of xi and xj in the low-dimensional representation space:(4)Sk=12∑i,j=1nfkxi−fkxj2WijV=∑i,j=1nfkxi2DiiV−∑i,j=1nfkxifkxjWijV=∑i,j=1nviTviDiiV−∑i,j=1nviTvjWijV=TrVTDVV−TrVTWVV=TrVTLVV,where fkxi=vi is the mapping function that maps the xi to its corresponding representation vi in the low-dimensional space. WV is the weight matrix to measure the pairwise similarity between the original data points, DV is a diagonal matrix defined by DiiV=∑jWijV, and LV=DV−WV [16] is the Laplacian matrix of data graph. Sk should be minimized to ensure the smoothness of data in the low-dimensional spatial relationship and make vi and vj as close as possible.The construction of weight matrixWV is an important factor affecting the performance of algorithms. In recent years, many scholars have proposed different methods to assign values to WV, such as robust adaptive neighbor [38], point product weighting [18], and kernel method [39, 40]. For simplicity, we use the 0-1 weighted scheme, but it is easy to be extended to other methods. Here, a k-nearest neighbor data graph with vertices corresponding to [x1,x2,…,xn] is constructed and the following data weight matrix is defined:(5)WijV=1,ifxi∈Nkxjorxj∈Nkxi,0,otherwise,where Nkxj represents the set of k-nearest neighbors of data point xj.Like the construction of the data graph, a k-nearest neighbor feature graph is also constructed by using the 0-1 weighted scheme and its vertices correspond toy1,y2,…,ym. Then, the feature weight matrix is defined:(6)WijU=1,ifyi∈Nkyjoryj∈Nkyi,0,otherwise,where DV is a diagonal matrix and DiiU=∑jWijU, and LU=DU−WU is the Laplacian matrix of the feature graph. Hence, TrUTLUU is obtained as the regularizer term.Assume thatn data points {xi,i=1,2,…,n} belong to c classes, in which the first l samples are labeled and the remaining n−l samples are unlabeled. We can construct an indicator matrix C∈R+l×c:(7)Cij=1,ifxiis labeled  withjthclass,0,otherwise.Then, a label hard constraints matrixA based on C is defined:(8)A=Cl×c00In−l,where A∈R+n×n−l+c and In−l is an n−l×n−l identity matrix.Combining the above two aspects with independent biorthogonal constraints, we can get the following objective function of CDONMTF:(9)DCDONMTFXUSVT=X−USVTF2+αTrVTLVV+TrUTLUU+λTrSTSs.t.U≥0,S≥0,V≥0,V=AZ,UTU=I,VTV=I,where U∈R+m×p, S∈R+p×q, and V∈R+n×q, α and λ are nonnegative regularized parameters, and TrSTS prevents all elements in U and V being too small.


### 2.3. Optimization and Update Rules

In order to solve the optimization problem, we first convert equation (9) to the following form:(10)DCDONMTFXUSVT=X−USZTATF2+λTrSTS+αTrZTATLVAZ+TrUTLUU+βTrAZEZTAT+TrUEUTs.t.U≥0,S≥0,Z≥0,where E=1¯−I(1¯ is a matrix of ones and I is an identity matrix), TrUEUT=∑i≠juiTuj and TrVEVT=∑i≠jviTvjV=AZ are penalty items to ensure the orthogonality of U and V, and β is a nonnegative orthogonal constraint balance parameter.The objective function equation (10) is not convex for U, S, and Z together. Therefore, it is unrealistic to find a global minimum. So, we try to deduce the multiplicative update rules to look for a local minimum.LetΨ=ψij, ϒ=υij, and Φ=ϕij be the Lagrange multipliers of Uij≥0, Sij≥0, and Zij≥0, respectively. By using matrix properties TrAB=TrBA and TrA=TrAT, the Lagrangian function ζ is(11)ζ=DCDONMTFXUSVT+TrΨUT+TrϒST+TrΦZT=TrXXT−2TrXAZSTUT+TrUSZTATAZSTUT+αTrZTATLVAZ+TrUTLUU+λTrSTS+βTrAZEZTAT+TrUEUT+TrΨUT+TrϒST+TrΦZT.The partial derivatives ofζ with respect to U, Z and S are(12)∂ζ∂U=−2XAZST+2USZTATAZST+2αLUU+2βUE+Ψ,∂ζ∂S=−2UTXAZ+2UTUSZTATAZ+2λS+ϒ,∂ζ∂Z=−2ATXTUS+2ATAZSTUTUS+2αATLvAZ+2βATAZE+Φ.Using the KKT conditions [41] ψijuij=0, υijsij=0, and ϕijzij=0 and substituting LV=DV−WV and LU=DU−WU, the following equations for uij, sij, and zij can be obtained:(13)XAZST+αWUU−USZTATAZST−αDUU−βUEijuij=0,UTXAZ−UTUSZTATAZ−λSijsij=0,ATXTUS+αATWVAZ−ATAZSTUTUS−αATDVAZ−βATAZEijzij=0.Finally, we get the following update rules:(14)U⟵U⋅XAZST+αWUUUSZTATAZST+αDUU+βUE,(15)S⟵S⋅UTXAZUTUSZTATAZ+λS,(16)Z⟵Z⋅ATXTUS+αATWVAZATAZSTUTUS+αATDVAZ+βATAZE,XAATXT,ATA,ATDVA, and ATWVA only need to be calculated before iteration and they are all fixed in the later updating procedure, which greatly reduces the running time. As a summary, we give the iterative update process of solving the optimization problem equation (10) in Algorithm 1.Algorithm 1: CDONMTF.
Input: data matrix X, constraints matrix A, dimension values p and q, the neighbor number k, parameters α, β, and λ, and number of iterations tmax.
Output: U,S, and Z.
(1)
randomly initialize three nonnegative matricesU,S,Z.
(2)
use the 0-1 weighted scheme to construct the k-nearest neighbor data graphWU and feature graph WV.
(3)
update matrixDiiU=∑jWijU,DiiV=∑jWijV,E=1¯−I.
(4)
fort=1,2,…,tmax, do
(a): updateU as
U⟵U⋅XAZST+αWUU/USZTATAZST+αDUU+βUE
(b): updateS as
S⟵S⋅UTXAZ/UTUSZTATAZ+λS
(c): updateZ as
Z⟵Z⋅ATXTUS+αATWVAZ/ATAZSTUTUS+αATDVAZ+βATAZE
End
(5)
end CDONMTF.


### 2.4. Convergence Analysis

We discuss the convergence of CDONMTF under the update rules equations (14)–(16) by the following theorem.Theorem 1.
The objective functionDCDONMTFXUSVT is nonincreasing under the update rules equations (14)–(16).In order to prove Theorem1, the following definition and four Lemmas are introduced first.Definition 1 (see [42]).
G
x
,
x
′ is an auxiliary function for Fx if the conditions(17)Gx,x′≥Fx,Gx,x=Fx,are satisfied.Lemma 1 (see [42]).
IfGx,x′ is an auxiliary function of Fx, then Fx is nonincreasing under the update:(18)xt+1=argminxGx,xt.Proof.
(19)
F
x
t
+
1
≤
G
x
t
+
1
,
x
t
≤
G
x
t
,
x
t
=
F
x
t
.
For any elementzij in Z, let Fijzij be the part about Z in the objective function equation (9).
Then, we can obtain(20)Fij′zij=−2ATXTUS+2ATAZUTUS+2αATLVAZ+2βATAZEij,Fij″zij=2ATAiiUTUSjj+2αATLVAii+2βATAiiEjj.Lemma 2.
The function(21)Gzij,zijt=Fijzijt+Fij′zij−zijt+ATAZUTUS+αATDVAZ+βATAZEijzijtzij−zijt2,is an auxiliary function for Fijzij.Proof.
The Taylor series expansion ofFijzij is defined as below:(22)Fijzij=Fijzijt+Fij′zij−zijt+ATAiiUTUSjj+αATLVAii+βATAiiEjjzij−zijt2.
Comparing equation (21) with equation (22), it is easy to find that Gzij,zijt≥Fijzij is equivalent to(23)ATAZUTUS+αATDVAZ+βATAZEijzijt≥ATAiiUTUSjj+αATLVAii+βATAiiEjj,since(24)ATAZUTU+βATAZEij=∑l=1rAATZilUTUSlj+β∑l=1rAATZilElj≥AATZijUTUSjj+βAATZijEjj≥∑l=1rATAilzijtUTUSlj+β∑l=1rATAilzijtElj≥zijtATAiiUTUSjj+βzijtATAiiEjj,αATDVAZij=α∑q=1rATDVAiqzqjt≥αATDVAiizijt≥αATDV−WVAiizijt=αATLVAiizijt.
Thus, equation (23) holds. Then, Gzij,zijt≥Fijzij and obviously Gzij,zij=Fijzij. According to the definition of auxiliary function, Lemma 2 is proved.
Denote the part aboutS in the objective function DCDONMTFXUSVT as Fijsij for any element sij in S, we can obtain(25)Fij′sij=−2UTXAZ+2UTUSZTATAZ+2λSij,Fij″sij=2UTUiiZTATAZjj+2λ.Lemma 3.
The function(26)Gsij,sijt=Fijsijt+Fij′sij−sijt+UTUSZTATAZ+λSijsijtsij−sijt2,is an auxiliary function for Fijsij.Proof.
ObviouslyGs,s=Fijs. According to the definition of auxiliary function, the only thing that needs to be proved is Gsij,sijt≥Fijsij. And Fijsij has a Taylor expansion as follows:(27)Fijsij=Fijsijt+Fij′sij−sijt+UTUiiZTATAZjj+λsij−sijt2.
Since(28)UTUSZTATAZ+λSij=∑l=1rUTUSilZTATAZlj+λSij≥UTUSijZTATAZjj+λSij≥∑l=1rUTUilsljtZTATAZjj+λsijt≥sijtUTUiiZTATAZjj+λsijt.
Then,Gsij,sijt≥Fijsij and Lemma 3 holds.
For any elementuij in U, the similar results can be obtained as zij in Z. They are(29)Fij′uij=−2XAZST+2USZTATAZST+2αLUU+2βUEij,Fij″uij=2SZTATAZSTjj+2αLUii+2βEjj.Lemma 4.
The function(30)Guij,uijt=Fijuijt+Fij′uij−uijt+USZTATAZST+αDUU+βUEijuijtuij−uijt2,is an auxiliary function for Fijuij.
The proof is similar to that of Lemma2.With these four Lemmas, the proof of Theorem1 can be carried out as follows.Proof of Theorem 1.
In order to get the minimum values ofGu,u′, Gs,s′, and Gz,z′, we let ∂Guij,uijt/∂uij=0, ∂Gsij,sijt/∂sij=0, and ∂Gzij,zijt/∂zij=0. Then, local minimum points can be obtained as follows:(31)uijt+1⟵uijt⋅XAZ+αWUUijUZTATAZ+αDUU+βUEij,(32)sijt+1⟵sijt⋅UTXAZijUTUSZTATAZ+λSij,(33)zijt+1⟵zijt⋅ATXTU+αATWVAZijATAZUTU+αATDVAZ+βATAZEij.
Equations (31)–(33) have the same form as equations (14)–(16) after rewriting while Guij,uijt, Gsij,sijt, and Gzij,zijt are the auxiliary functions of Fijuij, Fijsij, and Fijzij. It is easy to deduce that Fijuij, Fijsij, and Fijzij are not increased under the update rules of equations (14)–(16), so Theorem 1 ensures the convergence of our CDONMTF algorithm.


### 2.5. Connection to the Gradient Descent Method

Another common method of minimizing the objective functionDCDONMTFXUSVT is the gradient descent method [43]. Its update rules are(34)uij⟵uij+ηij∂D∂uij,sij⟵sij+θij∂D∂sij,zij⟵zij+φij∂D∂zij,where ηij, θij, and φij are the step parameters to control the gradient descent but usually difficult to set. However, due to the special form of partial derivatives, the step parameters can be set automatically by some techniques. Let(35)ηij=−uij2UZTATAZ+αDUU+βUEij,θij=−sij2UTUSZTATAZ+λSij,φij=−zij2ATAZUTU+αATDVAZ+βATAZEij.We can obtain(36)uij⟵uij+ηij∂D∂uij=uijXAZ+αWUUijUZTATAZ+αDUU+βUEij,sij⟵sij+θij∂D∂sij=sijUTXAZijUTUSZTATAZ+λSij,zij⟵zij+φij∂D∂zij=zij⋅ATXTU+αATWVAZijATAZUTU+αATDVAZ+βATAZEij.It can be seen that the multiplicative update rules in Equations (14)–(16) are special forms of the step parameters automatically selected by the gradient descent method, so the nonnegativity of U, S, and Z can be guaranteed simultaneously.


### 2.6. Computational Complexity

According to the iterative update rules, it is not difficult to calculate the arithmetic operations of DNMTF and CDONMTF in each iteration. The relative results of NMF and CNMF [27] are also listed here in Table 1. The complexity of each algorithm is Omnp. In addition to calculating the hard constraints matrix AV=AZ, the proposed CDONMTF framework is similar to the corresponding DNMTF structure [35]. The construction of k-nearest neighbor feature graph and data graph needs Omn2+nm2, and WU and WV are sparse. Both DU and DV are diagonal matrices, and A is a sparse matrix with only n nonzero elements. Therefore, it can be calculated effectively. When the number of iterations is t, the overall cost for CDONMTF is(37)Otmnp+mn2+nm2.Table 1
Arithmetic operation counts for each iteration.
Floating point operation
Addition (+)
Multiplication (×)
Division (÷)
Overall
NMF
2
m
n
p
+
2
m
+
n
p
2
2
m
n
p
+
2
m
+
n
p
2
+
m
+
n
p
(m + n)p
O (mnp)
CNMF [27]
2
m
+
2
n
−
l
+
c
n
−
l
+
c
p
+
2
m
+
n
p
2
2
m
+
2
n
−
l
+
c
n
−
l
+
c
p
+
2
m
+
n
p
2
+
m
+
n
−
l
+
c
p
m
+
n
−
l
+
c
p
O (mnp)
DNMTF [35]
5
m
n
p
+
7
m
p
2
+
4
n
p
2
+
m
+
n
k
+
3
p
+
p
2
2
p
+
1
5
m
n
p
+
7
m
p
2
+
4
n
p
2
+
m
+
n
k
+
2
p
+
p
2
2
p
+
1
m
+
n
+
p
p
O (mnp)
CDONMTF
m
n
2
p
+
q
+
m
p
3
p
+
4
q
+
k
+
4
+
n
q
4
p
+
3
q
+
k
+
4
+
2
m
q
2
+
n
p
2
+
12
n
n
−
l
+
c
q
+
p
q
m
n
2
p
+
q
+
m
p
3
p
+
4
q
+
k
+
2
+
n
q
4
p
+
3
q
+
k
+
2
+
2
m
q
2
+
n
p
2
+
12
n
n
−
l
+
c
q
+
p
q
m
p
+
p
+
n
−
l
+
c
q
O (mnp)


## 2.1. Related Works

Given a nonnegative data matrixX∈R+m×n, in which xi=X:,i∈R+m×1,i=1,2,…,n, are m-dimensional nonnegative column vectors, typically representing data points. yi=Xi,:T∈R+n×1,i=1,2,…,m, are n-dimensional nonnegative row vectors, typically representing features. The coclustering method can group the data space into p clusters, as well as the feature space into q clusters [7].A NMF algorithm can find two nonnegative matricesU∈R+m×r and V∈R+n×r, whose product UVT can be approximate to X. It can be further understood that individual original data point xi can be explicated as a linear combination of column vectors in U, which can be thought of as a base matrix, and matrix VT is an encoding matrix.The error between the matrixX and the product UVT is defined as a metric function DNMFXUVT, which should be minimized and usually called an objective function. A commonly used objective function is based on the square of the Frobenius norm [9]:(1)DNMFXUVT=X−UVTF2s.t.U≥0,V≥0.Shang et al. [35] considered the geometric structure of both data manifold and feature manifold by establishing two associated graphs embedded in the objective function of NMF and proposed a clustering method called DNMTF. The objective function is as follows:(2)DDNMTFXUSVT=X−USVTF2+αTrVTLVV+βTrUTLUUs.t.U≥0,S≥0,V≥0,where S is an auxiliary matrix, Tr⋅ denotes the trace of the matrix, LV and LU are graph Laplacians of the data graph and the feature graph, respectively, and α and β are regularized parameters.It is usually infeasible to obtain labels for all samples in a data set, but the cost of obtaining a small portion of data labels is relatively low. Under this consideration, Liu et al. [27] extended NMF to a semisupervised learning algorithm called constrained NMF (CNMF).The objective function of CNMF is(3)DCNMFXUVT=X−UZTATF2s.t.U≥0,Z≥0,where A represents a hard constraints matrix and Z is an auxiliary matrix.


## 2.2. Problem Formulation

DNMTF maintains the geometric structure of data manifold and feature manifold through two nearest neighbor graphs, and CNMF makes use of priori label information of original data effectively through hard constraints matrix. The clustering results of these two methods are much better than that of classical NMF methods. We try to integrate the idea of graph embedding and semisupervised learning and utilize biorthogonal constraints as an independent penalty to build the CDONMTF algorithm.Firstly, we construct a data graph and a feature graph to simulate the geometric structures of data space and feature space [37]. Following the graph embedding formula, the inner compactness of data is characterized by the nearest neighbor graph. Sk is used to measure the smoothness of xi and xj in the low-dimensional representation space:(4)Sk=12∑i,j=1nfkxi−fkxj2WijV=∑i,j=1nfkxi2DiiV−∑i,j=1nfkxifkxjWijV=∑i,j=1nviTviDiiV−∑i,j=1nviTvjWijV=TrVTDVV−TrVTWVV=TrVTLVV,where fkxi=vi is the mapping function that maps the xi to its corresponding representation vi in the low-dimensional space. WV is the weight matrix to measure the pairwise similarity between the original data points, DV is a diagonal matrix defined by DiiV=∑jWijV, and LV=DV−WV [16] is the Laplacian matrix of data graph. Sk should be minimized to ensure the smoothness of data in the low-dimensional spatial relationship and make vi and vj as close as possible.The construction of weight matrixWV is an important factor affecting the performance of algorithms. In recent years, many scholars have proposed different methods to assign values to WV, such as robust adaptive neighbor [38], point product weighting [18], and kernel method [39, 40]. For simplicity, we use the 0-1 weighted scheme, but it is easy to be extended to other methods. Here, a k-nearest neighbor data graph with vertices corresponding to [x1,x2,…,xn] is constructed and the following data weight matrix is defined:(5)WijV=1,ifxi∈Nkxjorxj∈Nkxi,0,otherwise,where Nkxj represents the set of k-nearest neighbors of data point xj.Like the construction of the data graph, a k-nearest neighbor feature graph is also constructed by using the 0-1 weighted scheme and its vertices correspond toy1,y2,…,ym. Then, the feature weight matrix is defined:(6)WijU=1,ifyi∈Nkyjoryj∈Nkyi,0,otherwise,where DV is a diagonal matrix and DiiU=∑jWijU, and LU=DU−WU is the Laplacian matrix of the feature graph. Hence, TrUTLUU is obtained as the regularizer term.Assume thatn data points {xi,i=1,2,…,n} belong to c classes, in which the first l samples are labeled and the remaining n−l samples are unlabeled. We can construct an indicator matrix C∈R+l×c:(7)Cij=1,ifxiis labeled  withjthclass,0,otherwise.Then, a label hard constraints matrixA based on C is defined:(8)A=Cl×c00In−l,where A∈R+n×n−l+c and In−l is an n−l×n−l identity matrix.Combining the above two aspects with independent biorthogonal constraints, we can get the following objective function of CDONMTF:(9)DCDONMTFXUSVT=X−USVTF2+αTrVTLVV+TrUTLUU+λTrSTSs.t.U≥0,S≥0,V≥0,V=AZ,UTU=I,VTV=I,where U∈R+m×p, S∈R+p×q, and V∈R+n×q, α and λ are nonnegative regularized parameters, and TrSTS prevents all elements in U and V being too small.


## 2.3. Optimization and Update Rules

In order to solve the optimization problem, we first convert equation (9) to the following form:(10)DCDONMTFXUSVT=X−USZTATF2+λTrSTS+αTrZTATLVAZ+TrUTLUU+βTrAZEZTAT+TrUEUTs.t.U≥0,S≥0,Z≥0,where E=1¯−I(1¯ is a matrix of ones and I is an identity matrix), TrUEUT=∑i≠juiTuj and TrVEVT=∑i≠jviTvjV=AZ are penalty items to ensure the orthogonality of U and V, and β is a nonnegative orthogonal constraint balance parameter.The objective function equation (10) is not convex for U, S, and Z together. Therefore, it is unrealistic to find a global minimum. So, we try to deduce the multiplicative update rules to look for a local minimum.LetΨ=ψij, ϒ=υij, and Φ=ϕij be the Lagrange multipliers of Uij≥0, Sij≥0, and Zij≥0, respectively. By using matrix properties TrAB=TrBA and TrA=TrAT, the Lagrangian function ζ is(11)ζ=DCDONMTFXUSVT+TrΨUT+TrϒST+TrΦZT=TrXXT−2TrXAZSTUT+TrUSZTATAZSTUT+αTrZTATLVAZ+TrUTLUU+λTrSTS+βTrAZEZTAT+TrUEUT+TrΨUT+TrϒST+TrΦZT.The partial derivatives ofζ with respect to U, Z and S are(12)∂ζ∂U=−2XAZST+2USZTATAZST+2αLUU+2βUE+Ψ,∂ζ∂S=−2UTXAZ+2UTUSZTATAZ+2λS+ϒ,∂ζ∂Z=−2ATXTUS+2ATAZSTUTUS+2αATLvAZ+2βATAZE+Φ.Using the KKT conditions [41] ψijuij=0, υijsij=0, and ϕijzij=0 and substituting LV=DV−WV and LU=DU−WU, the following equations for uij, sij, and zij can be obtained:(13)XAZST+αWUU−USZTATAZST−αDUU−βUEijuij=0,UTXAZ−UTUSZTATAZ−λSijsij=0,ATXTUS+αATWVAZ−ATAZSTUTUS−αATDVAZ−βATAZEijzij=0.Finally, we get the following update rules:(14)U⟵U⋅XAZST+αWUUUSZTATAZST+αDUU+βUE,(15)S⟵S⋅UTXAZUTUSZTATAZ+λS,(16)Z⟵Z⋅ATXTUS+αATWVAZATAZSTUTUS+αATDVAZ+βATAZE,XAATXT,ATA,ATDVA, and ATWVA only need to be calculated before iteration and they are all fixed in the later updating procedure, which greatly reduces the running time. As a summary, we give the iterative update process of solving the optimization problem equation (10) in Algorithm 1.Algorithm 1: CDONMTF.
Input: data matrix X, constraints matrix A, dimension values p and q, the neighbor number k, parameters α, β, and λ, and number of iterations tmax.
Output: U,S, and Z.
(1)
randomly initialize three nonnegative matricesU,S,Z.
(2)
use the 0-1 weighted scheme to construct the k-nearest neighbor data graphWU and feature graph WV.
(3)
update matrixDiiU=∑jWijU,DiiV=∑jWijV,E=1¯−I.
(4)
fort=1,2,…,tmax, do
(a): updateU as
U⟵U⋅XAZST+αWUU/USZTATAZST+αDUU+βUE
(b): updateS as
S⟵S⋅UTXAZ/UTUSZTATAZ+λS
(c): updateZ as
Z⟵Z⋅ATXTUS+αATWVAZ/ATAZSTUTUS+αATDVAZ+βATAZE
End
(5)
end CDONMTF.


## 2.4. Convergence Analysis

We discuss the convergence of CDONMTF under the update rules equations (14)–(16) by the following theorem.Theorem 1.
The objective functionDCDONMTFXUSVT is nonincreasing under the update rules equations (14)–(16).In order to prove Theorem1, the following definition and four Lemmas are introduced first.Definition 1 (see [42]).
G
x
,
x
′ is an auxiliary function for Fx if the conditions(17)Gx,x′≥Fx,Gx,x=Fx,are satisfied.Lemma 1 (see [42]).
IfGx,x′ is an auxiliary function of Fx, then Fx is nonincreasing under the update:(18)xt+1=argminxGx,xt.Proof.
(19)
F
x
t
+
1
≤
G
x
t
+
1
,
x
t
≤
G
x
t
,
x
t
=
F
x
t
.
For any elementzij in Z, let Fijzij be the part about Z in the objective function equation (9).
Then, we can obtain(20)Fij′zij=−2ATXTUS+2ATAZUTUS+2αATLVAZ+2βATAZEij,Fij″zij=2ATAiiUTUSjj+2αATLVAii+2βATAiiEjj.Lemma 2.
The function(21)Gzij,zijt=Fijzijt+Fij′zij−zijt+ATAZUTUS+αATDVAZ+βATAZEijzijtzij−zijt2,is an auxiliary function for Fijzij.Proof.
The Taylor series expansion ofFijzij is defined as below:(22)Fijzij=Fijzijt+Fij′zij−zijt+ATAiiUTUSjj+αATLVAii+βATAiiEjjzij−zijt2.
Comparing equation (21) with equation (22), it is easy to find that Gzij,zijt≥Fijzij is equivalent to(23)ATAZUTUS+αATDVAZ+βATAZEijzijt≥ATAiiUTUSjj+αATLVAii+βATAiiEjj,since(24)ATAZUTU+βATAZEij=∑l=1rAATZilUTUSlj+β∑l=1rAATZilElj≥AATZijUTUSjj+βAATZijEjj≥∑l=1rATAilzijtUTUSlj+β∑l=1rATAilzijtElj≥zijtATAiiUTUSjj+βzijtATAiiEjj,αATDVAZij=α∑q=1rATDVAiqzqjt≥αATDVAiizijt≥αATDV−WVAiizijt=αATLVAiizijt.
Thus, equation (23) holds. Then, Gzij,zijt≥Fijzij and obviously Gzij,zij=Fijzij. According to the definition of auxiliary function, Lemma 2 is proved.
Denote the part aboutS in the objective function DCDONMTFXUSVT as Fijsij for any element sij in S, we can obtain(25)Fij′sij=−2UTXAZ+2UTUSZTATAZ+2λSij,Fij″sij=2UTUiiZTATAZjj+2λ.Lemma 3.
The function(26)Gsij,sijt=Fijsijt+Fij′sij−sijt+UTUSZTATAZ+λSijsijtsij−sijt2,is an auxiliary function for Fijsij.Proof.
ObviouslyGs,s=Fijs. According to the definition of auxiliary function, the only thing that needs to be proved is Gsij,sijt≥Fijsij. And Fijsij has a Taylor expansion as follows:(27)Fijsij=Fijsijt+Fij′sij−sijt+UTUiiZTATAZjj+λsij−sijt2.
Since(28)UTUSZTATAZ+λSij=∑l=1rUTUSilZTATAZlj+λSij≥UTUSijZTATAZjj+λSij≥∑l=1rUTUilsljtZTATAZjj+λsijt≥sijtUTUiiZTATAZjj+λsijt.
Then,Gsij,sijt≥Fijsij and Lemma 3 holds.
For any elementuij in U, the similar results can be obtained as zij in Z. They are(29)Fij′uij=−2XAZST+2USZTATAZST+2αLUU+2βUEij,Fij″uij=2SZTATAZSTjj+2αLUii+2βEjj.Lemma 4.
The function(30)Guij,uijt=Fijuijt+Fij′uij−uijt+USZTATAZST+αDUU+βUEijuijtuij−uijt2,is an auxiliary function for Fijuij.
The proof is similar to that of Lemma2.With these four Lemmas, the proof of Theorem1 can be carried out as follows.Proof of Theorem 1.
In order to get the minimum values ofGu,u′, Gs,s′, and Gz,z′, we let ∂Guij,uijt/∂uij=0, ∂Gsij,sijt/∂sij=0, and ∂Gzij,zijt/∂zij=0. Then, local minimum points can be obtained as follows:(31)uijt+1⟵uijt⋅XAZ+αWUUijUZTATAZ+αDUU+βUEij,(32)sijt+1⟵sijt⋅UTXAZijUTUSZTATAZ+λSij,(33)zijt+1⟵zijt⋅ATXTU+αATWVAZijATAZUTU+αATDVAZ+βATAZEij.
Equations (31)–(33) have the same form as equations (14)–(16) after rewriting while Guij,uijt, Gsij,sijt, and Gzij,zijt are the auxiliary functions of Fijuij, Fijsij, and Fijzij. It is easy to deduce that Fijuij, Fijsij, and Fijzij are not increased under the update rules of equations (14)–(16), so Theorem 1 ensures the convergence of our CDONMTF algorithm.


## 2.5. Connection to the Gradient Descent Method

Another common method of minimizing the objective functionDCDONMTFXUSVT is the gradient descent method [43]. Its update rules are(34)uij⟵uij+ηij∂D∂uij,sij⟵sij+θij∂D∂sij,zij⟵zij+φij∂D∂zij,where ηij, θij, and φij are the step parameters to control the gradient descent but usually difficult to set. However, due to the special form of partial derivatives, the step parameters can be set automatically by some techniques. Let(35)ηij=−uij2UZTATAZ+αDUU+βUEij,θij=−sij2UTUSZTATAZ+λSij,φij=−zij2ATAZUTU+αATDVAZ+βATAZEij.We can obtain(36)uij⟵uij+ηij∂D∂uij=uijXAZ+αWUUijUZTATAZ+αDUU+βUEij,sij⟵sij+θij∂D∂sij=sijUTXAZijUTUSZTATAZ+λSij,zij⟵zij+φij∂D∂zij=zij⋅ATXTU+αATWVAZijATAZUTU+αATDVAZ+βATAZEij.It can be seen that the multiplicative update rules in Equations (14)–(16) are special forms of the step parameters automatically selected by the gradient descent method, so the nonnegativity of U, S, and Z can be guaranteed simultaneously.


## 2.6. Computational Complexity

According to the iterative update rules, it is not difficult to calculate the arithmetic operations of DNMTF and CDONMTF in each iteration. The relative results of NMF and CNMF [27] are also listed here in Table 1. The complexity of each algorithm is Omnp. In addition to calculating the hard constraints matrix AV=AZ, the proposed CDONMTF framework is similar to the corresponding DNMTF structure [35]. The construction of k-nearest neighbor feature graph and data graph needs Omn2+nm2, and WU and WV are sparse. Both DU and DV are diagonal matrices, and A is a sparse matrix with only n nonzero elements. Therefore, it can be calculated effectively. When the number of iterations is t, the overall cost for CDONMTF is(37)Otmnp+mn2+nm2.Table 1
Arithmetic operation counts for each iteration.
Floating point operation
Addition (+)
Multiplication (×)
Division (÷)
Overall
NMF
2
m
n
p
+
2
m
+
n
p
2
2
m
n
p
+
2
m
+
n
p
2
+
m
+
n
p
(m + n)p
O (mnp)
CNMF [27]
2
m
+
2
n
−
l
+
c
n
−
l
+
c
p
+
2
m
+
n
p
2
2
m
+
2
n
−
l
+
c
n
−
l
+
c
p
+
2
m
+
n
p
2
+
m
+
n
−
l
+
c
p
m
+
n
−
l
+
c
p
O (mnp)
DNMTF [35]
5
m
n
p
+
7
m
p
2
+
4
n
p
2
+
m
+
n
k
+
3
p
+
p
2
2
p
+
1
5
m
n
p
+
7
m
p
2
+
4
n
p
2
+
m
+
n
k
+
2
p
+
p
2
2
p
+
1
m
+
n
+
p
p
O (mnp)
CDONMTF
m
n
2
p
+
q
+
m
p
3
p
+
4
q
+
k
+
4
+
n
q
4
p
+
3
q
+
k
+
4
+
2
m
q
2
+
n
p
2
+
12
n
n
−
l
+
c
q
+
p
q
m
n
2
p
+
q
+
m
p
3
p
+
4
q
+
k
+
2
+
n
q
4
p
+
3
q
+
k
+
2
+
2
m
q
2
+
n
p
2
+
12
n
n
−
l
+
c
q
+
p
q
m
p
+
p
+
n
−
l
+
c
q
O (mnp)


## 3. Experiments and Result Analysis

The performance of the CDONMTF algorithm on some open data sets including 5 UCI machine learning data sets and 7 image data sets is tested. All the experiments are implemented on the Matlab R2014A platform, on an Intel Core i7-8700 (3.2 GHz).


### 3.1. Comparison of Algorithms and Evaluation Criteria

We compare our method to 8 existing representative methods that are similar, in some respects, to ours. These algorithms include GNMF [18], CNMF [27], DNMTF [35], PNMT [7], NMF2L2,0 [22], SODNMF [36], G-NBVD [19], and GSR-NBVD [19].For graph-based algorithms, GNMF, DNMF, G-NBVD,NMF2L2,0, GSR-NBVD, and our CDONMTF, are weighted with 0 and 1, and SODNMF is weighted with heat kernel. The nearest neighbor number is minm,n,5 as recommended in [18]. For CNMF, SODNMF, and CDONMTF, 20% of the samples in each data set are randomly selected as labeled samples. The number of iterations for all algorithms are set to 100.In all experiments, clustering accuracy (AC) [44], normalized mutual information (NMI) [44], and purity (PU) [45] are applied as criteria to evaluate the performance of clustering algorithms.


### 3.2. Parameter Analysis

The CDONMTF model has two types of parameters: decomposition dimension valuesp and q; regularization parameters α,β, and λ. Although many literature studies have studied how to choose the dimension [32, 46, 47] and regularization parameters [7, 48, 49], there is no unique way to select these parameters. For the choice of p and q, we use the following methods [19]:(38)p=θp×c,q=θq×c,where θp and θq are empirical, usually depending on the type of data sets, and c is the number of classes to cluster. In our experiments, θp=θq=3 are selected for UCI data sets and θp=3 and θq=2 for image data sets.For the regularization parametersα and β, their sensitivity is analyzed on 2 representative data sets, JAFFE and SEMEION. Values as α=0.01,0.1,1,10,100,1000 and β=0.01,0.05,0.1,0.5,1,5 are tested, while λ is fixed to 1. Figure 1 displays the clustering results of CDONMTF under different values of α and β. It can be seen that they are more robust for parameter α in a certain range, but more sensitive for parameter β. In the following experiments, we select α=100, β=0.5, and λ=1.Figure 1
Clustering performance of CDONMTF versus the value of parameters. (a) JAFFE. (b) SEMEION. (c) JAFFE. (d) SEMEION. (e) JAFFE. (f) SEMEION.
(a)
(b)
(c)
(d)
(e)
(f)


### 3.3. Experiments on UCI Data

In this section, we test clustering performance of different algorithms on 5 UCI nonnegative data sets (http://archive.ics.uci.edu/ml/index.phpl). Table 2 outlines the basic information of those data sets. The experimental results listed in Table 3 show that the proposed CDONMTF algorithm achieves the best clustering performance in the view of all 3 indexes: AC, NMI, and PU on three data sets (WDBC, lenses, and Chess), while the clustering results on Zoo and soybean data sets are not satisfying but acceptable.Table 2
Description of the five UCI data sets.
Data sets
No. of features (m)
No. of samples (n)
No. of classes (c)
WDBC
30
569
2
Lenses
4
24
3
Zoo
16
101
7
Chess
36
3196
2
Soybean
35
47
4Table 3
Comparison of clustering performance on UCI data sets.
Data sets
GNMF
CNMF
DNMTF
PNMT
NMF
2
L
2,0
SODNMF
G-NBVD
GSR-NBVD
CDONMTF
Accuracy (%)
WDBC
83.83
79.26
82.59
83.46
83.48
90.33
81.02
69.07
92.79
Lenses
41.67
56.17
50.13
57.33
41.67
37.51
41.67
45.83
57.83
Zoo
61.01
67.04
60.93
55.15
84.32
73.39
61.82
62.99
67.19
Chess
55.16
51.69
54.07
52.07
50.99
60.89
55.23
52.94
77.41
Soybean
74.47
72.34
78.72
55.28
89.36
81.02
76.61
78.09
87.23
Normalized mutual information (%)
WDBC
38.34
23.11
35.44
36.81
37.5
56.08
28.63
8.17
61.19
Lenses
2.01
18.52
3.71
13.89
2.02
4.99
2.58
2.52
21.73
Zoo
65.92
70.57
66.02
54.81
80.63
73.59
65.79
66.24
71.04
Chess
0.75
0.15
0.57
0.09
1.17
4.29
0.75
0.36
24.88
Soybean
71.38
71.56
72.44
36.74
81.49
79.47
72.12
72.31
77.85
Purity (%)
WDBC
83.83
79.26
82.61
83.46
83.48
90.33
81.02
69.07
92.79
Lenses
62.51
66.67
62.49
70.83
62.5
62.5
62.51
62.51
70.52
Zoo
82.18
85.39
82.28
74.12
87.05
83.29
82.04
82.36
85.91
Chess
55.16
52.22
54.07
52.22
52.22
60.89
55.23
52.94
77.41
Soybean
78.72
78.72
78.72
59.53
89.36
81.02
78.72
78.72
87.23


### 3.4. Experiments on Image Data

We continue to experiment on different types and sizes of some open image data sets: UMIST (http://images.ee.umist.ac.uk/danny/database.html), JAFFE [50], COIL20 [51], Optdigits (http://archive.ics.uci.edu/ml/machine-learning-datasets/optdigits/), MNIST1 (a subset of the first 10k samples) (http://www.cad.zju.edu.cn/home/dengcai/Data/MLData.html), and SEMEION (https://archive.ics.uci.edu/ml/machine-learning-databases/semeion/) in this part. Both types and sizes of the data sets are listed in Table 4.Table 4
Description of the six image data sets.
Data sets
No. of features (m)
No. of samples (n)
No. of classes (c)
Data types
UMIST
10304
575
20
Face image
JAFFE
4096
213
10
Facial expression
COIL20
1024
1440
20
Object image
Optdigit
64
3823
10
Handwritten digit
MNIST1
784
10000
10
Handwritten digit
SEMEION
256
1593
10
Handwritten digitWe use a subset of the firstc categories as the collection X for clustering. Different cluster numbers c are used to assess the efficiency of the algorithms under various sample scales of X in the experiments.Experiment results (Figures2–4) display the performance of different algorithms using different clustering numbers c. Table 5 shows the average performance of different algorithms under the same data set with different clustering numbers c. The optimal results are marked in bold in the table. Several interesting findings are noted as follows:(i)
SODNMF and the proposed CDONMTF achieve better clustering performance than PNMT, GNMF, DNMTF,NMF2L2,0, G-NBVD, and GSR-NBVD in most cases due to adding the priori label information of the data as an additional constraint(ii)
Compared with CNMF only considering priori label information, the proposed CDONMTF considers the geometric structure information contained in data and feature points and combines with additional biorthogonal constraints, which results in more efficient clustering performance(iii)
From the view of evaluation indexes (AC, NMI, and PU), the average clustering effect of CDONMTF is almost the best on all data sets except JAFFE, but it is also good enoughFigure 2
Comparisons of clustering AC (%) with different number of classes on six data sets. (a). UMIST. (b) JAFFE. (c) COIL20. (d) Optdigit. (e) MNIST1. (f) SEMEION.
(a)
(b)
(c)
(d)
(e)
(f)Figure 3
Comparisons of clustering NMI (%) with different number of classes on six data sets. (a). UMIST. (b) JAFFE. (c) COIL20. (d) Optdigit. (e) MNIST1. (f) SEMEION.
(a)
(b)
(c)
(d)
(e)
(f)Figure 4
Comparisons of clustering PU (%) with different number of classes on six data sets.(a). UMIST. (b) JAFFE. (c) COIL20. (d) Optdigit. (e) MNIST1. (f) SEMEION.
(a)
(b)
(c)
(d)
(e)
(f)Table 5
Comparison of average clustering performance on six image data sets.
Data sets
GNMF
CNMF
DNMTF
PNMT
NMF
2
L
2,0
SODNMF
G-NBVD
GSR-NBVD
CDONMTF
Accuracy (%)
UMIST
43.47
47.17
66.56
37.57
63.65
79.06
64.91
64.86
81.46
JAFFE
87.89
91.33
89.85
41.97
97.66
96.57
88.86
88.85
96.48
COIL20
59.71
62.36
74.71
43.91
72.54
83.29
77.18
77.21
86.11
Optdigit
89.53
88.23
95.46
54.04
92.64
95.06
93.91
95.13
97.77
MNIST1
56.52
64.14
65.41
40.71
72.64
90.21
72.76
67.83
90.88
SEMEION
70.56
73.63
79.16
42.74
76.73
82.69
77.45
79.45
89.27
Normalized mutual information (%)
UMIST
58.33
58.51
79.04
43.17
76.82
86.01
78.73
78.74
87.34
JAFFE
86.41
89.58
90.59
36.93
97.78
93.72
90.04
90.05
95.72
COIL20
64.98
65.42
83.12
49.44
80.22
88.95
85.61
85.56
90.54
Optdigit
86.16
77.18
92.72
41.09
89.88
92.01
93.18
93.55
96.61
MNIST1
45.73
54.31
66.61
29.41
71.34
84.34
71.04
69.27
84.62
SEMEION
68.94
60.12
73.96
33.94
71.83
78.76
73.71
74.15
81.89
Purity (%)
UMIST
50.04
52.25
74.33
42.49
71.05
82.91
73.21
73.19
84.47
JAFFE
88.66
92.15
91.31
44.98
98.06
94.53
90.58
90.57
96.76
COIL20
62.15
64.16
79.69
48.32
77.27
84.73
81.85
81.86
87.89
Optdigit
90.01
88.24
95.61
55.51
92.71
95.24
94.64
95.66
97.82
MNIST1
60.89
68.36
72.92
43.97
77.09
90.78
76.91
75.09
91.27
SEMEION
73.91
74.51
80.29
45.89
78.38
83.61
79.24
80.56
89.46


### 3.5. Experiments on Large-Scale Multiclass Data

PRID (Person Re-ID) 2011 data set [52] includes images extracted from multiple tracks recorded by two different static surveillance cameras. These images contain viewpoint changes and significant differences in lighting, background, and camera characteristics. 80 images are selected for each of 200 people (cam b),16,000 images in total and each image is reshaped to a 1024-dimensional column vector. Some representative images are illustrated in Figure 5.Figure 5
Image samples from PRID data set.The clustering results of different algorithms are demonstrated in Table6 when the clustering number c=50,80,110,140,170,200 is selected, respectively. Apparently CDONMTF and SODNMF perform much better than the other algorithms. And compared to SODNMF, CDONMTF shows slightly better performance in the view of the average value of AC, NMI, and PU. In addition, CDONMTF has a higher efficiency than SODNMF on the data sets shown in Table 7.Table 6
Comparison of clustering performance on PRID data set.
Algorithms
GNMF
CNMF
DNMTF
PNMT
NMF
2
L
2,0
SODNMF
G-NBVD
GSR-NBVD
CDONMTF
c = 50
AC
60.79
55.68
70.13
14.76
69.31
77.51
70.49
70.52
78.71
NMI
78.31
69.35
86.02
31.25
84.94
90.84
86.48
86.5
91.46
PU
65.14
58.84
76.21
16.16
75.44
81.38
76.62
76.62
81.99
c = 80
AC
52.72
49.76
65.14
11.58
63.48
74.74
64.63
64.65
76.25
NMI
74.78
68.29
83.37
31.18
81.61
89.76
83.08
83.07
90.23
PU
57.38
52.71
71.51
12.81
69.98
78.63
71.05
71.05
79.58
c = 110
AC
54.71
45.18
61.97
9.8
60.58
75.47
62.16
62.15
75.54
NMI
77.04
66.25
82.09
31.59
80.39
90.71
82.37
82.41
90.33
PU
59.32
48.35
68.89
10.88
67.61
79.32
69.03
69.05
79.11
c = 140
AC
53.02
44.35
59.59
9.05
59.42
74.81
59.51
59.71
73.78
NMI
76.71
66.51
80.79
32.51
79.93
90.31
80.91
81.03
88.96
PU
58.21
47.71
67.42
9.95
66.87
78.89
67.26
67.39
77.57
c = 170
AC
50.28
44.39
58.52
8.36
59.49
72.69
58.46
58.44
73.08
NMI
75.35
67.92
80.55
32.83
80.42
88.71
80.36
80.39
89.66
PU
55.33
47.74
66.32
9.19
66.95
76.42
66.23
66.22
77.33
c = 200
AC
51.31
42.12
58.24
7.73
58.42
71.21
58.47
58.38
73.08
NMI
76.64
66.59
79.97
32.98
79.82
87.02
79.93
79.81
89.65
PU
56.24
45.43
66.11
8.64
66.06
75.32
66.21
66.11
77.16
Average
AC
53.81
46.91
62.26
10.21
61.79
74.41
62.29
62.29
75.07
NMI
76.47
67.48
82.13
32.06
81.19
89.55
82.19
82.21
90.05
PU
58.59
50.13
69.41
11.27
68.82
78.33
69.41
69.41
78.79Table 7
Running time of different algorithms on seven image data sets (seconds).
Data sets
GNMF
CNMF
DNMTF
PNMT
NMF
2
L
2,0
SODNMF
G-NBVD
GSR-NBVD
CDONMTF
UMIST
0.75
6.29
1.48
4.81
22.32
81.8
4.41
16.31
9.26
JAFFE
0.08
0.64
0.19
0.57
2.91
12.06
0.54
1.63
1.02
COIL20
0.18
2.48
0.34
1.07
9.2
28.47
1.03
4.87
3.78
Optdigit
0.08
1.75
0.14
0.2
19.55
198.67
0.28
1.08
6.14
MNIST1
1.45
57.16
2.38
6.57
241.99
2738.62
8.7
23.36
87.99
SEMEION
0.06
0.8
0.1
0.26
5.06
25.86
0.28
1.41
1.58
PRID
15.69
223.73
28.11
45.59
672.37
9848.17
105.39
419.94
457.96


### 3.6. Clustering Visualization

In order to visualize the effectiveness of our CDONMTF algorithm, we use IRIS data set (https://archive.ics.uci.edu/ml/datasets/iris) for data analysis. The data set contains a total of 150 samples in 3 classes, and each sample contains 4 features. As shown in Figure 6(a), we use PCA to project sample data into two-dimensional space, display the distribution of the 3 classes by different colors and tag shapes, and randomly select 20% samples from each class as prior labeled samples and marked with red dots.Figure 6
IRIS data set clustering illustration based on different algorithms. (a) IRIS. (b) CNMF. (c) DNMTF. (d) CDONMTF.
(a)
(b)
(c)
(d)We use CDONMTF, DNMTF, and CNMF algorithms to cluster the original IRIS data and obtain a clustering label for each sample. The clustering label is marked by its class with the same color and a smaller tag shape. If the clustering label is the same as the originally displayed, then we get a data point with a correct classification (Figures6(b)–6(d)). For the red marked samples, DNMTF has many clustering errors, but CDONMTF and CNMF cluster correctly using the label information. For unlabeled samples, CDONMTF has the highest accuracy owing to its effective use of priori labels, geometric structures of data and feature manifolds, and independent biorthogonal constraints.


### 3.7. Analysis of Convergence

The update rules used to minimize the CDONMTF objective function are essentially iterative processes. We have proved in Theorem1 that the iterative rules are convergent. In this section, the convergence rate of these update rules is studied experimentally.Figure7 shows the convergence curves of CDONMTF on 6 image data sets. All samples of each data set are decomposed in the experiments. In each subfigure, the y−axis is the value of the objective function while the x−axis denotes the number of iterations. Although the sizes of the 6 data sets are different, the proposed algorithm usually converges within 100 iterations, which indicates that CDONMTF has relatively stable and fast convergence rate.Figure 7
Convergence curves of CDONMTF on six image data sets. (a) UMIST. (b) JAFFE. (c) COIL20. (d) Optdigit. (e) MNIST1. (f) SEMEION.
(a)
(b)
(c)
(d)
(e)
(f)


## 3.1. Comparison of Algorithms and Evaluation Criteria

We compare our method to 8 existing representative methods that are similar, in some respects, to ours. These algorithms include GNMF [18], CNMF [27], DNMTF [35], PNMT [7], NMF2L2,0 [22], SODNMF [36], G-NBVD [19], and GSR-NBVD [19].For graph-based algorithms, GNMF, DNMF, G-NBVD,NMF2L2,0, GSR-NBVD, and our CDONMTF, are weighted with 0 and 1, and SODNMF is weighted with heat kernel. The nearest neighbor number is minm,n,5 as recommended in [18]. For CNMF, SODNMF, and CDONMTF, 20% of the samples in each data set are randomly selected as labeled samples. The number of iterations for all algorithms are set to 100.In all experiments, clustering accuracy (AC) [44], normalized mutual information (NMI) [44], and purity (PU) [45] are applied as criteria to evaluate the performance of clustering algorithms.


## 3.2. Parameter Analysis

The CDONMTF model has two types of parameters: decomposition dimension valuesp and q; regularization parameters α,β, and λ. Although many literature studies have studied how to choose the dimension [32, 46, 47] and regularization parameters [7, 48, 49], there is no unique way to select these parameters. For the choice of p and q, we use the following methods [19]:(38)p=θp×c,q=θq×c,where θp and θq are empirical, usually depending on the type of data sets, and c is the number of classes to cluster. In our experiments, θp=θq=3 are selected for UCI data sets and θp=3 and θq=2 for image data sets.For the regularization parametersα and β, their sensitivity is analyzed on 2 representative data sets, JAFFE and SEMEION. Values as α=0.01,0.1,1,10,100,1000 and β=0.01,0.05,0.1,0.5,1,5 are tested, while λ is fixed to 1. Figure 1 displays the clustering results of CDONMTF under different values of α and β. It can be seen that they are more robust for parameter α in a certain range, but more sensitive for parameter β. In the following experiments, we select α=100, β=0.5, and λ=1.Figure 1
Clustering performance of CDONMTF versus the value of parameters. (a) JAFFE. (b) SEMEION. (c) JAFFE. (d) SEMEION. (e) JAFFE. (f) SEMEION.
(a)
(b)
(c)
(d)
(e)
(f)


## 3.3. Experiments on UCI Data

In this section, we test clustering performance of different algorithms on 5 UCI nonnegative data sets (http://archive.ics.uci.edu/ml/index.phpl). Table 2 outlines the basic information of those data sets. The experimental results listed in Table 3 show that the proposed CDONMTF algorithm achieves the best clustering performance in the view of all 3 indexes: AC, NMI, and PU on three data sets (WDBC, lenses, and Chess), while the clustering results on Zoo and soybean data sets are not satisfying but acceptable.Table 2
Description of the five UCI data sets.
Data sets
No. of features (m)
No. of samples (n)
No. of classes (c)
WDBC
30
569
2
Lenses
4
24
3
Zoo
16
101
7
Chess
36
3196
2
Soybean
35
47
4Table 3
Comparison of clustering performance on UCI data sets.
Data sets
GNMF
CNMF
DNMTF
PNMT
NMF
2
L
2,0
SODNMF
G-NBVD
GSR-NBVD
CDONMTF
Accuracy (%)
WDBC
83.83
79.26
82.59
83.46
83.48
90.33
81.02
69.07
92.79
Lenses
41.67
56.17
50.13
57.33
41.67
37.51
41.67
45.83
57.83
Zoo
61.01
67.04
60.93
55.15
84.32
73.39
61.82
62.99
67.19
Chess
55.16
51.69
54.07
52.07
50.99
60.89
55.23
52.94
77.41
Soybean
74.47
72.34
78.72
55.28
89.36
81.02
76.61
78.09
87.23
Normalized mutual information (%)
WDBC
38.34
23.11
35.44
36.81
37.5
56.08
28.63
8.17
61.19
Lenses
2.01
18.52
3.71
13.89
2.02
4.99
2.58
2.52
21.73
Zoo
65.92
70.57
66.02
54.81
80.63
73.59
65.79
66.24
71.04
Chess
0.75
0.15
0.57
0.09
1.17
4.29
0.75
0.36
24.88
Soybean
71.38
71.56
72.44
36.74
81.49
79.47
72.12
72.31
77.85
Purity (%)
WDBC
83.83
79.26
82.61
83.46
83.48
90.33
81.02
69.07
92.79
Lenses
62.51
66.67
62.49
70.83
62.5
62.5
62.51
62.51
70.52
Zoo
82.18
85.39
82.28
74.12
87.05
83.29
82.04
82.36
85.91
Chess
55.16
52.22
54.07
52.22
52.22
60.89
55.23
52.94
77.41
Soybean
78.72
78.72
78.72
59.53
89.36
81.02
78.72
78.72
87.23


## 3.4. Experiments on Image Data

We continue to experiment on different types and sizes of some open image data sets: UMIST (http://images.ee.umist.ac.uk/danny/database.html), JAFFE [50], COIL20 [51], Optdigits (http://archive.ics.uci.edu/ml/machine-learning-datasets/optdigits/), MNIST1 (a subset of the first 10k samples) (http://www.cad.zju.edu.cn/home/dengcai/Data/MLData.html), and SEMEION (https://archive.ics.uci.edu/ml/machine-learning-databases/semeion/) in this part. Both types and sizes of the data sets are listed in Table 4.Table 4
Description of the six image data sets.
Data sets
No. of features (m)
No. of samples (n)
No. of classes (c)
Data types
UMIST
10304
575
20
Face image
JAFFE
4096
213
10
Facial expression
COIL20
1024
1440
20
Object image
Optdigit
64
3823
10
Handwritten digit
MNIST1
784
10000
10
Handwritten digit
SEMEION
256
1593
10
Handwritten digitWe use a subset of the firstc categories as the collection X for clustering. Different cluster numbers c are used to assess the efficiency of the algorithms under various sample scales of X in the experiments.Experiment results (Figures2–4) display the performance of different algorithms using different clustering numbers c. Table 5 shows the average performance of different algorithms under the same data set with different clustering numbers c. The optimal results are marked in bold in the table. Several interesting findings are noted as follows:(i)
SODNMF and the proposed CDONMTF achieve better clustering performance than PNMT, GNMF, DNMTF,NMF2L2,0, G-NBVD, and GSR-NBVD in most cases due to adding the priori label information of the data as an additional constraint(ii)
Compared with CNMF only considering priori label information, the proposed CDONMTF considers the geometric structure information contained in data and feature points and combines with additional biorthogonal constraints, which results in more efficient clustering performance(iii)
From the view of evaluation indexes (AC, NMI, and PU), the average clustering effect of CDONMTF is almost the best on all data sets except JAFFE, but it is also good enoughFigure 2
Comparisons of clustering AC (%) with different number of classes on six data sets. (a). UMIST. (b) JAFFE. (c) COIL20. (d) Optdigit. (e) MNIST1. (f) SEMEION.
(a)
(b)
(c)
(d)
(e)
(f)Figure 3
Comparisons of clustering NMI (%) with different number of classes on six data sets. (a). UMIST. (b) JAFFE. (c) COIL20. (d) Optdigit. (e) MNIST1. (f) SEMEION.
(a)
(b)
(c)
(d)
(e)
(f)Figure 4
Comparisons of clustering PU (%) with different number of classes on six data sets.(a). UMIST. (b) JAFFE. (c) COIL20. (d) Optdigit. (e) MNIST1. (f) SEMEION.
(a)
(b)
(c)
(d)
(e)
(f)Table 5
Comparison of average clustering performance on six image data sets.
Data sets
GNMF
CNMF
DNMTF
PNMT
NMF
2
L
2,0
SODNMF
G-NBVD
GSR-NBVD
CDONMTF
Accuracy (%)
UMIST
43.47
47.17
66.56
37.57
63.65
79.06
64.91
64.86
81.46
JAFFE
87.89
91.33
89.85
41.97
97.66
96.57
88.86
88.85
96.48
COIL20
59.71
62.36
74.71
43.91
72.54
83.29
77.18
77.21
86.11
Optdigit
89.53
88.23
95.46
54.04
92.64
95.06
93.91
95.13
97.77
MNIST1
56.52
64.14
65.41
40.71
72.64
90.21
72.76
67.83
90.88
SEMEION
70.56
73.63
79.16
42.74
76.73
82.69
77.45
79.45
89.27
Normalized mutual information (%)
UMIST
58.33
58.51
79.04
43.17
76.82
86.01
78.73
78.74
87.34
JAFFE
86.41
89.58
90.59
36.93
97.78
93.72
90.04
90.05
95.72
COIL20
64.98
65.42
83.12
49.44
80.22
88.95
85.61
85.56
90.54
Optdigit
86.16
77.18
92.72
41.09
89.88
92.01
93.18
93.55
96.61
MNIST1
45.73
54.31
66.61
29.41
71.34
84.34
71.04
69.27
84.62
SEMEION
68.94
60.12
73.96
33.94
71.83
78.76
73.71
74.15
81.89
Purity (%)
UMIST
50.04
52.25
74.33
42.49
71.05
82.91
73.21
73.19
84.47
JAFFE
88.66
92.15
91.31
44.98
98.06
94.53
90.58
90.57
96.76
COIL20
62.15
64.16
79.69
48.32
77.27
84.73
81.85
81.86
87.89
Optdigit
90.01
88.24
95.61
55.51
92.71
95.24
94.64
95.66
97.82
MNIST1
60.89
68.36
72.92
43.97
77.09
90.78
76.91
75.09
91.27
SEMEION
73.91
74.51
80.29
45.89
78.38
83.61
79.24
80.56
89.46


## 3.5. Experiments on Large-Scale Multiclass Data

PRID (Person Re-ID) 2011 data set [52] includes images extracted from multiple tracks recorded by two different static surveillance cameras. These images contain viewpoint changes and significant differences in lighting, background, and camera characteristics. 80 images are selected for each of 200 people (cam b),16,000 images in total and each image is reshaped to a 1024-dimensional column vector. Some representative images are illustrated in Figure 5.Figure 5
Image samples from PRID data set.The clustering results of different algorithms are demonstrated in Table6 when the clustering number c=50,80,110,140,170,200 is selected, respectively. Apparently CDONMTF and SODNMF perform much better than the other algorithms. And compared to SODNMF, CDONMTF shows slightly better performance in the view of the average value of AC, NMI, and PU. In addition, CDONMTF has a higher efficiency than SODNMF on the data sets shown in Table 7.Table 6
Comparison of clustering performance on PRID data set.
Algorithms
GNMF
CNMF
DNMTF
PNMT
NMF
2
L
2,0
SODNMF
G-NBVD
GSR-NBVD
CDONMTF
c = 50
AC
60.79
55.68
70.13
14.76
69.31
77.51
70.49
70.52
78.71
NMI
78.31
69.35
86.02
31.25
84.94
90.84
86.48
86.5
91.46
PU
65.14
58.84
76.21
16.16
75.44
81.38
76.62
76.62
81.99
c = 80
AC
52.72
49.76
65.14
11.58
63.48
74.74
64.63
64.65
76.25
NMI
74.78
68.29
83.37
31.18
81.61
89.76
83.08
83.07
90.23
PU
57.38
52.71
71.51
12.81
69.98
78.63
71.05
71.05
79.58
c = 110
AC
54.71
45.18
61.97
9.8
60.58
75.47
62.16
62.15
75.54
NMI
77.04
66.25
82.09
31.59
80.39
90.71
82.37
82.41
90.33
PU
59.32
48.35
68.89
10.88
67.61
79.32
69.03
69.05
79.11
c = 140
AC
53.02
44.35
59.59
9.05
59.42
74.81
59.51
59.71
73.78
NMI
76.71
66.51
80.79
32.51
79.93
90.31
80.91
81.03
88.96
PU
58.21
47.71
67.42
9.95
66.87
78.89
67.26
67.39
77.57
c = 170
AC
50.28
44.39
58.52
8.36
59.49
72.69
58.46
58.44
73.08
NMI
75.35
67.92
80.55
32.83
80.42
88.71
80.36
80.39
89.66
PU
55.33
47.74
66.32
9.19
66.95
76.42
66.23
66.22
77.33
c = 200
AC
51.31
42.12
58.24
7.73
58.42
71.21
58.47
58.38
73.08
NMI
76.64
66.59
79.97
32.98
79.82
87.02
79.93
79.81
89.65
PU
56.24
45.43
66.11
8.64
66.06
75.32
66.21
66.11
77.16
Average
AC
53.81
46.91
62.26
10.21
61.79
74.41
62.29
62.29
75.07
NMI
76.47
67.48
82.13
32.06
81.19
89.55
82.19
82.21
90.05
PU
58.59
50.13
69.41
11.27
68.82
78.33
69.41
69.41
78.79Table 7
Running time of different algorithms on seven image data sets (seconds).
Data sets
GNMF
CNMF
DNMTF
PNMT
NMF
2
L
2,0
SODNMF
G-NBVD
GSR-NBVD
CDONMTF
UMIST
0.75
6.29
1.48
4.81
22.32
81.8
4.41
16.31
9.26
JAFFE
0.08
0.64
0.19
0.57
2.91
12.06
0.54
1.63
1.02
COIL20
0.18
2.48
0.34
1.07
9.2
28.47
1.03
4.87
3.78
Optdigit
0.08
1.75
0.14
0.2
19.55
198.67
0.28
1.08
6.14
MNIST1
1.45
57.16
2.38
6.57
241.99
2738.62
8.7
23.36
87.99
SEMEION
0.06
0.8
0.1
0.26
5.06
25.86
0.28
1.41
1.58
PRID
15.69
223.73
28.11
45.59
672.37
9848.17
105.39
419.94
457.96


## 3.6. Clustering Visualization

In order to visualize the effectiveness of our CDONMTF algorithm, we use IRIS data set (https://archive.ics.uci.edu/ml/datasets/iris) for data analysis. The data set contains a total of 150 samples in 3 classes, and each sample contains 4 features. As shown in Figure 6(a), we use PCA to project sample data into two-dimensional space, display the distribution of the 3 classes by different colors and tag shapes, and randomly select 20% samples from each class as prior labeled samples and marked with red dots.Figure 6
IRIS data set clustering illustration based on different algorithms. (a) IRIS. (b) CNMF. (c) DNMTF. (d) CDONMTF.
(a)
(b)
(c)
(d)We use CDONMTF, DNMTF, and CNMF algorithms to cluster the original IRIS data and obtain a clustering label for each sample. The clustering label is marked by its class with the same color and a smaller tag shape. If the clustering label is the same as the originally displayed, then we get a data point with a correct classification (Figures6(b)–6(d)). For the red marked samples, DNMTF has many clustering errors, but CDONMTF and CNMF cluster correctly using the label information. For unlabeled samples, CDONMTF has the highest accuracy owing to its effective use of priori labels, geometric structures of data and feature manifolds, and independent biorthogonal constraints.


## 3.7. Analysis of Convergence

The update rules used to minimize the CDONMTF objective function are essentially iterative processes. We have proved in Theorem1 that the iterative rules are convergent. In this section, the convergence rate of these update rules is studied experimentally.Figure7 shows the convergence curves of CDONMTF on 6 image data sets. All samples of each data set are decomposed in the experiments. In each subfigure, the y−axis is the value of the objective function while the x−axis denotes the number of iterations. Although the sizes of the 6 data sets are different, the proposed algorithm usually converges within 100 iterations, which indicates that CDONMTF has relatively stable and fast convergence rate.Figure 7
Convergence curves of CDONMTF on six image data sets. (a) UMIST. (b) JAFFE. (c) COIL20. (d) Optdigit. (e) MNIST1. (f) SEMEION.
(a)
(b)
(c)
(d)
(e)
(f)


## 4. Conclusions

A novel coclustering algorithm, namely, CDONMTF is proposed in this paper. The algorithm effectively utilizes the hidden geometric structure information in the data and feature spaces by constructing the associated data and feature graphs and incorporates the priori label information of the data and independent biorthogonal constraints to enhance the ability of data discrimination in the new representation space. The clustering experiments on multiple UCI machine learning data sets and image data sets show that the proposed algorithm has better performance over several existing clustering methods.



---

*Source: 7565640-2019-12-26.xml*