Skip to content

Multivariate Time Series Forecasting needs Cross Variable Loss

Conference: NeurIPS 2026 (Accepted, according to the conference list)
arXiv: 2608.05742
Code: https://github.com/Day333/CvLoss
Area: Time Series
Keywords: multivariate forecasting, cross-variable loss, residual consistency, graph total variation, direct forecasting

TL;DR

CvLoss leaves the forecasting backbone unchanged and supplements point-wise MSE with constraints on residual differences between forecast patches from different variables, achieving 111 wins, 2 ties, and 1 loss across the paper's 114 backbone-controlled comparison cells, with no loss computation at inference time.

Background & Motivation

Multivariate time series forecasting must preserve how variables evolve together, not merely predict each curve: electricity customers may respond to the same weather conditions, while traffic disturbances may propagate between locations. Methods such as iTransformer, Crossformer, and TimeFilter primarily organize cross-variable information in the historical input and then generate the entire future window through direct forecasting (DF). Even when the network shares representations, the usual point-wise MSE simply accumulates errors at individual time–variable positions without explicitly supervising relative relationships between variables.

This does not mean that an MSE-trained model cannot learn dependencies. Rather, a backbone's ability to exchange information is different from an objective's requirement to recover structure. Two predictors can have similar point-wise errors yet produce different correlation matrices for future variables. The paper illustrates this with an identical iTransformer on ECL: after adding CvLoss, 84.83% of cross-variable matrix entries have smaller structural errors. This percentage measures correlation-structure recovery, not the reduction in MSE.

Instead of estimating the full spatiotemporal covariance of future residuals, the authors divide the future window into patches and compare predicted and target relative differences on a cross-variable graph. The loss can therefore supervise both synchronous relationships within a patch index and lagged relationships across patch indices. Core idea: supplement point-wise supervision with graph-edge residual differences that constrain cross-variable structure, while retaining MSE as an absolute numerical anchor.

Method

Overall Architecture

The input is a historical window and the output is a forecast matrix covering all future variables; the forecasting backbone, input processing, and deployed forward path remain unchanged. During training, predictions and targets are partitioned identically into non-overlapping temporal patches, forming nodes indexed by patch and variable. A structural loss between nodes from different variables is then jointly optimized with MSE.

Each node stores the prediction residual vector for an entire patch, and graph edges check whether connected nodes have consistent residuals. This graph is a collection of relationships in the training objective, not an additional graph neural network: there is no message passing or inference-time correction. Prose and loss equations are therefore sufficient; a multi-module network diagram would be misleading.

The default is a predetermined complete cross-variable graph covering synchronous and asynchronous edges, not a graph learned from the correlation matrix of future labels. Training labels enter the default objective only through the usual supervised residuals. Testing requires neither the true future nor construction of this loss graph.

Key Designs

1. Relation-aware objective: distinguish point-wise error from joint error geometry

After flattening residuals across all future positions, MSE assigns the same uncoupled penalty to each coordinate. The authors assume zero-mean multivariate Gaussian residuals with a positive-definite covariance that does not depend on predictor parameters in the derivation. Under these conditions, the true negative log-likelihood weights errors through a precision matrix, whereas an isotropic Gaussian corresponds to ordinary MSE.

The paper's objective gap is:

\[ \Delta=\frac{1}{2TD}\mathbf{e}^{\top}\left(\mathbf{\Sigma}_{ST}^{-1}-\frac{1}{\sigma^{2}}\mathbf{I}\right)\mathbf{e}. \]

Here, the forecast horizon is \(T\), the number of variables is \(D\), and \(\mathbf{e}\) is the flattened forecast residual. Off-diagonal precision terms indicate that errors at different positions should be scored jointly. However, this derivation identifies a difference in objective form; it does not prove that a particular regularizer outperforms MSE on arbitrary data.

In particular, two quadratic forms being different functions does not imply a nonzero difference for every residual vector. The paper describes a nonspherical precision matrix as implying a strictly nonzero gap, but zero residuals or residuals in certain directions can still make the expression vanish. The expression also does not guarantee a positive gap or, by itself, a necessary reduction in prediction risk.

The authors additionally show that MSE plus squared graph-edge differences corresponds to a Gaussian precision family constructed from a graph Laplacian. This provides a valid structured family for relation-aware objectives, not recovery of the true precision matrix: uniform edge weights, fixed signs, and fixed topology remain substantial restrictions. A complete graph not excluding cross-variable edges should not be confused with the ability to represent arbitrary residual covariance.

2. Patch cross-variable graph: place synchronous and lagged relationships in one supervision set

The future window is partitioned into \(P\) non-overlapping patches of length \(L\), with \(T=PL\). Each patch of each variable becomes a node, giving \(N=PD\) nodes. The node residual is \(\mathbf{e}_v=\hat{\mathbf{z}}_v-\mathbf{z}_v\), the predicted patch vector minus its target; \(L=1\) recovers point nodes.

Synchronous edges connect different variables at the same patch index, while asynchronous edges connect different variables at different patch indices. The default graph is their union. Asynchronous comparisons still subtract corresponding within-patch positions, so they cover relative relationships induced by patch displacement rather than explicitly learning propagation direction, continuous delays, or a causal graph.

A complete cross-variable graph avoids prior knowledge of which variables interact and avoids estimating a high-dimensional precision matrix. For example, ECL with horizon 720 and 321 variables would require a full spatiotemporal precision matrix with more than 50 billion entries, making direct estimation impractical. Graph-edge constraints turn this into computable pairwise supervision, but the complete graph also includes irrelevant edges and weights all edges equally.

3. Residual graph total variation: preserve relative differences instead of making curves identical

For an edge, compare the difference between predicted patches with the difference between target patches. Their discrepancy is exactly the difference between the two node residuals, so CvLoss is mean graph total variation on the residual field:

\[ \mathcal{L}_{\mathrm{cv}}=\frac{1}{|\mathcal{E}|L}\sum_{(i,j)\in\mathcal{E}}\left\|(\hat{\mathbf{z}}_i-\hat{\mathbf{z}}_j)-(\mathbf{z}_i-\mathbf{z}_j)\right\|_1=\frac{1}{|\mathcal{E}|L}\sum_{(i,j)\in\mathcal{E}}\|\mathbf{e}_i-\mathbf{e}_j\|_1. \]

Here, \(\mathcal{E}\) is the graph-edge set, and the denominator normalizes by edge count and patch length. The constraint concerns how well predictions preserve true relative differences, not whether different variables have equal predictions. Level differences, amplitude differences, and opposing movements already present in the labels can be preserved.

This equivalence also exposes a blind spot: if every node has the same nonzero residual vector, every edge discrepancy is still zero. CvLoss can detect uncoordinated errors between variables but cannot independently identify a common offset. It must therefore accompany absolute-error supervision rather than be presented as a standalone guarantee of accurate forecasts.

The deployed loss uses the \(\ell_1\) norm because a few variable-specific shocks can create large edge discrepancies that dominate a squared penalty. The Gaussian proposition directly corresponds to a squared \(\ell_2\) graph penalty. Switching to \(\ell_1\) admits a pairwise Gibbs-field interpretation on the same graph, but no longer gives the objective derived by that Gaussian proposition. The norm choice has ablation support: it is an empirical decision, not a necessary consequence of the theorem.

Loss & Training

MSE remains the absolute anchor, and the final training objective is:

\[ \mathcal{L}_{\alpha}=(1-\alpha)\mathcal{L}_{\mathrm{df}}+\alpha\mathcal{L}_{\mathrm{cv}},\qquad \mathcal{L}_{\mathrm{df}}=\frac{1}{TD}\|\hat{Y}-Y\|_F^2. \]

In the main method, \(\alpha\in(0,1)\) is a scalar selected on validation and then fixed; patch length \(L\) is also selected on validation. The weight does not vary across batches, variables, or edges and is not learned adaptively during training. The learned unnormalized coefficients in Appendix D.2 belong only to an auxiliary diagnostic and are not this \(\alpha\).

Backbone settings, preprocessing, normalization, the Adam optimizer, and training protocols follow the public baselines; only the loss weight and patch length are tuned. Validation early stopping has a patience of 3 epochs. Test-time drop-last is disabled so that the final incomplete batch is retained. Sensitivity experiments additionally test \(\alpha=0\) and \(\alpha=1\) as boundary diagnostics, not as the definition of the main method's fixed-weight interval.

The full graph's edge count grows as \(O(P^2D^2)\) with patch and variable counts. Appendix C.10 separately studies an efficiency variant sampling at most 1000 random edges per batch; the main accuracy results do not come from this sampled variant. Structural-error top-K edge selection in Appendix D is likewise not the main-result configuration. Top-K requires correlation structure from training labels, and the statistical cost of selecting edges cannot be summarized solely by the cost of evaluating the retained K edges.

Deployment removes the entire structural-loss computation, adding no prediction parameters, buffers, or forward operations. “Zero inference overhead” refers to this unchanged prediction path, not to the absence of additional full-graph training cost. Backward overhead is more pronounced on high-dimensional ECL than on low-dimensional ETTh2.

Key Experimental Results

Main Results

Experiments cover 12 datasets: 4 ETT subsets, Weather, ECL, Traffic, Solar, and 4 PEMS subsets. The historical window is fixed at 96, and training, validation, and test splits are chronological. Long-horizon results average over horizons 96, 192, 336, and 720. The following excerpt from main-text Table 1 compares the identical TimeFilter backbone; lower MSE and MAE are better.

Dataset TimeFilter MSE +CvLoss MSE TimeFilter MAE +CvLoss MAE
ETTm1 0.377 0.372 0.393 0.381
ETTh1 0.420 0.419 0.428 0.427
Weather 0.240 0.236 0.270 0.260
ECL 0.159 0.157 0.256 0.252
Traffic 0.408 0.407 0.269 0.254
Solar 0.228 0.218 0.262 0.254
PEMS07 0.071 0.063 0.170 0.156

These pairs support attribution to the objective, whereas the other architectures in the same leaderboard cannot replace paired experiments. TQNet achieves 0.197 MSE on Solar, still below TimeFilter+CvLoss at 0.218, so the CvLoss column is not the best architecture on every dataset and metric.

Main-text Table 2 consolidates same-configuration comparisons across 7 backbones, defining a cell as a backbone–dataset–metric combination averaged over horizons. Among 114 cells, there are 111 wins, 2 ties, and 1 loss, with mean relative reductions of 3.39% in MSE and 3.61% in MAE. Backbones and datasets overlap between blocks, making this a descriptive evidence summary rather than a significance test over 114 independent tasks.

Ablation Study

The following excerpt from main-text Table 4 reports topology ablations on TQNet, averaged over 4 long horizons. Each cell displays MSE / MAE.

Graph configuration ETTm1 ETTh1 ECL Weather
DF, no structural edges 0.377 / 0.393 0.441 / 0.434 0.165 / 0.259 0.242 / 0.269
Synchronous edges only 0.375 / 0.389 0.444 / 0.435 0.164 / 0.255 0.244 / 0.269
Asynchronous edges only 0.374 / 0.387 0.444 / 0.435 0.164 / 0.256 0.241 / 0.266
Synchronous + asynchronous edges 0.372 / 0.387 0.438 / 0.430 0.162 / 0.253 0.240 / 0.264

Adding only one edge type does not guarantee gains: both single-support graphs perform worse than DF on ETTh1, while the synchronous-only graph does not improve Weather. The combined graph is more reliable overall, so the presence of dependencies does not imply that arbitrary relational constraints are beneficial.

Appendix Table 15 additionally compares edge-discrepancy norms under a fixed protocol, averaging over horizons 96, 192, 336, and 720. This table has its own baseline values and should not be merged with main-text Table 1 as if they were the same run.

Dataset / backbone Base MSE / MAE Squared L2 edge penalty L1 edge penalty
ECL / TimeFilter 0.158 / 0.256 0.157 / 0.254 0.156 / 0.252
ECL / TQNet 0.165 / 0.259 0.164 / 0.257 0.162 / 0.254
Weather / TimeFilter 0.241 / 0.271 0.241 / 0.270 0.238 / 0.262
Weather / TQNet 0.242 / 0.269 0.244 / 0.269 0.241 / 0.264

Key Findings

  • The L1 variant outperforms squared L2 in all 8 metric cells above, but L2 does not help everywhere: Weather/TQNet MSE changes from 0.242 to 0.244. This supports the norm preference under the evaluated protocol, not a universally optimal norm.
  • ECL's weight sensitivity demonstrates the practical value of the absolute anchor: in main-text Table 6, TimeFilter MSE is 0.159, 0.157, and 0.168 at weights 0, 0.5, and 1. Pure structural supervision is worse than the unregularized baseline in this example.
  • The appendix uses 7 paired random seeds, 2020–2026. Of the 13 cells with overlapping marginal intervals in Table 18, 11 have paired 95% intervals excluding zero. Weather/TQNet MSE improvements at horizons 336 and 720 remain inconclusive and should not uniformly be called significant wins.
  • Structural visualizations use the absolute Pearson correlation matrix of future sequences and compare element-wise absolute errors between predicted and true correlations. ECL's 84.83%, PEMS03's 61.98%, and Weather's 44% are proportions of entries with smaller structural errors, not forecasting-error reduction percentages; they also discard correlation signs.

The source contains numerical and protocol conflicts that must remain explicit. The prose states that Traffic MAE falls by 0.014, while Table 1's displayed values, 0.269 and 0.254, differ by 0.015. The ETTm1 baseline MAE is 0.393 in Table 1 but 0.394 in Appendix Table 8. Main-text Table 3 gives TQNet/ETTh1 CvLoss as MSE 0.430 and MAE 0.438, whereas the corresponding average row in Appendix Table 13 gives 0.438 and 0.430. This note does not use that conflicting pair for definitive win/loss attribution.

PEMS protocols are also inconsistent: Appendix B specifies horizons 12, 24, 36, and 48, while actual Table 9 lists only 12, 24, and 48 and its caption still calls them “four” horizons. The PEMS07 average above is retained from Table 1 without inventing horizon-36 results or recalculating the authors' average.

Highlights & Insights

  • Structural supervision can be independent of a structural backbone. Even with input-side graph filtering or cross-variable attention, the output objective can supply additional supervision. This makes backbone-controlled comparisons particularly important.
  • Smoothing residuals is more appropriate than directly smoothing forecasts. Different physical variables need not produce identical curves; predicted relative differences should instead match target relationships. Different variable levels therefore do not, by themselves, incur a penalty.
  • Relative constraints need an absolute anchor. The translation invariance of graph total variation directly explains the retention of MSE, rather than leaving the combined objective as an unexplained tuning heuristic.

Limitations & Future Work

  • Theoretical motivation is not a general optimality guarantee. Gaussian assumptions, fixed covariance, and the structured precision family are restrictive. The deployed L1 loss relies more heavily on empirical evidence, and graph total variation cannot arbitrarily represent signed, heterogeneous residual dependencies.
  • The complete graph includes irrelevant relationships. The authors acknowledge limited adaptivity in the predefined graph and fixed patches. Future work could investigate sparse or weighted graphs selected only from training history, reporting the full graph-construction and selection costs.
  • Evaluation remains deterministic forecasting with regular sampling. Probabilistic forecasts, irregular intervals, and missing values are untested. The deliberate exclusion of Exchange also limits extrapolation to weakly coupled, nearly random-walk settings.
  • Evidence for structural recovery is limited. Absolute Pearson correlation discards signs and does not directly validate the residual precision matrix. Signed correlations, joint-distribution calibration, and operational constraints would complement heatmaps.
  • Manuscript consistency affects reproducibility judgments. Small baseline discrepancies, metric reversals, and conflicting PEMS horizons require author code or clarification; a note should not silently reconcile them.
  • vs iTransformer / Crossformer / TimeFilter: These methods primarily change how historical variables exchange information, whereas CvLoss changes how future outputs receive structural supervision. The approaches can be combined, but strong input modeling does not mean every output regularizer will help.
  • vs FreDF / QDF / Time-o1: These objectives improve training through frequency-domain supervision, joint error modeling, or label transformations. CvLoss explicitly compares relative differences between cross-variable patches. Other objectives still win in specific settings, precluding a claim of universal replacement.
  • vs DBLoss / Patch-wise Structural Loss: Decomposition and patch-structure supervision emphasize sequence structure; this paper places relationships between nodes from different variables. Complementarity between temporal decomposition and cross-variable structure is worth testing, but requires identical backbones and tuning budgets; separate gains do not establish combined gains.

Rating

  • Novelty: 4/5. Organizes output-side cross-variable residual consistency into a plug-in graph loss with a simple mechanism and a clear problem focus.
  • Experimental Thoroughness: 4/5. Covers backbones, topology, norms, weights, and paired-seed analysis, but some protocols and tables conflict.
  • Writing Quality: 4/5. Version v3 distinguishes theory, L1 practice, and auxiliary diagnostics clearly, while the strictly nonzero-gap claim and table consistency need improvement.
  • Value: 4/5. Offers a minimally invasive training change for existing forecasters, chiefly providing consistent modest gains without extending the inference path.