Physics-Informed Long-Sequence Forecasting From Multi-Resolution Spatiotemporal Data

1 Introduction

Forecasting from spatiotemporal data has wide applications in domains such as transportation and energy. The input and output of such forecasting tasks are both sequences of graph signal frames, where each frame is a graph with multivariate features defined on nodes and edges. In these applications, long-run forecasting is usually required for planning and policy making. As illustrated by [Zhou et al., 2021], forecasting long sequences into the long-run future has higher requirements on the long-range alignment ability of models compared to short-run forecasting tasks.

One key aspect for forecasting long sequences is effectively modeling the complex patterns and dynamics of real-world spatiotemporal data aggregated in various spatial and temporal resolutions. Figure 1 shows an example from the taxi demand dataset. Data aggregated in different spatial (first 4 rows vs. the last row) and temporal (left/middle/right columns) resolutions demonstrate correlated but heterogeneous patterns, implying the necessity of refined modeling of multi-resolution data. In practical cases where data of high resolutions usually suffer from high missing rates and low signal-to-noise ratios due to the high cost of collection, correctly capturing the interaction among data in various resolutions is even more critical for forecasting.

While existing works have achieved great success in short-run predictions of uni-resolution spatiotemporal data within several hours into the future [Zhang et al., 2017; Geng et al., 2019; Yao et al., 2019; Wu et al., 2019; Zheng et al., 2020; Chen et al., 2020; Zhou et al., 2020; Cao et al., 2020], research on long-run forecasting of spatiotemporal data remains hardly developed. [Zhou et al., 2021] proposes an efficient Transformer-based [Vaswani et al., 2017] architecture for long-sequence time series forecasting, but still lacks modeling of spatial relations and ignores the interresolution dynamics. We include a more comprehensive summary of related works in the appendix 1 .

To bridge the gap between multi-resolution spatiotemporal data and long-sequence forecasting, the proposed model must address challenges in two aspects: (1) Inter-resolution modeling. The model must be able to capture interactions of data among resolutions to fully utilize them for forecasting. (2) Intra-resolution modeling. The model should effectively extract information representing the dynamics from either short or long sequences of each resolution.

In this paper, we propose Spatiotemporal Koopman MultiResolution Network (ST-KMRN) to address both challenges. For better inter-resolution modeling, we leverage the selfattention mechanism to fuse and communicate representations of input data in all temporal resolutions. The model is then jointly trained with forecasting losses over all spatial and temporal resolutions. We further construct downsampling and upsampling modules among forecasting outcomes of different temporal resolutions to enhance inter-resolution connections. To improve intra-resolution modeling, we combine the Koopman theory-based modeling of dynamic systems with deep learning-based prediction modules for more principled modeling of dynamics. The inferred Koopman matrix offers interpretations of heterogeneous dynamics in different resolutions at the same time.

1 https://github.com/mengcz13/mengcz13.github.io/raw/master/ pdf/ijcai2022-appendix.pdf

Our contributions are: (1) ST-KMRN captures interresolution dynamics of spatiotemporal data via connecting representations in different levels of resolutions via self-attention, providing an effective way of leveraging multi-resolution data. (2) ST-KMRN further improves inter-resolution modeling with upsampling and downsampling modules among predictions of various resolutions. (3) ST-KMRN improves the modeling of intra-resolution dynamics with the combination of Koopman theory-based modeling and deep learning-based forecasting models, which also provides interpretable information of different dynamics in multi-resolution data. (4) ST-KMRN achieves state-of-the-art performance on the longsequence forecasting tasks from real-world spatiotemporal datasets.

2 ProblemFormulation

In this work, we focus on the task of forecasting spatiotemporal multivariate time series with multiple resolutions of both space and time. We denote the set of regions by S = {s1, s2, . . . , sN}, the set of P historical time steps by TH = {t1, t2, . . . , tP}, the set of Q future time steps by TF = {tP+1, . . . , tP+Q}, and the multivariate variable of dimension D at region s and time step t by xt,s ∈ RD .

Definition 1 (Resolution). A resolution R(X) is a partition of the set X. In addition, for a temporal resolution R(T) (T = TH or T = TF ), we require |T| = |R(T)| × r, ∃r ∈ N, and R(T) =

{{t1, . . . , tr}, . . . , {t(|R(T)|-1)×r+1, . . . , t|R(T)|×r}}. We name r as the scale of the temporal resolution R(T).

Definition 2 (Aggregated Variable). The aggregated multivariate variable over a set of regions S′ and a set of time steps T′ is defined as xagg T′ ,S′ = agg({xt,s | t ∈ T′ , s ∈ S′ }), where agg is some aggregation function agg : RN×D → RD , such as summation or average. We omit agg in following notations since the aggregation function is usually the same for a given dataset and can be inferred from the context.

Definition 3 (Observation at Given Spatial and Temporal Resolutions). The observation at a given spatial resolution R(S) and a given temporal resolution R(T) is defined as a set of aggregated multivariate variables: OR(T),R(S) = {xrT ,rS |rT ∈ R(T), rS ∈ R(S)} ∈ R|R(T)|×|R(S)|×D .

We then formulate the input and output of the forecasting problem with multi-resolution spatiotemporal input data as follows:

Input (1) The set of regions S, the set of historical time steps TH, and the set of future time steps TF . (2) The list of available temporal resolutions TR = {RT1 (TH), . . . ,RTRT (TH)}, and the list of available spatial resolutions SR = {RS1 (S), . . . ,RSRS (S)}. RH1 (TH),RH1 (S) are the highest temporal and spatial resolution respectively. In other words, |rt| = 1, ∀ rt ∈ RT1 (TH), and |rs| = 1, ∀ rs ∈ RS1 (S). (3) The set of historical spatio-temporal resolution pairs will be HR = TR × SR. (4) The set of historical observations in multiple spatial and temporal resolutions O = {OTR,SR | (TR, SR) ∈ HR}.

Output The set of forecast values at the target spatial and temporal resolution ˆY = ORTout (TF ),RSout (S).

3 ProposedMethod: Spatiotemporal Koopman Multi-Resolution Network

In this section, we introduce our proposed model: Spatiotemporal Koopman Multi-Resolution Network (ST-KMRN). Figure 2a provides an overview of its architecture.

Input historical observations O in multiple spatial and temporal resolutions are first grouped by temporal resolutions into {OTR1 , . . . ,OTRRT }, where OTRi = concat({OTRi,SR1 , . . . ,OTRi,SRRS }) ∈ R|TRi|×NS×D ,NS = ∑RS k=1 |SRk| is the concatenated multivariate time series along the spatial dimension. For OTRi , we construct a spatial hierarchical graph Gi = (Vi, Ei), where Vi = ⋃RS k=1 SRk contains region sets at all spatial resolution levels and Ei is composed of 2 types of edges: (1) intra-spatial-resolution edges constructed with prior knowledge such as spatial information or connection strengths; (2) inter-spatial-resolution edges between all pairs of region sets (Sp, Sq) iff. Sp ⊂ Sq.

Each pair of (OTRi , Gi) is fed into the i-th spatiotemporal encoder to generate a temporal-resolution-specific embedding Ei ∈ RNS×H , where H is the embedding dimension. Embeddings of all temporal resolutions are then propagated and updated with the self-attention module and the updated embedding E′ i serves as the input of the i-th decoder to get the deep learning based forecasting result ˆYDL TRi . Meanwhile, a Koopman Theory based forecasting module takes OTRi as its input and gives its forecasting ˆYKTRi in parallel. ˆYDL T Ri and ˆYKTRi are then fused with a gating mechanism to get the first-stage forecasting result ˆY1TRi . We further design upsampling and downsampling modules for the second-stage forecasting: ˆY1TRi , ˆYds TRi = DownSamplingi-1(ˆY1TRi-1 ), ˆYups TRi = UpSamplingi+1(ˆY1TRi+1 ) is combined with the attention mechanism for the final forecast ˆYTRi (for i = 1 / N, ˆYds TRi-1 / ˆYups TRi+1 is omitted respectively).

(b) Neural network based module.

Figure 2: ST-KMRN's architecture.

3.1 Physics-InformedModeling of Intra-Resolution Dynamics

In this part, we introduce our proposed physics-informed modeling of intra-resolution dynamics, where a neural network based module and a Koopman theory based module are combined for encoding and forecasting in each temporal resolution.

Neural Network Based Module (ST-Encoder and ST-Decoder)

The architecture of the neural network based module is a modification of Graph-WaveNet [Wu et al., 2019]. Figure 2b shows the architecture of the encoder part, where multiple temporal convolutional layers (TCN) [Oord et al., 2016] and graph convolutional layers (GCN) [Kipf and Welling, 2017] are stacked alternately to encode spatial dependencies at various temporal scales. The output of the encoder is an embedding vector Ei summarizing the input sequence OTRi . Ei is then propagated with embeddings from other resolutions into E′ i via self-attention [Vaswani et al., 2017]. The decoder is composed of two 1D convolutional layers applied along the spatial dimension, which maps E′ i to the forecasting result ˆYDL T Ri .

Koopman Module

Koopman Theory Given a non-linear dynamic systemwith its state vector at time t denoted as st ∈ Rm. The system can be described as st+1 = F(st). As defined in [Koopman, 1931], the Koopman operator KF is a linear transformation defined on a function space F by KF g = g ◦ F for every g : Rm → R that belongs to the infinite-dimensional Hilbert space F. With this definition, we have KFg(st) = g ◦ F (st) = g(st+1).

The Koopman theory [Koopman, 1931] guarantees the existence of K, but in practice we often assume the existence of an invariant finite-dimensional subspace G of F spanned by k bases {g1, g2, . . . , gk}. Define gt = [g1(st), g2(st), . . . , gk(st)]T and gt+1 = [g1(st+1), g2(st+1), . . . , gk(st+1)]T , under the assumption we have gt, gt+1 ∈ G and there exists a Koopman matrix K ∈ Rk×k s.t. gt+1 =Kgt.

Koopman Theory Based Modeling of Temporal Dynamics Compared to directly modeling the state propagation function F with neural networks, the Koopman theory provides prior knowledge related to the system. In addition, the resulting Koopman matrix is more interpretable as the propagation of hidden states can be represented as a linear mapping.

The key problem is to find the pair of mappings between the state space Rm and the invariant subspace G ∈ Rk : g : Rm → Rk and g-1 : Rk → Rm. Here we model them with 2 multi-layer perceptrons (MLPs). g maps each frame of the input OTRi to the subspace, and generatesHTRi ∈ R|TRi|×k . Then, the Koopman matrixKTRi is estimated as

With HTRi [-1] as the starting state q0 for forecasting, the module forecasts Ti steps into the future in a recurrent way: ˆYKTRi = [g-1 (KTRi q0), g-1 (K2TRi q0), . . . , g-1 (KTi T Ri q0)].

Gated Fusion of ST-Encoder and Koopman Module

While the Koopman theory based module is more principled and can be incorporated with prior knowledge of the system, the neural network module is far more expressive and flexible. To combine the advantages of both modules, we design a fusion module that jointly leverages their forecasting. More specifically, we use a 2-layer 1D convolutional neural network (CNN) with the sigmoid activation function to generate a vector. The vector controls the ratio of each module's prediction in the combined forecasting results as follows:

where ⊗ is the element-wise product.

3.2 Inter-Resolution Dynamics Modeling

In this part, we introduce the self-attention module that fuses representations over all temporal resolutions as well as the upsampling and downsampling modules for the second stage forecasting.

Self-Attention Module

The self-attention mechanism, proposed in [Vaswani et al., 2017], captures dependencies of input elements without regard to their distance in the input sequence. With similar motivation, we adopt self-attention to model the interaction of representations frommultiple temporal resolutions. Given the list of resolution-specific embeddings [E1,E2, . . . ,ERT ], we concatenate all embeddings to construct the input representation as E ∈ RNS×RT×H . We then update the embeddings of all regions in parallel via multi-head self-attention across resolutions.

Upsampling and Downsampling Modules

Upsampling Module We set up a learnable upsampling module UpSamplingi for each ˆY1TRi , i > 1. Denote ri as the scale of the temporal resolution TRi. In our experimental settings we always have ri = kiri-1, ∃ki ∈ N, ∀i > 1. We implement the i-th upsampling module as a 1D CNN with D input channels and kiD output channels. Its output is then rearranged as the upsampled prediction in the same way as PixelShuffle [Shi et al., 2016] for the upsampled forecasting results ˆYups T Ri-1 .

Downsampling Module We also set up a downsampling module DownSamplingi for each ˆY1TRi , i < RT . In the same experimental settings as mentioned above, we implement the downsampling module as an aggregation operator, which is of the same type of aggregation for generating data in coarser temporal resolutions. For the temporal resolution TRi, i < RT , the operator aggregates the prediction in ˆY1TRi every ri+1/ri steps as the downsampled prediction ˆYds T Ri+1 .

Weighted Summation of Forecasting Results For the ith temporal resolution (1 < i < RT ), we use the weighted summation of the first-stage results, the upsampled results, and the downsampled results as the second-stage prediction. Weights for summation are the output of a learnable module implemented as CNN with softmax activation function defined for the i-th resolution:

where ⊗ is the element-wise product. For cases where i = 1[resp. RT ], the weighted summation is calculated in the same way but with the ˆYds TRi [resp. ˆYups TRi ] term removed.

3.3 Loss Function

Denote the forecasting loss function as L(ˆYpred, Ytrue), the total loss for training the model is

L(KM) TRi encourages the mapping from the state space to the invariant subspace in the Koopman module to preserve the distance.

4 Experiments

Datasets We evaluate the performance of ST-KMRN and all baselines on 3 datasets: (1) NewYorkYellowTaxi Trip Record Data (YellowCab) [NYCTLC, 2021] in 2017-2019; (2) New York Green Taxi Trip Record Data (GreenCab) [NYCTLC, 2021] in 2017-2019; and (3) Solar Energy Data (Solar Energy) [NREL, 2021] of Alabama in 2006. We use sliding windows to generate input/output sequence pairs ordered by starting time and divide all pairs into train/validation/test sets with the ratio 60%/20%/20%. Since all datasets are originally in a single resolution, we construct one extra spatial resolution (named as "agg regions") and two extra temporal resolutions on all datasets by aggregation to simulate the multi-resolution scenario.We list statistics and task settings of each dataset in the appendix.

Baselines We compare our model ST-KMRNwith the following baselines: (1) Historical Averaging (HA): We use the averaged value of historical frames with the 1-week period as the prediction. (2) Static: We use the value from the last available frame in the input sequence with the 1-week period as the prediction. (3) Gated Recurrent Unit (GRU) [Chung et al., 2014] (4) Informer [Zhou et al., 2021] (5) Graph WaveNet [Wu et al., 2019] (6) MTGNN [Wu et al., 2020] (7) KoopmanAE [Azencot et al., 2020]. For all baselines, we concatenate input in multiple temporal resolutions as features into one single temporal resolution. We have also considered other recent works on short-term spatiotemporal forcasting as baselines, such as AGCRN [Bai et al., 2020], GMAN [Zheng et al., 2020], DGCRN [Li et al., 2021]. However, their designs induce high computation and memory complexity on long sequential data and prohibit us from retrieving results practically.

4.1 Long-Sequence ForecastingWith Fully Observed Input

Set-up We train all models and evaluate their forecasting performance within certain horizons into the future. We repeat each experiment 3 times and report mean values and standard deviations of all metrics.

Discussion Table 1 shows the forecasting performance within the maximum possible horizon (10-day for YellowCab/GreenCab and 3-day for Solar Energy) of baselines and our proposed ST-KMRN measured in Mean Absolute Error (MAE) and Rooted Mean Squared Error (RMSE). Results within other horizons are in the appendix and we have similar observations. We observe that: (1) On all of 3 datasets, ST-KMRN outperforms all baselines, demonstrating the advantage of ST-KMRN in long-sequence forecasting with multiresolution data. (2) ST-KMRN achieves significantly decreased forecasting errors compared to the Graph-WaveNet baseline, with whom its encoders and decoders share similar neural network architectures. The results demonstrate that the gain of ST-KMRN comes from the enhanced modeling of both intra-resolution and inter-resolution dynamics. (3) On the YellowCab and Solar Energy datasets, our proposed ST-KMRN still holds the best performance compared to other baselines, but its advantages over HA and Static are not as large as on the GreenCab dataset. This is due to the strong periodicity within the fully observed input sequences. We further conduct experiments in a more challenging setting: input data is partially observed, where the periodicity in input is weaker.

4.2 Long-Sequence Forecasting Results With Partially Observed Input

Set-up We evaluate the effect of missing data in input sequences by training and evaluating all models with varying ratios of observable frames in the input of the finest spatial and temporal resolution (named as "obs ratio"). This setting aims to simulate the practical scenario when the model needs to forecast with low-quality input data of high resolutions.

Discussion Table 2 shows the 10-day forecasting performance of all models with various obs ratios on YellowCab. When the input data suffers from high missing ratios, the periodicity of data is less beneficial for forecasting, and the capability of capturing non-periodic patterns becomes more critical. Thus, trainable models start having advantages over the Static baseline. Under all obs ratios we select (0.8/0.6/0.4/0.2), our proposed ST-KMRN achieves the lowest forecasting errors, especially for higher ones (0.4 and 0.2), demonstrating its advantages with partially observed input. We also provide results for GreenCab and Solar Energy in the appendix, from which we have similar observations.

4.3 Ablation Study

Set-up We conduct the ablation study on the YellowCab dataset with 80% observation ratio to evaluate the contribution of proposed components. By removing each component from ST-KMRN, we have the following settings: (1) w/o SelfAttn: ST-KMRNwithout the self-attention module. (2) w/o Koopman: ST-KMRN without the Koopman module. The prediction will only be performed from the neural network module. (3) w/o ups/ds: ST-KMRN without upsampling and downsampling modules. Only the decoder of the target resolution will forecast. We repeat experiments of each setting for 3 times and report results in Table 3. Higher increase of errors after removing one module indicates larger contributions.

Discussion From Table 3, we observe that (1) the selfattention module applied on representations of different resolutions significantly improves performance. It allows information exchange among resolutions, and it contributes most on shorter horizons (6h). (2) The combination of both Koopman and neural network modules brings significant boosts to the forecasting performance for both short (6h) and long (10d) prediction horizons. It captures the correct patterns for each temporal resolution, which we will detail in Section "Interpretability of Koopman Module". (3) The upsampling and downsampling modules can improve performance on all forecasting horizons as (a) the downsampling module regularizes the output by minimizing the difference between its aggregated values with ground truth values in lower resolutions; (b) the upsampling module is trained to predict local temporal patterns and can fix the errors produced by the decoder, which forecasts values on all steps at once.

4.4 Interpretability of KoopmanModule: Revealing Dynamics in Each Resolution

Set-up Since the Koopman module models the temporal dynamics as a linear system in the hidden space, the derived Koopman matrix describing the system can provide rich interpretable information. With eigen decomposition applied to the Koopman matrix, we can decompose the system dynamics into components with different magnitudes and periods. An eigenvalue λejθ corresponds to a component f(t) = λt ejθt with magnitude λt and frequency θ (i.e. period T = 2π θ ).

Discussion Figure 3 displays the distribution of eigenvalues of Koopman matrices for the input data of each resolution. We observe that: (1) For data in the resolution of 30-min, 6hour, and 1-day, the maximum of magnitudes of eigenvalues reaches its peak value around the period of 24-hour, 1-day, and 7-day respectively (we ignore the eigenvalues with zero angles since they represent components with infinite periods and thus are non-periodic). This matches our observation that the taxi demand data shows strong daily and weekly patterns. (2) For data in the resolution of 30-min and 6-hour, we can still notice eigenvalues with large magnitudes around the 1week period, but it is hard to distinguish them from other components including non-periodic ones. The reason is that in data with high temporal resolutions, long periods correspond to large numbers of time steps T and frequency values θ close to 0. When we switch to data with a lower temporal resolution (1-day), the pattern of the same period can be represented by components with larger distances to non-periodic components in the spectral domain and are easier to identify. This again implies the necessity of utilizing multi-resolution data for better modeling temporal patterns with various periods.

5 Conclusion

We propose Spatiotemporal Koopman Multi-Resolution Network (ST-KMRN), which boosts the long-sequence forecasting with enhanced modeling of inter-resolution and intraresolution dynamics in multi-resolution data. To accomplish this, we introduce self-attention among representations of different resolutions together with upsampling and downsampling modules for better utilization of inter-resolution contextual information. Meanwhile, we improve the modeling of intra-resolution dynamics with the combination of Koopman theory based modeling. ST-KMRN achieves state-of-the-art performance on the long-sequence forecasting tasks frommultiple real-world spatiotemporal datasets. Limitations of our work include that the available resolutions used by the model are fixed and are usually dependent on the task and data availability, and the model only applies to data with regular time steps. Our future works include extending the framework to find optimal data resolutions to guide the data collection process and enabling it to model and forecast irregular spatiotemporal sequences.

Acknowledgements

This work is supported by KDDI Research, Inc., granted to co-author Yan Liu in her academic role at the University of Southern California. The views and conclusions are those of the authors and should not be interpreted as representing the official policies of the funding agency.