Graph Based Multi-Layer K-Means++ (G-MLKM) for Sensory Pattern Analysis in Constrained Spaces

In this paper, we focus on developing a novel unsupervised machine learning algorithm, named graph based multi-layer k-means++ (G-MLKM), to solve the data-target association problem when targets move on a constrained space and minimal information of the targets can be obtained by sensors. Instead of employing the traditional data-target association methods that are based on statistical probabilities, the G-MLKM solves the problem via data clustering. We first develop the multi-layer k-means++ (MLKM) method for data-target association at a local space given a simplified constrained space situation. Then a p-dual graph is proposed to represent the general constrained space when local spaces are interconnected. Based on the p-dual graph and graph theory, we then generalize MLKM to G-MLKM by first understanding local data-target association, extracting cross-local data-target association mathematically, and then analyzing the data association at intersections of that space. To exclude potential data-target association errors that disobey physical rules, we also develop error correction mechanisms to further improve the accuracy. Numerous simulation examples are conducted to demonstrate the performance of G-MLKM, which yields an average data-target association accuracy of 92.2%.


Introduction
Associating data with the right target in a multi-target environment is an important task in many research areas, such as object tracking [1], surveillance [2,3], and situational awareness [4]. Image sensors can be used to acquire rich information related to each target, which will significantly simplify the data-target association problem. For example, video cameras in a multi-target tracking mission can provide colors and shapes of targets as extra features in the association process [5]. However, considering the costs, security issues, and special environments (e.g., ocean tracking [6], military spying), a simple, reliable, and low-cost sensor network is often a preferred option [7]. Consequently, the data-target association problem needs to be further studied, especially in cases when the gathered data are cluttered and contains limited information related to the targets.
The existing approaches for data-target association, in general, consist of three procedures [8]: (i) Measurements collection-preparation before data association process, such as object identification in video frames, radar signals processing, or raw sensor data accumulation; (ii) measurements prediction-predict the potential future measurements based on history data, which yields an area (validation gate) that narrows down the search space; and (iii) optimal measurement selection-select the optimal measurement that matches history data according to a criterion (varies in different approaches) and update the history dataset. With the same procedures but different choices of the optimal measurement criteria, many data-target association techniques have already been developed. Among them, the wellknown techniques include the global nearest neighbor standard filter (Global NNSF) [9], joint probabilistic data association filter (JPDAF) [10][11][12][13], and multiple hypothesis tracking (MHT) [14].
The Global NNSF approach attempts to find the maximum likelihood estimate related to the possible measurements (non-Bayesian) at each scan (that measures the states of all targets simultaneously). For nearest neighbor correspondences, there is always a finite chance that the association is incorrect [15]. Besides that, the Global NNSF assumes a fixed number of targets and cannot adjust the target number during the data association process. A different well-known technique for data association is JPDAF, which computes association probabilities (weights) and updates the track with the weighted average of all validated measurements. Similar to Global NNSF, JPDAF cannot be applied in scenarios with targets birth and death [1]. The most successful algorithm based on this data-oriented view is the MHT [16], which takes a delayed decision strategy by maintaining and propagating a subset of hypotheses in the hope that future data will disambiguate decisions at present [1]. MHT is capable of associating noisy observations and is resistant to a dynamic number of targets during the association process. The main disadvantage of MHT is its computational complexity as the number of hypotheses increases exponentially over time.
There are other approaches available for data association. For example, the Markov chain Monte Carlo data association (MCMCDA) [5,17]. MCMCDA takes the data-oriented, combinatorial optimization approach to the data association problem but avoids the enumeration of tracks by applying a sampling method called Markov chain Monte Carlo (MCMC) [17], which implements statistical probabilities in the procedure of optimal measurement selection as well. In this paper, we assume an object generates at most a single detection in each sensor scan, namely, a point-target assumption. Hence, the approaches on multiple detections per object per time step, i.e., extended-target [18], are not discussed here. The data association in extended object tracking problems typically use data clustering techniques, such as k-means [19], to address the extended-target issue by specifying which measurements are from the same source. Then the corresponding association problems can be simplified as point-target tracking problems. For example, the authors in [20] proposed a clustering procedure and took into account the uncertainty and imprecision of similarity measures by using a geometric fuzzy representation, which shows the potential of applying clustering algorithms in the data association problem.
The main contribution of this paper is the development of an efficient unsupervised machine learning algorithm, called graph based multi-layer k-means++ (G-MLKM). The proposed G-MLKM differs from the existing data-target association methods in three aspects. First, in contrast to the previous developed data association approaches that estimate the potential measurement from history data for each target and select an optimal one from validated measurements based on statistical probabilities, G-MLKM solves the data-target association problem in the view of data clustering. Second, the previous approaches are mainly developed with respect to sensors that are capable of obtaining information from a multiple dimensional environment, such as radars, sonars, and video cameras. G-MLKM is proposed on sensors that only provide limited information. Interesting research on tracking targets with binary proximity sensors can be seen in [7], whose objective is only limited to target counting, while G-MLKM can associate data to targets. Third, G-MLKM can address the case that targets move in a constrained space, which requires dealing with data separation and merging.
The reminder of this paper is structured as follows. The data association problem in a constrained space and the corresponding tasks are described in Section 2. In Section 3, the multi-layer k-means++ (MLKM) method is developed for data-target association at local space given a simplified constrained space situation. The graph based multi-layer k-means++ (G-MLKM) algorithm is then developed in Section 4 for general constrained spaces. Simulation examples are then provided in Section 5. Section 6 provides a brief summary of the work presented in this paper.

Problem Formulation
In this paper, we consider the problem of data-target association when multiple targets move across a road network. Here, a road network is a set of connected road segments, along which low-cost sensors are spatially distributed. The sensors are used to collect information of targets, which, in particular, are the velocity of targets and the corresponding measured time. We assume (1) there is no false alarm in the sensor measurements, and (2) the target's velocity does not change rapidly within two adjacent sensors. The collected information about a target is normally disassociated with the target itself, meaning that the target from which the information was captured cannot be directly identified using the information. Hence, data-target associations is necessary. Figure 1 shows one road network example that consists of 6 road segments. Without loss of generality, let the total number of road segments in one road network be denoted as L. The road segments are denoted as R 1 , R 2 , · · · , R L , respectively. The length of road segment R i is denoted as D i for i = 1, 2, · · · , L. To simplify discussion, we assume the road segments are for one-way traffic, i.e., targets cannot change their moving directions within one road segment. However, when the road segment allows bidirectional traffic, we can separate it into two unidirectional road segments and the proposed approach in this paper directly applies. Let S i = {S i1 , S i2 , · · · , S iN i } be a set of N i ∈ R sensors placed along the direction of road segment R i . In other words, for sensor S ij ∈ S i , the larger the sub-notation j is, the further distance the sensor locates away from the starting point of road segment R i . We denote the corresponding distance between sensor S ij and the starting point of road segment R i as d ij . Hence, the position set for sensors in R i related to the starting point can be denoted as For each sensor S ij , its measurements are collected and stored in chronological order.
The collections are denoted as a column vector X ij , such that , where i ∈ {1, 2, · · · , L}, j ∈ {1, 2, · · · , N i }, the prime symbol represents the transpose operation for a vector, m ij is the total number of measurements in X ij , and x n ij , n ∈ {1, 2, · · · , m ij }, denotes an individual measurement in X ij . In particular, x n ij = [v n ij , t n ij ] stores the measured velocity v n ij when one target passed by sensor S ij at time t n ij . As the elements in X ij are stored in chronological order, the recorded time for each measurement satisfies t 1 ij < t 2 ij < · · · < t m ij ij , which can be distinguished based on the superscript n. All the measurement vectors stored by sensors that locate in the same road segment R i are stored into a matrix X i , such The added all-zero row vector inX ij is to unify the length of vectors in matrix X i considering that miss detection may happen or targets may remain (or stop) inside the road network for a given data collection period.
The road network collects X i , i = 1, · · · , L that only include information of target's velocity and the corresponding measurement time. In order to solve data-target association based on the L matrices, three tasks need to be accomplished. The first task (Task 1) is to cluster X i into m i groups for each road segment. Denote the data grouping result for each road segment as a new matrix T i , such that: whereT iz , z = 1, 2, · · · , m i , is a row vector consisting of N i measurements associated with the same target, defined as:T where τ u iz is an entry ofX iu for u = 1, 2, · · · , N i . Then a new row vector T iz is obtained fromT iz by excluding all zero elements.
The second task (Task 2) is to link the trajectories of targets at road intersections by pairing sensor S i1 /S iN i from multiple road segments that are connected geometrically. In particular, let O i nts T denote the index set of road segments that have outgoing targets related to one intersection i nts , and I i nts T denote the index set of road segments that have ingoing targets related to the same intersection. Since the road segments are unidirectional, the two index sets have no overlaps, i.e., O i nts T ∩ I i nts T = ∅. In particular, only the dataset that has a subscript of iN i (according to the unidirectional road segement setting) can be the candidate for O i nts T . Similarly, only the dataset that has a subscript notation of i1 can be the candidate for I i nts T . Therefore, datasets that belong to targets who move towards the intersection i nts are denoted as: while datasets that belong to targets who leave the intersection i nts are denoted as: where k O ∈ {1, 2, · · · , m iN i }, k I ∈ {1, 2, · · · , m i1 }, O i nts T ⊂ {1, 2, · · · , L}, and I i nts T ⊂ {1, 2, · · · , L}. Since targets may stop in the intersection or the data collection process terminates before targets exit the intersection, the total number of targets heading into an intersection i nts is always greater than or equal to the number of targets leaving the same intersection, i.e., |Q i nts I | ≥ |Q i nts O |. For simplicity of notation, denote |Q i nts I | and |Q i nts O | as n I and n O . Then we can calculate n I and n O via: The pairing task for intersection i nts can be denoted as a mapping function f , such that: where x k iN i ∈ Q i nts I and x l1 ∈ {Q i nts O , 0, 0, · · · , 0 n I −n O }. In particular, the function f for intersection i nts can be denoted as a permutation matrix G i nts ∈ R n I ×n I .
The last task (Task 3) is to merge data groups on the road network when loops may exist, i.e., targets may pass the same road segment several times. Hence, multiple data association groups may belong to the same target. The merged results can be denoted as L symmetric matrices G R i ∈ R m i ×m i for each road segment R i . If targets only pass the road segment R i once, G R i is an identity matrix.
In this paper, we are going to propose a new unsupervised machine learning algorithm to associate data-target for the collected L matrices. In particular, this algorithm first creates a new clustering structure for data grouping in each matrix (associated with each road segment), and then leverages graph theory and clustering algorithms to link the matrices from different road segments for each intersection. Finally, the entire dataset can be analyzed and associated properly to the targets. The output of this new algorithm will be a detail trajectory path for each target with the captured velocities along the road segments. In the next two sections, the new data-target associations algorithm will be explained in detail. We begin the discussion with a special case when the road network is consisted of a single road segment.

MLKM for a Single Road Segment
In this section, we consider the special case when L = 1, i.e., the road network only consists of one road segment, R 1 . In this special case, there are neither intersections nor loops in the road network. Therefore, the tasks in identifying data-target associations are simplified to cluster X 1 into m 1 groups (Task 1) only. One example of matrix X 1 ∈ R 10×9 is shown in Figure 2, which is the plot of measurements for 10 different targets that are captured by nine equally spaced sensors on road segment R 1 .

K-means++ Clustering and Deep Neural Network
K-means [19] and k-means++ [21] are perhaps the most common methods for data clustering. For a set of data points in a N-dimensional system, the two algorithms perform clustering by grouping the points that are closer to the optimally placed centroids. From the machine learning perspective, k-means learns where to optimally place a pre-defined number of centroids such that the cost function, defined as Φ Y (C) = ∑ y∈Y d 2 (y, C), is minimized, where d(y, C) = min y∈Y y − c i represents the distance between a sub-set of measurements Y and a centroid c i and C = {c 1 , ..., c k } represents the set of centroids. The associated cost function is the sum of the Euclidean distances from all data points to their closer centroid. The cost function and optimization algorithms are the same for k-means and k-means++ while the only difference between them is that k-means++ places the initial guesses for the centroids in places that have data concentration, and consequently improves the running time of Lloyd's algorithm and the quality of the final solution [21].
A much more complex boundary may exist between two data groups. Therefore, we also verify the potential performance of the deep neural network (DNN) algorithm [22] in the data association process, which is known for its capability of recognizing underlying patterns and defining better decision boundaries among data samples. For the purpose of evaluating the supervised DNN capabilities, a slight modification of the problem is considered. Instead of a complete unlabeled dataset X 1 , part of the measurements are prelabeled, i.e., data-target relations for part of the measurements are known. In addition, we extend the measurement's dimensions to further include vt, v 2 t, and vt 2 as extra features so that the inner structure of DNN can be simpler. Table 1 presents the detail settings of the DNN framework.

K-means++ with Data Preprocessing
While DNN can potentially provide better performance for the data association problem, it demands labeled datasets for training. In real scenarios, however, the training dataset may not be available. In contrast, k-mean++ can cluster data samples without the need for a labeled dataset. This unsupervised property of k-means++ enables a wider application domain. Hence, k-means++ is more practical for the task of clustering X 1 into m 1 groups. Moreover, when the dataset X 1 is small and sparse, k-means++ can perform well on the task of data-target association.
However, when the measurements are distributed along the time axis and velocity profiles are close, k-means++ tends to place the centroids in positions where data from different targets overlap and hence causes an inaccurate data-target pairing. This happens because k-means implements Euclidean distance to determine which centroid data sample (v, t) belongs, i.e., arg min where C is the set of centroids. When data samples distribute along time axis, the time difference becomes the determining factor for grouping results.
One natural way to balance the two components (time difference and velocity difference) in (7) is to process X 1 before applying k-means++. The idea of preprocessing is similar to the principal component analysis [23] that projects data into a main axis. The preprocessed data sample is denoted asx n 1j = [v n 1j ,t n 1j ], wheret n 1j is given by: where j ∈ {1, · · · , N 1 }, n ∈ {1, · · · , m 1j }, d 1j is the position of sensor S 1j with respect to the starting point of road segment R 1 , and d * is the reference point for projecting. In other words,t n 1j is the expected starting time for a constant velocity (v n 1j ) model given the current time t n 1j . Figure 3 is the preprocessed result for the dataset in Figure 2. In this example, the reference point d * is selected to be the starting point, and we can see clusters for each target have been formed after data preprocessing.

Multi-Layer K-means++
Through the preprocessing procedure, data can be roughly separated for different targets that provide dense and grouped subsets. The boundaries between two groups, however, maybe still too complex for k-means++ to define, especially, when X 1 is a large dataset and the grouped subsets are close to each other. Inspired by the DNN capability of defining classification boundaries via a multi-layer structure and a back-propagation philosophy, we propose a new multi-layer k-means++ (MLKM) method that integrates the DNN's multi-layer structure with the clustering capabilities of k-means++ to overcome the complex boundary challenge.  Target 1  Target 2  Target 3  Target 4   Target 5  Target 6  Target 7   Target 8  Target 9 Target 10 The proposed MLKM algorithm is performed via 3 layers: (i) Data segmentation and clustering-the dataset is sequentially partitioned into smaller groups for the purpose of creating sparse data samples for k-means++; (ii) error detection and correction-check the clustered data by searching for errors through predefined rules and re-cluster the data using nearest neighbor concepts [24] if an error is found. Note that the k-means++ associates the data closer to the optimally placed centroid based on the Euclidean distance between data point and centroid, which is a scalar quantity; and (iii) cluster matchingmatch the clusters of each segment by preprocessing the cluster centroids of all segments to the cluster centroid of the first segment and again grouping them based on k-means++. A detail explanation for these three layers are given as follows.

Layer 1 (Data Segmentation & Clustering)
Without loss of generality, we assume that there are K sensors per segment. The dataset X 1 ∈ R m 1 ×N 1 (m 1 and N 1 are the maximum number of measurements and the total sensor number in sensor set S 1 , respectively) is sequentially partitioned into E segments, such that: In other words, when N 1 %K = 0, the last segment will contain measurements from less than K sensors. In the following of the paper, we assume that N 1 %K = 0 in the following sections of this paper for the simplicity of presentation. When N 1 %K = 0, we can add some extra artificial sensors with all zero measurements. Then the data segment can be defined as X 1e = eK j=(e−1)K+1X 1j , where e = 1, 2, · · · , E. K-means++ algorithm is then applied to each X 1e by excluding all zero elements. By aggregating the clustering results, we can obtain a set of centroids for X 1e , e = 1, 2, · · · , E, defined as C 1e = {c e 11 , c e 12 , · · · , c e 1m 1 }, and the associated measurements with each c e 1k centroid are represented as T e 1k , where k ∈ {1, 2, · · · , m 1 }.

Layer 2 (Error Detection & Correction)
The first layer seeks to associate data for each data segment. Since the clustering standard used in k-means++ is a scalar quantity while the actual measurements are given by vectors, there are potential data association errors in T e 1k . Hence an additional layer to perform error detection and correction is needed. The error detection is to verify logic rules to determine if wrong data association appears in T e 1k . The error correction will conduct data re-association on the identified wrong associations. To avoid the same wrong reassociation again, the global nearest neighbor standard is chosen as the re-association technique instead of k-means++ given the assumption that the target's velocity does not change rapidly within two adjacent sensors.
We here proposed the following logic rules for error detection: 1k | indicates the cardinality of T e 1k . The first rule means that more than K measurements appear in T e 1k . The second rule means that more than one sensory measurements from the same sensor are associated with one target inT e 1k . The third rule means that target is recorded in a later time by a previous sensor. If one or more rules are satisfied, the corresponding T e 1k is then considered to be an erroneous data association and will be stored in Y * 1e , where Y * 1e refers to the wrong data associations in X 1e . The error correction is to re-associate data in Y * 1e for the purpose of breaking all the logic rules listed above. We propose to use the global nearest neighbor approach. Specifically, elements in Y * 1e that belongs to measurements of sensor S 1 are selected sequentially to be evaluated against with every measurement in Y * 1e that belongs to measurements of sensor S 1( +1) to obtain the best match. The evaluation is accomplished via the following optimization process: With this procedure, all T e 1k are updated with the corrected clusters and all c e 1k are re-calculated based on the updated T e 1k . The new corrected set of centroids C 1e is updated for all segments and grouped into C 1 = {C 11 C 12 ... C 1E }. The position of the centroid set C 1e is defined as:

Layer 3 (Cluster Matching)
Through the preceding two layers, data-target association can be accomplished for each data segment X 1e independently. However, the target associations are uncorrelated among each data segment. In particular, the unsupervised k-means++ only groups data samples that belong to the same target while the clusters of each target are anonymous. Hence, it is still unclear how to associate the clusters among different segments.
In Layer 3, we project C 1e , e = 1, · · · , E, using the preprocessing technique that is stated in Section 3.2. More precisely, the time component in c e 1k ∈ C 1e is preprocessed as: where c e 1k = [v e 1k , t e 1k ], and d 1e is the position of centroid set C 1e defined in (9). Then k-means++ is applied to the preprocessed C 1 to find the clusters that group cluster centroids in different data segments. Accordingly, the associated measurements T e 1k with respect to each centroid are merged together as T 1k and, hence, provides the complete data-target association result for the entire road segment.
Note that the proposed MLKM method may not be applied directly to the case when L > 1 (i.e., more than one road segments). Therefore, we propose a more general method, named G-MLKM, to solve the general data-target association problem for a general road network in the next section.

G-MLKM for a General Road Network
In this section, we consider the general case when the road network consists of multiple road segments. To solve the data-target association problem, we propose a new graph-based multi-layer k-means++ (G-MLKM) algorithm. In particular, G-MLKM uses graph theory to represent the road network as a graph, and then links data from different road segments at each intersection of the road network by analyzing the graph structure. The data-target association problem for a general road network is then solved by merging the clustering results at intersections with the MLKM results on each road segment.
We first briefly introduce graph theory and the representation of road networks using graphs as preliminaries. Then the procedures for G-MLKM are explained in detail. In particular, we begin with a new graph representation for the road network. Then the procedures for linking measurements at intersections (Task 2) are described. After that, we unify the results on road segments and intersections, and complete the data merging task (Task 3).

Graph Theory
For a system of L connected agents, its network topology can be modeled as a directed graph G = (V, E ), where V = {v 1 , v 2 , · · · , v L } and E ⊆ V × V are, respectively, the set of agents and the set of edges that connect the agents. An edge (v i , v j ) in set E means that the agent v j can access the state information of agent v i , but not necessarily vice versa [25].

Graph Representation of Road Networks
There are mainly two strategies to represent road networks using a graph, namely a primal graph and dual graph [26]. In a primal graph representation, road intersections or end points are represented by agents and road segments are represented by edges [27], while in a dual graph representation, road segments are represented by agents and an edge exists if two roads are intersected with each other [28]. Compared with a primal graph, dual graph concerns more on the topological relationship among road segments. As the data-target associations for each road segment can be solved by the MLKM method, the focus here is to cluster data at each intersection. As a consequence, the dual graph is a better option. However, the geometric properties such as road length are neglected by a dual graph. Hence, some further modification to the dual graph is needed.

G-MLKM Algorithm
In this subsection, we will provide the detail procedures for the G-MLKM algorithm that are composed of the following three steps.

Modified Graph Representation for Road Networks
Considering the cases when targets may stop in a road segment or data collection process may terminate before targets pass through a road segment, the total number of measurements collected by sensor S iN i (locates near the ending point of road segment R i ) may be less than the one collected by sensor S i1 (locates near the starting point of road segment R i ). If the entire road segment is abstracted as one single agent, the inequality of measurements in the road segment may create issues for the subsequent data-target associations process. Here, we modify the dual graph by incorporating the primal graph for the representation of the road segment. In other words, we propose to replace each road segment node in the dual graph by two agents with one directed edge connecting them and the direction of the edge is determined by the traffic direction. In particular, we use the sensor nodes S i1 and S iN i as the two agents. We may neglect the edge between S i1 and S iN i because we focus on data-target associations at intersections while the data-target associations within the road segment can be accomplished by the MLKM method without the need for the knowledge of the graph. Moreover, the connection between S i1 and S iN i is unidirectional when the traffic is unidirectional. We call the new graph the"p-dual graph", i.e., prime-based dual graph. An example of how to derive the p-dual graph is shown in Figure 4, where the original six agents in the dual graph are replaced by 12 agents and the edges between S i1 and S iN i are removed in the p-dual graph.  Figure 1 with the nodes and arrows representing, respectively, the agents and directed edges. (b) P-dual graph representation for the road network in Figure 1, where two sensor nodes represent one road segment and the edge within the two sensor nodes are ignored. In this example, there exist 3 subgraphs which are denoted as a, b, and c.
For a general road network with L edge segments, the edges of the new p-dual graph is given by V * = {S 11 , S 1N 1 , S 21 , · · · , S L1 , S LN L } with the corresponding adjacency matrix, A * ∈ R 2L×2L , given by:

Graph Analysis for Data Pairing at Intersections
From A * defined in (10), we can observe that the adjacency matrix A * has L columns and L rows that are all zeros. Hence, the sparse matrix A * can be further analyzed and decomposed to extract subgraphs related to different intersections. Then the task of linking the trajectories of targets at road intersections can be equivalently solved via pairing measurements of sensor S i1 /S iN i from road segments in the subgraphs, which is further decomposed into the following three procedures.

i. Subgraph Extraction
The first procedure is to extract subgraphs from A * . Let the letters in alphabet {a, b, c, ...} denote the names for different intersections. The subgraph extraction procedure begins with an intersection name as a, follows by b, c, and so on. For any intersection i nts , the subgraph extraction is conducted by cross-searching the non-zero entries of the matrix A * in a repeated row and column pattern. The corresponding indices of row and column containing non-zero entries, indicating the agents and edges that are included in that subgraph, are stored in the sets O i nts T and I i nts T , respectively. More precisely, O i nts T denotes the index set of road segments that have outgoing targets related to intersection i nts and I i nts T denotes the index set of road segments that have ingoing targets related to the same intersection. The index storing processes are defined as O i nts T = O i nts T ∪ {i}, and I i nts T = I i nts T ∪ {j}, where i, j are the corresponding row index and column index, respectively. The iterative search process will terminate and return (O i nts T , I i nts T ) when there is no more non-zero element in the recorded row and column indices. Algorithm 1 is the pseudo code for the subgraph extraction procedure. The extracted results are denoted as (O i nts T , I i nts T ), where i nts ∈ {a, b, c, · · · }. Algorithm 1 Subgraph Extraction 1: Input: ∀b ij ∈ A * ; 2: Output: (O i nts T , I i nts T ), i nts ∈ {a, b, c, · · · } 3: Idx row = Idx col = {1, 2, · · · , |A * |}; 4: i = 0; 5: for i nts in {a, b, c, · · · } do 6: O i nts T = I i nts T = ∅; 7: if |Idx row | ≥ 1 then 8: procedure INCREMENT(i) 9: i = i + 1; 10: if i ∈ Idx row then 11: return i; 12: else 13: INCREMENT(i); 14: procedure RECURSION(i) 15: if ∑ ∀j∈Idx col b ij ≥ 1 then 16: O i nts T = O i nts T ∪ {i}; 17: procedure EXTRACT(i) 18: for j in Idx col do 19: if b ij = 0 then 20: Idx col = Idx col \I i nts T ; 22: Idx row = Idx row \{i}; 23: for j ∈ I i nts T do 24: if ∑ ∀l∈Idx row b lj ≥ 1 then 25: for l in Idx row do 26: if b lj = 0 then 27: O i nts T = O i nts T ∪ {l}; 28: if Idx row ∩ O i nts T = ∅ then 29: EXTRACT(∃l ∈ (Idx row ∩ O i nts T ));

ii. Data Preprocessing at Intersections
Given that the subgraph that describes an intersection, i nts , is available from the preceding subgraph extraction procedure, datasets of X i1 /X iN i which are subjected to the pairing task for the corresponding intersection can be pinpointed. In particular, (3) and (4) define the dataset for the intersection i nts as an incoming dataset Q i nts I and an outgoing dataset Q i nts O , respectively. As we assume that (1) no false alarm in the measurements, and (2) the target's velocity does not change rapidly within two adjacent sensors, data pairing at intersections may interpret as data clustering. A potential machine learning technique for data clustering is the k-means++. However, the sensors S i1 /S iN i from different road segments are not guaranteed to locate near each other for a road intersection, which may contribute to a relatively large time difference in two sensors' measurements for one target. Hence, before applying k-means++, data preprocessing on Q i nts I and Q i nts O is necessary.
Based on the proposed preprocessing definition in (8), we here propose a new data preprocessing technique that first selects a virtual reference at the center of the intersection i nts and then recomputest k ij via projecting each element in I i nts T and O i nts T to the virtual reference as:t where k ∈ {1, 2, · · · , m ij } and r is the radius of the intersection circle centered at the virtual reference. An example of locating the virtual reference is shown in Figure 5, where the intersection consists of three road segments denoted as R i , R j , and R k . Figure 5. An intersection consists of three road segments denoted as R i , R j , and R k . The virtual reference for data preprocessing is in the center of the intersection with a radius of r to each road segment ending point.

iii. Data Pairing at Intersections and Error Correction
Denote the preprocessed datasets for Q i nts I and Q i nts O asQ i nts I andQ i nts O . Then k-means++ can be applied to the preprocessed intersection datasets {Q i nts I ,Q i nts O } for data pairing. Similar to the development of MLKM for the case of one road segment, errors may arise when conducting the data pairing/clustering. Error detection and correction are needed to further improve accuracy.
For an intersection i nts , the cardinalities of the preprocessedQ i nts I andQ i nts O remain the same as those of Q i nts I and Q i nts O . As defined in (5), |Q i nts I | = n I and |Q i nts O | = n O , where n I ≥ n O . The set of centroids is denoted as C i nts = {c i nts 1 , c i nts 2 , · · · , c i nts n I }, and the associated measurements with each centroid c i nts j , j ∈ {1, 2, · · · , n I }, are given as Y i nts j . The error correction is similar to Layer 2 in the MLKM method described in Section 3.3.2, and defines three logic rules for error detection: where |Y i nts j | is the cardinality of Y i nts j . The first rule means more than two measurements are associated in Y i nts j . Error can be determined in this case because each target has at most two measurements in one intersection. The second rule means either none or more than one sensory measurements can be found from the incoming datasetQ i nts I . The third rule means that the outgoing measurement in Y i nts j is recorded earlier than the incoming measurement. If one or more rules are satisfied, the corresponding Y i nts j is then considered to be an erroneous data association and will be stored in Y i nts . The error correction is to re-associate data in Y i nts for the purpose of breaking all three logic rules listed above. To achieve this goal, we separate Y i nts into two subsets denoted as Y I i nts and Y O i nts given by: where Y I i nts and Y O i nts store all measurements x iN i and x l1 in Y i nts , respectively. Re-associate data in Y i nts becomes a linear assignment problem [29] between Y I i nts and Y O i nts . The optimal pairing between Y I i nts and Y O i nts can be found when the matching score reaches to the minimum via solving the optimization problem of arg min are column vectors converted from subsets Y I i nts and Y O i nts , respectively. M ∈ R m O ×m I is a special binary matrix with the summation of each row being 1. After the error correction is accomplished, all Y i nts j will be updated to complete Task 2. Furthermore, a permutation matrix G i nts ∈ R n I ×n I can be created to record the pairing relationship between incoming dataset Q i nts I and outgoing dataset Q i nts O for each intersection.

Group Merging in the Road Network
K-means++ clustering on the preprocessed dataset at each intersection solves the task of linking the trajectories of targets at road intersections (Task 2) while the proposed MLKM method solves the task of data associations for each road segment (Task 1). If the clustering results at all intersections are combined with the MLKM results on all road segments, trajectory awareness for each target in the road network is achieved. This is valid for situations when targets only pass the same road segment once. However, when targets pass the same road segment and intersection multiple times, one target can be assigned to multiple associated data groups on the road segment. To determine the connections among all associated data groups, an extra task (Task 3) for merging data groups in the road network is needed. Given that the datasets at intersections are extracted from the L matrices collected from all road segments, clusters at the intersections can be classified based on the data groups for all road segments. Therefore, the task of determining the connections among the associated data groups in the road network can be focused on connections ofT iz defined in (2) for each road segment.
Let the symmetric matrix G R i ∈ R m i ×m i denote the connections among the m i association groups in road segment R i given by 1, ifT ip ,T iq belong to the same target, 0, otherwise.
To determine the entries in G R i , the depth-first search (DFS) [30] is implemented to detect cycles in the adjacency matrix A. If cycles do not exist, the non-diagonal entries are set to 0 and hence G R i is an identity matrix. Otherwise, further analysis on the connections among data groups at each road segment is operated sequentially in the following three steps.

i. Node Analysis on Dual Graph
The analysis starts with identifying road segments that have only outgoing flow, i.e., source nodes in the graph. The source nodes can be identified from the adjacency matrix A by checking the sum of each column. In particular, road segment R i is a source node when the sum of the ith column of A satisfies L ∑ l=1 a li = 0, where a li is the (l, i)th entry of the adjacency matrix, which represents the edge (R l , R i ).

ii. Trajectory Flow for Data Groups from Source Nodes
If the road segment R i is a source node, the m i data groups in R i resulting from the MLKM method are considered to be m i unique targets. Then the trajectories of these m i targets are traced in the road network. In particular, ifT iz ∩ X iN i = ∅, the target associated with data groupT iz does not contain any measurement from sensor S iN i , which corresponds to the case when target stops in the road segment or the data collection terminates before the target could approach to sensor S iN i . The trajectory tracking for this target is then completed. Otherwise, the permutation matrix G i of intersection i that is consisted of sensor S iN i is utilized to pinpoint the trajectory of the same target in the intersection, and its data groupT lz in the subsequent node or sink node R l where it is heading to. The trajectory tracking of the same target on the new road segments will keep on until the target stops or leaves the road network. The same process is used for tracing the flow of other targets.

iii. Matrix Description of Intermediate Nodes
After the trajectories of all targets from the road segments have been confirmed, data points for each target on different road segments can be merged. More precisely, the corresponding entry (p, q) in G R l that is assigned as 1 means that data groupsT l p and T lq belong to one target. Consequently, the corresponding matrix G R l can be determined.

Simulation
In this section, the performance of the proposed G-MLKM algorithm is evaluated. We first introduce the testing datasets generation process. Then the performance of the MLKM method on one road segment is evaluated and compared with k-means++ and DNN. Then the complete G-MLKM algorithm performance is evaluated. A detailed example presenting the output via using G-MLKM is given to show how matrices G i nts and G R i are created for data pairing at intersections and group merging.

Testing Data Generation
In order to obtain a quantitative performance evaluation of the data association techniques, labeled data is needed to obtain the percentage of true association between targets and their measurements. One convenient way to have accurate labeled dataset for data-target association is to generate it artificially. Let the generated testing dataset from the road network be M t = {T 1 , T 2 , · · · , T L }, where T i ∈ R m i ×m i has the same data structure as T i defined in (1). In particular, each element in T i is a data group that belongs to one target. Moreover, for any T i collected from road segments that have both incoming and outgoing flows, multiple rows may belong to the same target.
We utilize the road network structure shown in Figure 1 as a prototype for testing data generation. Moreover, N S sensors are assumed to be equally distributed on each road segment, where the length of the road segment is N S × d. The position set for sensors is selected as P i = {d, 2d, · · · , N S d} with respect to the starting point of road segment R i . The intersections are considered to have the same radius with the value of d/2. Hence, the distance between any two adjacency sensors is d. To further simplify the data generation process, we assume road segment R 1 is the only entrance of the road network during the data collection period with incoming targets number N A , and targets have equal possibilities of valid heading directions at each intersection. The targets are assumed to move with a constant velocity and the velocity is also discretely affected by Gaussian noise, such that, v ij = v 0 + N (µ, σ), where v ij is one velocity measurement at sensor S ij and v 0 is the velocity measurement at the previous sensor. The corresponding time measurement is calculated as t ij = t 0 + v ij /(j · d). The initial velocity and time for the N A targets are uniformly selected from the range (v min , v max ) and (t min , t max ), respectively (refer to Table 2). The testing dataset generating process stops when all targets move out of the road network.
With the generated testing datasets, we may evaluate the performance of the datatarget association techniques by calculating the data association accuracy, which is defined as the ratio between correctly classified number of data (M cr ) and the total number of data (M t ), such that, numel(M cr ) numel(M t ) × 100%, where numel(M) returns the number of elements in M. As multiple testing datasets are generated, the provided statistical information about performance includes the minimum (left -blue bar), average (middle -orange bar), and maximum (right -yellow bar) accuracies.

MLKM Performance and Comparisons
Before evaluating the entire accuracy of the proposed G-MLKM algorithm, the MLKM method is evaluated and compared with the other two common data clustering machine learning techniques, in particular, k-means++ and DNN, based on the collected dataset in road segment R 1 .

K-means++
The first set of simulations evaluate the performance of k-means++ based on two criteria: (i) Unprocessed vs. preprocessed data, and (ii) using different values of N A and N S . When the values of N A and N S increase, more data points are introduced into the dataset, leading to more overlapping among these data points. Figures 6 and 7 show the performance of K-means++ using the parameters listed in Table 2.  As can be observed, a higher accuracy is achieved using the preprocessed data than that using the unprocessed data. This can be seen by comparing the average, and maximum and minimum accuracy for the two methods that use the preprocessed data versus unprocessed data, as shown in Figure 6. Using the raw data, the measurements associated with a specific target are sparse along the time axis. However, the velocity measurements from the same sensor are closely grouped along the velocity axis. These conditions contribute to incorrect clustering of the data. The preprocessing technique reduces the distance between target related measurements, therefore reducing the effect of the velocity measurements on the clustering.
A low accuracy is obtained for large values of N A and N S . This can be observed by comparing average, maximum and minimum accuracy for different N A and N S , as shown in Figures 6 and 7. Similar to the unprocessed data, a large number of sensors/targets increases the density of measurement points. The concentration of measurements increases the probability that k-means/k-means++ clusters the data incorrectly (even with preprocessing).

DNN
The k-means++ fails to correctly cluster data when overlapping of measurements occurs. A deep neural networks (DNN) is used as an alternative approach because it has been shown to provide good results to uncover patterns for large dataset classification. One necessary condition for DNN is the availability of labeled datasets for training. To meet the requirements of DNN, it is assumed that labeled data is available for training.
The results for DNN are obtained using N A = 50 targets and N S = 50 sensors. Assuming that a portion of the data association has already been identified, the objective is to train a neural network to label the unidentified measurements. The number of 'training' sensors that provide labeled information and 'testing' sensors that provide unlabeled information are provided in Table 3. The accuracy is obtained for various proportions of 'training' sensors to 'testing' sensors. Table 3 also shows the accuracy obtained for different dataset configuration. It can be observed that the training (respectively, testing) accuracy is high (respectively, low), when the testing dataset is relatively small. However, when the testing dataset is relatively high, the testing performance increases significantly (up to 91%). A high training accuracy with a low testing accuracy means that DNN suffers from overfitting due to the small size of the training dataset. Given this comparison, DNN is applicable when a large portion of a training dataset is available to train the network for classifying a relatively small amount of measurements.

MLKM
K-means++ does not provide good accuracy for a high number of measurements but performs well when clustering small amounts of data. DNN can cluster large datasets but requires a large training dataset. MLKM combines the multi-layer back-propagation error correction from DNN and the clustering capabilities of k-means++. The DNN-inspired error correction significantly improves the performance of MLKM by preventing the clustering errors in layer 1 to propagate to the cluster association in layer 3.
The results for the MLKM method are obtained using N A = 50 number of targets and N S = 20 number of sensors. In addition, the time and velocity parameters are set to (t min , t max ) = U (−10, 30) and (v min , v max ) = N (50, 40), receptively. Figure 8 shows the performance of the MLKM method with and without error correction, as well as results using the standard k-means++ method with preprocessing. It can be observed that a higher accuracy is achieved using MLKM than that using k-means++. Figure 8 shows the average, and maximum and minimum accuracy for both methods. The error correction performed in layer 2 improves the average accuracy of MLKM by approximately 7% (MLKM w/ EC 91.65%; MLKM w/o EC 84.3%).

G-MLKM Overall Performance
The results for the G-MLKM method are obtained using N A = 20 number of targets and N S = 10 number of sensors. In addition, the time and velocity parameters are set to (t min , t max ) = U(0, 40) and (v min , v max ) = U(10, 50), respectively. Figure 9 shows the performance of the G-MLKM algorithm with and without error correction.  It can be observed that a higher accuracy is achieved using G-MLKM with error correction than the result without error correction. Figure 9 shows the average, and maximum and minimum accuracy for both methods. The second error correction performed in the algorithm improves the average accuracy of G-MLKM by approximately 11% (G-MLKM w/ EC 92.2%; G-MLKM w/o EC 81%).

Matrix Output of the G-MLKM Algorithm
The proposed G-MLKM algorithm implements multiple (determined by the structure of road networks) permutation matrices G i nts and L symmetric matrices G R i to represent the data cluster classification results at intersections and road segments, respectively. A detail example is illustrated to show the use of the proposed G-MLKM matrix output.
Take target N 1 as an example, it travels through road segment R 1 , R 2 , then heads to road segment R 5 . After that, it keeps on moving through road segment R 4 , R 2 and finally leaves the road network through road segment R 3 . The connections among associated data groups in each road segment that are related to target N 1 is represented as {1 1 , 1 2 , 6 2 , 1 5 , 3 4 , 5 3 }, which means data group 1 in road segment R 1 , data groups 1 and 6 in road segment R 2 , data group 1 in road segment R 5 , data group 3 in road segment R 4 , and data group 5 in road segment R 3 all belong to the measurements extracted from target N 1 . Table 4. Ground truth for five targets trajectories. The representation of A i denotes the associated data group A in road segment R i .
Let the intersection formed by road segments R 6 , R 5 , and R 4 be denoted as a. The incoming dataset Q a I ∈ R 3×1 can be stored in the sequence of {1 5 , 1 6 , 2 6 } and the outgoing dataset Q a O ∈ R 3×1 can be stored in the sequence of {1 4 , 2 4 , 3 4 }. Therefore, the permutation matrix G a may be determined as: Similarly, for the intersection formed by road segment R 1 , R 2 , and R 4 (named as intersection b), matrix G b may be determined as: O ∈ R 6×1 in the sequence of {1 2 , 2 2 , 3 2 , 4 2 , 5 2 , 6 2 }. For the intersection formed by road segment R 2 , R 3 , and R 5 (named as intersection c), G c may be determined as: , with Q c I ∈ R 6×1 stored in the sequence of {1 2 , 2 2 , 3 2 , 4 2 , 5 2 , 6 2 } and the outgoing dataset Q c O ∈ R 6×1 in the sequence of {1 3 , 2 3 , 3 3 , 4 3 , 5 3 ,1 5 }. With these matrices determined, the output result from G-MLKM can be clearly presented.

Conclusions and Future Work
This paper studied data pattern recognition for multi-targets in a constrained space, where the data were the minimal information provided by spatially distributed sensors. In contrast to the existing methods that rely on probabilistic hypothesis estimation, we proposed to utilize the machine learning approach for the data correlation analysis. Two common data clustering algorithms, namely, k-means++ and deep neural network, were first analyzed for data association given a simplified constrained space. Then the MLKM method was proposed via leveraging the structure advantage of DNN and the unsupervised clustering capability of k-means++. After that, graph theory was introduced in the purpose of extending the scope of MLKM for a general constrained space. In particular, we proposed a p-dual graph for data association at intersections and merged the results from local spaces and intersections through the dual graph of the constrained space. Simulation studies were provided to demonstrate the performance of the MLKM method and the proposed G-MLKM. Our future work will focus on releasing the assumptions in this paper to improve G-MLKM in the scenarios of false alarms.
Some interesting future work includes experimental verification of the proposed new approach in real-world environments and the consideration of constraints such as packet dropout, communication limitations, and other quality of service (QoS) parameters in sensor networks.