# Facial Expression Recognition Based on Discriminant Neighborhood Preserving Nonnegative Tensor Factorization and ELM

Authors: Gaoyun An; Shuai Liu; Yi Jin; Qiuqi Ruan; Shan Lu
Journal: Mathematical Problems in Engineering
 (2014)

---

## Abstract

A novel facial expression recognition algorithm based on discriminant neighborhood preserving nonnegative tensor factorization (DNPNTF) and extreme learning machine (ELM) is proposed. A discriminant constraint is adopted according to the manifold learning and graph embedding theory. The constraint is useful to exploit the spatial neighborhood structure and the prior defined discriminant properties. The obtained parts-based representations by our algorithm vary smoothly along the geodesics of the data manifold and have good discriminant property. To guarantee the convergence, the project gradient method is used for optimization. Then features extracted by DNPNTF are fed into ELM which is a training method for the single hidden layer feed-forward networks (SLFNs). Experimental results on JAFFE database and Cohn-Kanade database demonstrate that our proposed algorithm could extract effective features and have good performance in facial expression recognition.


---

## Body

## 1. Introduction

Facial expression recognition plays an important role in human-computer interaction, and 55% information is transferred by facial expression in face-to-face human communication [1]. Although many methods were proposed, recognizing facial expression is still challenging due to the complex, variable, and subtle facial expressions.One of the effective methods for facial expression recognition is the subspace-based algorithm [2–5]. It aims to project the samples into a lower dimensional space which preserves the needed information and discards the redundant information. There are many widely used subspace-based algorithms, for example, principle component analysis (PCA) [2], linear discriminant analysis (LDA) [3], neighborhood preserving embedding (NPE) [4], locality-preserving projection (LPP) [6], and singular value decomposition (SVD) [5], and so forth.Recently, the nonnegative matrix factorization (NMF) was introduced into facial expression recognition [7]. NMF decomposes the face samples into two nonnegative parts: the basis images and the corresponding weights. As many data in bases and weights were degenerated too close to zero, NMF derived the parts-based sparse representations. For facial expression recognition, the localized subtle features, such as the corners of mouth, upward or downward eyebrows, and change of eyes, are critical for the recognition performance. Since NMF yields the parts-based representations, it outperforms the subspace-based models. To further improve NMF, several variants have been presented by introducing different constraints to the objective function. Li et al. put forward a local NMF (LNMF) by adding a local constraint to the basis images [8], to learn a localized, parts-based representation. Hoyer gave a sparse constraint NMF (SNMF) by incorporating sparseness constraint into both the bases and the weights [9]. Cai et al. developed a graph constraint NMF (GNMF) by adding a graph preserving constraint to the weights [10]. Zafeiriou et al. used the discriminant NMF (DNMF) for frontal face verification [11].  Wang et al. extended NMF to PNMF with a PCA constraint and FNMF with a Fisher constraint [12].For facial expression recognition, NMF and its variants vectorize the samples before factorization, which may lose the local geometric structures. However, the spatial neighborhood relationships within pixels are critical for image representation, understanding, and recognition [13]. Another drawback of NMF is that it could not generate a unique decomposition result. Welling and Weber developed a positive tensor factorization (PTF) algorithm, which handled images as 2D matrices directly [14]. Shashua and Hazan proposed the nonnegative tensor factorization (NTF) which implemented the factorization in the rank-one tensor space [15]. The factorization in the tensor space could preserve the local structures and guarantee the uniqueness of the decomposition.On the other hand, the choice of classifier plays an important role for recognition. For facial expression recognition, nearest neighbor (NN) and support vector machine (SVM) are the commonly used methods [16]. The sparse representation classifier (SRC) was adopted in [17]. Recently, the extreme learning machine (ELM) was proposed for classification which is a training method for the single hidden layer feed-forward networks (SLFNs) [18]. The conventional methods need long time to converge or may lose the generalization property due to overfitting. However, ELM converges fast and provides good generalization performance. For ELM, the input weights and biases are randomly assigned, and the output weights can be simply calculated by the generalized inverse of the hidden layer output matrix. Therefore it converges extremely fast and obtains an excellent generalization capability. Many variants of ELM were proposed for different applications [19–26], including the Kernel-based ELM [21] and the incremental ELM (I-ELM) [23], which lead to the state-of-the-art results in different applications.In this paper, we propose a novel facial expression recognition algorithm based on discriminant neighborhood preserving nonnegative tensor factorization (DNPNTF) and ELM. It works well in the rank-one tensor space. The simple ELM is adopted to testify its effectiveness for facial expression recognition [18]. Our algorithm is composed of two stages: feature extraction and classification. Firstly, to extract the discriminant features, a neighborhood preserving constraint form of NTF is used. The constraint is derived according to the manifold learning and graph embedding theory [27–29]. Since the columns of the weighting matrix have a one-to-one correspondence with the columns of the original sample, the discriminant constraint is added to the weighting matrix. With the neighborhood preserving constraint, the obtained parts-based representations vary smoothly along the geodesics of the data manifold and are more discriminant. Secondly, the discriminant features extracted by DNPNTF are fed into ELM classifier to conduct the recognition task.The rest of this paper is organized as follows. The mathematical notations are given in Section2. In Section 3, we give the detailed analysis about DNPNTF and its optimization procedure. ELM is introduced in Section 4, and the experiments are given in Section 5. Finally, the conclusions are drawn in Section 6.


## 2. Basic Algebra and Notations

In this paper, a tensor is represented asA
∈
R
m
1
×
m
2
×
⋯
×
m
n, whose order is n. The element of A is denoted by A
i
1
i
2
⋯
i
n
>
0, where 1
≤
i
k
≤
m
k, 1
≤
k
≤
n.Definition 1 (inner product and tensor product [30]).
The inner product of two tensorsA
,
B
∈
R
m
1
×
m
2
×
⋯
×
m
n is defined as
(1)
〈
A
,
B
〉
=
∑
i
1
=
1
,
…
,
i
n
=
1
m
1
,
…
,
m
n
A
i
1
⋯
i
n
B
i
1
⋯
i
n
.
The tensor product of two tensorsA
∈
R
m
1
×
m
2
×
⋯
×
m
n and B
∈
R
p
1
×
p
2
×
⋯
×
p
q is
(2)
C
=
A
⊗
B
∈
R
m
1
×
m
2
×
⋯
×
m
n
×
p
1
×
⋯
×
p
q
.Definition 2 (rank-one tensor [30]).
Anth-order tensor A
∈
R
m
1
×
m
2
×
⋯
×
m
n could be represented as a tensor product of n tensors u
1
=
[
u
1
1
,
u
2
1
,
…
,
u
m
1
1
]
∈
R
m
1
,
…
,
u
n
=
[
u
1
n
,
u
2
n
,
…
,
u
m
n
n
]
∈
R
m
n as
(3)
A
=
u
1
⊗
u
2
⊗
⋯
⊗
u
n
.
HereA is called rank-one tensor, and A
i
1
i
2
⋯
i
n
=
u
i
1
1
u
i
2
2
⋯
u
i
n
n  (
1
≤
i
1
≤
m
1
,
…
,
1
≤
i
n
≤
m
n
).Definition 3 (mode product [30]).
The modek product of A
∈
R
m
1
×
m
2
×
⋯
×
m
n and U
∈
R
m
k
×
m
k
′ is
(4)
C
=
A
×
k
U
∈
R
m
1
×
m
k
-
1
×
m
k
′
×
m
k
+
1
×
⋯
×
m
n
,
where C
i
1
⋯
i
k
-
1
q
i
k
+
1
⋯
i
n
=
∑
p
=
1
m
k
A
i
1
⋯
i
k
-
1
p
i
k
+
1
⋯
i
n
·
U
p
,
q, 1
≤
i
j
≤
m
j, 1
≤
j
≤
n.


## 3. The DNPNTF Algorithm

In this section, we give a detailed description about the proposed DNPNTF algorithm. Instead of converting into vectors, it processes the samples in rank-one tensor space. The objective function of NTF is adopted, which could learn the parts-based representation and have the sparse property. To discover the spatial local geometric structure and the discriminant class-based information, a constraint is added in the objective function according to the manifold learning and graph embedding analysis. To guarantee the convergence, the project gradient method is used.


### 3.1. The Analysis of DNPNTF

Given am
×
N image database X, it contains N sample images X
=
[
x
1
,
…
,
x
n
]
∈
R
m
×
N. The dimension of each sample is m
(
m
=
m
1
×
n
2
). In NTF, the database is organized as a 3rd-order tensor A
∈
R
m
1
×
m
2
×
N, which is overlapped by {
x
i
∈
R
m
1
×
m
2
|
i
=
1
N
} sequentially. The objective function of NTF is
(5)
min
⁡
u
r
,
v
r
,
w
r
∥
A
-
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
∥
2
s
.
t
.
u
r
,
v
r
,
z
r
≥
0
,
where u
r
∈
R
m
1, v
r
∈
R
m
2, and z
r
∈
R
N
(
1
≤
r
≤
R
) describe the first, second, and third modules of A, respectively. Each sample x
i is approximated by
(6)
x
i
≈
∑
r
=
1
R
z
r
i
u
r
v
r
T
,
(
i
=
1
,
…
,
N
)
.By minimizing (5), the bases {
u
r
v
r
T
|
r
=
1
R
} and the corresponding weights {
z
r
|
r
=
1
R
} are conducted. The inner product of {
u
r
v
r
T
|
r
=
1
R
} and the sample image is calculated to derive the low-dimensional parts-based representation.To incorporate more properties into NTF, different constraints could be added into the objective function. The constraint form of objective function is(7)
min
⁡
u
r
,
v
r
,
z
r
∥
A
-
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
∥
2
+
α
F
(
U
,
V
)
+
β
J
(
Z
)
s
.
t
.
u
r
,
v
r
,
z
r
≥
0
,
∀
r
,
where U
=
[
u
1
,
u
2
,
…
,
u
R
]
∈
R
m
1
×
R, V
=
[
v
1
,
v
2
,
…
,
v
R
]
∈
R
m
2
×
R, Z
=
[
z
1
,
z
2
,
…
,
z
R
]
∈
R
N
×
R. F is the constraint function about {
U
,
V
} and J is the constraint about Z. α and β are the corresponding positive coefficients. To encode the spatial structure and discriminant class-based information into sparse representations, we propose a constraint function J according to the manifold learning and graph embedding analysis. In NTF, the columns of the weighting matrix have a one-to-one correspondence with the columns of the original image matrix. Therefore, we add the discriminant constraint to {
z
r
|
r
=
1
R
}, and J for z
r is defined as
(8)
J
(
Z
)
=
min
⁡
⁡
∑
∑
i
≠
j
∑
l
η
l
∥
z
i
-
z
j
∥
2
S
i
j
l
,
where S
i
j
l denote the graphs with different properties and η
l are the corresponding coefficients. By deriving different S
i
j, the graph embedding model could have different properties, such as the neighborhood preserving property and the discriminant property.Now, we discuss the selection ofS
i
j
l. The most commonly used graph is the Laplacian graph, and S
i
j
l is calculated in form of the heat kernel function as
(9)
S
i
j
=
e
-
∥
z
i
-
z
j
∥
2
/
t
,
t
∈
R
.Here,S
i
j measures the similarity between a pair of vertices and has neighborhood preserving property. To further incorporate the class-based discriminant information, we derive a universal penalty graph [27], where the similarity matrix S
i
j
p is defined as
(10)
S
i
j
p
=
{
1
or
1
k
,
if
(
i
,
j
)
∈
N
k
(
l
i
)
or
N
k
(
l
j
)
0
,
otherwise
,
where N
k
(
l
i
) represents the nearest k pairs of samples between class l
i and the other classes. The purpose of the penalty graph was to separate marginal samples between different classes.Now, the objective function of constrained NTF becomes(11)
min
⁡
u
r
,
v
r
,
z
r
∥
A
-
∑
r
=
1
R
u
r

⊗

v
r

⊗

z
r
∥
2

+

λ
∥
z
i

-

z
j
∥
2
S
i
j

-

σ
∥
z
i

-

z
j
∥
2
S
i
j
p
s
.
t
.
u
r
,
v
r
,
z
r
≥
0
,
∀
r
.By solving the generalized eigenvalue decomposition problem, the graph embedding criterion in (11) can be calculated as
(12)
1
2
λ
∥
z
i
-
z
j
∥
2
S
i
j
-
1
2
σ
∥
z
i
-
z
j
∥
2
S
i
j
p
=
1
2
λ
∑
i
,
j
t
r
[
(
z
i
-
z
j
)
(
z
i
-
z
j
)
T
]
S
i
j
-
1
2
σ
∑
i
,
j
t
r
[
(
z
i
-
z
j
)
(
z
i
-
z
j
)
T
]
S
i
j
p
=
1
2
t
r
[
λ
∑
i
,
j
[
(
z
i
-
z
j
)
(
z
i
-
z
j
)
T
]
S
i
j
-
σ
∑
i
,
j
[
(
z
i
-
z
j
)
(
z
i
-
z
j
)
T
]
S
i
j
p
]
=
t
r
[
λ
∑
i
,
j
(
z
i
S
i
j
(
z
i
)
T
-
z
i
S
i
j
(
z
j
)
T
)
-
σ
∑
i
,
j
(
z
i
S
i
j
p
(
z
i
)
T
-
z
i
S
i
j
p
(
z
j
)
T
)
]
=
t
r
[
λ
∑
i
,
j
(
z
i
E
i
j
(
z
i
)
T
-
z
i
E
i
j
(
z
j
)
T
)
-
σ
∑
i
,
j
(
z
i
E
i
j
p
(
z
i
)
T
-
z
i
E
i
j
p
(
z
j
)
T
)
]
=
t
r
[
λ
∑
i
,
j
z
i
L
i
j
(
z
j
)
T
-
σ
∑
i
,
j
z
i
L
i
j
p
(
z
j
)
T
]
=
λ
t
r
(
Z
T
L
Z
)
-
σ
t
r
(
Z
T
L
p
Z
)
.
And the final objective function of DNPNTF is
(13)
min
⁡
u
r
,
v
r
,
z
r
∥
A
-
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
∥
2

+

λ
t
r
(
Z
T
L
Z
)

-

σ
t
r
(
Z
T
L
p
Z
)
s
.
t
.
u
r
,
v
r
,
z
r
≥
0
,
∀
r
,
where λ
>
0 and σ
>
0 to make sure (13) should be nonnegative.


### 3.2. Projected Gradient Method of DNPNTF

The most popular approach to minimize NMF or NTF is the multiplicative update method. However, it cannot ensure the convergence of the constraint forms of NMF or NTF. In this paper, the projected gradient method is used to solve DNPNTF.The objective function of DNPNTF can be stated as(14)
D
obj
(
A
∥
U
V
Z
)
=
D
(
A
∥
U
V
Z
)
+
ς
D
(
Z
)
,
where ς is a positive constant. The goal of (14) is to find {
U
,
V
} and Z by solving the following problem:
(15)
min
⁡
U
,
V
,
Z
D
obj
(
A
∥
U
V
Z
)
s
.
t
.
u
i
,
k
≥
0
,
v
k
,
j
≥
0
,
z
i
,
j
≥
0
∀
k
,
i
,
j
.To find the optimal solution, (15) is divided into three subproblems: first, we fix V and Z and update U to arrive at the conditional optimal value of the subminimization problem; second, we fix U and Z and update V; last, we fix U and V and update Z. Three functions are defined as f
V
,
Z
(
U
)
=
D
obj
(
A
∥
U
V
Z
), f
U
,
Z
(
V
)
=
D
obj
(
A
∥
U
V
Z
), and f
U
,
V
(
Z
)
=
D
obj
(
A
∥
U
V
Z
). The update rules are defined as
(16)
U
(
t
+
1
)
=
[
U
(
t
)
-
α
t
∇
f
V
,
Z
(
U
(
t
)
)
]
,
V
(
t
+
1
)
=
[
V
(
t
)
-
α
t
∇
f
U
,
Z
(
V
(
t
)
)
]
,
Z
(
t
+
1
)
=
[
Z
(
t
)
-
α
t
∇
f
U
,
V
(
Z
(
t
)
)
]
.
Now the task is calculating ∇
f
V
,
Z
(
U
(
t
)
), ∇
f
U
,
Z
(
V
(
t
)
), and ∇
f
U
,
V
(
Z
(
t
)
).


#### 3.2.1. The Calculation ofU and V

Firstly, we discuss the calculation of∇
f
V
,
Z
(
U
(
t
)
) and ∇
f
U
,
Z
(
V
(
t
)
). The objective function could be written as
(17)
f
=
1
2
∥
A
-
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
∥
2
+
1
2
λ
t
r
(
Z
T
L
Z
)
-
1
2
σ
t
r
(
Z
T
L
p
Z
)
=
1
2
〈
A
-
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
,
A
-
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
〉
+
1
2
λ
∑
r
=
1
R
z
r
T
L
z
r
-
1
2
σ
∑
r
=
1
R
z
r
T
L
p
z
r
.The differential off is
(18)
d
(
f
u
s
)
=
〈
A
-
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
,
-
d
(
u
s
)
⊗
v
s
⊗
z
s
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
〉
+
0
-
0
=
〈
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
,
d
(
u
s
)
⊗
v
s
⊗
z
s
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
〉
-
〈
A
,
d
(
u
s
)
⊗
v
s
⊗
z
s
〉
.And the partial differential foru
s
p
(
1
≤
p
≤
m
1
) is
(19)
∂
f
∂
u
s
p
=
〈
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
,
e
p
⊗
v
s
⊗
z
s
〉
-
〈
A
,
e
p
⊗
v
s
⊗
z
s
〉
,
where the pth element in e
p
∈
R
m
1 is 1 and others are 0. That is, (
e
p
)
p
=
1 and (
e
p
)
k
≠
p
=
0. According to Definition 1, for any order tensors A
1
,
A
2
∈
R
a
1
×
a
2
×
⋯
×
a
n, B
1
,
B
2
∈
R
b
1
×
b
2
×
⋯
×
b
n, there is 〈
A
1
⊗
B
1
,
A
2
⊗
B
2
〉
=
〈
A
1
,
A
2
〉
〈
B
1
,
B
2
〉. Then (19) could be written as
(20)
∂
f
∂
u
s
p
=
∑
r
=
1
R
(
〈
u
r
,
e
p
〉
〈
v
r
⊗
z
r
,
v
s
⊗
z
s
〉
)
-
(
∑
k
=
1
,
q
=
1
,
i
=
1
m
1
,
m
2
,
N
A
k
q
i
(
e
p
)
k
v
s
q
z
s
i
)
=
∑
r
=
1
R
(
〈
u
r
,
e
p
〉
〈
v
r
⊗
z
r
〉
〈
v
s
⊗
z
s
〉
)
-
(
∑
q
=
1
,
i
=
1
m
2
,
N
A
p
q
i
(
e
p
)
p
v
s
q
z
s
i
+
∑
k
≠
p
A
k
q
i
(
e
p
)
k
v
s
q
z
s
i
)
=
∑
r
=
1
R
(
〈
u
r
,
e
p
〉
〈
v
r
,
v
s
〉
〈
z
r
,
z
s
〉
)
-
(
∑
q
=
1
,
i
=
1
m
2
,
N
A
p
q
i
v
s
q
z
s
i
+
0
)
=
∑
r
=
1
R
(
(
u
r
p
(
e
p
)
p
+
∑
k
≠
p
u
r
k
(
e
p
)
k
)
〈
v
r
,
v
s
〉
〈
z
r
,
z
s
〉
)
-
∑
q
=
1
,
i
=
1
m
2
,
N
A
p
q
i
v
s
q
z
s
i
=
∑
r
=
1
R
(
u
r
p
〈
v
r
,
v
s
〉
〈
z
r
,
z
s
〉
)
-
∑
q
=
1
,
i
=
1
m
2
,
N
A
p
q
i
v
s
q
z
s
i
.According to (16), the update rule for u
s
p is
(21)
u
s
p
=
u
s
p
-
μ
(
u
s
p
)
∂
f
∂
u
s
p
=
u
s
p
-
μ
(
u
s
p
)
(
∑
r
=
1
R
(
u
r
p
〈
v
r
,
v
s
〉
〈
z
r
,
z
s
〉
)
-
∑
q
=
1
,
i
=
1
m
2
,
N
A
p
q
i
v
s
q
z
s
i
)
.To confirm the nonnegative ofu
s
p, μ
(
u
s
p
) is set to be
(22)
μ
(
u
s
p
)
=
u
s
p
∑
r
=
1
R
(
u
r
p
〈
v
r
,
v
s
〉
〈
z
r
,
z
s
〉
)
.
Then the update function is
(23)
u
s
p
=
u
s
p
(
∑
q
=
1
,
i
=
1
m
2
,
N
A
p
q
i
v
s
q
z
s
i
)
∑
r
=
1
R
(
u
r
p
〈
v
r
,
v
s
〉
〈
z
r
,
z
s
〉
)
=
u
s
p
(
v
s
T
A
p
;
;
z
s
)
∑
r
=
1
R
(
u
r
p
(
v
r
T
v
s
)
(
z
r
T
z
s
)
)
=
u
s
p
(
v
s
T
A
p
;
;
z
s
)
(
U
p
;
(
(
V
T
v
s
)
⊙
(
Z
T
z
s
)
)
)
,
where U
p
;
∈
R
1
×
R represents the pth row of U
=
[
u
1
,
u
2
,
…
,
u
R
], ⊙ is the matrix Hadamard product, A
p
;
;
∈
R
m
2
×
N represents the matrix which fixes the first module of A, and traversals are the other two modules. It is defined as
(24)
(
A
p
;
;
)
q
i
=
A
p
q
i
,
1
≤
∀
q
≤
m
2
,
1
≤
∀
i
≤
N
.Similarly, the update rule of theqth element of v
s
(
v
s
q
,
1
≤
s
≤
R
,
1
≤
q
≤
m
2
) is
(25)
v
s
q
=
v
s
q
∑
p
=
1
,
i
=
1
m
1
,
N
(
A
p
q
i
u
s
p
z
s
i
)
∑
r
=
1
R
(
v
r
q
〈
u
r
,
u
s
〉
〈
z
r
,
z
s
〉
)
=
v
s
q
(
u
s
T
A
;
q
;
z
s
)
(
V
q
;
(
(
U
T
u
s
)
⊙
(
Z
T
z
s
)
)
)
,
where V
q
;
∈
R
1
×
R represents the pth row of V
=
[
v
1
,
v
2
,
…
,
v
R
]
∈
R
m
2
×
R, ⊙ is the matrix Hadamard product, A
;
q
;
∈
R
m
1
×
N represents the matrix which fixes the second module of A, and traversals are the other two modules. It is defined as
(26)
(
A
;
q
;
)
p
i
=
A
p
q
i
,
1
≤
∀
p
≤
m
1
,
1
≤
∀
i
≤
N
.Now,v
r and u
r are calculated.


#### 3.2.2. The Calculation ofZ

Then we discuss the calculation of∇
f
U
,
V
(
Z
(
t
)
). The differential of f along z
s
(
∀
s
,
1
≤
s
≤
R
) is
(27)
d
(
f
z
s
)
=
〈
A
-
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
,
-
d
(
z
s
)
⊗
u
s
⊗
v
s
〉
+
1
2
d
(
λ
z
s
T
L
z
s
)
-
1
2
d
(
σ
z
s
T
L
p
z
s
)
=
〈
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
,
d
(
z
s
)
⊗
u
s
⊗
v
s
〉
-
〈
A
,
d
(
z
s
)
⊗
u
s
⊗
v
s
〉
+
d
(
1
2
λ
z
s
T
(
D
-
S
)
z
s
)
-
d
(
1
2
σ
z
s
T
(
D
p
-
S
p
)
z
s
)
.For(
λ
z
s
T
(
D
-
S
)
z
s
)
/
2, the partial differential for z
s
i is
(28)
∂
(
(
1
/
2
)
λ
z
s
T
(
D
-
S
)
z
s
)
∂
z
s
i
=
∂
(
z
s
)
∂
z
s
i
∂
(
(
1
/
2
)
λ
z
s
T
(
D
-
S
)
z
s
)
∂
z
s
=
(
e
i
)
T
λ
(
D
-
S
)
z
s
,
where the ith element in e
i
∈
R
N is 1 and others are 0. That is, (
e
i
)
i
=
1 and (
e
i
)
k
≠
i
=
0. Then the partial differential for z
s
i is
(29)
∂
f
∂
z
s
i
=
〈
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
,
u
s
⊗
v
s
⊗
e
i
〉
-
〈
A
,
u
s
⊗
v
s
⊗
e
i
〉
+
λ
(
e
i
)
T
(
D
-
S
)
z
s
-
σ
(
e
i
)
T
(
D
p
-
S
p
)
z
s
=
∑
r
=
1
R
〈
u
r
,
u
s
〉
〈
v
r
,
v
s
〉
z
r
i
-
∑
p
=
1
,
q
=
1
m
1
,
m
2
A
p
q
i
u
s
p
v
s
q
+
(
λ
D
i
i
z
s
i
-
λ
∑
k
=
1
N
S
i
k
z
s
k
)
-
(
σ
D
i
i
p
z
s
i
-
σ
∑
k
=
1
N
S
i
k
p
z
s
k
)
.According to (16), the update rule for z
s
i is
(30)
z
s
i
=
z
s
i
-
μ
(
z
s
i
)
∂
f
∂
z
s
i
=
z
s
i
-
μ
(
z
s
i
)
(
∑
r
=
1
R
〈
u
r
,
u
s
〉
〈
v
r
,
v
s
〉
z
r
i
+
λ
D
i
i
z
s
i
+
σ
∑
k
=
1
N
S
i
k
p
z
s
k
-
∑
p
=
1
,
q
=
1
m
1
,
m
2
A
p
q
i
u
s
p
v
s
q
-
λ
∑
k
=
1
N
S
i
k
z
s
k
-
σ
D
i
i
p
z
s
i
)
.The update stepμ
(
z
s
i
) is set as
(31)
μ
(
z
s
i
)
=
z
s
i
(
∑
r
=
1
R
〈
u
r
,
u
s
〉
〈
v
r
,
v
s
〉
z
r
i
+
λ
D
i
i
z
s
i
+
σ
∑
k
=
1
N
S
i
k
p
z
s
k
)
.And the final update rule ofz
s
i is
(32)
z
s
i
=
z
s
i
∑
p
=
1
,
q
=
1
m
1
,
m
2
A
p
q
i
u
s
p
v
s
q
+
λ
∑
k
=
1
N
S
i
k
z
s
k
+
σ
D
i
i
p
z
s
i
∑
r
=
1
R
〈
u
r
,
u
s
〉
〈
v
r
,
v
s
〉
z
r
i
+
λ
D
i
i
z
s
i
+
σ
∑
k
=
1
N
S
i
k
p
z
s
k
=
z
s
i
u
s
T
A
;
;
i
v
s
+
λ
S
i
;
z
s
+
σ
D
i
i
p
z
s
i
Z
i
;
(
(
U
T
u
s
)
⊙
(
V
T
v
s
)
)
+
λ
D
i
i
z
s
i
+
σ
λ
S
i
;
p
z
s
,
where Z
i
;
∈
R
1
×
R represents the ith row of Z
=
[
z
1
,
z
2
,
…
,
z
R
]; S
i
;
∈
R
1
×
N and S
i
;
p
∈
R
1
×
N represent the ith row of S and S
~
p, respectively; ⊙ is Hadamard product, A
;
;
i
∈
R
m
1
×
m
2, and is defined as
(33)
(
A
;
;
i
)
p
q
=
A
p
q
i
,
1
≤
∀
p
≤
m
1
,
1
≤
∀
q
≤
m
2
.Nowu
r, v
r, and z
r in the objective function are all calculated.


## 3.1. The Analysis of DNPNTF

Given am
×
N image database X, it contains N sample images X
=
[
x
1
,
…
,
x
n
]
∈
R
m
×
N. The dimension of each sample is m
(
m
=
m
1
×
n
2
). In NTF, the database is organized as a 3rd-order tensor A
∈
R
m
1
×
m
2
×
N, which is overlapped by {
x
i
∈
R
m
1
×
m
2
|
i
=
1
N
} sequentially. The objective function of NTF is
(5)
min
⁡
u
r
,
v
r
,
w
r
∥
A
-
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
∥
2
s
.
t
.
u
r
,
v
r
,
z
r
≥
0
,
where u
r
∈
R
m
1, v
r
∈
R
m
2, and z
r
∈
R
N
(
1
≤
r
≤
R
) describe the first, second, and third modules of A, respectively. Each sample x
i is approximated by
(6)
x
i
≈
∑
r
=
1
R
z
r
i
u
r
v
r
T
,
(
i
=
1
,
…
,
N
)
.By minimizing (5), the bases {
u
r
v
r
T
|
r
=
1
R
} and the corresponding weights {
z
r
|
r
=
1
R
} are conducted. The inner product of {
u
r
v
r
T
|
r
=
1
R
} and the sample image is calculated to derive the low-dimensional parts-based representation.To incorporate more properties into NTF, different constraints could be added into the objective function. The constraint form of objective function is(7)
min
⁡
u
r
,
v
r
,
z
r
∥
A
-
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
∥
2
+
α
F
(
U
,
V
)
+
β
J
(
Z
)
s
.
t
.
u
r
,
v
r
,
z
r
≥
0
,
∀
r
,
where U
=
[
u
1
,
u
2
,
…
,
u
R
]
∈
R
m
1
×
R, V
=
[
v
1
,
v
2
,
…
,
v
R
]
∈
R
m
2
×
R, Z
=
[
z
1
,
z
2
,
…
,
z
R
]
∈
R
N
×
R. F is the constraint function about {
U
,
V
} and J is the constraint about Z. α and β are the corresponding positive coefficients. To encode the spatial structure and discriminant class-based information into sparse representations, we propose a constraint function J according to the manifold learning and graph embedding analysis. In NTF, the columns of the weighting matrix have a one-to-one correspondence with the columns of the original image matrix. Therefore, we add the discriminant constraint to {
z
r
|
r
=
1
R
}, and J for z
r is defined as
(8)
J
(
Z
)
=
min
⁡
⁡
∑
∑
i
≠
j
∑
l
η
l
∥
z
i
-
z
j
∥
2
S
i
j
l
,
where S
i
j
l denote the graphs with different properties and η
l are the corresponding coefficients. By deriving different S
i
j, the graph embedding model could have different properties, such as the neighborhood preserving property and the discriminant property.Now, we discuss the selection ofS
i
j
l. The most commonly used graph is the Laplacian graph, and S
i
j
l is calculated in form of the heat kernel function as
(9)
S
i
j
=
e
-
∥
z
i
-
z
j
∥
2
/
t
,
t
∈
R
.Here,S
i
j measures the similarity between a pair of vertices and has neighborhood preserving property. To further incorporate the class-based discriminant information, we derive a universal penalty graph [27], where the similarity matrix S
i
j
p is defined as
(10)
S
i
j
p
=
{
1
or
1
k
,
if
(
i
,
j
)
∈
N
k
(
l
i
)
or
N
k
(
l
j
)
0
,
otherwise
,
where N
k
(
l
i
) represents the nearest k pairs of samples between class l
i and the other classes. The purpose of the penalty graph was to separate marginal samples between different classes.Now, the objective function of constrained NTF becomes(11)
min
⁡
u
r
,
v
r
,
z
r
∥
A
-
∑
r
=
1
R
u
r

⊗

v
r

⊗

z
r
∥
2

+

λ
∥
z
i

-

z
j
∥
2
S
i
j

-

σ
∥
z
i

-

z
j
∥
2
S
i
j
p
s
.
t
.
u
r
,
v
r
,
z
r
≥
0
,
∀
r
.By solving the generalized eigenvalue decomposition problem, the graph embedding criterion in (11) can be calculated as
(12)
1
2
λ
∥
z
i
-
z
j
∥
2
S
i
j
-
1
2
σ
∥
z
i
-
z
j
∥
2
S
i
j
p
=
1
2
λ
∑
i
,
j
t
r
[
(
z
i
-
z
j
)
(
z
i
-
z
j
)
T
]
S
i
j
-
1
2
σ
∑
i
,
j
t
r
[
(
z
i
-
z
j
)
(
z
i
-
z
j
)
T
]
S
i
j
p
=
1
2
t
r
[
λ
∑
i
,
j
[
(
z
i
-
z
j
)
(
z
i
-
z
j
)
T
]
S
i
j
-
σ
∑
i
,
j
[
(
z
i
-
z
j
)
(
z
i
-
z
j
)
T
]
S
i
j
p
]
=
t
r
[
λ
∑
i
,
j
(
z
i
S
i
j
(
z
i
)
T
-
z
i
S
i
j
(
z
j
)
T
)
-
σ
∑
i
,
j
(
z
i
S
i
j
p
(
z
i
)
T
-
z
i
S
i
j
p
(
z
j
)
T
)
]
=
t
r
[
λ
∑
i
,
j
(
z
i
E
i
j
(
z
i
)
T
-
z
i
E
i
j
(
z
j
)
T
)
-
σ
∑
i
,
j
(
z
i
E
i
j
p
(
z
i
)
T
-
z
i
E
i
j
p
(
z
j
)
T
)
]
=
t
r
[
λ
∑
i
,
j
z
i
L
i
j
(
z
j
)
T
-
σ
∑
i
,
j
z
i
L
i
j
p
(
z
j
)
T
]
=
λ
t
r
(
Z
T
L
Z
)
-
σ
t
r
(
Z
T
L
p
Z
)
.
And the final objective function of DNPNTF is
(13)
min
⁡
u
r
,
v
r
,
z
r
∥
A
-
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
∥
2

+

λ
t
r
(
Z
T
L
Z
)

-

σ
t
r
(
Z
T
L
p
Z
)
s
.
t
.
u
r
,
v
r
,
z
r
≥
0
,
∀
r
,
where λ
>
0 and σ
>
0 to make sure (13) should be nonnegative.


## 3.2. Projected Gradient Method of DNPNTF

The most popular approach to minimize NMF or NTF is the multiplicative update method. However, it cannot ensure the convergence of the constraint forms of NMF or NTF. In this paper, the projected gradient method is used to solve DNPNTF.The objective function of DNPNTF can be stated as(14)
D
obj
(
A
∥
U
V
Z
)
=
D
(
A
∥
U
V
Z
)
+
ς
D
(
Z
)
,
where ς is a positive constant. The goal of (14) is to find {
U
,
V
} and Z by solving the following problem:
(15)
min
⁡
U
,
V
,
Z
D
obj
(
A
∥
U
V
Z
)
s
.
t
.
u
i
,
k
≥
0
,
v
k
,
j
≥
0
,
z
i
,
j
≥
0
∀
k
,
i
,
j
.To find the optimal solution, (15) is divided into three subproblems: first, we fix V and Z and update U to arrive at the conditional optimal value of the subminimization problem; second, we fix U and Z and update V; last, we fix U and V and update Z. Three functions are defined as f
V
,
Z
(
U
)
=
D
obj
(
A
∥
U
V
Z
), f
U
,
Z
(
V
)
=
D
obj
(
A
∥
U
V
Z
), and f
U
,
V
(
Z
)
=
D
obj
(
A
∥
U
V
Z
). The update rules are defined as
(16)
U
(
t
+
1
)
=
[
U
(
t
)
-
α
t
∇
f
V
,
Z
(
U
(
t
)
)
]
,
V
(
t
+
1
)
=
[
V
(
t
)
-
α
t
∇
f
U
,
Z
(
V
(
t
)
)
]
,
Z
(
t
+
1
)
=
[
Z
(
t
)
-
α
t
∇
f
U
,
V
(
Z
(
t
)
)
]
.
Now the task is calculating ∇
f
V
,
Z
(
U
(
t
)
), ∇
f
U
,
Z
(
V
(
t
)
), and ∇
f
U
,
V
(
Z
(
t
)
).


### 3.2.1. The Calculation ofU and V

Firstly, we discuss the calculation of∇
f
V
,
Z
(
U
(
t
)
) and ∇
f
U
,
Z
(
V
(
t
)
). The objective function could be written as
(17)
f
=
1
2
∥
A
-
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
∥
2
+
1
2
λ
t
r
(
Z
T
L
Z
)
-
1
2
σ
t
r
(
Z
T
L
p
Z
)
=
1
2
〈
A
-
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
,
A
-
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
〉
+
1
2
λ
∑
r
=
1
R
z
r
T
L
z
r
-
1
2
σ
∑
r
=
1
R
z
r
T
L
p
z
r
.The differential off is
(18)
d
(
f
u
s
)
=
〈
A
-
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
,
-
d
(
u
s
)
⊗
v
s
⊗
z
s
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
〉
+
0
-
0
=
〈
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
,
d
(
u
s
)
⊗
v
s
⊗
z
s
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
〉
-
〈
A
,
d
(
u
s
)
⊗
v
s
⊗
z
s
〉
.And the partial differential foru
s
p
(
1
≤
p
≤
m
1
) is
(19)
∂
f
∂
u
s
p
=
〈
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
,
e
p
⊗
v
s
⊗
z
s
〉
-
〈
A
,
e
p
⊗
v
s
⊗
z
s
〉
,
where the pth element in e
p
∈
R
m
1 is 1 and others are 0. That is, (
e
p
)
p
=
1 and (
e
p
)
k
≠
p
=
0. According to Definition 1, for any order tensors A
1
,
A
2
∈
R
a
1
×
a
2
×
⋯
×
a
n, B
1
,
B
2
∈
R
b
1
×
b
2
×
⋯
×
b
n, there is 〈
A
1
⊗
B
1
,
A
2
⊗
B
2
〉
=
〈
A
1
,
A
2
〉
〈
B
1
,
B
2
〉. Then (19) could be written as
(20)
∂
f
∂
u
s
p
=
∑
r
=
1
R
(
〈
u
r
,
e
p
〉
〈
v
r
⊗
z
r
,
v
s
⊗
z
s
〉
)
-
(
∑
k
=
1
,
q
=
1
,
i
=
1
m
1
,
m
2
,
N
A
k
q
i
(
e
p
)
k
v
s
q
z
s
i
)
=
∑
r
=
1
R
(
〈
u
r
,
e
p
〉
〈
v
r
⊗
z
r
〉
〈
v
s
⊗
z
s
〉
)
-
(
∑
q
=
1
,
i
=
1
m
2
,
N
A
p
q
i
(
e
p
)
p
v
s
q
z
s
i
+
∑
k
≠
p
A
k
q
i
(
e
p
)
k
v
s
q
z
s
i
)
=
∑
r
=
1
R
(
〈
u
r
,
e
p
〉
〈
v
r
,
v
s
〉
〈
z
r
,
z
s
〉
)
-
(
∑
q
=
1
,
i
=
1
m
2
,
N
A
p
q
i
v
s
q
z
s
i
+
0
)
=
∑
r
=
1
R
(
(
u
r
p
(
e
p
)
p
+
∑
k
≠
p
u
r
k
(
e
p
)
k
)
〈
v
r
,
v
s
〉
〈
z
r
,
z
s
〉
)
-
∑
q
=
1
,
i
=
1
m
2
,
N
A
p
q
i
v
s
q
z
s
i
=
∑
r
=
1
R
(
u
r
p
〈
v
r
,
v
s
〉
〈
z
r
,
z
s
〉
)
-
∑
q
=
1
,
i
=
1
m
2
,
N
A
p
q
i
v
s
q
z
s
i
.According to (16), the update rule for u
s
p is
(21)
u
s
p
=
u
s
p
-
μ
(
u
s
p
)
∂
f
∂
u
s
p
=
u
s
p
-
μ
(
u
s
p
)
(
∑
r
=
1
R
(
u
r
p
〈
v
r
,
v
s
〉
〈
z
r
,
z
s
〉
)
-
∑
q
=
1
,
i
=
1
m
2
,
N
A
p
q
i
v
s
q
z
s
i
)
.To confirm the nonnegative ofu
s
p, μ
(
u
s
p
) is set to be
(22)
μ
(
u
s
p
)
=
u
s
p
∑
r
=
1
R
(
u
r
p
〈
v
r
,
v
s
〉
〈
z
r
,
z
s
〉
)
.
Then the update function is
(23)
u
s
p
=
u
s
p
(
∑
q
=
1
,
i
=
1
m
2
,
N
A
p
q
i
v
s
q
z
s
i
)
∑
r
=
1
R
(
u
r
p
〈
v
r
,
v
s
〉
〈
z
r
,
z
s
〉
)
=
u
s
p
(
v
s
T
A
p
;
;
z
s
)
∑
r
=
1
R
(
u
r
p
(
v
r
T
v
s
)
(
z
r
T
z
s
)
)
=
u
s
p
(
v
s
T
A
p
;
;
z
s
)
(
U
p
;
(
(
V
T
v
s
)
⊙
(
Z
T
z
s
)
)
)
,
where U
p
;
∈
R
1
×
R represents the pth row of U
=
[
u
1
,
u
2
,
…
,
u
R
], ⊙ is the matrix Hadamard product, A
p
;
;
∈
R
m
2
×
N represents the matrix which fixes the first module of A, and traversals are the other two modules. It is defined as
(24)
(
A
p
;
;
)
q
i
=
A
p
q
i
,
1
≤
∀
q
≤
m
2
,
1
≤
∀
i
≤
N
.Similarly, the update rule of theqth element of v
s
(
v
s
q
,
1
≤
s
≤
R
,
1
≤
q
≤
m
2
) is
(25)
v
s
q
=
v
s
q
∑
p
=
1
,
i
=
1
m
1
,
N
(
A
p
q
i
u
s
p
z
s
i
)
∑
r
=
1
R
(
v
r
q
〈
u
r
,
u
s
〉
〈
z
r
,
z
s
〉
)
=
v
s
q
(
u
s
T
A
;
q
;
z
s
)
(
V
q
;
(
(
U
T
u
s
)
⊙
(
Z
T
z
s
)
)
)
,
where V
q
;
∈
R
1
×
R represents the pth row of V
=
[
v
1
,
v
2
,
…
,
v
R
]
∈
R
m
2
×
R, ⊙ is the matrix Hadamard product, A
;
q
;
∈
R
m
1
×
N represents the matrix which fixes the second module of A, and traversals are the other two modules. It is defined as
(26)
(
A
;
q
;
)
p
i
=
A
p
q
i
,
1
≤
∀
p
≤
m
1
,
1
≤
∀
i
≤
N
.Now,v
r and u
r are calculated.


### 3.2.2. The Calculation ofZ

Then we discuss the calculation of∇
f
U
,
V
(
Z
(
t
)
). The differential of f along z
s
(
∀
s
,
1
≤
s
≤
R
) is
(27)
d
(
f
z
s
)
=
〈
A
-
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
,
-
d
(
z
s
)
⊗
u
s
⊗
v
s
〉
+
1
2
d
(
λ
z
s
T
L
z
s
)
-
1
2
d
(
σ
z
s
T
L
p
z
s
)
=
〈
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
,
d
(
z
s
)
⊗
u
s
⊗
v
s
〉
-
〈
A
,
d
(
z
s
)
⊗
u
s
⊗
v
s
〉
+
d
(
1
2
λ
z
s
T
(
D
-
S
)
z
s
)
-
d
(
1
2
σ
z
s
T
(
D
p
-
S
p
)
z
s
)
.For(
λ
z
s
T
(
D
-
S
)
z
s
)
/
2, the partial differential for z
s
i is
(28)
∂
(
(
1
/
2
)
λ
z
s
T
(
D
-
S
)
z
s
)
∂
z
s
i
=
∂
(
z
s
)
∂
z
s
i
∂
(
(
1
/
2
)
λ
z
s
T
(
D
-
S
)
z
s
)
∂
z
s
=
(
e
i
)
T
λ
(
D
-
S
)
z
s
,
where the ith element in e
i
∈
R
N is 1 and others are 0. That is, (
e
i
)
i
=
1 and (
e
i
)
k
≠
i
=
0. Then the partial differential for z
s
i is
(29)
∂
f
∂
z
s
i
=
〈
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
,
u
s
⊗
v
s
⊗
e
i
〉
-
〈
A
,
u
s
⊗
v
s
⊗
e
i
〉
+
λ
(
e
i
)
T
(
D
-
S
)
z
s
-
σ
(
e
i
)
T
(
D
p
-
S
p
)
z
s
=
∑
r
=
1
R
〈
u
r
,
u
s
〉
〈
v
r
,
v
s
〉
z
r
i
-
∑
p
=
1
,
q
=
1
m
1
,
m
2
A
p
q
i
u
s
p
v
s
q
+
(
λ
D
i
i
z
s
i
-
λ
∑
k
=
1
N
S
i
k
z
s
k
)
-
(
σ
D
i
i
p
z
s
i
-
σ
∑
k
=
1
N
S
i
k
p
z
s
k
)
.According to (16), the update rule for z
s
i is
(30)
z
s
i
=
z
s
i
-
μ
(
z
s
i
)
∂
f
∂
z
s
i
=
z
s
i
-
μ
(
z
s
i
)
(
∑
r
=
1
R
〈
u
r
,
u
s
〉
〈
v
r
,
v
s
〉
z
r
i
+
λ
D
i
i
z
s
i
+
σ
∑
k
=
1
N
S
i
k
p
z
s
k
-
∑
p
=
1
,
q
=
1
m
1
,
m
2
A
p
q
i
u
s
p
v
s
q
-
λ
∑
k
=
1
N
S
i
k
z
s
k
-
σ
D
i
i
p
z
s
i
)
.The update stepμ
(
z
s
i
) is set as
(31)
μ
(
z
s
i
)
=
z
s
i
(
∑
r
=
1
R
〈
u
r
,
u
s
〉
〈
v
r
,
v
s
〉
z
r
i
+
λ
D
i
i
z
s
i
+
σ
∑
k
=
1
N
S
i
k
p
z
s
k
)
.And the final update rule ofz
s
i is
(32)
z
s
i
=
z
s
i
∑
p
=
1
,
q
=
1
m
1
,
m
2
A
p
q
i
u
s
p
v
s
q
+
λ
∑
k
=
1
N
S
i
k
z
s
k
+
σ
D
i
i
p
z
s
i
∑
r
=
1
R
〈
u
r
,
u
s
〉
〈
v
r
,
v
s
〉
z
r
i
+
λ
D
i
i
z
s
i
+
σ
∑
k
=
1
N
S
i
k
p
z
s
k
=
z
s
i
u
s
T
A
;
;
i
v
s
+
λ
S
i
;
z
s
+
σ
D
i
i
p
z
s
i
Z
i
;
(
(
U
T
u
s
)
⊙
(
V
T
v
s
)
)
+
λ
D
i
i
z
s
i
+
σ
λ
S
i
;
p
z
s
,
where Z
i
;
∈
R
1
×
R represents the ith row of Z
=
[
z
1
,
z
2
,
…
,
z
R
]; S
i
;
∈
R
1
×
N and S
i
;
p
∈
R
1
×
N represent the ith row of S and S
~
p, respectively; ⊙ is Hadamard product, A
;
;
i
∈
R
m
1
×
m
2, and is defined as
(33)
(
A
;
;
i
)
p
q
=
A
p
q
i
,
1
≤
∀
p
≤
m
1
,
1
≤
∀
q
≤
m
2
.Nowu
r, v
r, and z
r in the objective function are all calculated.


## 3.2.1. The Calculation ofU and V

Firstly, we discuss the calculation of∇
f
V
,
Z
(
U
(
t
)
) and ∇
f
U
,
Z
(
V
(
t
)
). The objective function could be written as
(17)
f
=
1
2
∥
A
-
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
∥
2
+
1
2
λ
t
r
(
Z
T
L
Z
)
-
1
2
σ
t
r
(
Z
T
L
p
Z
)
=
1
2
〈
A
-
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
,
A
-
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
〉
+
1
2
λ
∑
r
=
1
R
z
r
T
L
z
r
-
1
2
σ
∑
r
=
1
R
z
r
T
L
p
z
r
.The differential off is
(18)
d
(
f
u
s
)
=
〈
A
-
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
,
-
d
(
u
s
)
⊗
v
s
⊗
z
s
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
〉
+
0
-
0
=
〈
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
,
d
(
u
s
)
⊗
v
s
⊗
z
s
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
〉
-
〈
A
,
d
(
u
s
)
⊗
v
s
⊗
z
s
〉
.And the partial differential foru
s
p
(
1
≤
p
≤
m
1
) is
(19)
∂
f
∂
u
s
p
=
〈
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
,
e
p
⊗
v
s
⊗
z
s
〉
-
〈
A
,
e
p
⊗
v
s
⊗
z
s
〉
,
where the pth element in e
p
∈
R
m
1 is 1 and others are 0. That is, (
e
p
)
p
=
1 and (
e
p
)
k
≠
p
=
0. According to Definition 1, for any order tensors A
1
,
A
2
∈
R
a
1
×
a
2
×
⋯
×
a
n, B
1
,
B
2
∈
R
b
1
×
b
2
×
⋯
×
b
n, there is 〈
A
1
⊗
B
1
,
A
2
⊗
B
2
〉
=
〈
A
1
,
A
2
〉
〈
B
1
,
B
2
〉. Then (19) could be written as
(20)
∂
f
∂
u
s
p
=
∑
r
=
1
R
(
〈
u
r
,
e
p
〉
〈
v
r
⊗
z
r
,
v
s
⊗
z
s
〉
)
-
(
∑
k
=
1
,
q
=
1
,
i
=
1
m
1
,
m
2
,
N
A
k
q
i
(
e
p
)
k
v
s
q
z
s
i
)
=
∑
r
=
1
R
(
〈
u
r
,
e
p
〉
〈
v
r
⊗
z
r
〉
〈
v
s
⊗
z
s
〉
)
-
(
∑
q
=
1
,
i
=
1
m
2
,
N
A
p
q
i
(
e
p
)
p
v
s
q
z
s
i
+
∑
k
≠
p
A
k
q
i
(
e
p
)
k
v
s
q
z
s
i
)
=
∑
r
=
1
R
(
〈
u
r
,
e
p
〉
〈
v
r
,
v
s
〉
〈
z
r
,
z
s
〉
)
-
(
∑
q
=
1
,
i
=
1
m
2
,
N
A
p
q
i
v
s
q
z
s
i
+
0
)
=
∑
r
=
1
R
(
(
u
r
p
(
e
p
)
p
+
∑
k
≠
p
u
r
k
(
e
p
)
k
)
〈
v
r
,
v
s
〉
〈
z
r
,
z
s
〉
)
-
∑
q
=
1
,
i
=
1
m
2
,
N
A
p
q
i
v
s
q
z
s
i
=
∑
r
=
1
R
(
u
r
p
〈
v
r
,
v
s
〉
〈
z
r
,
z
s
〉
)
-
∑
q
=
1
,
i
=
1
m
2
,
N
A
p
q
i
v
s
q
z
s
i
.According to (16), the update rule for u
s
p is
(21)
u
s
p
=
u
s
p
-
μ
(
u
s
p
)
∂
f
∂
u
s
p
=
u
s
p
-
μ
(
u
s
p
)
(
∑
r
=
1
R
(
u
r
p
〈
v
r
,
v
s
〉
〈
z
r
,
z
s
〉
)
-
∑
q
=
1
,
i
=
1
m
2
,
N
A
p
q
i
v
s
q
z
s
i
)
.To confirm the nonnegative ofu
s
p, μ
(
u
s
p
) is set to be
(22)
μ
(
u
s
p
)
=
u
s
p
∑
r
=
1
R
(
u
r
p
〈
v
r
,
v
s
〉
〈
z
r
,
z
s
〉
)
.
Then the update function is
(23)
u
s
p
=
u
s
p
(
∑
q
=
1
,
i
=
1
m
2
,
N
A
p
q
i
v
s
q
z
s
i
)
∑
r
=
1
R
(
u
r
p
〈
v
r
,
v
s
〉
〈
z
r
,
z
s
〉
)
=
u
s
p
(
v
s
T
A
p
;
;
z
s
)
∑
r
=
1
R
(
u
r
p
(
v
r
T
v
s
)
(
z
r
T
z
s
)
)
=
u
s
p
(
v
s
T
A
p
;
;
z
s
)
(
U
p
;
(
(
V
T
v
s
)
⊙
(
Z
T
z
s
)
)
)
,
where U
p
;
∈
R
1
×
R represents the pth row of U
=
[
u
1
,
u
2
,
…
,
u
R
], ⊙ is the matrix Hadamard product, A
p
;
;
∈
R
m
2
×
N represents the matrix which fixes the first module of A, and traversals are the other two modules. It is defined as
(24)
(
A
p
;
;
)
q
i
=
A
p
q
i
,
1
≤
∀
q
≤
m
2
,
1
≤
∀
i
≤
N
.Similarly, the update rule of theqth element of v
s
(
v
s
q
,
1
≤
s
≤
R
,
1
≤
q
≤
m
2
) is
(25)
v
s
q
=
v
s
q
∑
p
=
1
,
i
=
1
m
1
,
N
(
A
p
q
i
u
s
p
z
s
i
)
∑
r
=
1
R
(
v
r
q
〈
u
r
,
u
s
〉
〈
z
r
,
z
s
〉
)
=
v
s
q
(
u
s
T
A
;
q
;
z
s
)
(
V
q
;
(
(
U
T
u
s
)
⊙
(
Z
T
z
s
)
)
)
,
where V
q
;
∈
R
1
×
R represents the pth row of V
=
[
v
1
,
v
2
,
…
,
v
R
]
∈
R
m
2
×
R, ⊙ is the matrix Hadamard product, A
;
q
;
∈
R
m
1
×
N represents the matrix which fixes the second module of A, and traversals are the other two modules. It is defined as
(26)
(
A
;
q
;
)
p
i
=
A
p
q
i
,
1
≤
∀
p
≤
m
1
,
1
≤
∀
i
≤
N
.Now,v
r and u
r are calculated.


## 3.2.2. The Calculation ofZ

Then we discuss the calculation of∇
f
U
,
V
(
Z
(
t
)
). The differential of f along z
s
(
∀
s
,
1
≤
s
≤
R
) is
(27)
d
(
f
z
s
)
=
〈
A
-
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
,
-
d
(
z
s
)
⊗
u
s
⊗
v
s
〉
+
1
2
d
(
λ
z
s
T
L
z
s
)
-
1
2
d
(
σ
z
s
T
L
p
z
s
)
=
〈
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
,
d
(
z
s
)
⊗
u
s
⊗
v
s
〉
-
〈
A
,
d
(
z
s
)
⊗
u
s
⊗
v
s
〉
+
d
(
1
2
λ
z
s
T
(
D
-
S
)
z
s
)
-
d
(
1
2
σ
z
s
T
(
D
p
-
S
p
)
z
s
)
.For(
λ
z
s
T
(
D
-
S
)
z
s
)
/
2, the partial differential for z
s
i is
(28)
∂
(
(
1
/
2
)
λ
z
s
T
(
D
-
S
)
z
s
)
∂
z
s
i
=
∂
(
z
s
)
∂
z
s
i
∂
(
(
1
/
2
)
λ
z
s
T
(
D
-
S
)
z
s
)
∂
z
s
=
(
e
i
)
T
λ
(
D
-
S
)
z
s
,
where the ith element in e
i
∈
R
N is 1 and others are 0. That is, (
e
i
)
i
=
1 and (
e
i
)
k
≠
i
=
0. Then the partial differential for z
s
i is
(29)
∂
f
∂
z
s
i
=
〈
∑
r
=
1
R
u
r
⊗
v
r
⊗
z
r
,
u
s
⊗
v
s
⊗
e
i
〉
-
〈
A
,
u
s
⊗
v
s
⊗
e
i
〉
+
λ
(
e
i
)
T
(
D
-
S
)
z
s
-
σ
(
e
i
)
T
(
D
p
-
S
p
)
z
s
=
∑
r
=
1
R
〈
u
r
,
u
s
〉
〈
v
r
,
v
s
〉
z
r
i
-
∑
p
=
1
,
q
=
1
m
1
,
m
2
A
p
q
i
u
s
p
v
s
q
+
(
λ
D
i
i
z
s
i
-
λ
∑
k
=
1
N
S
i
k
z
s
k
)
-
(
σ
D
i
i
p
z
s
i
-
σ
∑
k
=
1
N
S
i
k
p
z
s
k
)
.According to (16), the update rule for z
s
i is
(30)
z
s
i
=
z
s
i
-
μ
(
z
s
i
)
∂
f
∂
z
s
i
=
z
s
i
-
μ
(
z
s
i
)
(
∑
r
=
1
R
〈
u
r
,
u
s
〉
〈
v
r
,
v
s
〉
z
r
i
+
λ
D
i
i
z
s
i
+
σ
∑
k
=
1
N
S
i
k
p
z
s
k
-
∑
p
=
1
,
q
=
1
m
1
,
m
2
A
p
q
i
u
s
p
v
s
q
-
λ
∑
k
=
1
N
S
i
k
z
s
k
-
σ
D
i
i
p
z
s
i
)
.The update stepμ
(
z
s
i
) is set as
(31)
μ
(
z
s
i
)
=
z
s
i
(
∑
r
=
1
R
〈
u
r
,
u
s
〉
〈
v
r
,
v
s
〉
z
r
i
+
λ
D
i
i
z
s
i
+
σ
∑
k
=
1
N
S
i
k
p
z
s
k
)
.And the final update rule ofz
s
i is
(32)
z
s
i
=
z
s
i
∑
p
=
1
,
q
=
1
m
1
,
m
2
A
p
q
i
u
s
p
v
s
q
+
λ
∑
k
=
1
N
S
i
k
z
s
k
+
σ
D
i
i
p
z
s
i
∑
r
=
1
R
〈
u
r
,
u
s
〉
〈
v
r
,
v
s
〉
z
r
i
+
λ
D
i
i
z
s
i
+
σ
∑
k
=
1
N
S
i
k
p
z
s
k
=
z
s
i
u
s
T
A
;
;
i
v
s
+
λ
S
i
;
z
s
+
σ
D
i
i
p
z
s
i
Z
i
;
(
(
U
T
u
s
)
⊙
(
V
T
v
s
)
)
+
λ
D
i
i
z
s
i
+
σ
λ
S
i
;
p
z
s
,
where Z
i
;
∈
R
1
×
R represents the ith row of Z
=
[
z
1
,
z
2
,
…
,
z
R
]; S
i
;
∈
R
1
×
N and S
i
;
p
∈
R
1
×
N represent the ith row of S and S
~
p, respectively; ⊙ is Hadamard product, A
;
;
i
∈
R
m
1
×
m
2, and is defined as
(33)
(
A
;
;
i
)
p
q
=
A
p
q
i
,
1
≤
∀
p
≤
m
1
,
1
≤
∀
q
≤
m
2
.Nowu
r, v
r, and z
r in the objective function are all calculated.


## 4. Extreme Learning Machine

ELM is proposed by Huang et al. [18] for SLFNs. Unlike the traditional feedforward neural network training methods, such as the gradient-descent method, the standard optimization method, and the least-square based method, ELM need not tune the hidden layer of SLFNs which may cause learning complicated and inefficient. It could reach the smallest training error and have better generalization performance. The learning speed of ELM is fast, and the parameters have not to be tuned manually. In our proposed algorithm, the extracted features by DNPNTF are fed into ELM for classification.Given a training setX
=
{
(
x
i
,
t
i
)
∣
x
i
∈
R
n
,
t
i
∈
R
n
,
i
=
1,2
,
…
,
N
}, where x
i is the n
×
1 input feature vector and t
i is a m
×
1 target vector. ELM with K hidden nodes and activation function g
(
x
) is modeled as
(34)
∑
i
=
1
K
β
i
g
(
x
j
)
=
∑
i
=
1
K
β
i
g
(
w
i
·
x
j
+
b
i
)
=
o
j
,
j
=
1
,
…
,
N
,
where w
i
=
(
w
i
1
,
w
i
2
,
…
,
w
i
d
)
T represents the input weight, which is the ith neuron in the hidden layer and the input layer; β
i
=
(
β
i
1
,
β
i
2
,
…
,
β
i
m
)
T is the weight vector between the ith hidden neuron and the output layer; o
j is the target vector of the jth input data. In training step, ELM aims to approximate N training samples with zero error, which means ∑
j
=
1
N
∥
o
j
-
t
j
∥
=
0. Then there exist β
i, w
j, and b
j satisfying that
(35)
∑
i
=
1
K
β
i
g
(
w
i
·
x
j
+
b
i
)
=
t
j
,
j
=
1
,
…
,
N
.Equation (35) can be reformulated compactly as
(36)
H
β
=
T
,
where
(37)
H
(
w
1
,
…
,
w
K
,
b
1
,
…
,
b
K
,
x
1
,
…
,
x
N
)
=
(
g
(
w
1
·
x
1
+
b
1
)
⋯
g
(
w
K
·
x
1
+
b
K
)
⋮
⋱
⋮
g
(
w
1
·
x
N
+
b
1
)
⋯
g
(
w
K
·
x
N
+
b
K
)
)
N
×
K
,
β
=
(
β
1
T
⋮
β
K
T
)
K
×
m
,
T
=
(
t
1
T
⋮
t
N
T
)
N
×
m
.
H is called the hidden layer output matrix of the neural network, and the ith column of H is the ith hidden neuron output with respect to inputs x
1
,
x
2
,
…
,
x
N. It is proved by Huang et al. [18] that weights and biases need not be adjusted and can be arbitrarily given. Therefore, the output weights could be determined by finding the least-square solution β
^ as
(38)
H
β
=
T
⟹
β
^
=
H
†
T
,
where H
† is the Moore-Penrose generalized inverse of matrix H. Furthermore, the smallest training error can be obtained by β
^ as
(39)
∥
H
β
^
-
T
∥
=
∥
H
H
†
T
-
T
∥
=
min
⁡
β
⁡
∥
H
β
-
T
∥
.
As analyzed by Huang, ELM could obtain a good generalization performance with a dramatically increased learning speed by solving (39).


## 5. Experiments

In this section, we apply DNPNTF via ELM to facial expression recognition. We compare DNPNTF with NMF [7], DNMF [11], and NTF [9] and give the experimental results by employing ELM, NN, SVM [16], and SRC [17]. Two facial expression databases are used: the JAFFE database [31] and the Cohn-Kanade database [32]. Raw facial images are cropped according to the position of eyes and normalized to 32 × 32 pixels. Figure 1 shows an example of the original face image and the corresponding cropped image. According to the rank-one tensor theory, gray level images are encoded in tensor space.Figure 1
Illustration of the preprocess for original images.Since the results of ELM may vary during each different execution, we repeat the execution for 5 times and take the average value as the final result. It is proved by theory analysis and experiments that the classification performance of ELM is affected by the hidden activation function and the number of hidden nodes [23]. However, in this paper we just focus on the application of ELM to facial expression recognition. The activation function used in our algorithm is a simple sigmoidal function. The number of hidden nodes is set to be the same as the number of facial expression classes (e.g., 7 for the JAFFE database and 6 for the Cohn-Kanade database). For SVM, the radial basis function (RBF) is used, which is exp
⁡
⁡
(
-
γ
∥
x
i
-
x
j
∥
2
)
·
γ that is set to be 3 as an empirical value. For SRC, “Homotopy” algorithm is used to solve the minimization of l
1 norm constraint.


### 5.1. Experiments on JAFFE Database

The JAFFE database [31] is an expression database which contains 213 static facial images captured from 10 Japanese females. Each person poses 2 to 4 examples for each of the 6 prototypic expressions (anger, disgust, fear, happiness, sadness, and surprise) plus the natural face. To evaluate the algorithms, we randomly partition all images into 10 groups, with roughly 70 samples in each group. We take any 9 groups for training and calculate the recognition rates with the remaining one. We repeat it for all the 10 possible choices. Finally, the average result over 10 times’ testing was taken.The average recognition rates of different feature extraction algorithms are shown in Figure2, where the vertical axis represents the correct recognition rate in percentage and the horizontal axis represents the corresponding dimensions (from 1 to 120). Here, only the NN classifier is used. In the lower range of dimensions, the recognition rates of DNPNTF are similar to other algorithms. This is because DNPNTF extracts the parts-based sparse representations, and only a few features could be generated for recognition in the low range of dimensions. In the higher range of dimensions, DNPNTF outperforms the others. With the increase of the extracted parts-based features, DNPNTF could achieve good recognition performance. Since different constraints were added, the improved versions of NMF, including DNMF and NTF, outperform the conventional NMF.Figure 2
The recognition rate versus dimension curves achieved on the JAFFE database.The top recognition rates of different algorithms with corresponding dimensions are illustrated in Table1. NMF achieves the highest rate at a low dimension, while DNMF achieves the highest rate at a high dimension. Although more dimensions are needed, DNPNTF achieves the highest recognition rate compared with others. This is because the constraints about manifold structure and discriminant information are considered, which are critical for classification.Table 1
The top recognition rate and the corresponding dimension achieved on the JAFFE database.
DNPNTF
NMF
DNMF
NTF
Recognition rates (%)
90.71
87.43
89.33
88.91
Dimensions
98
16
66
101Figure3 shows the basis images obtained on the JAFFE database by NMF, NTF, and DNPNTF. Based on the principle of NMF, the face images are represented by combining multiple basis images with addition only, and the basis images are expected to represent facial parts. In this database, the basis images calculated by NMF are not sparse. NTF and DNPNTF which execute in the tensor space could generate parts-based sparse representations. Since more constraints were adopted, DNPNTF generate sparser basis images which reflect distinct features for recognition.Basis images obtained on the JAFFE database.
(a)
NMF
(b)
NTF
(c)
DNPNTFThen we give the experiments to prove the effectiveness of DNPNTF via ELM. The average recognition rates of DNPNTF with ELM, NN, SVM, and SRC are given in Figure4, where the vertical axis represents the correct recognition rate in percentage and the horizontal axis represents the corresponding dimensions (from 1 to 120). ELM and SRC achieve better recognition performance compared with NN and SVM, and ELM achieves the highest recognition rate. The top recognition rates with the corresponding dimensions are given in Table 2.Table 2
The top recognition rate of different classifiers achieved on the JAFFE database.
ELM
NN
SRC
SVM
Recognition rates (%)
94.39
90.71
93.54
91.03
Dimensions
55
98
63
101Figure 4
The recognition rate versus dimension curves of DNPNTF by using different classifiers on the JAFFE database.


### 5.2. Experiments on Cohn-Kanade Database

The Cohn-Kanade database [32] consists of a large amount of image sequences starting from natural face and ending with the peak of the corresponding expression. 104 subjects with different ages, genders, and races are instructed to pose a series of 23 facial displays, including the 6 prototypic expressions. In our experiments, for every image sequence, we take 2 to 8 continuous frames near the peak expression as the static samples. We use the face images of all subjects. We partition the subject to 3 exclusive groups, and in each group, for each of the prototypic expression, we select 100 samples; that is, there are 600 samples in each group and the size of the total set is 1800. During the experiment, we adopt the leave-one-group-out strategy and 3-fold cross-validation: each time two groups are taken as training set and the remaining group is left for testing. This procedure is repeated for 3 times.The average recognition rates of different algorithms on the Cohn-Kanade database are shown in Figure5, where the vertical axis represents the correct recognition rate in percentage and the horizontal axis represents the corresponding dimensions (from 1 to 120). Here, only the NN classifier is used. Table 3 shows the top recognition rates with the corresponding dimensions. The recognition rates obtained on the Cohn-Kanade database are lower than those obtained on the JAFFE database. It can be explained that the experiments on the Cohn-Kanade database are person-independent, which are more difficult than the person-dependent experiments on the JAFFE database. From Figure 5, we can see that the performance of DNPNTF is superior to others with nearly all dimensions. Its recognition rates improve with the increase of dimensions.Table 3
The top recognition rate and the corresponding dimensions achieved on the Cohn-Kanade database.
DNPNTF
NMF
DNMF
NTF
Recognition rates (%)
63.61
58.17
58.83
59.72
Dimensions
101
76
26
120Figure 5
The recognition rate versus dimension curves achieved on the Cohn-Kanade database.Lastly, we give the experiments about different classifiers on the Cohn-Kanade database. The average recognition rates of DNPNTF via ELM, NN, SVM, and SRC are shown in Figure6, and the top recognition rates are given in Table 4. SVM and SRC achieve better performance compared with NN. ELM achieves the best recognition accuracy among all tested algorithms on almost all dimensions. It means ELM could use the information contained in extracted features better than other classifiers.Table 4
The top recognition rate and the corresponding dimensions achieved on the Cohn-Kanade database.
ELM
NN
SRC
SVM
Recognition rates (%)
86.00
63.61
73.00
72.06
Dimensions
51
101
116
106Figure 6
The recognition rate versus dimension curves of DNPNTF by using different classifiers on the Cohn-Kanade database.


## 5.1. Experiments on JAFFE Database

The JAFFE database [31] is an expression database which contains 213 static facial images captured from 10 Japanese females. Each person poses 2 to 4 examples for each of the 6 prototypic expressions (anger, disgust, fear, happiness, sadness, and surprise) plus the natural face. To evaluate the algorithms, we randomly partition all images into 10 groups, with roughly 70 samples in each group. We take any 9 groups for training and calculate the recognition rates with the remaining one. We repeat it for all the 10 possible choices. Finally, the average result over 10 times’ testing was taken.The average recognition rates of different feature extraction algorithms are shown in Figure2, where the vertical axis represents the correct recognition rate in percentage and the horizontal axis represents the corresponding dimensions (from 1 to 120). Here, only the NN classifier is used. In the lower range of dimensions, the recognition rates of DNPNTF are similar to other algorithms. This is because DNPNTF extracts the parts-based sparse representations, and only a few features could be generated for recognition in the low range of dimensions. In the higher range of dimensions, DNPNTF outperforms the others. With the increase of the extracted parts-based features, DNPNTF could achieve good recognition performance. Since different constraints were added, the improved versions of NMF, including DNMF and NTF, outperform the conventional NMF.Figure 2
The recognition rate versus dimension curves achieved on the JAFFE database.The top recognition rates of different algorithms with corresponding dimensions are illustrated in Table1. NMF achieves the highest rate at a low dimension, while DNMF achieves the highest rate at a high dimension. Although more dimensions are needed, DNPNTF achieves the highest recognition rate compared with others. This is because the constraints about manifold structure and discriminant information are considered, which are critical for classification.Table 1
The top recognition rate and the corresponding dimension achieved on the JAFFE database.
DNPNTF
NMF
DNMF
NTF
Recognition rates (%)
90.71
87.43
89.33
88.91
Dimensions
98
16
66
101Figure3 shows the basis images obtained on the JAFFE database by NMF, NTF, and DNPNTF. Based on the principle of NMF, the face images are represented by combining multiple basis images with addition only, and the basis images are expected to represent facial parts. In this database, the basis images calculated by NMF are not sparse. NTF and DNPNTF which execute in the tensor space could generate parts-based sparse representations. Since more constraints were adopted, DNPNTF generate sparser basis images which reflect distinct features for recognition.Basis images obtained on the JAFFE database.
(a)
NMF
(b)
NTF
(c)
DNPNTFThen we give the experiments to prove the effectiveness of DNPNTF via ELM. The average recognition rates of DNPNTF with ELM, NN, SVM, and SRC are given in Figure4, where the vertical axis represents the correct recognition rate in percentage and the horizontal axis represents the corresponding dimensions (from 1 to 120). ELM and SRC achieve better recognition performance compared with NN and SVM, and ELM achieves the highest recognition rate. The top recognition rates with the corresponding dimensions are given in Table 2.Table 2
The top recognition rate of different classifiers achieved on the JAFFE database.
ELM
NN
SRC
SVM
Recognition rates (%)
94.39
90.71
93.54
91.03
Dimensions
55
98
63
101Figure 4
The recognition rate versus dimension curves of DNPNTF by using different classifiers on the JAFFE database.


## 5.2. Experiments on Cohn-Kanade Database

The Cohn-Kanade database [32] consists of a large amount of image sequences starting from natural face and ending with the peak of the corresponding expression. 104 subjects with different ages, genders, and races are instructed to pose a series of 23 facial displays, including the 6 prototypic expressions. In our experiments, for every image sequence, we take 2 to 8 continuous frames near the peak expression as the static samples. We use the face images of all subjects. We partition the subject to 3 exclusive groups, and in each group, for each of the prototypic expression, we select 100 samples; that is, there are 600 samples in each group and the size of the total set is 1800. During the experiment, we adopt the leave-one-group-out strategy and 3-fold cross-validation: each time two groups are taken as training set and the remaining group is left for testing. This procedure is repeated for 3 times.The average recognition rates of different algorithms on the Cohn-Kanade database are shown in Figure5, where the vertical axis represents the correct recognition rate in percentage and the horizontal axis represents the corresponding dimensions (from 1 to 120). Here, only the NN classifier is used. Table 3 shows the top recognition rates with the corresponding dimensions. The recognition rates obtained on the Cohn-Kanade database are lower than those obtained on the JAFFE database. It can be explained that the experiments on the Cohn-Kanade database are person-independent, which are more difficult than the person-dependent experiments on the JAFFE database. From Figure 5, we can see that the performance of DNPNTF is superior to others with nearly all dimensions. Its recognition rates improve with the increase of dimensions.Table 3
The top recognition rate and the corresponding dimensions achieved on the Cohn-Kanade database.
DNPNTF
NMF
DNMF
NTF
Recognition rates (%)
63.61
58.17
58.83
59.72
Dimensions
101
76
26
120Figure 5
The recognition rate versus dimension curves achieved on the Cohn-Kanade database.Lastly, we give the experiments about different classifiers on the Cohn-Kanade database. The average recognition rates of DNPNTF via ELM, NN, SVM, and SRC are shown in Figure6, and the top recognition rates are given in Table 4. SVM and SRC achieve better performance compared with NN. ELM achieves the best recognition accuracy among all tested algorithms on almost all dimensions. It means ELM could use the information contained in extracted features better than other classifiers.Table 4
The top recognition rate and the corresponding dimensions achieved on the Cohn-Kanade database.
ELM
NN
SRC
SVM
Recognition rates (%)
86.00
63.61
73.00
72.06
Dimensions
51
101
116
106Figure 6
The recognition rate versus dimension curves of DNPNTF by using different classifiers on the Cohn-Kanade database.


## 6. Conclusions

In this paper, a novel DNPNTF algorithm with the application to facial expression recognition was proposed, which adopts ELM as the classifier. To incorporate the spatial information and the discriminant class information, a discriminant constraint is added to the objective function according to the manifold learning and graph embedding theory. To guarantee the convergence, the project gradient method is used for optimization. Theoretical analysis and experimental results demonstrate that DNPNTF could achieve better performance compared with NTF, NMF, and its variant. Then the discriminant features generated by DNPNTF are fed into ELM to learn an optimal model for recognition. In our experiments, DNPNTF via ELM achieves higher recognition rate compared with NN, SVM, and SRC.



---

*Source: 390328-2014-10-19.xml*