# Nonnegative Matrix Factorization-Based Spatial-Temporal Clustering for Multiple Sensor Data Streams

Authors: Di-Hua Sun; Chun-Yan Sang
Journal: Journal of Sensors
 (2014)

---

## Abstract

Cyber physical systems have grown exponentially and have been attracting a lot of attention over the last few years. To retrieve and mine the useful information from massive amounts of sensor data streams with spatial, temporal, and other multidimensional information has become an active research area. Moreover, recent research has shown that clusters of streams change with a comprehensive spatial-temporal viewpoint in real applications. In this paper, we propose a spatial-temporal clustering algorithm (STClu) based on nonnegative matrix trifactorization by utilizing time-series observational data streams and geospatial relationship for clustering multiple sensor data streams. Instead of directly clustering multiple data streams periodically, STClu incorporates the spatial relationship between two sensors in proximity and integrates the historical information into consideration. Furthermore, we develop an iterative updating optimization algorithm STClu. The effectiveness and efficiency of the algorithm STClu are both demonstrated in experiments on real and synthetic data sets. The results show that the proposed STClu algorithm outperforms existing methods for clustering sensor data streams.


---

## Body

## 1. Introduction

Cyber physical systems (CPS) have grown exponentially and have been attracting a lot of attention over the last few years [1]. Emerging applications include transportation, energy, manufacturing, environment monitoring, and many parts of our social infrastructure. It consists of a large number of sensors to monitor physical or environmental conditions [2, 3]. The incoming information in CPS needs to be considered as a stream and not so much as an ever-growing data set [4, 5].To obtain interesting relationships and useful information from multiple data streams, a variety of methods have developed over the past decades. Beringer and Hüllermeier [6] used a discrete Fourier transformation (DFT) to summarize the data streams and presented an online version of the classic K-means clustering algorithm to cluster parallel data streams. To avoid recalculating the DFT coefficients and minimize processing time, an incremental update mechanism for clustering large number of data streams was proposed in [7]. To discover cross relationships among multiple data streams, a clustering-on-demand framework was presented in [8]. Yeh et al. [9] proposed a framework for clustering over multiple evolving streams by correlations and events, which monitored the distribution of clusters over multiple data streams based on their correlation. To find communities in a large set of interacting entities, Aggarwal and Yu [10] proposed an online analytical processing framework for community detection of data streams. Evolutionary clustering was a relatively new topic and was formulated by Chakrabarti et al. [11] in 2006. To discover communities and capture their evolution with temporal smoothness, incremental and evolutionary clustering technologies were designed to handle dynamic data streams [10, 12–14]. Wang et al. [15] proposed an ECKF framework for evolutionary clustering large-scale data based on low-rank kernel matrix factorization. Chi et al. [16] proposed two algorithms to incorporate temporal smoothness in evolutionary spectral clustering. Ning et al. [14] presented an incremental approach for spectral clustering to handle not only insertion/deletion of data points but also similarity changes between existing points.The values of a sensor in the temporal domain usually have spatial correlation features, meaning that the sensors’ readings are influenced by near sensors. Such as in transportation systems, the congestion swiftly expands along the street and influences nearby sensors. A serious congestion usually lasts for a few hours and covers hundreds of sensors when reaching the full size. As time passes by, it shrinks slowly, eventually reduces the coverage, and finally disappears. As shown in Figure1, we find that those records in an atypical event are spatially close and timely relevant. To retrieve the useful information from massive sensor data streams, Tang et al. [5] proposed a model of atypical cluster to describe the atypical events and retrieve them with spatial and temporal attribute. To explore the spatial correlation between data sampled by different sensors for sensor networks, α-local spatial clustering algorithm was proposed [17]. To explore the temporal locality features of sensor data, Wei and Peng [18] proposed an incremental clustering algorithm that clusters objects by inducing clustering results from previous time windows.Figure 1
Time-series of upstream and downstream traffic flow.In this paper, we propose a spatial-temporal clustering algorithm (STClu) based on nonnegative matrix trifactorization by utilizing time-series observational data streams and geospatial relationship for clustering multiple sensor data streams. Our proposal was motivated by recent progress in nonnegative matrix factorization and its applications [19–25], with the aim being to discover clusters of objects with similar behavior and, consequently, discover potential anomalies and atypical events which may be revealed by the relationship evolving over time. We assume that the spatial feature is the summary of the atypical event in temporal dimension, and the temporal feature is the summary of the event in spatial dimension. STClu simultaneously incorporates the spatial correlation between two sensors in proximity and integrates the historical information into consideration. To maintain a consistent clustering result, there is a trade-off between the spatially close and timely relevant cost embedded in STClu and the benefit from current observations and the spatially close historical results. Thus, the clustering results at a given time step are determined by current snapshot, prior knowledge, and spatial correlation between two sensors in proximity.The rest of this paper is organized as follows. Related research work is introduced in Section2. Then, the preliminaries of STClu are given in Section 3. Section 4 presents the proposed algorithm in detail, while convergence and correctness proofs of the STClu are provided in the Appendix. Experimental results using synthetic and real world data sets are presented in Section 5. Finally, Section 6 gives the conclusion.


## 2. A Brief Review of Matrix Factorization and Its Applications

Graphs are used in a wide range of applications, such as transportation networks, sensor networks, and social networks. So, various problems are studied under graph mining. To find clusters in graphs is a new challenge if the graph is evolving over time. Matrix is a usual representation of a graph; the relationships between matrix factorization andK-means clustering have been explored in [26].Singular value decomposition (SVD) has served as a building block for many important applications, such as PCA and LSI [27]. SVD factorizes a matrix with the general form of A
≈
SV
D
T, where S is a unitary basis consisting of left-singular vectors of A, D is a unitary basis consisting of right-singular vectors of A, and V is a diagonal matrix with singular values on the diagonal. In SVD, matrices S and D are allowed to have negative entries; the projected data might have negative values in spite of the original data being positive. That can prevent from the clustering results to be intuitive for applications such as documents or images that have a positive data input.Nonnegative matrix factorization (NMF) [28, 29] is a linear, nonnegative approximate data representation technique. NMF focuses on the analysis of data matrices whose elements are nonnegative, a common occurrence in data sets derived from text and images. The nonnegative data matrix A is factorized into matrixes F and G as A
≈
F
G
T, with the constraints that F
∈
R
d
×
k and G
∈
R
n
×
k are nonnegative.In particularly, the result of aK-means clustering run can be written as a matrix factorization A
=
F
G
T, where A is a data matrix, F contains the cluster centroids, and G contains the cluster membership indicators. Nonnegative matrix factorization focuses on the analysis of data matrices whose elements are nonnegative [28].To extend the applicable range of NMF methods, when the data matrix is unconstrained, Semi-NMF is motivated from the perspective of clustering [29], in which it restricts G to be nonnegative while placing no restriction on the signs of F. For reasons of interpretability, Convex-NMF constrains the basis vectors F
=
(
f
1
,
…
,
f
k
). The vectors defining F lie within the column space of A: f
l
=
w
1
l
a
1
+
⋯
+
w
n
l
a
n
=
A
w
l, or F
=
A
W. It can be applied to both nonnegative and mixed-sign data matrices. This constraint could interpret the columns f
l as weighted sums of certain data points and these columns would capture a notion of centroids.To expand the range of application of NMF, Semi-NMF and Convex-NMF algorithm have been proposed in [29]. Semi-NMF offers a low-dimensional representation of data points which lends itself to a convenient clustering interpretation. Convex-NMF restricts the columns of F to be convex combinations of data points in X. The advantage of this constraint is that the columns as weighted sums of certain data points could interpret. Based on NMF, a clustering algorithm has been proposed in [29].However, the traditional NMF, Semi-NMF, and Convex-NMF are linear model and they may fail to discover the nonlinearities of data streams. In real world, the data streams have potential nonlinear structure. Kernel method is a powerful technique for extracting useful information form nonlinear correlations. To solve the problem of nonlinearities, the kernel method is to map the data nonlinearly into a kernel feature space. Then, the Convex-NMF method can be accomplished in the kernel feature space to process the nonlinear data.Indeed, several important applications can be modeled as large sparse graphs, such as transportation network analysis and social network analysis. Low-rank approximation for the matrix of a graph is essential in finding patterns and detecting anomalies. And it can extract correlations and remove noise from matrix structured data. This has led to proposal of methods such as CUR [30], CMD [31], and the family of Colibri [32].


## 3. Preliminaries

To monitor the traffic status along a freeway, sensors are deployed and utilized to collect readings, such as vehicle speed and volume of traffic. It is assumed that data records arrive synchronously, which means that all data streams will be updated simultaneously. The sensors are the points in the graph and an edge is formed between each pair of sensors; these data records can be used to construct graphs of transportation networks periodically.Given two sets of spatial objectsO
(
t
-
1
)
=
{
o
1
(
t
-
1
)
,
…
,
o
m
(
t
-
1
)
} and O
(
t
)
=
{
o
1
(
t
)
,
…
,
o
n
(
t
)
}, the spatial relationships are depicted as Figure 2. Let G
(
t
-
1
)
=
(
S
(
t
)
,
E
(
t
)
) be a graph associated with the sensor network of transportation system at timestamp t, where s
i
(
t
)
∈
S
(
t
) represents the set of sensors associated with graph G
(
t
) and e
i
j
(
t
)
∈
E
(
t
) represents the set of edges between s
i
(
t
) and s
j
(
t
). And similar definitions are depicted for G
(
t
-
1
)
=
(
S
(
t
-
1
)
,
E
(
t
-
1
)
) associated with the sensor network of transportation system at timestamp t
-
1.Figure 2
Spatial-temporal cluster: graphG
(
t
-
1
)
=
(
S
(
t
-
1
)
,
E
(
t
-
1
)
), graph G
(
t
)
=
(
S
(
t
)
,
E
(
t
)
) and the spatial relationship between timestamps t
-
1 and t.Without loss of generality, we use the adjacency matricesX
(
t
-
1
)
∈
R
m
×
m and X
(
t
)
∈
R
n
×
n to represent two graphs G
(
t
-
1
)
=
(
S
(
t
-
1
)
,
E
(
t
-
1
)
) and G
(
t
)
=
(
S
(
t
)
,
E
(
t
)
), respectively. X
(
i
)
(
i
,
j
) is the element at the ith row and jth column of the matrix X
(
i
); X
(
i
)
(
:
,
j
) is the jth column of X
(
i
). Every row or column in X
(
i
) corresponds to a node in S
(
i
). If there is an edge from node s
i
(
i
)
∈
S
(
i
) to node s
j
(
i
)
∈
S
(
i
) with similarity ρ
i
j
(
i
), we set the value to X
(
i
)
(
i
,
j
). Otherwise, we set it to zero. The similarity matrix is computed as follows:
(1)
ρ
i
j
(
t
)
=
exp
⁡
(
-
∥
s
i
(
t
)
-
s
j
(
t
)
∥
2
2
σ
2
)
,
where σ is a scaling parameter controlling how rapidly similarity ρ
i
j
(
t
) reduces with the distance between s
i
(
t
) and s
j
(
t
).The spatial relationships between two graphs are denoted as a set of edges connecting the nodes betweenG
(
t
) and G
(
t
-
1
). We use the adjacency matrix X
(
t
-
1
,
t
)
∈
R
n
×
m to represent the spatial relationships. Given two objects s
i
(
t
-
1
)
∈
S
(
t
-
1
) and s
j
(
t
)
∈
S
(
t
), the locations of s
i
(
t
-
1
) and s
j
(
t
) are denoted as L
·
s
i
(
t
-
1
) and L
·
s
j
(
t
), respectively. The spatial relationship between L
·
s
i
(
t
-
1
) and L
·
s
j
(
t
) is defined as follows:
(2)
D
(
L
·
s
i
(
t
)
,
L
·
s
j
(
t
-
1
)
)
=
(
(
L
·
s
i
1
(
t
)
-
L
·
s
j
1
(
t
-
1
)
)
2
+
(
L
·
s
i
2
(
t
)
-
L
·
s
j
2
(
t
-
1
)
)
2
+
⋯
+
(
L
·
s
i
n
(
t
)
-
L
·
s
j
n
(
t
-
1
)
)
2
)
1
/
2
.At timestampt, the objective of STClu is to create sets of partitions {
Δ
i
(
t
)
}
i
=
1
k. ∀
i
∈
{
1
,
…
,
k
}, Δ
i
(
t
) should meet the following conditions: ⋃
i
=
1
k
Δ
i
(
t
)
=
{
s
1
(
t
)
,
s
2
(
t
)
,
…
,
s
n
(
t
)
} and ⋂
i
=
1
k
Δ
i
(
t
)
=
⌀. The community Δ
i
(
t
) is a set of similar data streams such that Δ
i
(
t
)
=
{
s
1
(
t
)
,
…
,
s
|
Δ
i
(
t
)
|
(
t
)
}, where |
Δ
i
(
t
)
| is the total number of streams in Δ
i
(
t
).


## 4. Our Approach

### 4.1. Objective Function of STClu

By utilizing time-series observational data streams and geospatial relationship for clustering multiple sensor data streams, the proposed methodSTClu takes the spatial correlation between two sensors in proximity and the historical information into consideration. The cost of spatially close and timely relevant is embedded in the objective function. The cluster partitions are obtained by decomposing the adjacency matrix representation of the graph into a product of its latent factors. The square of the Euclidean distance and the generalized Kullback-Leibler divergence are common metrics for the approximation error [24]. In STClu, we seek to minimize the distance function ∥
A
-
B
∥
F
2 between two nonnegative matrices A and B, defined as
(3)
∥
A
-
B
∥
F
2
=
∑
i
j
(
A
i
j
-
B
i
j
)
2
,
where ∥
·
∥
F
2 is the Frobenius norm. This is lower bounded by zero and clearly vanishes if, and only if, A
=
B.The objective function can be written as(4)
J
STClu
=
∥
X
(
t
)
-
C
(
t
)
U
(
t
)
R
(
t
)
T
∥
F
2
+
∥
X
(
t
-
1
)
-
C
(
t
-
1
)
U
(
t
-
1
)
R
(
t
)
T
∥
F
2
+
∥
X
(
t
-
1
,
t
)
-
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
∥
F
2
,
where C
(
t
) and C
(
t
-
1
) are the representative column matrices of X
(
t
) and X
(
t
-
1
), respectively, U
(
t
) and U
(
t
-
1
) are the weight matrices, respectively, R
(
t
)
∈
R
n
×
k
1 and R
(
t
-
1
)
∈
R
m
×
k
2 are the corresponding cluster membership matrices for graphs G
(
t
) and G
(
t
-
1
), respectively, and U
(
t
-
1
,
t
) reflects the correspondence between the subgroups derived from two time steps.Specifically, assume cluster numbersk
1 and k
2 at time steps t and t
-
1, respectively. If k
2
<
k
1, then the extra clusters at time t must be computed using only the current data and we let U
(
t
-
1
)
=
[
U
(
t
-
1
)
,
0
c
×
(
k
2
+
1
:
k
2
)
]. Otherwise, k
2
>
k
1 and we simply remove the deleted cluster weight vectors in U
(
t
-
1
).


### 4.2. Optimization of STClu

As can be seen, the objective function in (4) is minimized with respect to C
(
t
), U
(
t
), R
(
t
), and U
(
t
-
1
,
t
), and it is unrealistic to expect an algorithm to find the global minimum. In the following, we introduce an alternating schema to optimize the objective function, which can achieve a local minimum.We optimize the objective function with respect to one variable while fixing the other variables. This iterative procedure repeats until convergence. Using the matrix propertiesTr
⁡
(
A
B
)
=
Tr
⁡
(
B
A
) and Tr
⁡
(
A
)
=
Tr
⁡
(
A
T
), the objective function J in (4) can be rewritten as
(5)
J
STClu
=
Tr
⁡
(
(
X
(
t
)
-
C
(
t
)
U
(
t
)
R
(
t
)
T
)
(
X
(
t
)
-
C
(
t
)
U
(
t
)
R
(
t
)
T
)
T
)
+
Tr
⁡
(
(
X
(
t
-
1
)
-
C
(
t
-
1
)
U
(
t
-
1
)
R
(
t
)
T
)
T
)

(
X
(
t
-
1
)
-
C
(
t
-
1
)
U
(
t
-
1
)
R
(
t
)
T
)
h
h
h
h
h
h
h
h
h
×
(
X
(
t
-
1
)
-
C
(
t
-
1
)
U
(
t
-
1
)
R
(
t
)
T
)
T
)
+
Tr
⁡
(
(
X
(
t
-
1
,
t
)
-
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
)
T
)

(
X
(
t
-
1
,
t
)
-
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
)
h
h
h
h
h
h
h
h
h
h
×
(
X
(
t
-
1
,
t
)
-
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
)
T
)
=
Tr
⁡
(
X
(
t
)
X
(
t
)
T
)
-
2
Tr
⁡
(
X
(
t
)
R
(
t
)
U
(
t
)
T
C
(
t
)
T
)
+
Tr
⁡
(
C
(
t
)
U
(
t
)
R
(
t
)
T
R
(
t
)
U
(
t
)
T
C
(
t
)
T
)
+
Tr
⁡
(
X
(
t
-
1
)
X
(
t
-
1
)
T
)
-
2
Tr
⁡
(
X
(
t
-
1
)
R
(
t
)
U
(
t
-
1
)
T
C
(
t
-
1
)
T
)
+
Tr
⁡
(
C
(
t
-
1
)
U
(
t
-
1
)
R
(
t
)
T
R
(
t
)
U
(
t
-
1
)
T
C
(
t
-
1
)
T
)
+
Tr
⁡
(
X
(
t
-
1
,
t
)
X
(
t
-
1
,
t
)
T
)
-
2
Tr
⁡
(
X
(
t
-
1
,
t
)
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
R
(
t
)
T
)
+
Tr
⁡
(
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
R
(
t
)
T
)
.To alternately update the entries ofC
(
t
), U
(
t
), R
(
t
), and U
(
t
-
1
,
t
), we resort to a Lagrangian function. Let ζ
i
j, ξ
j
l, ψ
k
j, and ϑ
i
j be the Lagrange multiplier for constraints C
i
j
(
t
)
≥
0, U
j
l
(
t
)
≥
0, R
k
j
(
t
)
≥
0, and U
i
j
(
t
-
1
,
t
)
≥
0, respectively. Then the Lagrange function L is expressed as
(6)
L
=
Tr
⁡
(
X
(
t
)
X
(
t
)
T
)
-
2
Tr
⁡
(
X
(
t
)
R
(
t
)
U
(
t
)
T
C
(
t
)
T
)
+
Tr
⁡
(
C
(
t
)
U
(
t
)
R
(
t
)
T
R
(
t
)
U
(
t
)
T
C
(
t
)
T
)
+
Tr
⁡
(
X
(
t
-
1
)
X
(
t
-
1
)
T
)
-
2
Tr
⁡
(
X
(
t
-
1
)
R
(
t
)
U
(
t
-
1
)
T
C
(
t
-
1
)
T
)
+
Tr
⁡
(
C
(
t
-
1
)
U
(
t
-
1
)
R
(
t
)
T
R
(
t
)
U
(
t
-
1
)
T
C
(
t
-
1
)
T
)
+
Tr
⁡
(
X
(
t
-
1
,
t
)
X
(
t
-
1
,
t
)
T
)
-
2
Tr
⁡
(
X
(
t
-
1
,
t
)
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
R
(
t
)
T
)
+
Tr
⁡
(
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
R
(
t
)
T
)
+
Tr
⁡
(
ζ
C
(
t
)
T
)
+
Tr
⁡
(
ξ
U
(
t
)
T
)
+
Tr
⁡
(
ψ
R
(
t
)
T
)
+
Tr
⁡
(
ϑ
U
(
t
-
1
,
t
)
T
)
.Since the derivatives ofL with respect to U
(
t
), C
(
t
), R
(
t
), and U
(
t
-
1
,
t
) must be zero, we have
(7)
∂
L
∂
U
(
t
)
=
-
2
C
(
t
)
T
X
(
t
)
R
(
t
)
+
2
C
(
t
)
T
C
(
t
)
U
(
t
)
R
(
t
)
T
R
(
t
)
+
ξ
,
∂
L
∂
C
(
t
)
=
-
2
X
(
t
)
R
(
t
)
U
(
t
)
T
+
2
C
(
t
)
U
(
t
)
R
(
t
)
T
R
(
t
)
U
(
t
)
T
+
ζ
,
∂
L
∂
R
(
t
)
=
-
2
X
(
t
)
T
C
(
t
)
U
(
t
)
+
2
R
(
t
)
U
(
t
)
T
C
(
t
)
T
C
(
t
)
U
(
t
)
-
2
X
(
t
-
1
)
T
C
(
t
-
1
)
U
(
t
-
1
)
+
2
R
(
t
)
U
(
t
-
1
)
T
C
(
t
-
1
)
T
C
(
t
-
1
)
U
(
t
-
1
)
-
2
X
(
t
-
1
,
t
)
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
+
2
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
+
ψ
,
∂
L
∂
U
(
t
-
1
,
t
)
=
-
2
R
(
t
)
T
X
(
t
-
1
,
t
)
R
(
t
-
1
)
+
2
R
(
t
)
T
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
R
(
t
-
1
)
+
ϑ
.Using the Karush-Kuhn-Tucker conditions [29, 33] and letting ξ
j
l
U
j
l
(
t
)
=
0, ζ
i
j
C
i
j
(
t
)
=
0, ψ
k
j
R
k
j
(
t
)
=
0, and ϑ
i
j
U
i
j
(
t
-
1
,
t
)
=
0, we obtain the following equations for U
(
t
), C
(
t
), R
(
t
), and U
(
t
-
1
,
t
):
(8)
[
-
C
(
t
)
T
X
(
t
)
R
(
t
)
+
C
(
t
)
T
C
(
t
)
U
(
t
)
R
(
t
)
T
R
(
t
)
]
j
l
U
j
l
(
t
)
=
0
,
[
-
X
(
t
)
R
(
t
)
U
(
t
)
T
+
C
(
t
)
U
(
t
)
R
(
t
)
T
R
(
t
)
U
(
t
)
T
]
i
j
C
i
j
(
t
)
=
0
,
[
-
X
(
t
)
T
C
(
t
)
U
(
t
)
+
R
(
t
)
U
(
t
)
T
C
(
t
)
T
C
(
t
)
U
(
t
)
-
X
(
t
-
1
)
T
C
(
t
-
1
)
U
(
t
-
1
)
+
R
(
t
)
U
(
t
-
1
)
T
C
(
t
-
1
)
T
C
(
t
-
1
)
U
(
t
-
1
)
-
X
(
t
-
1
,
t
)
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
+
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
]
k
j
R
k
j
(
t
)
=
0
,
[
-
R
(
t
)
T
X
(
t
-
1
,
t
)
R
(
t
-
1
)
+
R
(
t
)
T
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
R
(
t
-
1
)
]
i
j
U
i
j
(
t
-
1
,
t
)
=
0
.According to (8), we present the following updating rules:
(9)
U
j
l
(
t
)
⟵
U
j
l
(
t
)
[
C
(
t
)
T
X
(
t
)
R
(
t
)
]
j
l
[
C
(
t
)
T
C
(
t
)
U
(
t
)
R
(
t
)
T
R
(
t
)
]
j
l
,
(10)
C
i
j
(
t
)
⟵
C
i
j
(
t
)
[
X
(
t
)
R
(
t
)
U
(
t
)
T
]
i
j
[
C
(
t
)
U
(
t
)
R
(
t
)
T
R
(
t
)
U
(
t
)
T
]
i
j
,
(11)
R
k
j
(
t
)
⟵
R
k
j
(
t
)
(
+
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
]
k
j
)
-
1

(
+
X
(
t
-
1
,
t
)
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
]
k
j
[
X
(
t
)
T
C
(
t
)
U
(
t
)
+
X
(
t
-
1
)
T
C
(
t
-
1
)
U
(
t
-
1
)
h
h
h
h
h
h
h
h
h
h
h
h
h
h
+
X
(
t
-
1
,
t
)
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
]
k
j
)
h
h
h
h
h
h
h
h
h
h
h
×
(
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
]
k
j

[
R
(
t
)
U
(
t
)
T
C
(
t
)
T
C
(
t
)
U
(
t
)
h
h
h
h
h
h
h
h
h
h
h
h
h
h
h
+
R
(
t
)
U
(
t
-
1
)
T
C
(
t
-
1
)
T
C
(
t
-
1
)
U
(
t
-
1
)
h
h
h
h
h
h
h
h
h
h
h
h
h
h
h
+
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
h
h
h
h
h
h
h
h
h
h
h
h
h
h
h
×
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
]
k
j
)
-
1
)
,
(12)
U
i
j
(
t
-
1
,
t
)
⟵
U
i
j
(
t
-
1
,
t
)
[
R
(
t
)
T
X
(
t
-
1
,
t
)
R
(
t
-
1
)
]
i
j
[
R
(
t
)
T
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
R
(
t
-
1
)
]
i
j
.Proofs of the convergence of update rules (9), (10), (11), and (12) are given in the Appendix.


### 4.3. Algorithm Description

To eliminate linearly dependent columns and construct the subspace used for low-rank approximation from data matrixX, a unique and independent subspace C is constructed using Colibri-S and Colibri-D [32]. Algorithm 1 summarizes the details of the iterative updating algorithm for (4). It is assumed that the adjacency matrices  X
(
t
), X
(
t
-
1
), X
(
t
-
1
,
t
), R
(
t
-
1
), and U
(
t
-
1
) are given. To improve the clustering efficiency, the variable factors are initialized by the previous clustering results instead of random values [29]. Applying the update rules in (9), (10), (11), and (12) until convergence, we obtain a new C
(
t
), U
(
t
), R
(
t
), and U
(
t
-
t
,
t
).Algorithm 1:STClu.
Input: MatricesX
(
t
), X
(
t
-
1
), X
(
t
-
1
,
t
), R
(
t
-
1
), and U
(
t
-
1
), cluster number k
(
t
), and maximum iteration T
Output: MatricesC
(
t
), U
(
t
), R
(
t
), and U
(
t
-
1
,
t
)
(1)   use the family of Colibri methods to getC
(
t
)
(2)  determine cluster numberk
(
t
)
(3)  initializeU
(
t
), R
(
t
), and U
(
t
-
1
,
t
)  using R
(
t
-
1
), U
(
t
-
1
), and U
(
t
-
1
,
t
), respectively
(4)  if we need to form a new partition
(5)    go to (8)
(6)  else
(7)    letR
(
t
)
=
R
(
t
-
1
) and U
(
t
)
=
U
(
t
-
1
), and return
(8)  while not converging andt
≤
T do
(9)    updateU
(
t
) by (9)
(10)  update C
(
t
) by (10)
(11)   updateR
(
t
) by (11)
(12)update U
(
t
-
1
,
t
) by (12)
(13) end while
(14) returnU
(
t
), C
(
t
), R
(
t
), and U
(
t
-
1
,
t
)The computational cost of the proposedSTClu algorithm is discussed as follows. To compute the complexity of the STClu algorithm, we need to inspect its main operations, namely, the decomposition matrix and multiplicative updates. In the decomposition step, the computation is actually done using the three small matrices, U, C, and R. Each of  U, C, and R is multiplied by an update factor, which involves several matrix multiplications. We assume that k is the cluster number at time step t, c
1 and c
0 are the sizes of the data subspaces at time steps t and t
-
1, respectively, and n
1 and n
0 are the sizes of the input matrices at time steps t and t
-
1, respectively; the multiplicative update stops after T iterations, and the cost for STClu is O
(
T
(
c
2
k
+
c
n
k
+
c
n
2
+
c
2
n
)
), where n
=
max
⁡
(
n
0
,
n
1
) and c
=
max
⁡
(
c
0
,
c
1
).


## 4.1. Objective Function of STClu

By utilizing time-series observational data streams and geospatial relationship for clustering multiple sensor data streams, the proposed methodSTClu takes the spatial correlation between two sensors in proximity and the historical information into consideration. The cost of spatially close and timely relevant is embedded in the objective function. The cluster partitions are obtained by decomposing the adjacency matrix representation of the graph into a product of its latent factors. The square of the Euclidean distance and the generalized Kullback-Leibler divergence are common metrics for the approximation error [24]. In STClu, we seek to minimize the distance function ∥
A
-
B
∥
F
2 between two nonnegative matrices A and B, defined as
(3)
∥
A
-
B
∥
F
2
=
∑
i
j
(
A
i
j
-
B
i
j
)
2
,
where ∥
·
∥
F
2 is the Frobenius norm. This is lower bounded by zero and clearly vanishes if, and only if, A
=
B.The objective function can be written as(4)
J
STClu
=
∥
X
(
t
)
-
C
(
t
)
U
(
t
)
R
(
t
)
T
∥
F
2
+
∥
X
(
t
-
1
)
-
C
(
t
-
1
)
U
(
t
-
1
)
R
(
t
)
T
∥
F
2
+
∥
X
(
t
-
1
,
t
)
-
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
∥
F
2
,
where C
(
t
) and C
(
t
-
1
) are the representative column matrices of X
(
t
) and X
(
t
-
1
), respectively, U
(
t
) and U
(
t
-
1
) are the weight matrices, respectively, R
(
t
)
∈
R
n
×
k
1 and R
(
t
-
1
)
∈
R
m
×
k
2 are the corresponding cluster membership matrices for graphs G
(
t
) and G
(
t
-
1
), respectively, and U
(
t
-
1
,
t
) reflects the correspondence between the subgroups derived from two time steps.Specifically, assume cluster numbersk
1 and k
2 at time steps t and t
-
1, respectively. If k
2
<
k
1, then the extra clusters at time t must be computed using only the current data and we let U
(
t
-
1
)
=
[
U
(
t
-
1
)
,
0
c
×
(
k
2
+
1
:
k
2
)
]. Otherwise, k
2
>
k
1 and we simply remove the deleted cluster weight vectors in U
(
t
-
1
).


## 4.2. Optimization of STClu

As can be seen, the objective function in (4) is minimized with respect to C
(
t
), U
(
t
), R
(
t
), and U
(
t
-
1
,
t
), and it is unrealistic to expect an algorithm to find the global minimum. In the following, we introduce an alternating schema to optimize the objective function, which can achieve a local minimum.We optimize the objective function with respect to one variable while fixing the other variables. This iterative procedure repeats until convergence. Using the matrix propertiesTr
⁡
(
A
B
)
=
Tr
⁡
(
B
A
) and Tr
⁡
(
A
)
=
Tr
⁡
(
A
T
), the objective function J in (4) can be rewritten as
(5)
J
STClu
=
Tr
⁡
(
(
X
(
t
)
-
C
(
t
)
U
(
t
)
R
(
t
)
T
)
(
X
(
t
)
-
C
(
t
)
U
(
t
)
R
(
t
)
T
)
T
)
+
Tr
⁡
(
(
X
(
t
-
1
)
-
C
(
t
-
1
)
U
(
t
-
1
)
R
(
t
)
T
)
T
)

(
X
(
t
-
1
)
-
C
(
t
-
1
)
U
(
t
-
1
)
R
(
t
)
T
)
h
h
h
h
h
h
h
h
h
×
(
X
(
t
-
1
)
-
C
(
t
-
1
)
U
(
t
-
1
)
R
(
t
)
T
)
T
)
+
Tr
⁡
(
(
X
(
t
-
1
,
t
)
-
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
)
T
)

(
X
(
t
-
1
,
t
)
-
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
)
h
h
h
h
h
h
h
h
h
h
×
(
X
(
t
-
1
,
t
)
-
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
)
T
)
=
Tr
⁡
(
X
(
t
)
X
(
t
)
T
)
-
2
Tr
⁡
(
X
(
t
)
R
(
t
)
U
(
t
)
T
C
(
t
)
T
)
+
Tr
⁡
(
C
(
t
)
U
(
t
)
R
(
t
)
T
R
(
t
)
U
(
t
)
T
C
(
t
)
T
)
+
Tr
⁡
(
X
(
t
-
1
)
X
(
t
-
1
)
T
)
-
2
Tr
⁡
(
X
(
t
-
1
)
R
(
t
)
U
(
t
-
1
)
T
C
(
t
-
1
)
T
)
+
Tr
⁡
(
C
(
t
-
1
)
U
(
t
-
1
)
R
(
t
)
T
R
(
t
)
U
(
t
-
1
)
T
C
(
t
-
1
)
T
)
+
Tr
⁡
(
X
(
t
-
1
,
t
)
X
(
t
-
1
,
t
)
T
)
-
2
Tr
⁡
(
X
(
t
-
1
,
t
)
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
R
(
t
)
T
)
+
Tr
⁡
(
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
R
(
t
)
T
)
.To alternately update the entries ofC
(
t
), U
(
t
), R
(
t
), and U
(
t
-
1
,
t
), we resort to a Lagrangian function. Let ζ
i
j, ξ
j
l, ψ
k
j, and ϑ
i
j be the Lagrange multiplier for constraints C
i
j
(
t
)
≥
0, U
j
l
(
t
)
≥
0, R
k
j
(
t
)
≥
0, and U
i
j
(
t
-
1
,
t
)
≥
0, respectively. Then the Lagrange function L is expressed as
(6)
L
=
Tr
⁡
(
X
(
t
)
X
(
t
)
T
)
-
2
Tr
⁡
(
X
(
t
)
R
(
t
)
U
(
t
)
T
C
(
t
)
T
)
+
Tr
⁡
(
C
(
t
)
U
(
t
)
R
(
t
)
T
R
(
t
)
U
(
t
)
T
C
(
t
)
T
)
+
Tr
⁡
(
X
(
t
-
1
)
X
(
t
-
1
)
T
)
-
2
Tr
⁡
(
X
(
t
-
1
)
R
(
t
)
U
(
t
-
1
)
T
C
(
t
-
1
)
T
)
+
Tr
⁡
(
C
(
t
-
1
)
U
(
t
-
1
)
R
(
t
)
T
R
(
t
)
U
(
t
-
1
)
T
C
(
t
-
1
)
T
)
+
Tr
⁡
(
X
(
t
-
1
,
t
)
X
(
t
-
1
,
t
)
T
)
-
2
Tr
⁡
(
X
(
t
-
1
,
t
)
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
R
(
t
)
T
)
+
Tr
⁡
(
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
R
(
t
)
T
)
+
Tr
⁡
(
ζ
C
(
t
)
T
)
+
Tr
⁡
(
ξ
U
(
t
)
T
)
+
Tr
⁡
(
ψ
R
(
t
)
T
)
+
Tr
⁡
(
ϑ
U
(
t
-
1
,
t
)
T
)
.Since the derivatives ofL with respect to U
(
t
), C
(
t
), R
(
t
), and U
(
t
-
1
,
t
) must be zero, we have
(7)
∂
L
∂
U
(
t
)
=
-
2
C
(
t
)
T
X
(
t
)
R
(
t
)
+
2
C
(
t
)
T
C
(
t
)
U
(
t
)
R
(
t
)
T
R
(
t
)
+
ξ
,
∂
L
∂
C
(
t
)
=
-
2
X
(
t
)
R
(
t
)
U
(
t
)
T
+
2
C
(
t
)
U
(
t
)
R
(
t
)
T
R
(
t
)
U
(
t
)
T
+
ζ
,
∂
L
∂
R
(
t
)
=
-
2
X
(
t
)
T
C
(
t
)
U
(
t
)
+
2
R
(
t
)
U
(
t
)
T
C
(
t
)
T
C
(
t
)
U
(
t
)
-
2
X
(
t
-
1
)
T
C
(
t
-
1
)
U
(
t
-
1
)
+
2
R
(
t
)
U
(
t
-
1
)
T
C
(
t
-
1
)
T
C
(
t
-
1
)
U
(
t
-
1
)
-
2
X
(
t
-
1
,
t
)
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
+
2
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
+
ψ
,
∂
L
∂
U
(
t
-
1
,
t
)
=
-
2
R
(
t
)
T
X
(
t
-
1
,
t
)
R
(
t
-
1
)
+
2
R
(
t
)
T
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
R
(
t
-
1
)
+
ϑ
.Using the Karush-Kuhn-Tucker conditions [29, 33] and letting ξ
j
l
U
j
l
(
t
)
=
0, ζ
i
j
C
i
j
(
t
)
=
0, ψ
k
j
R
k
j
(
t
)
=
0, and ϑ
i
j
U
i
j
(
t
-
1
,
t
)
=
0, we obtain the following equations for U
(
t
), C
(
t
), R
(
t
), and U
(
t
-
1
,
t
):
(8)
[
-
C
(
t
)
T
X
(
t
)
R
(
t
)
+
C
(
t
)
T
C
(
t
)
U
(
t
)
R
(
t
)
T
R
(
t
)
]
j
l
U
j
l
(
t
)
=
0
,
[
-
X
(
t
)
R
(
t
)
U
(
t
)
T
+
C
(
t
)
U
(
t
)
R
(
t
)
T
R
(
t
)
U
(
t
)
T
]
i
j
C
i
j
(
t
)
=
0
,
[
-
X
(
t
)
T
C
(
t
)
U
(
t
)
+
R
(
t
)
U
(
t
)
T
C
(
t
)
T
C
(
t
)
U
(
t
)
-
X
(
t
-
1
)
T
C
(
t
-
1
)
U
(
t
-
1
)
+
R
(
t
)
U
(
t
-
1
)
T
C
(
t
-
1
)
T
C
(
t
-
1
)
U
(
t
-
1
)
-
X
(
t
-
1
,
t
)
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
+
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
]
k
j
R
k
j
(
t
)
=
0
,
[
-
R
(
t
)
T
X
(
t
-
1
,
t
)
R
(
t
-
1
)
+
R
(
t
)
T
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
R
(
t
-
1
)
]
i
j
U
i
j
(
t
-
1
,
t
)
=
0
.According to (8), we present the following updating rules:
(9)
U
j
l
(
t
)
⟵
U
j
l
(
t
)
[
C
(
t
)
T
X
(
t
)
R
(
t
)
]
j
l
[
C
(
t
)
T
C
(
t
)
U
(
t
)
R
(
t
)
T
R
(
t
)
]
j
l
,
(10)
C
i
j
(
t
)
⟵
C
i
j
(
t
)
[
X
(
t
)
R
(
t
)
U
(
t
)
T
]
i
j
[
C
(
t
)
U
(
t
)
R
(
t
)
T
R
(
t
)
U
(
t
)
T
]
i
j
,
(11)
R
k
j
(
t
)
⟵
R
k
j
(
t
)
(
+
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
]
k
j
)
-
1

(
+
X
(
t
-
1
,
t
)
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
]
k
j
[
X
(
t
)
T
C
(
t
)
U
(
t
)
+
X
(
t
-
1
)
T
C
(
t
-
1
)
U
(
t
-
1
)
h
h
h
h
h
h
h
h
h
h
h
h
h
h
+
X
(
t
-
1
,
t
)
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
]
k
j
)
h
h
h
h
h
h
h
h
h
h
h
×
(
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
]
k
j

[
R
(
t
)
U
(
t
)
T
C
(
t
)
T
C
(
t
)
U
(
t
)
h
h
h
h
h
h
h
h
h
h
h
h
h
h
h
+
R
(
t
)
U
(
t
-
1
)
T
C
(
t
-
1
)
T
C
(
t
-
1
)
U
(
t
-
1
)
h
h
h
h
h
h
h
h
h
h
h
h
h
h
h
+
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
h
h
h
h
h
h
h
h
h
h
h
h
h
h
h
×
R
(
t
-
1
)
U
(
t
-
1
,
t
)
T
]
k
j
)
-
1
)
,
(12)
U
i
j
(
t
-
1
,
t
)
⟵
U
i
j
(
t
-
1
,
t
)
[
R
(
t
)
T
X
(
t
-
1
,
t
)
R
(
t
-
1
)
]
i
j
[
R
(
t
)
T
R
(
t
)
U
(
t
-
1
,
t
)
R
(
t
-
1
)
T
R
(
t
-
1
)
]
i
j
.Proofs of the convergence of update rules (9), (10), (11), and (12) are given in the Appendix.


## 4.3. Algorithm Description

To eliminate linearly dependent columns and construct the subspace used for low-rank approximation from data matrixX, a unique and independent subspace C is constructed using Colibri-S and Colibri-D [32]. Algorithm 1 summarizes the details of the iterative updating algorithm for (4). It is assumed that the adjacency matrices  X
(
t
), X
(
t
-
1
), X
(
t
-
1
,
t
), R
(
t
-
1
), and U
(
t
-
1
) are given. To improve the clustering efficiency, the variable factors are initialized by the previous clustering results instead of random values [29]. Applying the update rules in (9), (10), (11), and (12) until convergence, we obtain a new C
(
t
), U
(
t
), R
(
t
), and U
(
t
-
t
,
t
).Algorithm 1:STClu.
Input: MatricesX
(
t
), X
(
t
-
1
), X
(
t
-
1
,
t
), R
(
t
-
1
), and U
(
t
-
1
), cluster number k
(
t
), and maximum iteration T
Output: MatricesC
(
t
), U
(
t
), R
(
t
), and U
(
t
-
1
,
t
)
(1)   use the family of Colibri methods to getC
(
t
)
(2)  determine cluster numberk
(
t
)
(3)  initializeU
(
t
), R
(
t
), and U
(
t
-
1
,
t
)  using R
(
t
-
1
), U
(
t
-
1
), and U
(
t
-
1
,
t
), respectively
(4)  if we need to form a new partition
(5)    go to (8)
(6)  else
(7)    letR
(
t
)
=
R
(
t
-
1
) and U
(
t
)
=
U
(
t
-
1
), and return
(8)  while not converging andt
≤
T do
(9)    updateU
(
t
) by (9)
(10)  update C
(
t
) by (10)
(11)   updateR
(
t
) by (11)
(12)update U
(
t
-
1
,
t
) by (12)
(13) end while
(14) returnU
(
t
), C
(
t
), R
(
t
), and U
(
t
-
1
,
t
)The computational cost of the proposedSTClu algorithm is discussed as follows. To compute the complexity of the STClu algorithm, we need to inspect its main operations, namely, the decomposition matrix and multiplicative updates. In the decomposition step, the computation is actually done using the three small matrices, U, C, and R. Each of  U, C, and R is multiplied by an update factor, which involves several matrix multiplications. We assume that k is the cluster number at time step t, c
1 and c
0 are the sizes of the data subspaces at time steps t and t
-
1, respectively, and n
1 and n
0 are the sizes of the input matrices at time steps t and t
-
1, respectively; the multiplicative update stops after T iterations, and the cost for STClu is O
(
T
(
c
2
k
+
c
n
k
+
c
n
2
+
c
2
n
)
), where n
=
max
⁡
(
n
0
,
n
1
) and c
=
max
⁡
(
c
0
,
c
1
).


## 5. Experiments and Results

In this section, we use several synthetic and real world data sets to evaluate the effectiveness and efficiency of theSTClu algorithm.


### 5.1. Baseline Algorithms and Evaluation Methods

To demonstrate the efficiency ofSTClu, we compare it with the following three popular clustering algorithms. The first is K-means clustering in the principle component analysis (PCA) subspace. The other two algorithms are the typical spectral clustering algorithm normalized cut (Ncut) [34] and NMF-based clustering [19].To evaluate clustering quality, all our comparisons are based on clustering accuracy (ACC) and normalized mutual information (NMI) measurements.Clustering accuracy (ACC) discovers the one-to-one relationships between clusters and classes and measures the extent to which each cluster contains data points from the corresponding class. It is defined as [18, 19]
(13)
ACC
=
∑
i
=
1
n
δ
(
map
(
r
i
)
,
l
i
)
n
,
where r
i denotes the cluster label, l
i denotes the true class label, n is the total number of documents, δ
(
x
,
y
) is the delta function that is equal to one if x
=
y and zero otherwise, and map
(
r
i
) is the permutation mapping function that maps each cluster label r
i to the equivalent label from the data set.Between two random variablesi (category label) and j (cluster label), NMI is defined as [35]
(14)
NMI
=
∑
i
=
1
k
(
a
)
∑
j
=
1
k
(
b
)
n
i
,
j
log
⁡
(
(
n
·
n
i
,
j
)
/
(
n
i
·
n
j
)
)
(
∑
i
=
1
k
(
a
)
n
i
log
⁡
(
n
i
/
n
)
)
(
∑
j
=
1
k
(
b
)
n
j
log
⁡
(
n
j
/
n
)
)
,
where n is the number of documents, n
i and n
j denote the number of documents in category i and cluster j, respectively, n
i
,
j denotes the number of documents in category i as well as in clusterj, and k
(
a
) and k
(
b
) are the cluster numbers in the true and predicted clusters. The NMI score is 1 if the clustering results perfectly match the category labels whereas the score is 0 if the data are randomly partitioned. The higher the NMI score is, the better the clustering quality is.Since there are no predefined categories in our data, we have to design an alternative way to carry out the evaluation. In this paper, we use theK-means clustering result of the combined feature data as the target labels for feature clustering. Similarly, the K-means clustering result of the combined user data is used as the target labels for user clustering. To randomize the experiments, we conduct the evaluation using different cluster numbers. Note that all the experiments ran on the same machine. At each time step, for each given cluster number, 30 test runs are carried out on different randomly chosen clusters. For each experiment, we run it 10 times and report the average.


### 5.2. Synthetic Data Sets

In this section, we use synthetic datasets to demonstrate thatSTClu can significantly cluster data for multiple evolving data streams over time. The number of nodes and number of clusters in each time step are given. The synthetic datasets are generated by two steps [18]. The first step generates the objects’ geographic information, and the second step generates the values of the objects’ nongeographic attributes. The objects’ nongeographic attributes are generated by applying prototype systems used in [6]. It can be defined as follows: f
(
t
+
Δ
t
)
=
f
(
t
)
+
f
′
(
t
+
Δ
t
), f
′
(
t
+
Δ
t
)
=
f
′
(
t
)
+
u
(
t
), in which f
(
·
) is a stochastic process, t
=
0, and Δ
t
,
2
Δ
t
,
… are independent random variables, uniformly distributed in an interval [
-
a
,
a
]. So, the data stream S
(
·
) can be defined by S
(
t
)
=
f
(
t
+
h
(
t
)
)
+
g
(
t
), where h
(
·
) and g
(
·
) are stochastic processes that are generated in the same way as the prototype f
(
·
). The constant a determines the smoothness of a data stream, which can be different for p
(
·
), h
(
·
), and g
(
·
), such as 0.05, 0.1, and 0.2, respectively. To examine the quality of STClu for clustering multiple data streams, we generate random datasets with the number of streams ranging from 100 to 1000. Each stream has a length of 2,000 points.We show the performance of four methods while the time window evolving over time and the number of clusters in the data varies from 2 to 20. The stream number is fixed at 400, and each of them contains 1,000 points. As can be seen in Table1, the average ACC and NMI among four methods with the different numbers on synthetic dataset are shown. It can be clearly found that the STClu outperforms all other approaches and having highest NMIs at evolving time steps. The reason is that it incorporates the previous clustering results of the upstream sensors and the spatial relationship between the upstream and downstream sensors.Table 1
Comparison of ACC and NMI results for the four methods based on synthetic data sets.
k
ACC
NMI
K-means
Ncut
NMF
STClu
K-means
Ncut
NMF
STClu
2
0.3959
0.4954
0.5449
0.6542
0.3480
0.5130
0.5895
0.6749
4
0.3855
0.4740
0.5120
0.6436
0.4148
0.5233
0.4416
0.6967
5
0.3773
0.4685
0.5205
0.6757
0.3810
0.5016
0.5781
0.6122
6
0.3710
0.4866
0.5425
0.5603
0.4162
0.5430
0.6041
0.6059
8
0.3864
0.4636
0.4893
0.6038
0.4257
0.5247
0.5098
0.7022
10
0.3711
0.4540
0.4919
0.6231
0.3960
0.5154
0.4082
0.7238
12
0.3896
0.4802
0.5231
0.5868
0.4051
0.5088
0.5916
0.5421
15
0.3819
0.4728
0.5102
0.4834
0.4180
0.5207
0.5252
0.6281
18
0.4121
0.4827
0.5082
0.5293
0.4315
0.5786
0.5850
0.5823
20
0.4019
0.5275
0.6110
0.6516
0.4171
0.5823
0.4390
0.6299
Average
0.3873
0.4805
0.5254
0.6012
0.4053
0.5311
0.5272
0.6398The average ACC and NMI value for the proposed algorithm varies with the varied number of clusters or data streams shown in Figures3(a) and 3(b), respectively. As shown in Figure 3, we can observe that the performance of the proposed algorithm STClu is very robust to the number of clusters and data streams.The performance of the proposed algorithm varying with the size of (a) clusters, (b) data streams, respectively.
(a)
(b)Finally, we evaluate the average processing time for a round of clustering multiple data streams. There are different factors affecting the execution time. In this experiment, we evaluate the effect of the window size and the number of data streams on the response time for the compared algorithms.The following set of experiments evaluate the effect of window size on the execution time of these algorithms. Figure4(a) shows that as the window size increases, the processing time of these algorithms increases. The y-axis shows the execution time, while the x-axis represents the window size w
(
t
) from 60 to 600. In Figure 4(b), the y-axis shows the execution time, while the x-axis represents the number of data streams being varied from 100 to 1000. The other parameters are fixed as follows: w
(
t
)
=
300 s, k
=
4. Note that as the number of streams increases, the execution time of these methods increases. However, the time of the STClu and NMF increases much more slowly than others. This is because of the execution efficiency of STClu and NMF for clustering multiple data streams in the approximate subspace. Although STClu requires more time than NMF, as can be seen from the experiments, it is able to process extremely large sensor data streams with spatial, temporal, and other multidimensional information.Average processing time: (a) while the window sizew
(
t
) varies from 60 to 600, (b) while the number of streams varies from 100 to 1000, respectively.
(a)
(b)


### 5.3. Real World Data Sets

The real world data sets were obtained from the PeMS (http://pems.dot.ca.gov/) database. With sensor devices installed on road networks, a system monitors traffic flows on major US highways 24 hours a day, 7 days a week, thus acquiring huge volumes of data. We collected data from 963 sensors positioned along the freeway. We selected 15 months of daily data records giving the occupancy rates (between 0 and 1) of different car lanes on freeways in the San Francisco bay area over time. The data set covered records sampled every 10 min from 1 January 2008 to 30 March 2009. Data for public holidays were removed from the data set, as well as data collected on two days with anomalies (8 and 9 March 2008), where all sensors were muted between 2:00 and 3:00 am.First, we discuss the experiments on the PeMS data set with the number of clusters varying from 2 to 20. For each given cluster numberk, the clustering results of all these algorithms on the PeMS data set are shown in Table 2; Ncut, NMF, and the proposed STClu algorithm achieve better performance than the K-means algorithm, since these three algorithms consider the geometrical structure information contained in the data.Table 2
Comparison of the ACC and NMI results for the four methods on PeMS with the number of clusters varying from 2 to 20.
k
ACC
NMI
K-means
Ncut
NMF
STClu
K-means
Ncut
NMF
STClu
2
0.2716
0.3014
0.3163
0.5060
0.2808
0.3027
0.3260
0.5127
3
0.2685
0.2950
0.3064
0.5065
0.2762
0.2947
0.3197
0.5167
5
0.2660
0.2933
0.3089
0.5107
0.2741
0.2952
0.3256
0.5540
8
0.2641
0.2988
0.3155
0.5123
0.2759
0.3011
0.3119
0.5679
10
0.2687
0.2919
0.2996
0.5078
0.2748
0.2899
0.3105
0.5282
12
0.2642
0.2890
0.3004
0.5115
0.2711
0.2889
0.3137
0.5611
15
0.2697
0.2969
0.3097
0.5057
0.2777
0.2973
0.3129
0.5098
19
0.2674
0.2946
0.3059
0.5104
0.2754
0.2943
0.2959
0.5510
20
0.2764
0.2976
0.3052
0.5164
0.2813
0.2955
0.3023
0.6043
Average
0.2685
0.2954
0.3075
0.5097
0.2764
0.2955
0.3132
0.5451Next, we discuss the experiments on the PeMS data set with different time steps. For each time step, the clustering results of all these algorithms on the PeMS data set are shown in Table3. We observe that the STClu algorithm performs better than the other three algorithms for most time steps, except t
=
1. The reason for this is that STClu integrates more historical knowledge and the spatial relationship between two sensors in proximity.Table 3
Comparison of the ACC and NMI results for the four methods on PeMS with different time steps.
Time step
ACC
NMI
K-means
Ncut
NMF
STClu
K-means
Ncut
NMF
STClu
1
0.2671
0.3033
0.3262
0.5052
0.2796
0.3085
0.3248
0.5052
2
0.2777
0.3157
0.3340
0.5053
0.2908
0.3184
0.3276
0.5065
3
0.2805
0.3102
0.3057
0.5100
0.2895
0.3019
0.3279
0.5477
4
0.2716
0.3074
0.2753
0.5083
0.2838
0.2856
0.3162
0.5323
5
0.2744
0.3054
0.3303
0.5140
0.2842
0.3115
0.3165
0.5836
6
0.2782
0.3090
0.3104
0.5028
0.2878
0.3035
0.3193
0.4839
7
0.2823
0.3263
0.3283
0.5077
0.2983
0.3208
0.3214
0.5277
8
0.2779
0.3275
0.2845
0.5091
0.2967
0.2999
0.3069
0.5395
9
0.2576
0.3409
0.3810
0.5027
0.2881
0.3796
0.3346
0.4834
10
0.2480
0.3256
0.3627
0.5086
0.2779
0.3759
0.3423
0.5350
11
0.2390
0.3118
0.3471
0.5072
0.2698
0.3678
0.3245
0.5227
12
0.2305
0.2989
0.3315
0.5120
0.2637
0.3656
0.3445
0.5653
Average
0.2654
0.3152
0.3264
0.5077
0.2842
0.3283
0.3255
0.5277Finally, we report the experimental results on PeMS data set over a different number of time steps with a varying number of clusters. The clustering results of all these algorithms on the PeMS data set are shown in Table4. STClu outperforms the other three algorithms in terms of both ACC and NMI, because it simultaneously considers the geometric structure and incorporates the prior knowledge and the spatial correlation between sensors in proximity.Table 4
Comparison of the ACC and NMI results for the four methods on PeMS with the number of clusters varying from 2 to 20 at different time steps.
Time step
ACC
NMI
K-means
Ncut
NMF
STClu
K-means
Ncut
NMF
STClu
1
0.2694
0.3024
0.3213
0.5056
0.2802
0.3056
0.3254
0.5090
2
0.2731
0.3054
0.3202
0.5059
0.2835
0.3066
0.3237
0.5116
3
0.2733
0.3018
0.3073
0.5104
0.2818
0.2986
0.3268
0.5509
4
0.2679
0.3031
0.2954
0.5103
0.2799
0.2934
0.3141
0.5501
5
0.2716
0.2987
0.3150
0.5109
0.2795
0.3007
0.3135
0.5559
6
0.2712
0.2990
0.3054
0.5072
0.2795
0.2962
0.3165
0.5225
7
0.2760
0.3116
0.3190
0.5067
0.2880
0.3091
0.3172
0.5188
8
0.2727
0.3111
0.2952
0.5098
0.2861
0.2971
0.3014
0.5453
9
0.2670
0.3193
0.3431
0.5096
0.2847
0.3376
0.3185
0.5439
10
0.2607
0.3184
0.3494
0.5107
0.2822
0.3466
0.3389
0.5536
11
0.2481
0.3093
0.3384
0.5085
0.2731
0.3398
0.3301
0.5344
12
0.2539
0.3044
0.3084
0.5129
0.2757
0.3286
0.3308
0.5731
Average
0.2654
0.3152
0.3264
0.5077
0.2842
0.3283
0.3255
0.5277


## 5.1. Baseline Algorithms and Evaluation Methods

To demonstrate the efficiency ofSTClu, we compare it with the following three popular clustering algorithms. The first is K-means clustering in the principle component analysis (PCA) subspace. The other two algorithms are the typical spectral clustering algorithm normalized cut (Ncut) [34] and NMF-based clustering [19].To evaluate clustering quality, all our comparisons are based on clustering accuracy (ACC) and normalized mutual information (NMI) measurements.Clustering accuracy (ACC) discovers the one-to-one relationships between clusters and classes and measures the extent to which each cluster contains data points from the corresponding class. It is defined as [18, 19]
(13)
ACC
=
∑
i
=
1
n
δ
(
map
(
r
i
)
,
l
i
)
n
,
where r
i denotes the cluster label, l
i denotes the true class label, n is the total number of documents, δ
(
x
,
y
) is the delta function that is equal to one if x
=
y and zero otherwise, and map
(
r
i
) is the permutation mapping function that maps each cluster label r
i to the equivalent label from the data set.Between two random variablesi (category label) and j (cluster label), NMI is defined as [35]
(14)
NMI
=
∑
i
=
1
k
(
a
)
∑
j
=
1
k
(
b
)
n
i
,
j
log
⁡
(
(
n
·
n
i
,
j
)
/
(
n
i
·
n
j
)
)
(
∑
i
=
1
k
(
a
)
n
i
log
⁡
(
n
i
/
n
)
)
(
∑
j
=
1
k
(
b
)
n
j
log
⁡
(
n
j
/
n
)
)
,
where n is the number of documents, n
i and n
j denote the number of documents in category i and cluster j, respectively, n
i
,
j denotes the number of documents in category i as well as in clusterj, and k
(
a
) and k
(
b
) are the cluster numbers in the true and predicted clusters. The NMI score is 1 if the clustering results perfectly match the category labels whereas the score is 0 if the data are randomly partitioned. The higher the NMI score is, the better the clustering quality is.Since there are no predefined categories in our data, we have to design an alternative way to carry out the evaluation. In this paper, we use theK-means clustering result of the combined feature data as the target labels for feature clustering. Similarly, the K-means clustering result of the combined user data is used as the target labels for user clustering. To randomize the experiments, we conduct the evaluation using different cluster numbers. Note that all the experiments ran on the same machine. At each time step, for each given cluster number, 30 test runs are carried out on different randomly chosen clusters. For each experiment, we run it 10 times and report the average.


## 5.2. Synthetic Data Sets

In this section, we use synthetic datasets to demonstrate thatSTClu can significantly cluster data for multiple evolving data streams over time. The number of nodes and number of clusters in each time step are given. The synthetic datasets are generated by two steps [18]. The first step generates the objects’ geographic information, and the second step generates the values of the objects’ nongeographic attributes. The objects’ nongeographic attributes are generated by applying prototype systems used in [6]. It can be defined as follows: f
(
t
+
Δ
t
)
=
f
(
t
)
+
f
′
(
t
+
Δ
t
), f
′
(
t
+
Δ
t
)
=
f
′
(
t
)
+
u
(
t
), in which f
(
·
) is a stochastic process, t
=
0, and Δ
t
,
2
Δ
t
,
… are independent random variables, uniformly distributed in an interval [
-
a
,
a
]. So, the data stream S
(
·
) can be defined by S
(
t
)
=
f
(
t
+
h
(
t
)
)
+
g
(
t
), where h
(
·
) and g
(
·
) are stochastic processes that are generated in the same way as the prototype f
(
·
). The constant a determines the smoothness of a data stream, which can be different for p
(
·
), h
(
·
), and g
(
·
), such as 0.05, 0.1, and 0.2, respectively. To examine the quality of STClu for clustering multiple data streams, we generate random datasets with the number of streams ranging from 100 to 1000. Each stream has a length of 2,000 points.We show the performance of four methods while the time window evolving over time and the number of clusters in the data varies from 2 to 20. The stream number is fixed at 400, and each of them contains 1,000 points. As can be seen in Table1, the average ACC and NMI among four methods with the different numbers on synthetic dataset are shown. It can be clearly found that the STClu outperforms all other approaches and having highest NMIs at evolving time steps. The reason is that it incorporates the previous clustering results of the upstream sensors and the spatial relationship between the upstream and downstream sensors.Table 1
Comparison of ACC and NMI results for the four methods based on synthetic data sets.
k
ACC
NMI
K-means
Ncut
NMF
STClu
K-means
Ncut
NMF
STClu
2
0.3959
0.4954
0.5449
0.6542
0.3480
0.5130
0.5895
0.6749
4
0.3855
0.4740
0.5120
0.6436
0.4148
0.5233
0.4416
0.6967
5
0.3773
0.4685
0.5205
0.6757
0.3810
0.5016
0.5781
0.6122
6
0.3710
0.4866
0.5425
0.5603
0.4162
0.5430
0.6041
0.6059
8
0.3864
0.4636
0.4893
0.6038
0.4257
0.5247
0.5098
0.7022
10
0.3711
0.4540
0.4919
0.6231
0.3960
0.5154
0.4082
0.7238
12
0.3896
0.4802
0.5231
0.5868
0.4051
0.5088
0.5916
0.5421
15
0.3819
0.4728
0.5102
0.4834
0.4180
0.5207
0.5252
0.6281
18
0.4121
0.4827
0.5082
0.5293
0.4315
0.5786
0.5850
0.5823
20
0.4019
0.5275
0.6110
0.6516
0.4171
0.5823
0.4390
0.6299
Average
0.3873
0.4805
0.5254
0.6012
0.4053
0.5311
0.5272
0.6398The average ACC and NMI value for the proposed algorithm varies with the varied number of clusters or data streams shown in Figures3(a) and 3(b), respectively. As shown in Figure 3, we can observe that the performance of the proposed algorithm STClu is very robust to the number of clusters and data streams.The performance of the proposed algorithm varying with the size of (a) clusters, (b) data streams, respectively.
(a)
(b)Finally, we evaluate the average processing time for a round of clustering multiple data streams. There are different factors affecting the execution time. In this experiment, we evaluate the effect of the window size and the number of data streams on the response time for the compared algorithms.The following set of experiments evaluate the effect of window size on the execution time of these algorithms. Figure4(a) shows that as the window size increases, the processing time of these algorithms increases. The y-axis shows the execution time, while the x-axis represents the window size w
(
t
) from 60 to 600. In Figure 4(b), the y-axis shows the execution time, while the x-axis represents the number of data streams being varied from 100 to 1000. The other parameters are fixed as follows: w
(
t
)
=
300 s, k
=
4. Note that as the number of streams increases, the execution time of these methods increases. However, the time of the STClu and NMF increases much more slowly than others. This is because of the execution efficiency of STClu and NMF for clustering multiple data streams in the approximate subspace. Although STClu requires more time than NMF, as can be seen from the experiments, it is able to process extremely large sensor data streams with spatial, temporal, and other multidimensional information.Average processing time: (a) while the window sizew
(
t
) varies from 60 to 600, (b) while the number of streams varies from 100 to 1000, respectively.
(a)
(b)


## 5.3. Real World Data Sets

The real world data sets were obtained from the PeMS (http://pems.dot.ca.gov/) database. With sensor devices installed on road networks, a system monitors traffic flows on major US highways 24 hours a day, 7 days a week, thus acquiring huge volumes of data. We collected data from 963 sensors positioned along the freeway. We selected 15 months of daily data records giving the occupancy rates (between 0 and 1) of different car lanes on freeways in the San Francisco bay area over time. The data set covered records sampled every 10 min from 1 January 2008 to 30 March 2009. Data for public holidays were removed from the data set, as well as data collected on two days with anomalies (8 and 9 March 2008), where all sensors were muted between 2:00 and 3:00 am.First, we discuss the experiments on the PeMS data set with the number of clusters varying from 2 to 20. For each given cluster numberk, the clustering results of all these algorithms on the PeMS data set are shown in Table 2; Ncut, NMF, and the proposed STClu algorithm achieve better performance than the K-means algorithm, since these three algorithms consider the geometrical structure information contained in the data.Table 2
Comparison of the ACC and NMI results for the four methods on PeMS with the number of clusters varying from 2 to 20.
k
ACC
NMI
K-means
Ncut
NMF
STClu
K-means
Ncut
NMF
STClu
2
0.2716
0.3014
0.3163
0.5060
0.2808
0.3027
0.3260
0.5127
3
0.2685
0.2950
0.3064
0.5065
0.2762
0.2947
0.3197
0.5167
5
0.2660
0.2933
0.3089
0.5107
0.2741
0.2952
0.3256
0.5540
8
0.2641
0.2988
0.3155
0.5123
0.2759
0.3011
0.3119
0.5679
10
0.2687
0.2919
0.2996
0.5078
0.2748
0.2899
0.3105
0.5282
12
0.2642
0.2890
0.3004
0.5115
0.2711
0.2889
0.3137
0.5611
15
0.2697
0.2969
0.3097
0.5057
0.2777
0.2973
0.3129
0.5098
19
0.2674
0.2946
0.3059
0.5104
0.2754
0.2943
0.2959
0.5510
20
0.2764
0.2976
0.3052
0.5164
0.2813
0.2955
0.3023
0.6043
Average
0.2685
0.2954
0.3075
0.5097
0.2764
0.2955
0.3132
0.5451Next, we discuss the experiments on the PeMS data set with different time steps. For each time step, the clustering results of all these algorithms on the PeMS data set are shown in Table3. We observe that the STClu algorithm performs better than the other three algorithms for most time steps, except t
=
1. The reason for this is that STClu integrates more historical knowledge and the spatial relationship between two sensors in proximity.Table 3
Comparison of the ACC and NMI results for the four methods on PeMS with different time steps.
Time step
ACC
NMI
K-means
Ncut
NMF
STClu
K-means
Ncut
NMF
STClu
1
0.2671
0.3033
0.3262
0.5052
0.2796
0.3085
0.3248
0.5052
2
0.2777
0.3157
0.3340
0.5053
0.2908
0.3184
0.3276
0.5065
3
0.2805
0.3102
0.3057
0.5100
0.2895
0.3019
0.3279
0.5477
4
0.2716
0.3074
0.2753
0.5083
0.2838
0.2856
0.3162
0.5323
5
0.2744
0.3054
0.3303
0.5140
0.2842
0.3115
0.3165
0.5836
6
0.2782
0.3090
0.3104
0.5028
0.2878
0.3035
0.3193
0.4839
7
0.2823
0.3263
0.3283
0.5077
0.2983
0.3208
0.3214
0.5277
8
0.2779
0.3275
0.2845
0.5091
0.2967
0.2999
0.3069
0.5395
9
0.2576
0.3409
0.3810
0.5027
0.2881
0.3796
0.3346
0.4834
10
0.2480
0.3256
0.3627
0.5086
0.2779
0.3759
0.3423
0.5350
11
0.2390
0.3118
0.3471
0.5072
0.2698
0.3678
0.3245
0.5227
12
0.2305
0.2989
0.3315
0.5120
0.2637
0.3656
0.3445
0.5653
Average
0.2654
0.3152
0.3264
0.5077
0.2842
0.3283
0.3255
0.5277Finally, we report the experimental results on PeMS data set over a different number of time steps with a varying number of clusters. The clustering results of all these algorithms on the PeMS data set are shown in Table4. STClu outperforms the other three algorithms in terms of both ACC and NMI, because it simultaneously considers the geometric structure and incorporates the prior knowledge and the spatial correlation between sensors in proximity.Table 4
Comparison of the ACC and NMI results for the four methods on PeMS with the number of clusters varying from 2 to 20 at different time steps.
Time step
ACC
NMI
K-means
Ncut
NMF
STClu
K-means
Ncut
NMF
STClu
1
0.2694
0.3024
0.3213
0.5056
0.2802
0.3056
0.3254
0.5090
2
0.2731
0.3054
0.3202
0.5059
0.2835
0.3066
0.3237
0.5116
3
0.2733
0.3018
0.3073
0.5104
0.2818
0.2986
0.3268
0.5509
4
0.2679
0.3031
0.2954
0.5103
0.2799
0.2934
0.3141
0.5501
5
0.2716
0.2987
0.3150
0.5109
0.2795
0.3007
0.3135
0.5559
6
0.2712
0.2990
0.3054
0.5072
0.2795
0.2962
0.3165
0.5225
7
0.2760
0.3116
0.3190
0.5067
0.2880
0.3091
0.3172
0.5188
8
0.2727
0.3111
0.2952
0.5098
0.2861
0.2971
0.3014
0.5453
9
0.2670
0.3193
0.3431
0.5096
0.2847
0.3376
0.3185
0.5439
10
0.2607
0.3184
0.3494
0.5107
0.2822
0.3466
0.3389
0.5536
11
0.2481
0.3093
0.3384
0.5085
0.2731
0.3398
0.3301
0.5344
12
0.2539
0.3044
0.3084
0.5129
0.2757
0.3286
0.3308
0.5731
Average
0.2654
0.3152
0.3264
0.5077
0.2842
0.3283
0.3255
0.5277


## 6. Conclusions

Clustering multiple sensor data streams in CPS has been extensively studied in various applications, including transportation systems, sensor networks, and social networks. To extract and retain meaningful information from multiple sensor data streams, we assume that the spatial feature is the summary of the atypical event in temporal dimension, and the temporal feature is the summary of the event in spatial dimension. In this work, we have proposed a spatial-temporal clustering algorithm, calledSTClu. It utilizes time-series observational data streams and geospatial relationship for clustering multiple sensor data streams. It aims to discover clusters of objects with similar behavior and discover potential anomalies which may be revealed by the relationship evolving over time. STClu simultaneously incorporates the spatial relationship between two sensors in proximity and integrates the historical information into consideration. Experimental results on synthetic and real world data sets show that STClu achieves a better performance for clustering multiple sensor data streams evolving over time.



---

*Source: 824904-2014-07-17.xml*