Beyond Attention: Hierarchical Mamba Models for Scalable Spatiotemporal Traffic Forecasting

1. Introduction

Accurate traffic forecasting is critical to predictive network resource allocation, network planning, and network optimization, directly influencing the operational expenditure of telecom providers [1]. Traditional time series forecasting methods such as ARIMA [2] and machine learning (ML) based models such as LSTM [3] treat each base station independently and fail to account for spatial dependencies among neighboring cells. To address this challenge, recent work has focused on employing ML-based models jointly with data pre-processing and feature engineering techniques to incorporate spatial knowledge in the training procedure [4,5].

A key challenge in spatiotemporal modeling is the trade-off between using a single global model and multiple cell-specific models. While cell-specific models can capture local dynamics with higher fidelity, they are costly to train, validate, deploy, and maintain, especially at the scale of modern cellular networks. In contrast, a global model is trained to cover the entire spatial grid. A global model incorporates spatial dependencies as features and formulates the forecasting task as a multivariate time series problem. Although global models can exploit shared patterns across cells, they often underperform due to distributional heterogeneity [1]. In 5G networks, user mobility and overlapping coverage introduce more complex spatial correlations [6], rendering purely temporal modeling insufficient, which leads to large forecasting errors.

By explicitly modeling spatial context alongside temporal dynamics, spatiotemporal architectures have the potential to mitigate the limitations of standalone temporal approaches. Indeed, a single spatiotemporal model could exploit correlations among neighboring cells, leading to higher forecasting accuracy. Recent advancements using graph neural networks (GNNs) and attention-based models have shown promise in learning complex spatiotemporal dependencies [7]. However, these architectures suffer from high computational cost, which limits the feasibility of their real-time large-scale deployment.

To address these challenges, we propose the Hierarchical SpatioTemporal Mamba (HiSTM), a novel spatiotemporal forecasting model based on the Mamba architecture [8] that employs a hierarchical design to capture both local variability and global patterns. HiSTM formulates cellular traffic forecasting as a multivariate time series over spatial frames, where each cell is treated as a pixel and its traffic measurements evolve over time. Frame-wise spatial encoding applies spatial convolutions independently at each time step to capture local neighborhood context, which is subsequently modeled using Mamba-based state-space dynamics for temporal forecasting. Unlike graph-based or attention-heavy models, the HiSTM achieves strong predictive accuracy while maintaining low computational overhead. It provides a favorable trade-off between performance and scalability.

Our main contributions are summarized as follows:

• We propose HiSTM, a hierarchical spatiotemporal forecasting architecture that performs frame-wise spatial encoding using CNNs, models temporal dynamics via Mamba state-space blocks, and aggregates temporal information through attention-based mechanisms.
• We show that HiSTM achieves state-of-the-art performance on the Milan cellular traffic dataset, reducing MAE by 29.4% compared to the STN baseline and achieving the lowest RMSE among all evaluated models. The model exhibits a strong cross-dataset generalization, achieving a 47.3% MAE improvement on the unseen Trentino dataset and successfully adapting to correlation-based network topologies (Liverpool 5G RU).
• We introduce cluster-aware variants (Global-Aware and ClusterFiLM) to explicitly model spatial heterogeneity. Our results show that a single Global-Aware HiSTM model outperforms 10,000 independently-trained LSTMs, demonstrating the advantage of learning shared representations of traffic regimes over localized training.
• We evaluate robustness under missing and noisy data conditions and show that HiSTM maintains stable performance with up to 30% data loss in comparison to baselines.
• We show that HiSTM provides improved long-horizon forecasting stability, exhibiting a 58% slower error accumulation rate compared to STN in autoregressive multi-step prediction.

This article substantially extends our previous conference work [9]. The description of the proposed model architecture and a section of the Milan evaluation are based on the conference version. However, this study goes significantly further by introducing: a bigger set of strong contemporary baselines, cluster-based models/experiments, benchmarking against 10,000 local LSTMs, correlation-based datasets, ablation/granularity studies, and explainability analyses.

The rest of this paper is structured as follows: Section II reviews 5G traffic prediction methods. Section III formalizes the forecasting problem and introduces HiSTM along with baseline models. Section IV describes the datasets and experimental setup. Section V presents results and analysis, including ablation studies and large-scale comparisons. Section VI provides an explainability analysis of the model’s internal dynamics. Section VII concludes the study and outlines directions for future work.

2. Related Work

Traditionally, 5G traffic forecasting is treated as a purely temporal forecasting task, relying on recurrent neural networks. AI-based models such as vanilla RNNs, LSTMs [3], and GRUs [10] were applied to learn sequential traffic patterns, often achieving gains over classical time-series methods such as ARIMA [2]. However, these approaches inherently ignore spatial correlations. Prior studies note that while RNN-based methods can capture long-term dependencies, they neglect spatial context, reducing their effectiveness [11,12]. For example, Chen et al. proposed a multivariate LSTM (MuLSTM) using dual streams for traffic and handover sequences, but the model remains fundamentally temporal [13].

Temporal models employ hybrid or preprocessing strategies for improved accuracy. For instance, Hachemi et al. [4] applied an FFT filter before LSTM to separate periodic signals, while Wang and Zhang [5] used Gaussian Process Regression alongside LSTM to manage bursty patterns. These methods often operate at hourly aggregation levels and improve accuracy in long-term forecasts. However, they still lack spatial awareness and suffer from high computational and memory costs, limiting their applicability in large-scale or real-time cellular deployments.

More recently, transformer-based architectures have gained traction for temporal sequence modeling due to their ability to capture global dependencies and adapt to irregular patterns. Models such as Informer [14], Autoformer [15], TimesNet [16], and PatchTST [17] achieve state-of-the-art results in long-sequence time-series forecasting. These methods introduce efficient attention mechanisms, decomposition blocks, and frequency-domain representations to handle non-stationary patterns and long-range dependencies. However, despite their performance gain in long-range temporal forecasting, they still treat traffic data as independent sequences and overlook spatial correlations among base stations.

Recently, spatial models have become popular as they can exploit spatial correlations within the data. For example, grid-based methods apply ConvLSTM [18] or 3D-CNNs to traffic maps [19], though these require regular cell layouts. Graph-based architectures, such as STCNet and its multi-component variant A-MCSTCNet, combine CNNs or GRUs with attention mechanisms [20]. STGCN-HO uses graph convolutions based on handover-derived adjacency to model spatial links across base stations, with Gated Linear Units (GLUs) for time dependencies [21]. Many of these models supplement deep networks with auxiliary features, such as weekday labels or slice-specific parameters, to improve generalization. More advanced spatiotemporal architectures, such as DSTL, adopt a dual-step transfer learning scheme that clusters gNodeBs and fine-tunes shared RNN models per group [22]. It forecasts at a 10-minute interval while reducing training overhead. Spatiotemporal models that incorporate spatial structure can capture richer dependencies, but this often comes with increased architectural complexity and computational cost. In prior work, we investigated an enhanced Spatiotemporal Network (STN) for cellular traffic forecasting that combines a temporal sLSTM branch (an xLSTM implementation [23]) with a Conv3D spatial branch and fuses both streams via attention-based integration [24]. Although the model achieved strong spatiotemporal forecasting accuracy and improved generalization to unseen regions, these gains came with a substantial computational footprint, including higher memory usage, increased MACs, larger model size, and elevated inference latency (especially at the grid level). Consequently, the approach can be incompatible with deployment scenarios that impose strict latency or resource constraints. This motivates our investigation of structured state-space and Mamba-based models for spatiotemporal forecasting, which aim to retain predictive quality while reducing computational overhead.

To the best of our knowledge, Mehrabian et al. [25] published the first model to incorporate the Mamba framework into a spatiotemporal graph-based predictor for 5G traffic. It adapts a dynamic graph structure within a bidirectional Mamba block to model complex spatial and temporal dependencies.

In contrast to graph-based Mamba predictors that operate on an explicit graph structure, HiSTM models spatiotemporal dynamics using kernelized neighborhood tensors. For grid datasets, a $K\times K$ neighborhood centered kernel is extracted on each target cell (c.f. 4.3.1), while non-grid datasets, the $K\times K$ neighborhoods are constructed by ranking other cells using Pearson correlation computed on the training interval (c.f. 4.3.2). This design avoids reliance on a predefined graph while retaining spatial context in a consistent tensor form.

Despite recent advances, existing spatiotemporal forecasting models face a fundamental trade-off between modeling capacity and scalability. Attention- and graph-based architectures capture complex dependencies but incur high computational overhead, while simpler models fail to account for spatial heterogeneity. This gap motivates the need for a scalable spatiotemporal architecture that jointly models spatial and temporal dependencies without relying on expensive attention mechanisms or explicit graph construction.

In contrast, our model, HiSTM, uses the Mamba framework and adapts data input and output processing to work directly on spatial grids with temporal attention layers. This design improves accuracy and efficiency, especially for long-term forecasts and generalization to unseen data. Accordingly, we evaluate HiSTM against strong representatives of temporal Transformers, recurrent spatiotemporal predictors, and modern sequence backbones (c.f.  4.2) under a unified experimental protocol.

3. System Design

This section formalizes the spatiotemporal forecasting problem and describes the proposed model architectures.

3.1. Problem Formulation

Cellular network traffic forecasting aims to predict future traffic volume based on historical observations. We formulate this task as a spatiotemporal sequence prediction problem. Given a sequence of traffic measurements $M=\{M_1,M_2,...,M_T\}$ where each $M_t\in {\mathbb{R}}^{H\times W}$ represents a spatial grid of $H\times W$ cells at time step t, our goal is to predict the next traffic volume grid $M_{T+1}$.

Beyond grid-structured maps, the same formulation can be applied by constructing semantic neighborhood tensors; in this paper we instantiate this by defining neighborhoods using Pearson correlation computed on the training interval and arranging the top-correlated cells into a $K\times K$ kernel. Each grid cell contains a scalar value representing the traffic volume. This volume correlates with resource demands in the cellular network.

We define our input tensor as $X_t^{(i,j)}\in {\mathbb{R}}^{T\times K\times K}$, where T is the sequence length (window of observation time steps), $K\times K$ represents the spatial dimensions of the input kernel, and $(i,j)$ are the spatial coordinate in the grid $M_t$. The prediction target is $x_{t+1}^{(i,j)}\in \mathbb{R}$,

representing the traffic volume at the center cell of the kernel for the next time step. Formally, we aim to learn a function $f:{\mathbb{R}}^{T\times K\times K}\to \mathbb{R}$ that minimizes the prediction error ${min}_{\theta }L(f_{\theta }\left(X_t^{(i,j)}\right),x_{t+1}^{(i,j)})$, where L is the loss function (i.e., Mean Absolute Error), and $\theta$ represents the learnable parameters of the model. The spatiotemporal nature of the data introduces a unique challenge: capturing both spatial correlations between neighboring cells and temporal dependencies across the sequence.

3.2. Proposed Architecture: HiSTM

We propose the Hierarchical SpatioTemporal Mamba (HiSTM), an architecture that combines hierarchical spatiotemporal processing with attention-based temporal aggregation (Figure 1). The term "hierarchical" here refers to the progressive abstraction of features: starting from raw local spatial interactions (Convolution), to temporal dynamics (Mamba), and finally aggregating into a high-level global context (Attention). Given an input tensor $\mathbf{X}\in {\mathbb{R}}^{T\times K\times K}$ (where T is the number of time steps and $K\times K$ is the spatial kernel), the model predicts target values through three key components.

3.2.1. Hierarchical Spatiotemporal Encoding

The input $\mathbf{X}$ is first augmented with an initial channel dimension (i.e., $D_{in}=1$), then passes through N stacked Encoder Layers. Each Encoder Layer l transforms its input ${\mathbf{X}}^{(l-1)}$ (where ${\mathbf{X}}^{\left(0\right)}$ is the initial augmented input) into an output ${\mathbf{X}}^{\left(l\right)}\in {\mathbb{R}}^{T\times K\times K\times C}$. The operations within each layer are:

• Spatial Convolution: Input features are first reshaped appropriately (e.g., to $T\times D_{in}^{\prime }\times K\times K$). A 2D convolution followed by a ReLU activation is applied. The first layer up-projects $D_{in}^{\prime }$ to C channels. Subsequent layers take $D_{in}^{\prime }=C$ channels and output C channels.
• Temporal Mamba Processing: The C-channel output from the convolution is reshaped to ${\mathbf{X}}_{\mathrm{flat}}\in {\mathbb{R}}^{\left(K^2\right)\times T\times C}$. A Mamba SSM [8] then models the temporal dependencies for each of the $K^2$ spatial locations treated as sequences of length T, with $d_{\mathrm{mamba}}=C$.

We use Mamba for its state space foundation, which models sequences through continuous-time dynamics rather than attention. This enables precise control over temporal structure and inductive bias, making it suitable for forecasting tasks where long-range temporal dependencies interact with fine-grained spatial patterns. Its ability to selectively retain and propagate information aligns with the demands of spatiotemporal modeling.

The Mamba output is reshaped back to ${\mathbb{R}}^{T\times K\times K\times C}$, forming ${\mathbf{X}}^{\left(l\right)}$. The output of this encoding stage is ${\mathbf{X}}_{\mathrm{encoded}}={\mathbf{X}}^{\left(N\right)}$.

3.2.2. Temporal Attention-Based Aggregation

From the encoded features ${\mathbf{X}}_{\mathrm{encoded}}\in {\mathbb{R}}^{T\times K\times K\times C}$, features corresponding to the center spatial cell are extracted across all T time steps. This yields a sequence ${\mathbf{X}}_{\mathrm{center}}\in {\mathbb{R}}^{T\times C}$. An attention mechanism then computes an aggregated context vector $\mathbf{c}\in {\mathbb{R}}^{B\times C}$:

$$e_t={\mathrm{Linear}}_{\mathrm{att}}\left({\mathbf{h}}_t\right),{\alpha }_t={\displaystyle \frac{exp\left(e_t\right)}{{\sum }_{j=1}^{T}exp\left(e_j\right)}},\mathbf{c}=\sum _{t=1}^T{\alpha }_t{\mathbf{h}}_t$$

(1)

where ${\mathbf{h}}_t\in {\mathbb{R}}^C$ is the feature vector from ${\mathbf{X}}_{\mathrm{center}}$ (for a given batch instance) at time step t. The ${\mathrm{Linear}}_{\mathrm{att}}$ layer maps from ${\mathbb{R}}^C\to \mathbb{R}$, producing a scalar energy $e_t$. The softmax function normalizes these energies across all T time steps to obtain attention weights ${\alpha }_t$.

3.2.3. Prediction Head

Finally, the aggregated context vector $\mathbf{c}\in {\mathbb{R}}^C$ is fed into a Multilayer Perceptron (MLP) head. This MLP consists of two linear layers with a ReLU activation in between, passing the dimensionality from C to $ML{P}_{in}$, and then to 1, to produce the final prediction $\widehat{y}\in {\mathbb{R}}^1$.

3.2.4. Cluster-Aware Variants

  To investigate the impact of spatial heterogeneity on model performance, we introduce two cluster-aware extensions to the base HiSTM architecture. These variants leverage external knowledge about traffic regimes (e.g., Urban vs. Rural) during training:

Global-Blind (GB): 

A single HiSTM model trained on all regions, but without any explicit cluster identifier. This setup ignores spatial heterogeneity and serves as a global baseline.

Global-Aware (GA): 

This variant incorporates a learnable cluster embedding $E_c\in {\mathbb{R}}^C$ associated with the cluster label c of the target cell. The embedding is added to the feature representation before the final MLP head, enabling the model to condition its predictions on region-specific characteristics.

Cluster-Specific (CS): 

Four separate HiSTM models are trained, one for each spatial cluster, using only data belonging to that cluster. This configuration represents cluster-localized learning without any information sharing across regions.

ClusterFiLM (CF): 

This variant employs Feature-wise Linear Modulation (FiLM) [26]. A FiLM generator network takes the cluster label c as input and produces scaling (${\gamma }_c$) and shifting (${\beta }_c$) parameters, which modulate intermediate feature maps $\mathbf{F}$ in the encoder layers via an affine transformation:

$$\mathrm{FiLM}(\mathbf{F}\mid c)={\gamma }_c·\mathbf{F}+{\beta }_c.$$

This mechanism allows the model to dynamically adapt its internal representations based on spatial context.

In addition, we consider CF-Unbalanced, a variant of ClusterFiLM trained with a cluster-weighted loss to compensate for data imbalance across regions.

4. Datasets & Experimental Setup

This section outlines the datasets, preprocessing steps, training configuration, and evaluation process used to assess model performance in spatiotemporal data forecasting.

4.1. Datasets

We utilize two distinct datasets to evaluate HiSTM’s performance across different spatial topologies: geospatial grid and correlation-based (relational network) data.

4.1.1. Geospatial Data (Milan and Trentino)

  Introduced by Barlacchi et al. [27], this dataset represents a classic Euclidean topology where traffic is aggregated into a regular spatial grid ($100\times 100$ for Milan, $117\times 117$ for Trentino). Spatial dependencies are inherently geometric; the "neighbors" of a cell are its physically adjacent grid cells. The data exhibits significant spatial heterogeneity (Figure 2), with activity concentrated in urban centers and dispersed in rural peripheries. Traffic is sampled at 10-minute intervals (144 measurements per day).

Lag plot analysis (Figure 3) reveals that these datasets maintain correlation at higher lag values. Individual cell traffic demonstrates higher volatility (Approximate Entropy [28] of 1.386 for a single cell) compared to spatially aggregated traffic (0.196), with the latter exhibiting enhanced cyclical patterns and passing the Augmented Dickey Fuller test [29] for stationarity. It is evident that the aggregated series is more correlated with itself for different lags; hence, it is easier to predict. Consequently, incorporating the spatial element in the prediction can help the model to capture the distribution more effectively, reduce the influence of individual events, and magnify predictable cyclical patterns.

To characterize spatial heterogeneity in traffic patterns, we partition the study area into four spatial clusters using mean traffic volume computed on the training set. This avoids information leakage from the evaluation period. The resulting clusters correspond to Downtown, Urban, Suburban, and Rural regions, exhibiting distinct traffic distributions and temporal dynamics (Table 1).

We characterize the temporal behavior within these spatial regions by comparing traffic distributions under different temporal conditions. We define nighttime as 22:00–06:00 and daytime as 06:00–22:00. Figure 4 shows a diurnal separation across all spatial types, with reduced volumes at night. The effect is strongest in Downtown, while Suburban and Rural areas show comparatively tighter distributions.

Figure 5 highlights a secondary weekly-cycle effect: activity decreases on weekends across all clusters, with the strongest shift in Downtown and comparatively stable behavior in Rural regions.

These spatial clusters are used in the cluster-aware HiSTM evaluation in 5.4.

4.1.2. Correlation-Based Relational Data (Liverpool 5G RU)

The Liverpool City Region High Density Demand (LCR HDD) project [30] tracks performance metrics at the Radio Unit (RU) level for high-density event venues (Salt & Tar and ACC Arena). Traffic is sampled at 1 second in Salt & Tar for 3,000 users and every 3 seconds in ACC Arena for 12,000 users, for a total of 10,000 timestamps. Unlike the grid-based Milan dataset, the "spatial" dependencies here are not defined by geometric proximity on a map, but by functional relationships (e.g., user load balancing, interference, and shared user mobility). This dataset challenges the model to capture dependencies in a graph-like structure where the "neighbors" are RUs with correlated traffic patterns, regardless of their physical distance.

4.2. Baseline Models

To assess the performance of our approach, we compare HiSTM with selected baselines that span a representative range of spatiotemporal forecasting paradigms, including convolutional models, recurrent architectures, Transformer-based methods, and state-space–inspired hybrids. These baselines cover models that emphasize spatial–temporal dependency modeling, long-range temporal dynamics, or computational efficiency in large-scale forecasting settings.

To assess the impact of cluster-aware conditioning, we compare HiSTM against two additional baseline configurations. Global-Blind (GB) refers to the standard HiSTM trained on the full dataset without cluster information, while Cluster-Specific (CS) consists of separate HiSTM models trained independently on data from each spatial cluster. These baselines allow us to isolate the effect of shared representation learning and explicit spatial conditioning.

Unless a baseline’s formulation restricts inputs, we train all models on the same forecasting task with identical temporal and spatial inputs: for each target cell, the input is a length-T history of $K\times K$ neighborhoods and the target is the next-step value of the center cell. For the RU dataset, neighborhoods are constructed by correlation and provided in the same tensorized $T\times K\times K$ form. For purely temporal models (e.g., xLSTM), the spatial dimension is flattened to match the model’s expected input shape. Baseline architectures are configured following their original implementations, with hyperparameters selected to balance performance and computational efficiency. All models are trained under the same optimization and early-stopping criteria.

• Informer [14]: A Transformer variant utilizing the ProbSparse self-attention mechanism to efficiently handle long sequence forecasting with reduced computational overhead.
• AutoFormer [15]: An auto-correlation-based Transformer that decomposes time series into trend and seasonal components to improve long-term forecasting accuracy.
• TimesNet [16]: A unified temporal convolutional architecture that models multi-periodic patterns through 2D time–frequency representations, achieving strong results across various forecasting tasks.
• PatchTST [17]: A patch-based Transformer for time series forecasting that tokenizes continuous temporal segments into patches, enabling efficient long-sequence modeling.
• STN [19]: A deep neural network capturing spatiotemporal correlations for long-term mobile traffic forecasting.
• xLSTM [23]: A scalable LSTM variant with exponential gating and novel memory structures (sLSTM/mLSTM) designed as Transformer alternatives. For fair comparison, our implementation uses only one mLSTM layer (denoted as xLSTM[1:0] in their paper) to prioritize computational efficiency while retaining its parallelizable architecture.
• STTRE [31]: A Transformer-based architecture leveraging relative embeddings to model dependencies in multivariate time series.
• VMRNN-B & VMRNN-D [32]: Vision Mamba-LSTM hybrids addressing CNN’s limited receptive fields and ViT’s computational costs. VMRNN-B (basic) and VMRNN-D (deep) use Mamba’s selective state-space mechanisms for compact yet competitive performance.
• PredRNN++ [33]: An improved version of PredRNN that introduces gradient highway units and spatiotemporal memory flow to better capture long-term dependencies in sequence forecasting.

4.3. Data Preprocessing

To ensure robust spatiotemporal feature extraction, we apply different preprocessing strategies tailored to the topology of each dataset. For each prediction target $x_{t+1}^{(i,j)}$, we construct the input by stacking the T most recent observations ($t-T-1,\dots ,t$) of the cell’s $K\times K$ neighborhood, i.e., we use a length-T sequence of past frames to forecast the next step. Depending on the dataset, the $K\times K$ neighborhood is defined geospatially or semantically. Unless stated otherwise, all the model are trained as a spatiotemporal univariate timeseries, where the model has visibility to the target variable for T length and a $K\times K$ frame.

4.3.1. Geospatial Grid Processing

For the grid-based Milan dataset, we preserve the inherent grid structure. We set a spatial kernel size of $K=11$, exploiting the physical adjacency of cells. Boundary effects are mitigated by cropping the $100\times 100$ grid to $96\times 96$, ensuring divisibility by the stride. Temporal sequences are constructed by concatenating six consecutive time steps ($T=6$).

4.3.2. Correlation-Based Spatial Construction

For the RU dataset, which lacks a grid topology, we reconstruct a semantic spatial neighborhood to bridge the gap between correlation-based relationships and our kernel-based architecture. We define "spatial" proximity using Pearson correlation. For each target RU, we compute its temporal correlation with all other RUs and then construct a virtual spatial kernel where:

• The center pixel is the target RU itself.
• The surrounding pixels are populated by the top-correlated RUs, sorted by correlation strength (filling the inner ring first, then the outer ring).

This transformation allows the layers of spatiotemporal models to operate on "functional neighbors" (units that behave similarly) just as they would on physical neighbors.

4.3.3. Normalization and Splitting

For all datasets, inputs and targets are normalized to the [0, 1] range via Min-Max scaling fitted on the training set to prevent leakage. The data is chronologically partitioned into training (70%), validation (15%), and test (15%) sets.

4.4. Implementation and Training Configuration

We implement the HiSTM model in PyTorch, leveraging GPU-optimized operations for its Mamba and Conv2D modules to ensure computational efficiency. For reproducability, we release the complete source code for our HiSTM implementation in a public repository [34].

For training and inference, we use an AI-server with a single NVIDIA A100 80GB GPU with 64 CPU cores and 512 GB RAM, using CUDA 12.4 and PyTorch 2.6.0+cu124. All models are trained under the same optimization setup (Adam, learning rate ${10}^{-4}$, learning-rate scheduler, batch size 128) and the same early-stopping criterion (patience 15 with best checkpoint selected by validation loss).

4.5. Evaluation Metrics

We evaluate prediction accuracy using five metrics. Mean Absolute Error (MAE) measures the average absolute difference between predictions and ground truth, providing a straightforward assessment of overall accuracy. Root Mean Squared Error (RMSE) penalizes larger deviations more heavily, emphasizing sensitivity to outliers. The coefficient of determination (R²) quantifies the proportion of variance explained by the model, reflecting its explanatory power. Structural Similarity Index (SSIM) evaluates the spatial coherence and perceptual similarity of predicted traffic maps. In addition, Mean Absolute Percentage Error (MAPE) and Mean Average Scaled Error (MASE) are employed in the cell-level analysis to quantify relative accuracy across regions with varying traffic volumes. their percentage- and scaled-based interpretation is used for comparing model behavior between urban, suburban, and rural cells, where absolute traffic magnitudes differ in magnitude.

All metrics are computed after reversing the normalization to the original scale.

5. Results and Analysis

This section evaluates forecasting accuracy, cross-dataset generalization, robustness under missing observations, transfer to non-grid network topologies, and computational efficiency, comparing our proposed HiSTM model against the baselines.

To make the evaluation easier to follow, we organize the results from core predictive performance to increasingly deployment-relevant considerations. We first report single-step and autoregressive multi-step accuracy. We then evaluate cross-dataset generalization and robustness under controlled missing-data scenarios. Next, we test whether HiSTM transfers beyond grid-based traffic forecasting by using correlation-defined network neighborhoods on a 5G Radio Unit (RU) dataset. Finally, we analyze spatial heterogeneity via cluster-aware variants, contrast HiSTM with a large collection of independent local models, study architectural granularity (spatial kernel size and temporal look-back), and report computational efficiency.

Evaluation protocol. Unless otherwise stated, metrics are computed on the test set and averaged across all spatial locations (for sequence metrics: across test time indices). For multi-step forecasting, we use autoregressive rollout and report per-horizon averages. Best values within each table are highlighted in bold.

5.1. Prediction Accuracy

This subsection quantifies predictive performance under standard (non-cluster-aware) training, focusing on (i) single-step accuracy on Milan, (ii) multi-step autoregressive stability (error accumulation), and (iii) cross-dataset generalization to Trentino (unseen, yet similar, dataset).

5.1.1. Single-Step Prediction Results

On the Milan dataset, under the standard setting without cluster conditioning, HiSTM achieves the best single-step forecasting performance in terms of MAE and provides competitive improvements across complementary metrics (Table 2). The model attains the lowest RMSE and the highest $R^2$ score among the compared methods. The combination of frame-wise spatial encoding with Mamba-based temporal modeling yields strong next-step accuracy while preserving spatial structure.

5.1.2. Multi-Step Autoregressive Forecasting

Having established single-step accuracy, we evaluate the stability of HiSTM under autoregressive forecasting (multi-step on Milan dataset), where error accumulation poses a greater challenge. Due to the increased computational cost of autoregressive rollout, these experiments are conducted on a randomly sampled 20% subset of the test data. As this subset differs from the full test set used for single-step prediction, the first-step autoregressive results are not directly comparable to single-step scores. The same subset is used across all models, ensuring a fair and consistent comparison between the models in this evaluation.

HiSTM demonstrates improved stability over extended forecasting horizons (Table 3). While all models accumulate error as the forecast horizon increases, HiSTM maintains the lowest MAE and RMSE at each of the six temporal steps indicating improved temporal stability under autoregressive rollout. At step 6, HiSTM achieves an MAE that is 36.8% lower than STN and 11.3% lower than STTRE. The increasing performance gap under autoregressive rollout suggests that HiSTM accumulates less compounding error over multiple steps.

5.1.3. Cross-Dataset Generalization on Trentino Dataset

To assess whether the observed performance extends beyond the training distribution, we evaluate all models on the unseen Trentino dataset.

On the Trentino dataset, HiSTM achieves strong out-of-distribution performance (Table 4). Relative to STN, HiSTM reduces MAE by 47.3% and reduces RMSE by 36.9%. HiSTM also attains the highest SSIM and the highest $R^2$ score. In terms of RMSE, TimesNet scores lower (4.7210 vs. 4.8134), while HiSTM remains best on MAE and ties or leads on the other reported metrics. Overall, HiSTM exhibits a strong generalization capability to new spatial environments with different activity distributions and scales.

5.1.4. Cell-Specific Modeling and Spatially-Aware Accuracy

While aggregate metrics quantify average performance, they obfuscate local behavior; therefore, we examine cell-level forecasts to assess spatially-aware accuracy across representative traffic regimes.

We select a 7-day window (1008 time steps) from the Milan test set and generate one-step-ahead forecasts using a 6-step temporal memory window and an $11\times 11$ spatial kernel on the $100\times 100$ grid. To highlight spatial heterogeneity, we inspect four representative cells corresponding to (a) urban, (b) suburban, (c) rural, and (d) the cell with the maximum temporal variance (Figure 6). Predicted versus actual trajectories are shown in Figure 7.

The urban cell exhibits the lowest MAPE (8.49%), with close alignment between predicted and true trajectories even during peak periods. Performance remains stable in the suburban cell (11.62% MAPE). As expected, relative error increases in the rural cell (24.86% MAPE), where low absolute volumes can inflate percentage-based metrics. The high-variance cell remains challenging but is tracked reasonably well (15.30% MAPE), suggesting that HiSTM can handle non-stationary dynamics, with residual underestimation around sharp peaks indicating headroom for improving extreme-event modeling.

5.2. Missing Data

Next, We evaluate robustness to incomplete telemetry by injecting controlled missing values and imputing with last observation carried forward (LOCF) prior to forecasting. Table 5 reports MAE/RMSE across combinations of missingness in training and testing, isolating how performance degrades under increasing observation loss.

We synthetically introduce missing values into the Milan traffic frames ${\left\{M_t\right\}}_{t=1}^{T_{\mathrm{all}}}$, where $M_t(i,j)$ denotes the scalar traffic value at spatial location $(i,j)$ at time t. Missingness is injected at the cell–time granularity, i.e., individual spatiotemporal entries $M_t(i,j)$. For each dataset split S with temporal length $L_S$, we treat every spatial location independently as a univariate time series

$$x^{(i,j)}=\left(M_1(i,j),\dots ,M_{L_S}(i,j)\right).$$

Given a target missing-data rate $m\in [0,1]$ (0%, 10%, 20%, and 30%), we partition the missing budget into two components intended to reflect common operational failure modes: random noise dropouts and contiguous outages. Specifically, $0.8m$ of the missing values are assigned to random dropouts and $0.2m$ to outages, yielding per-series counts $M_n=\lfloor 0.8m{L}_S\rfloor$ and $M_o=\lfloor 0.2m{L}_S\rfloor$, respectively.

For the random-dropout component, we uniformly sample without replacement an index set $I_n\subseteq \{2,\dots ,L_S\}$ with $|I_n|=M_n$ and set $M_t(i,j)=\mathrm{NaN}$ for all $t\in I_n$. For the outage component, we simulate contiguous telemetry losses by masking disjoint time intervals. Each interval is generated by sampling a start time $s\ge 2$ and a duration

$$\ell \sim \mathrm{Uniform}\left(\right\{1,\dots ,2T\left\}\right),$$

where T denotes the model lookback (sequence) length used to construct inputs. This choice enforces that any single outage spans at most two consecutive input windows. Intervals that overlap previously masked indices are rejected, and the procedure is repeated until exactly $M_o$ additional time steps are removed. The resulting outage index set can be written as

$$I_o=\bigcup _k\{s_k,\dots ,s_k+{\ell }_k-1\},\left|I_o\right|=M_o,I_o\cap I_n=\varnothing .$$

Missing patterns are generated independently for each spatial location $(i,j)$, i.e., missingness is not correlated across space, and they are generated separately within each dataset split to avoid temporal leakage. The validation split remains uncorrupted.

After masking, we apply LOCF imputation independently along time for each $(i,j)$, replacing every missing value with the most recent previously observed value in the same series. Because the first index of each split ($t=1$) is never masked, LOCF is always well-defined. Model inputs $X_t^{(i,j)}\in {\mathbb{R}}^{T\times K\times K}$ are then constructed from the imputed frames $\left\{{\tilde{M}}_t\right\}$ using the same temporal window length T and spatial kernel size K as in the main experimental setup.

We report MAE and RMSE over combinations of training and test missing-data conditions. Models are trained with 0% and 10% missing data and evaluated on test sets with 0–30% missing data. As shown in Table 5, HiSTM consistently achieves the best performance across all evaluated train/test combinations. When trained on complete data, HiSTM yields the lowest MAE and RMSE across test missing-data rates (MAE 5.1462–5.1451; RMSE 10.9949–10.9915), outperforming the next-best baseline (Informer). When trained with 10% missing data, HiSTM maintains its advantage (MAE 5.1762–5.1754) under all tested test-set missing-data levels.

Finally, we note that LOCF imputation can smooth short-term fluctuations; consequently, increasing the proportion of missing data may occasionally reduce apparent error for some methods by attenuating high-frequency variance in the inputs. The consistent ranking across all configurations nevertheless indicates that HiSTM is comparatively robust when operating on imputed streams, which reflects realistic deployment scenarios involving telemetry dropouts and temporary backhaul failures.

5.3. Performance on 5G High Density Demand (RU) Dataset

To test transfer beyond gridded maps, we evaluate all models on a non-grid 5G RU dataset where spatial neighborhoods are defined by correlation rather than physical adjacency. We report MAE per Key Performance Indicator (KPI) to assess whether HiSTM retains advantages under heterogeneous targets and graph-like dependencies.

We evaluate HiSTM on the Liverpool 5G High Density Demand dataset to assess its capability in forecasting complex network metrics beyond traffic volume. We predict five KPIs: Throughput (ThPut), Physical Resource Blocks (PRB), and Block Error Rate (BLER). For suffixes: S is sum, A is average, and C is count.

In this experiment, we use a lookback (sequence) length of $T=20$ past observations to forecast the next time step, and we construct each input neighborhood with kernel size $K=5$ (i.e., a $5\times 5$ local neighborhood in the correlation-defined adjacency). This dataset contains 33 cells, which motivates the choice $K=5$ as increasing the kernel to $K=6$ would require padding to reach 36 cells, introducing additional complexity beyond the correlation-based structuring. Spatial neighborhoods are constructed using correlation-based kernels rather than geographic grids, enabling the models to exploit functional dependencies between Radio Units. To avoid temporal leakage, correlation statistics are computed using only the training split of the time series (i.e., without access to the evaluation interval).

Table 6 presents the MAE for each target metric. HiSTM achieves the lowest error across all five KPIs. For PRB Average (PRB-A), HiSTM achieves the lowest MAE, improving upon the next best model, Autoformer. A similar trend is observed for BLER-C, where HiSTM again records the best performance, outperforming both TimesNet and STTRE. These results show that correlation-defined neighborhoods provide a useful proxy for spatial structure in this setting, and that HiSTM can transfer to non-grid topologies.

5.4. Cluster-Aware HiSTM Variants

We evaluate whether explicitly encoding region identity improves HiSTM under heterogeneous traffic regimes. Using the spatial clusters (defined in 4.1.1), we compare global and cluster-conditioned HiSTM variants (defined in 3.2.4), reporting cluster-wise and aggregate errors to quantify trade-offs between shared representation learning and local specialization. Table 7 summarizes the evaluation between clusters using MAE, RMSE, and MASE.

Several patterns emerge:

• The Downtown (High volume, high variance) region exhibits the largest absolute errors in all models (MAE ≈ 32–42), consistent with its high magnitude and volatile traffic conditions. Both GB and GA achieve the lowest MAE values (32.45 and 32.57). Forecasting in Downtown regions benefits most from shared global learning, as strong variability and frequent regime shifts are better captured by Global-Aware HiSTM than by cluster-specific models.
• In the Urban (Medium–high volume, structured dynamics) and Suburban (Moderate volume, stable dynamics) clusters, GA consistently shows the best performance, with MAE 21.24 and 8.91, respectively. Urban regions exhibit consistent improvements with cluster-aware conditioning, indicating that combining global patterns with region-specific embeddings effectively balances generalization and local adaptation. In Suburban regions, performance differences across model variants are small, suggesting that moderate spatial heterogeneity can be handled effectively by both global and cluster-aware models.
• For the Rural (Low volume, homogeneous traffic) region, the locally trained CS model performs best (MAE 2.49, RMSE 4.28). Rural regions favor localized modeling, as their low variance and homogeneous traffic patterns reduce the benefit of global context, making cluster-specific training competitive.
• The ClusterFiLM models perform competitively, particularly CF-U, which reduces aggregate MAE from 5.33 (CF) to 5.26 and achieves the lowest overall MASE (0.091) among all FiLM variants. The re-weighted loss helps preserve accuracy in small clusters without degrading performance on large ones.

The Global-Aware HiSTM attains the best aggregate accuracy (MAE 5.13, RMSE 11.02, MASE 0.089). Our results indicate that the Global-Aware HiSTM provides the best trade-off across heterogeneous regions, leveraging shared representations where beneficial while retaining robustness in homogeneous areas.

5.5. 10k LSTMs vs HiSTM

To isolate the value of joint spatiotemporal learning, we compare HiSTM against a detailed baseline consisting of 10,000 independently trained single-cell LSTM models, trained and evaluated on the Milan dataset.

To assess the benefit of spatial contexts, we conduct a large-scale comparative experiment between HiSTM and a collection of 10,000 independent LSTM models. Each LSTM is trained separately for a single cellular tower, using only its local temporal data without any spatial information. This set of 10k LSTMs is evaluated against one global HiSTM.

Figure 8 presents a heatmap of the spatial distribution of model improvement across the Milan traffic grid, where each cell indicates whether HiSTM or the corresponding single-tower LSTM achieved a lower MAE/RMSE. HiSTM outperforms independent LSTMs in most spatial locations, demonstrating the consistent advantage of joint spatiotemporal learning.

The distribution of RMSE differences (Figure 9) further supports this observation. The histogram shows that most cells exhibit positive RMSE differences (i.e., ${\mathrm{RMSE}}_{\mathrm{LSTM}}-{\mathrm{RMSE}}_{\mathrm{HiSTM}}>0$), indicating that HiSTM achieves lower prediction error in nearly all regions. Only a small fraction of cells show negative differences. The spatial context and shared regime learning are more valuable than increasing the number of local models, even at extreme scale, as one HiSTM model is more effective than ten thousand independently trained LSTMs.

5.6. Model Granularity

To examine the sensitivity of HiSTM to architectural and temporal design choices, we conduct an ablation study varying both the spatial kernel size and the temporal look-back window. The goal is to assess how the model’s forecasting accuracy changes with different spatial receptive fields and temporal dependencies.

5.6.1. Spatial Granularity (Variation in Kernel Sizes)

Table 8 reports HiSTM performance under different spatial kernel sizes ranging from $3\times 3$ to $13\times 13$. Kernel size affects both MAE and RMSE. MAE improves up to $9\times 9$ (best MAE 5.1303), after which changes are marginal. RMSE continues to improve slightly with larger kernels, reaching its minimum at $13\times 13$ (10.9666). Performance saturates beyond medium-sized kernels, indicating that HiSTM effectively captures relevant spatial context without requiring large receptive fields.

5.6.2. Temporal Granularity (Variation in Look-Back Window Size)

To investigate the impact of temporal receptive field length, we vary the look-back window from 10 minutes (1 step) to 6 hours (36 steps). We additionally incorporate two forms of periodic temporal context: the last-day (LD) kernel and the last-week (LW) kernel. As summarized in Table 9, performance improves as the window length increases up to 1.5 hours, after which it stabilizes.

Incorporating periodic historical kernels (+LD and +LD+LW) consistently improves predictive accuracy across all temporal window configurations, highlighting the importance of explicit seasonal context. Performance peaks at a moderate look-back window of 1.5 hours when both daily and weekly kernels are included. This indicates that combining short-term dynamics with periodic information yields the most effective temporal representation.

5.7. Computational Efficiency Analysis

Table 10 presents a comparison of the computational characteristics across eleven models, including both classical baselines and recent transformer-based architectures ($T=6$ and $K=11\times 11$).

With only 33.8K parameters (0.13 MB), HiSTM is the most lightweight model in the benchmark—approximately $45\times$ smaller than PredRNNpp, $44\times$ smaller than VMRNN-D, and $18\times$ smaller than xLSTM. Despite its compact size, HiSTM achieves a competitive inference latency of 1.19ms, outperforming larger recurrent and transformer-based models such as PredRNNpp (6.38ms), VMRNN-D (18.58ms), and Autoformer (2.06ms). Only xLSTM (1.01ms) and PatchTST (0.46ms) deliver marginally faster inference, although at a significantly higher model or memory cost.

In terms of computational workload, HiSTM performs at $1.36\times {10}^7$ MACs—lower than transformer-heavy models like PatchTST ($2.54\times {10}^8$) and PredRNNpp ($1.57\times {10}^9$). This balance between model compactness, low latency, and moderate MACs underscores HiSTM’s suitability for resource-constrained or real-time forecasting environments, where both efficiency and responsiveness are critical.

Overall, HiSTM delivers a favorable efficiency profile compared to the state-of-the-art temporal models in this benchmark, supporting the architectural design as a practical option for spatiotemporal prediction tasks in resource-constrained or real-time settings.

6. Explainability Analysis

This section analyzes the internal behavior of HiSTM to understand how spatial and temporal information contribute to its predictive performance. We analyze the internal dynamics of its two core temporal components: the Mamba State Space Model (SSM) and the Temporal Attention mechanism. This analysis reveals distinct roles for each component in capturing traffic patterns. While the following analyses do not provide formal explainability guarantees, they offer qualitative evidence of how HiSTM leverages spatial and temporal structure.

6.1. Holistic Context vs. Selective Attention

We examined the attention weights ${\alpha }_t$ assigned to the input sequence of length $T=6$. Surprisingly, the attention distribution is highly uniform, with a normalized entropy of $1.0000$. As shown in Figure 10, the mean weights range narrowly from $0.1661$ at $t-6$ to $0.1690$ at $t-1$, with negligible standard deviation ($<0.0002$). This confirms that the model pays nearly equal attention to all time steps in the lookback window.

This uniformity indicates that HiSTM does not rely on "picking" a single informative time step (e.g., simply copying the last observation). Instead, the Temporal Attention layer acts as a global context aggregator, treating the entire observation window as a single, cohesive feature set. This places the burden of temporal dynamics modeling largely on the Mamba layer: Mamba compresses the sequential history into a rich hidden state, while attention summarizes that history rather than performing hard selection.

To test whether the Temporal Attention is merely acting as mean pooling, we performed an inference-time ablation where we forced ${\alpha }_t=1/T$ (uniform weights), keeping the same trained checkpoint and all other components fixed. This replacement degraded performance: MAE increased by $+2.3\%$, RMSE increased by $+15.3\%$, and $R^2$ dropped by $0.66\%$ relative. These results suggest that even small, sample-specific deviations from perfect uniformity encode useful temporal emphasis, and the attention layer is not fully redundant.

We further investigate whether this behavior holds during dynamic traffic shifts by analyzing two distinct scenarios.

First, we examine a traffic sequence with high variance and fluctuation (Figure 11). Despite the erratic input pattern, the attention mechanism maintains a remarkably stable distribution. This suggests that the model’s robustness to noise is partly due to its refusal to overfit to short-term volatility.

Second, we analyze a "ramp-up" phase where traffic volume consistently increases over the window (Figure 12). Even in the presence of a strong positive trend, the attention weights do not skew towards the most recent time step ($t-1$). This confirms that HiSTM aggregates the entire temporal context to determine the trend’s trajectory, rather than reacting solely to the latest observation. In both cases, the uniformity of the attention weights (red dashed lines) validates the role of the Mamba layer as the primary driver of temporal dynamics modeling.

Although the learned attention weights appear close to uniform, forcing uniform weights at inference degrades performance, indicating that the attention module contributes meaningfully beyond simple averaging.

This behavior explains why HiSTM exhibits stable multi-step forecasts, as temporal emphasis is distributed across the entire context window rather than being dominated by short-term fluctuations.

6.2. Latent State Dynamics

We visualized the evolution of the Mamba hidden state by projecting the feature trajectories into a 2D latent space using PCA (Figure 13). The latent states exhibit a low variance ($\sigma \approx 0.02$), suggesting a robust representation that is resistant to input noise.

Analysis of the state transition magnitudes reveals a non-linear processing of time. While most step-to-step transitions are small ($\approx 0.010$), there is a distinct peak in magnitude ($0.0165$) during the transition from $t-4$ to $t-3$. This suggests a "critical update" phase where the SSM integrates sufficient history to reduce uncertainty in the evolving traffic trajectory, effectively locking in a prediction trajectory midway through the sequence.

Such controlled latent state transitions are consistent with the observed reduction in error accumulation during autoregressive forecasting.

6.3. Spatial Dependency Analysis

To quantify the contribution of spatial context to prediction accuracy, we performed a spatial ablation study across the four identified clusters. We measured the relative increase in Mean Absolute Error (MAE) as the input spatial kernel was progressively restricted from the full $11\times 11$ context down to a single $1\times 1$ pixel (Figure 14).

The results align with the traffic characteristics identified in Sub 5.4. In the Downtown cluster, removing spatial context causes a sharp performance degradation, with MAE increasing by $+228.9\%$ (from 32.57 to 107.13) when only the center pixel is visible. This heavy reliance on neighbors is consistent with the high volume and variance observed in this region, where complex dynamics and user mobility increase sensitivity to neighboring cells, making spatial context critical for accurate forecasting.

In contrast, the Rural cluster is remarkably insensitive to spatial context, showing only a $+12.4\%$ increase in error (from 2.51 to 2.82) when stripped of its surroundings. This reflects the sparse and homogeneous nature of rural traffic, where patterns are driven by local, static behavior rather than complex spatial interactions. The Urban ($+46.1\%$) and Suburban ($+45.6\%$) clusters occupy a middle ground, indicating a moderate dependency proportional to their density.

This cluster-specific sensitivity validates the efficacy of our Global-Aware approach ( 5.4): by learning cluster embeddings, HiSTM can dynamically attend to spatial features in Downtown areas while avoiding noise overfitting in Rural regions. This adaptive sensitivity to spatial context explains why HiSTM achieves larger accuracy gains in high-density regions while maintaining stable performance in sparse areas.

6.4. Data Distribution and Model Behavior

To explain the systematic advantages observed for HiSTM in the previous experiments, we analyze the underlying data distribution and its implications for spatial modeling.

The performance gap between the 10,000 independent LSTM models and the unified HiSTM approach is consistent with the dataset’s distributional structure (from 5.5). Analysis of the Milan traffic dataset reveals a hierarchical structure where both the overall distribution and individual cluster distributions exhibit bimodality, suggesting that shared representations can be beneficial.

The overall traffic distribution is strongly bimodal (${\mathrm{BIC}}_{2-component}$ = 947,771,715 vs ${\mathrm{BIC}}_{1-component}$ = 1,108,616,079), with two distinct regimes: a low-traffic component (mean=27.48 units, weight=81.7%) and a high-traffic component (mean=217.83 units, weight=18.3%). This bimodality persists at the cluster level, with each of the four spatial clusters showing clear two-component structure:

Rural: 

Low-mode dominated (71.8% in lower component) with means at 14.78 and 48.82 units

Suburban: 

Strong low-mode preference (81.9%) with means at 81.56 and 190.08 units

Urban: 

Mixed distribution (74.8% low-mode) with means at 229.82 and 459.76 units

Downtown: 

High-variance bimodal (78.2% low-mode) with means at 462.41 and 1151.31 units

This hierarchical bimodality is a main contributor to the observed gains of HiSTM over independent local models. While an individual LSTM can learn the local regime structure of a single cell, it cannot leverage correlated regime transitions across space. In contrast, a shared spatiotemporal model can transfer information about regime shifts across locations.

In contrast, HiSTM can leverage spatial dependencies to build shared representations that transfer knowledge about regime transitions. This is particularly relevant when regime changes co-occur across neighboring or functionally related cells.

The improvement metrics are consistent with this interpretation. HiSTM shows larger RMSE improvements in high-variance regions (e.g., Downtown and Urban) and smaller but still positive gains in more homogeneous areas (Rural), indicating that spatial context can help across diverse regimes.

In 5.4, a seemingly paradoxical finding emerges when comparing cluster-specific versus global approaches. The Global-Blind model (trained on all data without explicit cluster identifiers) outperforms the Cluster-Specific models in aggregate metrics (MAE: 5.15 vs. 5.31). This suggests that learning from the full dataset can improve representation quality when dynamics (e.g., regime structure) are shared across regions, even if local distributions differ in magnitude.

However, the way cluster information is incorporated matters. Comparing ClusterFiLM variants with the Global-Aware model shows that learned cluster embeddings can be an effective mechanism for capturing spatial heterogeneity, yielding the best overall performance in this set of variants (MAE: 5.13). Overall, the Global-Aware variant provides a favorable balance between sharing information across the full grid and adapting to region-specific characteristics.

Together, these findings clarify why a unified spatiotemporal model can generalize more effectively than isolated local models when regime structure is shared across space.

7. Conclusion

In this work, we presented HiSTM, a hierarchical spatiotemporal forecasting model that integrates local spatial convolutions, Mamba-based temporal modeling, and global attention aggregation. Our extensive evaluation on real-world cellular traffic datasets shows that HiSTM exhibits better performance in comparison to the baselines, including Transformer-, RNN-, and Mamba-based architectures, in both accuracy and computational efficiency.

A key finding of this study is the scalability of global spatiotemporal modeling. We showed that a single HiSTM model outperforms 10,000 independent, cell-specific LSTMs. This result challenges the conventional reliance on localized modeling for precision, revealing that the ability to learn shared representations of traffic regimes (e.g., bimodal distributions) across the entire network outweighs the specificity of local training. By leveraging spatial correlations, HiSTM generalizes well to unseen environments (e.g., Trentino) and correlation-based topologies (e.g., Liverpool 5G RU), maintaining better accuracy where isolated models fail.

The explainability analysis shows that HiSTM’s architectural components align with known traffic dynamics, providing qualitative support for the observed gains in accuracy, robustness, and scalability. The broadly distributed temporal attention weights indicate that the model aggregates contextual information across the entire look-back window, rather than relying on heuristic selection of specific time steps. In parallel, the Mamba-based state-space layer summarizes complex temporal dynamics into a compact latent representation, which is consistent with the model’s observed stability under high-variance traffic conditions and robustness to up to 30% missing data.

Future work will focus on extending Mamba-based spatiotemporal models to graph-based implementations to handle irregular network topologies without kernel approximations. Additionally, we aim to explore the integration of exogenous features (e.g., weather, social events) to enhance predictive capability in anomaly-heavy scenarios. HiSTM stands as a scalable, lightweight, and robust solution for the next generation of predictive network management.