A Dynamic Stiefel Graph Neural Network for Efficient Spatio-Temporal Time Series Forecasting

1 Introduction

In recent years, spatio-temporal time series (STTS) have become increasingly important in fields such as traffic flow forecasting, financial market analysis, and weather forecasting [Saad et al., 2024; Li et al., 2024b]. These data encompass changes and relationships in both time and space dimensions. The spatial relationships, in a broader sense, not only represent geographical locations but also denote the relationships between other types of entities, such as companies or sensors. Accurately forecasting STTS is crucial for decision support and resource optimization in these areas.

In contrast to traditional time series forecasting, spatiotemporal forecasting faces the challenge of analyzing complex and dynamic spatio-temporal correlations among STTS data. The key problem is how to effectively model both the temporal dependencies and the spatial interactions between variables. Existing representative methods for multivariate time series forecasting, such as TimeMixer [Wang et al., 2024a] and CycleNet [Lin et al., 2024], have demonstrated excellent performance in capturing periodicity and trend changes in time series data. These methods utilize techniques such as decomposition to model complex temporal patterns effectively. However, they primarily focus on temporal dimensions and often overlook the spatial correlations in STTS, limiting their optimal performance in spatiotemporal forecasting. Graph Neural Networks (GNNs) [Kipf and Welling, 2016], on the other hand, have gained popularity in spatio-temporal analysis due to their ability to model spatial relationships through graph structures. Recent advances in this area include the development of extended message-passing mechanisms for simultaneous spatial and temporal feature extraction [Kong et al., 2024b], integrating spatio-temporal graph convolution networks for traffic prediction [Gupta et al., 2023], and geometric masking techniques for interpreting GNNs in spatio-temporal contexts [Yu et al., 2024]. These approaches fully utilize the advantages of GNNs in handling spatial interactions between nodes, thereby providing more accurate and comprehensive results for spatio-temporal forecasting.

∗Corresponding Author: Liang Xie

However, although GNNs have great potential in spatiotemporal forecasting, their application faces two major challenges. Firstly, existing graph-based methods are computationally expensive, limiting their performance on large-scale datasets. Secondly, many graph methods struggle to adapt to the dynamic changes in spatio-temporal data, which further increases the computational costs [Bastos et al., 2024]. Although some graph Fourier methods, such as FourierGNN [Yi et al., 2024a], have alleviated these issues to some extent, they still fall short in analyzing dynamic spatio-temporal data.

To address these challenges, this paper proposes the Dynamic Spatio-temporal Stiefel Graph Neural Network (DSTSGNN), which introduces the Stiefel Graph Spectral Convolution (SGSC) to effectively handle the dynamic changes and high-dimensional characteristics of spatio-temporal data. DST-SGNN first reduces the dimensionality and complexity of spatio-temporal data through patch-based sequence decomposition, breaking it into smaller, manageable patches to reduce redundancy and simplify the structure. Building on this, the model introduces Multi-layer Stiefel Graph Spectral Convolution (MSGSC) to effectively capture complex spatiotemporal correlations while maintaining high efficiency. To efficiently exploit dynamic graphs, Linear Dynamic Graph Optimization on Stiefel Manifold (LDGOSM) is proposed to efficiently learn the transformation matrix for SGSC.

The main advantages of DST-SGNN are as follows:

• DST-SGNN can dynamically model the evolving spatial correlations over time, thereby providing more accurate forecasting results.
• DST-SGNN achieves significant optimization in terms of time and space complexity compared to traditional graph neural networks, making it more efficient for processing large-scale spatio-temporal data.
• Extensive experiments on seven spatio-temporal datasets demonstrate DST-SGNN's superiority in accuracy and efficiency in comparison with state-of-the-art spatio-temporal and time series forecasting methods.

2 RelatedWork

2.1 Spatio-Temporal Time Series Forecasting

Multivariate time series forecasting is a key application in time series analysis, aiming to predict future values by modeling multiple correlated variables. Patch-based methods (such as PatchTST [Nie et al., 2022] and MTST [Zhang et al., 2024] ) have made significant progress by dividing time series into local segments and using the Transformer architecture to capture both local and global temporal dependencies. Meanwhile, methods based on time series decomposition, such as TimesNet [Wu et al., 2022] and TimeMixer [Wang et al., 2024a], have further improved the performance of longsequence forecasting by optimizing the decomposition process, and they better capture the intrinsic temporal structure of time series.

In the field of STTS, researchers have developed several innovative methods. STPGNN [Kong et al., 2024a] identifies pivotal nodes and captures their complex spatio-temporal dependencies to enhance traffic prediction. MiTSformer [Chen and Zhao, 2024] addresses mixed time series by recovering latent continuous variables and leveraging multi-scale context for balanced modeling. ModWaveMLP [Sun et al., 2024] offers an efficient and simple solution for traffic forecasting using MLPs combined with mode decomposition and wavelet denoising. Lastly, GPT-ST [Li et al., 2024b] introduces a spatio-temporal pre-training framework that enhances downstream models' learning capabilities through masked autoencoders and adaptive masking strategies. These methods, each with unique strengths, provide diverse solutions for spatiotemporal time series forecasting.

2.2 Graph Neural Networks for Spatio-Temporal Analysis

GNNs are widely applied in spatio-temporal prediction. Their powerful modeling capabilities can effectively capture the complex relationships within the data. In the context of spatio-temporal prediction, the graph structure needs to be dynamically updated according to the data to adapt to spatio-temporal changes [Li et al., 2024a; Ben et al., 2024]. Meanwhile, as demonstrated by LLGformer [Jin et al., 2025] and MPFGSN [Zheng et al., 2025], full spatiotemporal graphs can significantly improve the prediction performance and help uncover potential patterns. However, the dynamic update and the construction of full spatio-temporal graphs will incur significant computational overhead. Especially when dealing with large-scale dynamic data, this limits the widespread application of GNNs.

Graph spectral convolution methods perform convolution in the spectral domain based on the graph Fourier transform, which can effectively capture the complex relationships in STTS [Cao et al., 2020]. Nevertheless, traditional methods require eigenvalue decomposition of the adjacency matrix, resulting in high computational costs. Methods such as GCN [Kipf and Welling, 2016; Peng and Zhang, 2023] and ChebNet [Defferrard et al., 2016; Lei et al., 2024] replace eigenvalue decomposition with a spatial-like convolution that aggregates information via the adjacency matrix. However, the construction and dynamic update of the adjacency matrix remain computationally expensive. The FourierGNN [Yi et al., 2024a] uses the standard Fourier transform to improve computational efficiency, but it is only applicable to fixed graph structures and is not suitable for dynamic spatio-temporal relationships.

3 ProblemDefinition

The problem of STTS forecasting involves forecasting future values of multiple variables based on their past observations across both time and space dimensions. STTS is similar to multivariate time series, and the main difference between them is that there exist spatial correlations among variables of STTS. Specifically, given a time series dataset with N variables observed over T consecutive time points, the goal is to predict the values of these variables at future time points t + 1 to t + τ, where τ is the prediction length. The prediction function Fθ takes historical data Xt as input and outputs the predicted values Yt for the next τ time points, capturing complex patterns and dependencies in the STTS data.

Formally, the prediction problem can be defined as finding a function Fθ such that Yt = Fθ(Xt), where Yt represents the predicted multivariate values for future time points t + 1 to t + τ .

4 Methodology

This section details the DST-SGNN model methodology. Section 4.1 presents the overall framework, while Sections 4.2-4.5 cover module implementation, with Section 4.4 focusing on the core theory of the Stiefel graph neural network. Section 4.6 describes training and inference procedures.

4.1 Overall Framework

Figure 1 illustrates the overall framework of the DST-SGNN model, which mainly includes four modules: Patching, Decomposition, LDGOSM, and MSGSC. First, the Patching module decomposes the original high-dimensional spatiotemporal data into multiple smaller sub-sequences (patches)

to reduce data complexity. Next, the Decomposition module further decomposes each patch to extract seasonal and trend components, which are used as nodes to construct a hyperpatch graph. After constructing the hyperpatch graph, an initial transformation matrix is formed from the graph's data. Subsequently, the LDGOSM module dynamically optimizes the transformation matrix on the Stiefel manifold. Finally, the MSGSC module performs multi-step graph spectral convolutions, where HP denotes the element-wise Hadamard product operation in the spectral convoluation process. The seasonal and trend components are concatenated through the Concat module, and the resulting features are fed into the Linear Layers module to generate the final predictions. In the following, we will introduce each module in detail.

4.2 Patching and Decomposition

Patching

Given a multivariate time series window, the input is Xt ∈ RT×N of N variables at timestamp t. Let p denote the patch size, and s denote the stride (the length of the nonoverlapping region between two consecutive patches). Then, Xt is fed into the Patching module T : RT×N → RJ×p×N . This module transforms Xt into J = ⌊(T -p)/s⌋ + 2 patches. The specific calculation for the Patching module ˆXt = T (Xt) is as follows:

where ˆXt ∈ RJ×p×N represents the output of the Patching module, and ˆXt,j is its j-th column. To ensure that J is an integer, we pad the last s values of the sequence to the end of the original sequence.

Intra-patch Sequence Decomposition for Constructing the Hyperpatch Graph

First, the time series ˆXt is decomposed into seasonal component ˆX(1) t and trend component ˆX(2) t using the series decomposition block of Autoformer [Wu et al., 2021].

In the hyperpatch graph structure, each patch in the seasonal component ˆX(1) t and trend component ˆX(2) t is defined as a node in the graph. For the sliding window of the time series, each element within the window is considered to have a connection with every other element, thus it is represented in the form of a spatio-temporal fully connected graph. Based on ˆXm t ∈ RJ×N×p ,m = (1), (2), a hyperpatch graph structure Gmt = (XG,m t , AG,m t ) is further constructed. This structure is initialized as a fully connected graph containing (J ×N)×p nodes, where XG,m t ∈ R(J×N)×p represents the node features, and AG,m t ∈ R((J×N)×p)×((J×N)×p) is the adjacency matrix, used to describe the relationships between nodes.

4.3 Embedding and Normalization

In the DST-SGNN model, the Embedding module maps the data XG,m t ∈ R(J×N)×p after Decomposition into the hidden space ˜Xm t ∈ R(J×N)×K , where K is the dimension of the hidden layer. This process reduces the data dimensionality and enhances the feature representation capability. The Normalization module normalizes each time series instance to have zero mean and unit variance, thereby reducing the distribution differences between training and testing data [Kim et al., 2021]. This design improves the model's generalization ability and adaptability to different data distributions.

4.4 Stiefel Graph Spectral Convolution

Basic Concepts

Definition 1. The Stiefel Graph Spectral Convolution of a signal x ∈ Rn×1 with a kernel g ∈ Rn×1 is defined as:

where F satisfies the following optimization problem:

In Definition 1, ⊙ denotes the Hadamard product, St(n, d) ≜ {F ∈ Rn×d : F T F = Id} is the Stiefel manifold, L = In - D-1/2AD-1/2 is the Laplacian matrix, A ∈ Rn×n is the graph adjacency matrix, and D is the diagonal matrix, and its diagonal element is the sum of all elements in each row of A, which represents the degree of the corresponding node in an undirected graph (i.e., the number of edges connected to the node).

Theorem 1. The matrix F can be obtained by solving the eigenvalue decomposition of D-1/2AD-1/2 and selecting the eigenvectors corresponding to the d largest eigenvalues.

The proof of Theorem 1 is provided in Appendix A.1. 1

The F T x in Definition 1 can be considered as a pseudograph Fourier transform, which we refer to as the Stiefel Graph Fourier Transform (SGFT), defined as:

We can also define the Inverse Stiefel Graph Fourier Transform (ISGFT) as:

The SGSC of x and g can be obtained by first applying the SGFT to x and g separately, then performing element-wise multiplication, and finally applying the ISGFT. Since F is not orthogonal, it does not satisfy the definition of the standard graph Fourier transform. However, we can prove that SGSC is still a specific type of standard graph spectral convolution [Defferrard et al., 2016].

Theorem 2. The Stiefel Graph Spectral Convolution can be viewed as the following specific filtered graph spectral convolution:

where Λ = diag(λ1, λ2, · · · , λn), and λ1, λ2, · · · , λn are the eigenvalues of A sorted in descending order. P = (p1, p2, · · · , pn), where pi ∈ Rn×1 is the eigenvector corresponding to λi. gθ(Λ) = diag(θ), and θ = (θ1, θ2, · · · , θn) satisfies:

The proof of Theorem 2 is provided in Appendix A.2.

According to Theorem 2, it can be seen that the Stiefel Graph Spectral Convolution satisfies all the properties of graph spectral convolution. Moreover, setting the coefficients θi corresponding to smaller eigenvalues to zero reduces the time and space complexity of the graph spectral convolution. And it also serves as a filtering and noise-reduction mechanism.

The SGSC we define can be easily extended to the case of multi-dimensional features. For X ∈ Rn×k , if we choose the convolution kernel G ∈ Rn×k , then the SGSC is given by:

Algorithm 1 Implementation of LDGOSM

Input: X ∈ Rn×d , E ∈ Rn×d Output: W ∈ Rd×d

• 1: Perform eigenvalue decomposition on B = XT X to obtain eigenvector matrix D and eigenvalue diagonal matrix Λ;
• 2: Let M =DΛ-1/2DT ;
• 3: Compute C = ET X;
• 4: Compute H = B +C T C;
• 5: Performeigenvalue decomposition onM T HM to obtain eigenvector matrix U ;
• 6: returnW =MU;

Linear Dynamic Graph Optimization on Stiefel Manifold

In the SGSC, the introduction of the Stiefel manifold simplifies the convolution computation, but the eigenvalue decomposition of L still results in a time complexity of O(N 3 ), where N represents the total number of nodes across all variables and all patches in our paper. Moreover, for spatiotemporal data, the spatial relationships between variables change over time, requiring graph adjacent matrix A to be recalculated frequently. Each computation of the graph requires solving for eigenvalues, leading to a more complex computation. To simplify the dynamic computation of A, inspired by [Yang et al., 2024], we introduce the following calculation of the dynamic graph adjacency matrix:

Furthermore, to reduce the complexity of the SGSC, we improve the computation of F in Theorem 1. We propose an efficient Linear Dynamic Graph Optimization on Stiefel Manifold (LDGOSM). Inspired by [Huang et al., 2018], we introduce a linear transformation: F = XW, where W ∈ Rd×d . By substituting equation (10) into the equivalent form of (4) (as proven in Theorem 1 in Appendix A.1), the optimization problem can be transformed into the following objective function:

The optimal solution W for the objective function (11) can be obtained through Algorithm 1 (see Appendix A.3 for detailed derivation).

In Algorithm 1, the complexity matrix B and the eigenvalue decomposition are O(nd2 ) and O(d3 ), respectively. Computing matrices M, C, and H are O(d3 ), O(nd2 ), and O(d2 ). Computing M T HM and the subsequent eigenvalue decomposition are both O(d3 ). Finally, computing matrix W is O(d3 ). The total complexity is O(nd2 + kd3 ), where k is a small constant. If d ≪ n, the complexity is dominated by O(nd2 ), making it approximately linear with nodes number n.

1 The appendix can be accessed at:

https://github.com/komorebi424/DST-SGNN.

5 Experiments

Multi-layer SGSC

To enhance the effectiveness of dynamic spatio-temporal graph learning, we employ the following Multi-layer SGSC (MSGSC):

where Gi represents the convolutional kernel of each layer in the MSGSC.

Theorem 3. The MSGSC has the following equivalent computational form:

The proof of Theorem 3 is provided in Appendix A.4.

Through Theorem 3, when computing the MSGSC, multiple SGFT and ISGFT operations are avoided, thereby improving computational efficiency while retaining the effects of the original multi-step convolution.

In our implementation, S(Gj) = F T Gj can be learned via a linear fully connected layer (FC), i.e., ˆGj = FC(F T ). Therefore, the computation of MSGSC can be expressed as:

where, F is optimized using the Algorithm 1.

4.5 Concat Module and Linear Layers Module

In DST-SGNN, the seasonal and trend components y(1) t and y(2) t output from the MSGSC module are processed through inverse RevIN, resulting in ˆ Y (1) t and ˆ Y (2) t . These features are then concatenated via the Concat module to form ˆ Yt, which is subsequently fed into the Linear Layers module to generate the final prediction Yt.

4.6 The Inference and Training of DST-SGNN

In DST-SGNN, spatio-temporal data is decomposed into patches and split into seasonal and trend components. The LDGOSM module optimizes the transformation matrix F on the Stiefel manifold. The MSGSC module captures complex patterns via graph spectral convolutions. Predictions are generated by concatenating features and passing them through linear layers. The training process involves only standard eigenvalue decomposition, which can be efficiently handled using PyTorch's linalg module, and uses MAE loss for optimization.

5.1 Datasets and Settings

Datasets

To evaluate the performance of DST-SGNN 2 , this study conducted experiments on seven spatio-temporal benchmark datasets: PEMS03, PEMS07, Electricity, CSI300, exchange rate, Solar, and METR-LA. Table 1 provides a summary of the dataset statistics, with detailed information on the seven public benchmark datasets as follows:

• PEMS03 and PEMS07: These two datasets consist of California highway traffic flow data, collected in realtime every 30 seconds by the California Performance Measurement System (PEMS) [Song et al., 2020]. The original traffic flow data is aggregated at 5-minute intervals. The datasets record the geographical information of the sensor stations. PEMS03 and PEMS07 are datasets from two specific regions.
• Electricity: Electricity contains hourly time series of the electricity consumption of 321 customers from 2012 to 2014 [Asuncion et al., 2007].
• CSI300: This is a stock dataset, compiling the closing prices of one hundred stocks that have been continuously listed in the CSI 300 from 2012 to 2024 [Zheng et al., 2025].
• exchange rate: This dataset records the daily exchange rates of eight different countries from 1990 to 2016 [Lai et al., 2018].
• Solar 3 : The dataset selects power plant data points in Florida as the dataset, which includes 593 points. The data collection period spans from January 6, 2001, to December 31, 2006, with samples taken every hour, totaling 8,760 in length.
• METR-LA: The dataset includes traffic information collected by loop detectors on highways in Los Angeles County from March 1, 2012, to June 30, 2012. It contains data from 207 sensors, with samples taken every 5 minutes. [Yi et al., 2024a].

We divided all the datasets chronologically into a training set (70%), a validation set (10%), and a test set (20%).

Experimental Settings

The DST-SGNN model is implemented in Python using PyTorch 2.2.0 and is trained on a GPU (NVIDIA GeForce RTX 4090). MSE and MAE are used as evaluation metrics. We set the input window size T to 96, 192, or 336, and the learning rate to 0.0001. All of the models follow the same experimental setup with forecasting window sizes H ∈ {96, 192, 336, 720} for the PEMS03, Electricity, Sola and METR-LA datasets, H ∈ {48, 96, 192, 336} for the exchange rate datasets, H ∈ {36, 48, 60, 96} for the CSI300 dataset, and H ∈ {60, 96, 192, 336} for the PEMS07 dataset.

2 Our code at https://github.com/komorebi424/DST-SGNN.

3 https://www.nrel.gov/grid/solar-power-data.html.

5.2 Methods for Comparison

Our experiments comprehensively compared the forecasting performance of our model with several representative state-of-the-art (SOTA) models on seven datasets, including the advanced long-term forecasting models TimeMixer [Wang et al., 2024a] and TimeXer [Wang et al., 2024b], the Transformer-based model iTransformer [Liu et al., 2023], the spatio-temporal forecasting model MiTSformer [Chen and Zhao, 2024], and the Fourier transform based models FourierGNN [Yi et al., 2024a], FreTS [Yi et al., 2024b], Fredformer [Piao et al., 2024], and FITS [Xu et al., 2023]. We used the code released in the original papers, and all models followed the same experimental settings as described in the original papers.

5.3 Main Results

Table 2 presents the experimental results of DST-SGNN compared with other methods on the seven spatio-temporal datasets, leading to the following explanations:

• The DST-SGNNmodel shows excellent performance on various datasets, with an average rank of 2.7857 for MSE and 1.5000 for MAE. In particular, on the CSI300 and exchange rate datasets, the DST-SGNN significantly outperforms other comparative methods across all forecasting window lengths. This result demonstrates that DST-SGNN has a significant advantage in handling
• financial data and enables more accurate forecasting of stock prices and exchange rates, which is crucial in financial decision-making.
• FreTS based on Fourier transform performs well on the traffic datasets PEMS03, PEMS07, and METR-LA, but underperforms on other datasets. In contrast, iTransformer, which is based on the Transformer architecture, shows better performance on non-traffic datasets and falls short on traffic datasets. The reason is that traffic data exhibit clear periodicity and trends, which can be easily captured by Fourier transform-based methods, leading to better results on traffic datasets. Moreover, transformer-based models excel at modeling long-term dependencies and complex nonlinear patterns on nonperiodic data. DST-SGNN is superior on all datasets by incorporating the strengths of MSGSC, which effectively extracts complex relationships in both temporal and spatial dimensions. Whether dealing with periodic traffic data or non-periodic financial and power data, DST-SGNN maintains high forecasting accuracy and efficiency.
• TimeXer also performs well on all datasets, but as the forecasting length increases, the advantages of DSTSGNN become increasingly evident. This is attributed to the optimization of SGSC on the Stiefel manifold, which enables the model to more effectively handle the dynamic changes in spatio-temporal data. Additionally, DST-SGNN uses patch-based sequence decomposition to break down the original data into multiple smaller subsequences. This process further reduces the high di-

mensionality and complexity of the data. As a result, our model can maintain high efficiency and accuracy even with long forecasting windows. These characteristics make the advantages of DST-SGNN more prominent in long forecasting sequences and provide a more efficient and accurate solution for STTS forecasting.

5.4 Ablation Studies

To further investigate the specific role of the core module SGSC in DST-SGNN, we conducted ablation studies with the following three variants:

• Rep-StdGSC: Variant replacing SGSC with standard graph spectral convolution, using an orthogonal matrix F ∈ Rn×n for Graph Fourier Transform to perform spectral convolutions.
• Rep-SpatConv: Variant replacing SGSC with spatial graph convolution, where the convolution operation is defined as ∑i AiXiW, aggregating neighboring node features directly in the spatial domain.
• w/o-Stiefel: Variant removing the Stiefel manifold constraint and replacing the original LDGOSM with a multilayer perceptron (MLP).

The experimental results shown in Table 3 indicate that the three variants-Rep-StdGSC, Rep-SpatConv, and w/oStiefel-perform worse than DST-SGNN, although they still demonstrate the effectiveness of GNNs in spatio-temporal analysis. DST-SGNN's superior performance can be attributed to its usage of the Stiefel manifold, which provides filtering capabilities by discarding low-information eigenvalues. This mechanism not only reduces noise but also enhances the model's ability to capture complex spatiotemporal patterns efficiently. Compared to Rep-StdGSC and Rep-SpatConv, DST-SGNN's other advantage lies in its computational efficiency, which will be further analyzed in the subsequent sections regarding time complexity.

5.5 Time Consumption Analysis

Table 4 presents the experimental time consumption and parameter count for DST-SGNN, Rep-StdGSC, Rep-SpatConv, MiTSformer, Fredformer, and iTransformer on the PEMS07 dataset. The results show that DST-SGNN achieves lower time and space complexity than most methods due to its optimized SGSC. Although DST-SGNN has a longer training time than iTransformer due to the eigenvalue decomposition in the Stiefel manifold optimization, its inference time and parameter count are lower. This is attributed to the reduced complexity of eigenvalue decomposition to O(d3 ) with d ≪ n, which enables faster inference despite the additional training overhead. In addition, DST-SGNN has the fewest parameters. This highlights DST-SGNN's advantage in spatial complexity, making it more memory-efficient while also reducing inference overhead. For more experimental results, please refer to Appendix C.

6 Conclusions

The DST-SGNN model proposed in this paper effectively handles the dynamic changes and high-dimensional characteristics of spatio-temporal data by combining the advantages of graph neural networks and SGSC. DST-SGNN reduces data dimensionality and complexity through patchbased sequence decomposition and optimizes its core module, the SGSC, on the Stiefel manifold, reducing computational complexity while maintaining sensitivity to dynamic spatio-temporal relationships. Experiments show that DSTSGNN outperforms existing methods in forecasting accuracy and computational efficiency on seven public benchmark datasets, particularly in handling financial time series data and long forecasting windows. Our experiments have covered a variety of spatio-temporal data. In the future, we will consider more spatio-temporal application scenarios and further utilize and improve our model to apply it to these scenarios.

Acknowledgments

This research was supported by the PCL-CMCC Foundation for Science and Innovation (Grant No. 2024ZY2B0050).