Dynamic multiple-graph spatial-temporal synchronous aggregation framework for traffic prediction in intelligent transportation systems

Graph neural network

Compared to CNN, graph neural network (GNN) can handle non-Euclidean data, which has led to its wide application in many fields such as feature extraction, node classification, etc. The classifications of GNN are categorized into two types: spectral domain and spatial domain. The former, also known as GCN, has undergone three important developments. Bruna et al. (2014) implemented graph convolution operations by replacing the convolution kernel with a learnable diagonal matrix. To overcome the problem of excessive computation, Defferrard, Bresson & Vandergheynst (2016) introduced Chebyshev polynomials to approximate the convolutional kernel. Kipf & Welling (2017) further reduced Chebyshev polynomials to the first order, which greatly simplified the computation and obtained the most common GCN expressions. The latter aims to define GNN by iteratively updating the representation of nodes from spatial neighbor aggregation. In order to avoid excessive nodes participating in the computation, graph sample and aggregate (GraphSAGE) et al. (Hamilton, Ying & Leskovec, 2017) limited the number of neighboring nodes by sampling and then achieved information aggregation through pooling operations. Atwood & Towsley (2016) defined the weights of neighbors using the K-hop transfer probabilities obtained after a random walk. Graph attention network (GAT) (Velivčkovič et al., 2018) employed the attention mechanism to define the weights for various neighbors, which allowed for a more flexible characterization of the aggregation in different scenarios.

Traffic prediction

In recent years, neural network has been extensively utilized for traffic prediction, which offers superior modeling performance compared to traditional methods (Guo et al., 2021). Most studies utilized graph convolution network to capture spatial feature (Li et al., 2022; Wang et al., 2022; Zhu et al., 2022; Yao et al., 2023; Yang et al., 2022). GraphSAGE has also been used to model spatial dependency for inductive learning (Liu et al., 2023; Liu, Ong & Chen, 2022). Temporal modules based on RNN and its LSTM, as well as GRU, have been introduced to learn temporal dependence (Pan et al., 2022; Subramaniyan et al., 2023; Bao et al., 2022; Shu, Cai & Xiong, 2022; Wan et al., 2022). To improve computational efficiency, some studies employed CNN instead of RNN to model temporal correlation (Ji, Yu & Lei, 2023; Zhang et al., 2022). Li et al. (2018) designed an encoder-decoder architecture that employed a diffusion process characterized by a bidirectional walk of a graph to learn spatial dependency and proposed the diffusion convolutional gated recurrent unit to model temporal dynamics. Based on this, Wu et al. (2019) further introduced an adaptive adjacency matrix in diffusion convolution to discover hidden spatial features and then employed dilated stacked 1D convolutions with larger receptive fields to capture temporal trends. Yu, Yin & Zhu (2018) defined the problem on a graph and then modeled spatial correlation by utilizing graph convolution with Chebyshev polynomials approximation, and took advantage of gated CNN to extract temporal features. On this foundation, Guo et al. (2019) incorporated GCN and attention mechanisms to enhance the representation of features, then deployed three parallel sets of components to learn different temporal trends. T-GCN (Zhao et al., 2020) integrated GCN with GRU as a way to learn the complex dynamics of traffic features. Graph multi-attention network (GMAN) (Zheng et al., 2020) adopted exclusively attention mechanisms rather than convolutional network to transform input features into predictions, which had a high computational complexity but improved the prediction accuracy.

To achieve synchronous modeling of spatial-temporal correlations, STSGCN (Song et al., 2020) constructed localized synchronous graphs, which aggregated information from neighbor nodes at the current time step and from themselves at adjacent time steps. Based on this, STFGNN (Li & Zhu, 2021) further fused spatial graphs and temporal graphs to learn the hidden correlations simultaneously. Auto-DSTSGN (Jin et al., 2022) designed an automated dilated spatial-temporal synchronous graph module to extract the short-range and long-range correlations by stacked layers with dilated factors. STGSA (Wei et al., 2023) proposed a specialized graph aggregation to capture spatial-temporal dynamics.

To summarize, most of the studies adopted separate components to learn spatial and temporal correlations, failing to achieve synchronous modeling. A few works implemented simultaneous aggregation of spatial and temporal features, but their models only took into account the hidden dependencies of traffic data itself, failing to model the impacts of external factors. In addition, their spatial-temporal synchronous graphs represented only a single correlation in terms of space and time, with a failure to model spatial and temporal dynamics from multiple perspectives, resulting in the overlook of certain dependencies. To address these shortcomings in spatial-temporal synchronous modeling, this article proposes a novel approach that learns the influence of external factors by fusing them with traffic features through attention mechanisms, and constructs diverse spatial-temporal synchronous graphs to extract spatial-temporal features in multiple ways.

Feature augmentation module

Traffic features are affected by external factors such as weather, temperature, wind, etc. In order to integrate external factors with traffic features for accurate traffic prediction, we design FAW, whose structure is shown in Fig. 2. The comfort of the temperature is closely related to people’s willingness to travel, and the traffic flow in the road network shows varying intensities under different weather conditions (e.g., sunny, rainy). In addition, wind also affects people’s travel plans. Therefore, we take account of the effects of four external factors that have the most impact on traffic, namely, maximum temperature, minimum temperature, weather, and wind. The external factors matrix is denoted as $E\in {\mathbb{R}}^{T\times N\times 4}$. The historical traffic feature matrix is denoted as $X_H\in {\mathbb{R}}^{T\times N\times C}$, where C is the dimension of traffic features. To better illustrate the point in spatial and temporal attention, we define a fully connected layer as:

(4) $$f\left(x\right)=ReLU\left(xW+b\right),$$where W and $b$ are learnable parameters, and $ReLU$ is the activation function.

$S\in {\mathbb{R}}^{T\times N\times N}$ represents the spatial attention matrix of N nodes in T time steps. We consider both traffic features and external factors to measure the spatial attention of different nodes. To be specific, we concatenate traffic features with external factors, and calculate the spatial correlation coefficient between node $v_i$ and $v$ at time step $t_j$ based on the scaled dot-product, which can be formulated as:

(5) $${\alpha }_{v_{i,v}}^{t_j}=\frac{⟨f_{s,1}(x_{v_i,t_j}\parallel x_{v_i,t_j}),f_{s,2}(x_{v_i,t_j}\parallel x_{v_i,t_j})⟩}{\sqrt{C+4}}$$where $⟨\cdot ,\cdot ⟩$ represents the inner product, $\parallel$ denotes the concatenation, $f_{s,1}(\cdot )$ and $f_{s,2}(\cdot )$ represent two different fully connected layers respectively, $x_{v_i,t_j},x_{v,t_j}\in X_H$, and $e_{v_i,t_j},e_{v,t_j}\in E$. The spatial attention coefficient is then obtained by normalizing ${\alpha }_{v_i,v}^{t_j}$ via softmax:

(6) $$s_{v_i,v}^{t_j}=\frac{\mathrm{e}\mathrm{x}\mathrm{p}\left({\alpha }_{{\mathrm{v}}_{\mathrm{i}},\mathrm{v}}^{{\mathrm{t}}_{\mathrm{j}}}\right)}{{\sum }_{\mathrm{v}\in \mathrm{V}}exp\left({\alpha }_{v_i,v}^{t_j}\right)}.$$

$U\in {\mathbb{R}}^{N\times T\times T}$ denotes the temporal attention matrix of T time steps among N nodes. We take into consideration both traffic features and external factors to learn the temporal attention of various time steps. Specifically, we first concatenate traffic features with external factors, and then employ the scaled dot-product to compute the temporal correlation coefficient. After that, the temporal attention coefficient is obtained by softmax, which can be described as:

(7) $$\begin{gathered} {\beta }_{t_j,t}^{v_i}=\frac{⟨f_{u,1}(x_{v_i,t_j}\parallel e_{v_i,t_j}),f_{u,2}(x_{v_i,t}\parallel e_{v_i,t})⟩}{\sqrt{C+4}}, \\ {u}_{t_j,t}^{v_i}=\frac{\exp ({\beta }_{t_j,t}^{v_i})}{{\sum }_{t\in T}\exp ({\beta }_{t_j,t}^{v_i})}, \end{gathered}$$where ${\beta }_{t_j,t}^{v_i}$ denotes the temporal correlation coefficient of node $v_i$ between time step $t_j$ and $t$, $u_{t_j,t}^{v_i}$ is the temporal attention coefficient, $f_{u,1}(\cdot ),f_{u,2}(\cdot )$ represent two independent fully connected layers respectively.

After obtaining the spatial attention matrix S and the temporal attention matrix U, traffic features and external factors are fused and mapped to the high-dimensional space, and the mathematical expression can be defined as:

(8) $$X=f_x\left({(U\cdot (S\cdot (X_H\parallel E){)}^T)}^T\right)$$where $X\in {\mathbb{R}}^{T\times N\times D}$ is the fused feature.

Dynamic spatial-temporal synchronous graph construction

Traffic features are related in multiple ways over time and space. On the one hand, nodes in different locations have adjacency and distance relationships. On the other hand, different traffic features may have the same temporal pattern or trend. A single spatial or temporal graph can only focus on one aspect of spatial-temporal dependencies while ignoring others. Aiming to model spatial-temporal correlations from multiple perspectives, we first propose two spatial graphs, i.e., the neighbor graph $G^{Ne}$ and the distance graph $G^{Di}$, as well as two temporal graphs, i.e., the trend graph $G^{Tr}$ and the pattern graph $G^{Pa}$, and then design two dynamic spatial-temporal synchronous graphs based on them, i.e., the trend spatial-temporal synchronous graph $G^{TN}$ and the pattern spatial-temporal synchronous graph $G^{PD}$.

Neighbor graph. $G^{Ne}$ denotes the spatial adjacency of nodes, and the elements of its adjacency matrix $A^{Ne}$ are defined as:

(9) $$A_{ij}^{Ne}=\{\begin{array}{cc}1,\hfill & \mathrm{i}\mathrm{f}{v}_i\mathrm{c}\mathrm{o}\mathrm{n}\mathrm{n}\mathrm{e}\mathrm{c}\mathrm{t}\mathrm{s}\mathrm{t}\mathrm{o}{v}_j,\hfill \\ 0,\hfill & \mathrm{e}\mathrm{l}\mathrm{s}\mathrm{e}.\hfill \end{array}$$

Distance graph. $G^{Di}$ represents the distance relationship between nodes, the elements of $A^{Di}$ are fomulated as:

(10) $$A_{ij}^{Di}=\{\begin{array}{cc}\mathrm{e}\mathrm{x}\mathrm{p}\left(-\frac{d_{ij}^2}{{\sigma }^2}\right),\hfill & \mathrm{i}\mathrm{f}{d}_{ij}\le {\epsilon }_{\mathrm{D}},\hfill \\ 0,\hfill & \mathrm{e}\mathrm{l}\mathrm{s}\mathrm{e},\hfill \end{array}$$where $d_{ij}$ denotes the Euclidean distance between node $v_i$ and $v_j$, $\sigma$ denotes the standard deviation of Euclidean distance, ${\epsilon }_D$ is hyperparameters used to control the sparsity of $A^{Di}$.

Trend graph. In order to capture the trend similarity between traffic features, we propose the trend graph $G^{Tr}$ based on dynamic time warping (DTW) algorithm. Given two sequences $P=\left[p_1,p_2,\dots ,p_m\right]$ and $Q=\left[q_1,q_2,\dots ,q_n\right]$, the element $M\left(i,j\right)$ of distance matrix $M_{m\times n}$ can be calculated by $|p_i-q_j|$ (Li & Zhu, 2021), and the DTW distance between P and Q can be obtained by the iteration of the following equations:

(11) $$\begin{gathered} M_D\left(i,j\right)=M\left(i,j\right)+M_{min}\left(i,j\right), \\ {M}_{min}\left(i,j\right)=min\left(M_D\left(i-1,j\right),M_D\left(i,j-1\right),M_D\left(i-1,j-1\right)\right). \end{gathered}$$

The adjacency matrix $A^{Tr}$ can be formulated as:

(12) $$A_{ij}^{Tr}=\{\begin{array}{cc}1,\hfill & \mathrm{i}\mathrm{f}{M}_D\left(i,j\right)\le {\epsilon }_T,\hfill \\ 0,\hfill & \mathrm{e}\mathrm{l}\mathrm{s}\mathrm{e},\hfill \end{array}$$where $M_D\left(i,j\right)$ come from (11), ${\epsilon }_T$ denotes the hyperparameter to control the sparsity of $A^{Tr}$.

Pattern graph. Some traffic nodes exhibit strong linear correlation in their features because they are located in the same functional area (e.g., residential area, commercial area, etc.), as shown in Fig. 1B. For the purpose of extracting pattern similarity, we design the pattern graph $G^{Pa}$ by Pearson correlation coefficient. The traffic feature series of nodes $v_i$ and $v_j$ are denoted as $X_i=\left[x_i^1,x_i^2,\dots ,x_i^T\right],X_j=\left[x_j^1,x_j^2,\dots ,x_j^T\right]$, respectively. The Pearson correlation coefficient between $v_i$ and $v_j$ can be defined as:

(13) $${\rho }_{ij}=\frac{{\sum }_{t=1}^T\left(x_i^t-\overline{x_i}\right)\left(x_j^t-\overline{x_j}\right)}{\sqrt{{\sum }_{t=1}^T{\left(x_i^t-\overline{x_i}\right)}^2}\sqrt{{\sum }_{t=1}^T{\left(x_j^t-\overline{x_j}\right)}^2}},$$and the adjacency matrix $A^{Pa}$ of $G^{Pa}$ can be formulated as:

(14) $$A_{ij}^{Pa}=\{\begin{array}{l}{\rho }_{ij},\text{\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}if\hspace{0.17em}\hspace{0.17em}}{\rho }_{ij}\ge {\epsilon }_P,\\ 0,\text{\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}else},\end{array}.$$where ${\epsilon }_P$ denotes the hyperparameter to manipulate the sparsity of $A^{Pa}$.

Dynamic spatial-temporal synchronous graphs. The spatial and temporal correlations of traffic features exist simultaneously, and to model this intricate spatial-temporal dependencies synchronously from multiple perspectives, we propose the trend spatial-temporal synchronous graphs $G^{TN}$ and the pattern spatial-temporal synchronous graphs $G^{PD}$ inspired by STSGCN (Song et al., 2020) and STFGNN (Li & Zhu, 2021). First, we design two predefined spatial-temporal synchronous graphs $G_{Pre}^{TN}$ as well as $G_{Pre}^{PD}$, each of which contains two time steps, and the adjacency matrices are denoted as $A_{Pre}^{TN}\in {\mathbb{R}}^{2N\times 2N}$ and $A_{Pre}^{PD}\in {\mathbb{R}}^{2N\times 2N}$, respectively. The structure of $A_{Pre}^{TN}$ and $A_{Pre}^{PD}$ is illustrated in Fig. 3A. The main diagonals are $\left[A^{Tr},A^{Tr}\right]$ and $\left[A^{Pa},A^{Pa}\right]$, separately, denoting that each node has connectivity to nodes with the same trend or pattern at the same time step, while the counter-diagonal are $\left[A^{Ne},A^{Ne}\right]$ and $\left[A^{Di},A^{Di}\right]$, respectively, indicating that nodes are connected to their neighboring nodes or proximity nodes at the adjacent time step. Then, the two predefined adjacency matrices are multiplied by the learnable parameters of the same shape to obtain the adjacency matrices of the dynamic spatial-temporal synchronous graphs as follows:

(15) $$A^{TN}=W^{TN}\otimes A_{Pre}^{TN}\in {\mathbb{R}}^{2N\times 2N},A^{PD}=W^{PD}\otimes A_{Pre}^{PD}\in {\mathbb{R}}^{2N\times 2N},$$where $A^{TN},A^{PD}$ denote the adjacency matrices of $G^{TN}$ and $G^{PD}$ respectively, $\otimes$ is element-wise product, $W^{TN},W^{PD}$ are learnable parameters.

Spatial-temporal synchronous aggregation module

Aiming to learn the hidden spatial and temporal correlations synchronously, we build spatial-temporal synchronous aggregation modules corresponding to $A^{TN}$ and $A^{PD}$, that is, the trend spatial-temporal synchronous aggregation module (T-STSAM) and the pattern spatial-temporal synchronous aggregation module (P-STSAM), which have the same architecture but with different adjacency matrices, as shown in Fig. 2. Multiple gated graph convolutions are stacked in each STSAM to learn spatial-temporal correlations simultaneously, and each gated graph convolution can be formulated as:

(16) $$GGC\left(h\right)=tanh\left(h{W}_1^h+b_1^h\right)\otimes \sigma \left(h{W}_2^h+b_2^h\right),$$where $h$ denotes the input, $tanh(\cdot )$ denotes the tanh activation function, $\sigma (\cdot )$ denotes the sigmoid activation function, and $W_1^h,W_2^h\in {\mathbb{R}}^{D\times D},b_1^h,b_2^h\in {\mathbb{R}}^D$ are learnable parameters. The output of the previous gated graph convolution is used as the input of the next one. In addition, max pooling and a cropping operation are also included in a STSAM.

The gated graph convolution in T-STSAM allows simultaneous aggregation of information from nodes with similar temporal trends and nodes with spatial adjacencies, the mathematical equation can be formulated as:

(17) $$h_m^{TN}=GGC\left(A^{TN}{h}_{m-1}^{TN}\right),$$where $h_{m-1}^{TN},h_m^{TN}\in {\mathbb{R}}^{2N\times D}$ denote the input and output of the $m$-th gated graph convolution in T-STSAM, respectively.

The gated graph convolution in P-STSAM implements synchronous modeling of dependencies from nodes with the same temporal pattern and nodes with short spatial distances, and the computational formula is defined as:

(18) $$h_m^{PD}=GGC\left(A^{PD}{h}_{m-1}^{PD}\right),$$where $h_{m-1}^{PD},h_m^{PD}\in {\mathbb{R}}^{2N\times D}$ denote the input and output of the $m$-th gated graph convolution in P-STSAM, respectively.

In each STSAM, we take advantage of jump knowledge network (JK-Net) to aggregate the outputs of various gated graph convolutions and then retain the most potent representation by max pooling, which can be described as:

(19) $$\eqalign{ & h_{max}^{TN} = MaxPooling\left( {h_1^{TN},h_2^{TN}, \ldots ,h_M^{TN}} \right), \cr & h_{max}^{PD} = MaxPooling\left( {h_1^{PD},h_2^{PD}, \ldots ,h_M^{PD}} \right), \cr} $$$$\begin{array}{l}{h}_{max}^{TN}=MaxPooling\left(h_1^{TN},h_2^{TN},\dots ,h_M^{TN}\right),\\ h_{max}^{PD}=MaxPooling\left(h_1^{PD},h_2^{PD},\dots ,h_M^{PD}\right),\end{array}$$where $h_{max}^{TN},h_{max}^{PD}\in {\mathbb{R}}^{2N\times D}$ denote the outputs of max pooling in T-STSAM and P-STSAM, respectively, M denotes the number of gated graph convolutions.

The cropping operation is performed after the max pooling to reduce the hidden features from $2N$ dimensions to N dimensions, leading to the output of STSAM, which can be expressed as:

(20) $$\begin{array}{l}{O}^{TN}=h_{max}^{TN}\left[N:2N,:\right]\in {ℝ}^{N\times D},\\ O^{PD}=h_{max}^{PD}\left[N:2N,:\right]\in {ℝ}^{N\times D},\end{array}$$where $O^{TN},O^{PD}$ denote the outputs of T-STSAM and P-STSAM, respectively.

Dilated gated convolution module

Regarding time dimension, parallel STSAMs are conducive to extracting short-term dependencies due to their independent parameters while still lacking in modeling long-term correlations. To address this limitation, gated 1D convolution with shared parameters is utilized to extract long-term temporal features, which can be formulated as:

(21) $$O^{GC}=tanh({\mathrm{\Theta }}_1\ast X^{l-1}+b_1)\otimes \sigma ({\mathrm{\Theta }}_2\ast X^{l-1}+b_2),$$where $O^{GC}$ denotes the output of gated 1D convolution, $X^{l-1}$ denotes the input of $l$-th layer, ${\mathrm{\Theta }}_1,{\mathrm{\Theta }}_2$ denote two 1D convolution respectively, and $b_1,b_2$ are learnable parameters.

Inspired by Graph WaveNet (Wu et al., 2019), dilated instead of fixed step sizes are introduced to ${\mathrm{\Theta }}_1$ and ${\mathrm{\Theta }}_2$ for expanding the receptive fields of 1D convolutions. However, the three-layer convolutions with step sizes of [1,2,4] in Graph WaveNet only cover eight historical time steps, while the length of the input sequences in our model is 12. Therefore, we deploy four layers in which the step sizes of 1D convolutions are set as [1,2,4,4] respectively, while the kernel size of 1D convolution in each layer is fixed as 2, as shown in Fig. 3B. With this set of mechanisms, the receptive field of four-layer 1D convolutions can be expanded to the length of the input time series in our model.

Dynamic multiple-graph spatial-temporal synchronous aggregation layer

T-STSAMs and P-STSAMs with gated fusion, as well as a DGCM, compose a DMSTSAL, four DMSTSALs are stacked in a DMSTSAF, and the output of the previous layer is used as the input to the next layer, which is presented in Fig. 2. In each DMSTSAL, we first slide from the input sequence to obtain the time step pairs, and then construct two kinds of dynamic spatial-temporal synchronous graphs for every time step pair. In order to expand the receptive field and reduce the number of layers, the distances between two time steps in four layers are set to be dilated with the same step sizes as in DGCMs, which are [1,2,4,4], as shown in Fig. 3B.

Denoting the input of the $l$-th DMSTSAL as $X^{l-1}=\left[x_1,x_2,\dots ,x_{T^l}\right]$ and the distance of the two time steps as $d^l\in \left[1,2,4,4\right]$, the time step pairs generated by its sliding in the $l$-th layer can be described as $\left[\left(x_1,x_{1+d^l}\right),\left(x_2,x_{2+d^l}\right),\dots ,\left(x_{T^l-d^l},x_{T^l}\right)\right]$, then the number of time step pairs can be formulated as:

(22) $$M^l=T^l-d^l.$$

To achieve simultaneous learning of spatial-temporal dependencies from multiple perspectives, we construct a trend spatial-temporal synchronous graph $G^{TN}$ and a pattern spatial-temporal synchronization graph $G^{PD}$ for each time step pair. The number of each type of dynamic spatial-temporal graph is the same as the number of time step pairs, which can be obtained from (22).

Traffic features are heterogeneous across time, and to capture hidden features more accurately, we design parallel rather than shared STSMs to extract spatial-temporal dependencies. Specifically, a T-STSAM is deployed for each $G^{TN}$, and a P-STSAM is allocated for each $G^{PD}$. Thus, $M^l$ T-STSAMs and $M^l$ P-STSAMs are laid out in the $l$-th DMSTSL.

Both T-STSAMs and P-STSAMs are able to represent spatial-temporal correlations in some way, so it is necessary to fuse the outputs of two types of STSAMs. To achieve this goal, we start by aggregating the outputs of each kind of STSAM into a sequence, which can be defined as:

(23) $$\begin{gathered} O_{AG}^{TN}=\left[O_1^{TN},O_2^{TN},\dots ,O_{M^l}^{TN}\right]\in {\mathbb{R}}^{M^l\times N\times D}, \\ {O}_{AG}^{PD}=\left[O_1^{PD},O_2^{PD},\dots ,O_{M^l}^{PD}\right]\in {\mathbb{R}}^{M^l\times N\times D}, \end{gathered}$$where $O_{AG}^{TN},O_{AG}^{PD}$ denote the outputs of T-STSAMs and P-STSAMs, respectively, they represent the results of spatial-temporal synchronous aggregation, which is the core component that forms the prediction of our model. Their elements $O_m^{TN},O_m^{PD},1\le m\le M^l$ can be obtained from (20). Then, the gating mechanism is applied to incorporate the outputs of two concatenation operations, which can be formulated as:

(24) $$O^F=z\otimes O_{AG}^{TN}+\left(1-z\right)\otimes O_{AG}^{PD}\in {\mathbb{R}}^{M^l\times N\times D},$$where $O^F$ is one of the two parts that compose the output of the current layer, $z$ is employed to manipulate the proportion of information in the fusion, whose mathematical equation can be defined as:

(25) $$z=\sigma \left(O_{AG}^{TN}{W}_1^z+O_{AG}^{PD}{W}_2^z+b^z\right),$$where $W_1^z,W_2^z,b^z$ are learnable parameters. Another part of the current layer’s outcome is the result $O^{GC}$ of the DGCM, which can be obtained from (21). The outcome $O^F$ of the fusion operation and the output $O^{GC}$ of the DGCM are added together to form the result of the $l$-th DMSTSAL, which can be formulated as:

(26) $$X^l=O^F+O^{GC}\in {\mathbb{R}}^{M^l\times N\times D},$$

The results of the four DMSTSALs are the inputs to the output module.

Output module

The output module is responsible for generating the final predictions of our model, in which a concatenation that aggregates the outputs of four layers is first employed to capture comprehensive spatial-temporal correlations, which can be defined as:

(27) $$X_{CAT}=(X^1\parallel X^2\parallel X^3\parallel X^4)\in {\mathbb{R}}^{(M^1,M^2,M^3,M^4)\times N\times D},$$where $X^l,1\le l\le 4$ can be obtained from (26). A series of fully connected layers are then utilized to produce the final predictions of $T^{\mathrm{\prime }}$ time steps. Specifically, we design two-fully-connected-layers to produce the prediction of time step $t$, which can be formulated as:

(28) $${\hat{y}}_t=ReLU\left(\left(X_{CAT}{W}_1^t+b_1^t\right)W_2^t+b_2^t\right)\in {\mathbb{R}}^{1\times N\times C},$$where $W_1^t,b_1^t,W_2^t,b_2^t$ are learnable parameters. Finally, the results of (28) repeated $T^{\mathrm{\prime }}$ times compose the final predictions of our model, which can be described as:

(29) $$\hat{Y}=\left[{\hat{y}}_1,{\hat{y}}_2,\dots ,{\hat{y}}_{T^{\mathrm{\prime }}}\right]\in R^{T^{\mathrm{\prime }}\times N\times C}.$$

Smooth L1 loss rather than L1 loss is chosen as the loss function, which deals with the unsmooth disadvantage at the zero-point and can be defined as:

(30) $$L(Y,\hat{Y})=\{\begin{array}{cc}\frac{1}{2}{\left(Y-\hat{Y}\right)}^2/\delta ,\hfill & if\left|Y-\hat{Y}\right|\le \delta ,\hfill \\ \left|Y-\hat{Y}\right|-\frac{1}{2}\delta ,\hfill & otherwise\hfill \end{array},$$where $Y,\hat{Y}$ denote the true value and the predicted value, respectively, and $\delta$ denotes a threshold parameter to determine the sensitivity.

Datasets

Four public datasets are chosen for evaluating the prediction performance of DMSTSAF: PEMS03, PEMS04, PEMS07, and PEMS08. All datasets are generated by Caltrans Performance Measurement System (PeMS). Concretely, traffic features in these datasets are collected by sensors located on California highways at 5-min intervals, which means there are 12 time steps in an hour. Spatial graphs are constructed from the distribution of sensors. The distinctions among the datasets are the geographic locations of sensors and the temporal ranges of data. The detailed information of four datasets can be found in Table 2.

Experiment settings and evaluation metrics

In our experimental implementation, each dataset is partitioned for training, validation, and testing in a ratio of 60%, 20%, and 20%. We employ the traffic flow of 12 historical time steps (1 h) to predict the traffic flow of 12 future time steps. The traffic flow in each dataset is standardized using Z-score normalization.

DMSTSAF is performed by Pytorch using a PC with NVIDIA RTX 3080. We chose hyperparameters of our model experimentally to ensure superior performance. The sparsity of the adjacency matrices for the distance graph, the trend graph, and the pattern graph are set as 0.01. Each STSAM consists of two graph convolutions, and the output dimension of the fully connected layer in FAW and the hidden dimensions of graph convolutions in STSAM are set as $D=64$. The hidden representations of two fully-connected layers to generate predictions are tuned as [128,1], respectively. The optimizer in experiments is set as Adam. The batch sizes for four datasets are tuned as [32,32,8,64], respectively. The initial learning rate is set as 0.003 and scaled down to 0.3 times every 30 epochs to shorten the training time. Enroll up to 200 epochs in each training.

We take three metrics to evaluate the performance of models, which are mean absolute error (MAE), mean absolute percentage error (MAPE), and root mean square error (RMSE). Smaller values indicate better performance for these metrics, and their mathematical equation can be formulated as:

(31) $$\begin{gathered} MAE\left(y_i,\widehat{y_i}\right)=\frac{1}{N\times T^{\mathrm{\prime }}}{\sum }_{i=1}^{N\times T^{\mathrm{\prime }}}\left|y_i-\widehat{y_i}\right|, \\ MAPE\left(y_i,\widehat{y_i}\right)=\frac{1}{N\times T^{\mathrm{\prime }}}{\sum }_{i=1}^{N\times T^{\mathrm{\prime }}}\frac{\left|y_i-\widehat{y_i}\right|}{y_i}, \\ RMSE\left(y_i,\widehat{y_i}\right)=\sqrt{\frac{1}{N\times T^{\mathrm{\prime }}}{\sum }_{i=1}^{N\times T^{\mathrm{\prime }}}{\left(y_i-\widehat{y_i}\right)}^2}, \end{gathered}$$where $y_i$, $\widehat{y_i}$ denote the true value and the predicted value, respectively.

Baseline methods

We compare DMSTSAF with the following seven models for traffic flow prediction:

• FC-LSTM (Sutskever, Vinyals & Le, 2014): Long short-term memory network with fully connected layers, a variant of the recurrent neural network, which consists of the forget gate, the input gate, and the output gate.

• GRU (Fu, Zhang & Li, 2016): Gate recurrent unit, a variant of LSTM, which simplifies the three gates in LSTM to two, i.e., the reset gate and the update gate.

• T-GCN (Zhao et al., 2020): Temporal graph convolutional network, which incorporates GCN with GRU to extract spatial-temporal features.

• DCRNN (Li et al., 2018): Diffusion convolution recurrent neural network, in which random walks of a graph are utilized to model spatial dependence, GRU is employed to learn temporal correlation.

• STGCN (Yu, Yin & Zhu, 2018): Spatio-temporal graph convolution network, which takes advantage of GCN to extract spatial features, and makes use of 1D convolution rather than recurrent neural network to model temporal dependency.

• STSGCN (Song et al., 2020): Spatial-temporal graph convolutional network, which constructs localized spatial-temporal graphs based on spatial and temporal adjacencies, and designs spatial-temporal synchronous graph convolutional modules to learn spatial-temporal correlations synchronously.

• STFGNN (Li & Zhu, 2021): Spatial-temporal fusion graph neural network, an improved study based on STSGCN, which designs temporal graphs with dynamic time warping algorithm and proposes spatial-temporal fusion graphs to model localized spatial-temporal dependencies, and then designs gated CNN to extract long-range dependencies.

Experimental results

The experiments of traffic flow prediction are performed on four datasets, and a comparison of the results for all models is shown in Table 3. DMSTSAF obtains the smallest metrics on all four datasets, demonstrating the consistent superiority of DMSTSAF over the baselines. Compared to the state-of-the-art baseline STFGNN, DMSTSAF improves 4.27%, 3.72%, 8.54%, and 7.89% in terms of MAE on PEMS03, PEMS04, PEMS07, and PEMS08 respectively, while the improvements of MAPE are 6.63%, 7.49%, 7.95%, and 7.85%. In addition, our model also achieves 3.68%, 4.23%, 8.38%, and 7.65% improvements in terms of RMSE.

Moreover, for the purpose of evaluating the ability of models on multi-step prediction, we conduct prediction experiments for 12 future time steps on PEMS03. In comparison to STFGNN, our model shows 1.11–8.89% improvement in terms of MAE, 1.81–13.19% improvement for MAPE, and 1.38–7.48% improvement for RMSE. To achieve a more intuitive comparison, we illustrate the results of each model with a line plot, as shown in Fig. 4, which verifies that DMSTSAF overwhelmingly outperforms all the baselines.

To further demonstrate the prediction performances of DMSTSAF and the best baselines, we visualize the ground-truth and predictions of STSGCN, STFGNN, and DMSTSAF on PEMS03 test set for both sunny and rainy days, as presented in Fig. 5. It can be seen that compared to the best baselines, our proposed DMSTSAF can fit the ground-truth more accurately in different weather conditions.

Ablation study

To evaluate the validity of different modules in DMSTSAF, ablation studies are implemented on PEMS04 and PEMS08. We design six variants, a short illustration is introduced as follows:

• -FAW, which uses only traffic features as the input to the model, and utilizes an FC rather than a feature augmentation module to transform the input from C dimensions to D dimensions.

• - $G^{TN}$, which removes trend spatial-temporal synchronous graphs $G^{TN}$ and the corresponding T-STSAM, and retains only pattern spatial-temporal synchronous graphs $G^{PD}$ as well as P-STSAM.

• - $G^{PD}$, which is the opposite of - $G^{TN}$, removing $G^{PD}$ and P-STSAM, while retaining $G^{TN}$ as well as T-STSAM.

• -Dilation, the step sizes of 1D convolutions in DGCM and the distance between two time steps to construct dynamic spatial-temporal synchronous graphs are fixed to 1.

• -Concatenation, which removes the concatenation on the outputs of four DMSTSALs, and treats the output of the $4$-th DMSTSAL as the input to the output module.

• -DGCM, which removes dilated gated convolution modules from four DMSTSALs.

Table 4 demonstrates the results of ablation studies, which illustrates that DMSTSAF outperforms variants on PEMS04 as well as PEMS08. Compared to the three indicators of -FAW, our model achieves 1.92%, 2.60% as well as 1.67% improvement on PEMS04. Furthermore, it improves by 3.04%, 2.79%, and 1.39% on PEMS08. These improvements indicate the effectiveness of fusing traffic features and external factors in learning spatial-temporal correlations. In contrast with - $G^{TN}$, DMSTSAF achieves 1.02%, 1.53%, 1.01% improvements on PEMS04, as well as 1.39%, 2.23%, 0.88% improvements on PEMS08. In comparison to - $G^{PD}$, improvements of 0.6%, 0.94%, 0.54% on PEMS04 are obtained by our model, while improvements of 0.70%, 1.56%, 0.52% are achieved on PEMS08. These achievements verify the validity of dynamic multiple-graph in modeling spatial-temporal dependencies. Compared to -DGCM, our model achieves 1.87%, 2.68%, as well as 2.74% improvement on PEMS04. Furthermore, it improves 2.49%, 4.54%, and 1.67% on PEMS08. These improvements illustrate the positive role of DGCM in capturing long-term correlations. In addition, there are also different levels of improvement in our model compared to other variants, which demonstrates the efficiencies of the corresponding components.

Parameters study

In order to further validate our model, parameters study are implemented on PEMS04 and PEMS08. The number of graph convolutions and the hidden dimensions D in STSAM are under consideration, Fig. 6 presents the results of experiments. Our model achieves minimal errors as each STSAM consists of two graph convolutions. As for the hidden dimensions, it can be found in Fig. 6 that our model obtains optimal metrics when D equals to 64.